Search in sources :

Example 31 with Variable

use of cbit.vcell.math.Variable in project vcell by virtualcell.

the class SimulationData method getVarAndFunctionDataIdentifiers.

/**
 * This method was created in VisualAge.
 * @return java.lang.String[]
 */
public synchronized DataIdentifier[] getVarAndFunctionDataIdentifiers(OutputContext outputContext) throws IOException, DataAccessException {
    // Is this zip format?
    boolean bIsChombo = false;
    try {
        bIsChombo = isChombo();
    } catch (FileNotFoundException e) {
        e.printStackTrace(System.out);
    }
    File zipFile1 = getZipFile(bIsChombo, null);
    File zipFile2 = getZipFile(bIsChombo, 0);
    bZipFormat1 = false;
    bZipFormat2 = false;
    if (zipFile1.exists()) {
        bZipFormat1 = true;
    } else if (zipFile2.exists()) {
        bZipFormat2 = true;
    }
    refreshLogFile();
    boolean bIsComsol = false;
    try {
        bIsComsol = isComsol();
    } catch (FileNotFoundException e) {
        e.printStackTrace(System.out);
    }
    if (!bIsComsol) {
        try {
            refreshMeshFile();
        } catch (MathException e) {
            e.printStackTrace(System.out);
            throw new DataAccessException(e.getMessage());
        }
    }
    if (!isRulesData && !getIsODEData() && !bIsComsol && dataFilenames != null) {
        // read variables only when I have never read the file since variables don't change
        if (dataSetIdentifierList.size() == 0) {
            File file = getPDEDataFile(0.0);
            DataSet dataSet = getPDEDataSet(file, 0.0);
            String[] varNames = dataSet.getDataNames();
            int[] varTypeInts = dataSet.getVariableTypeIntegers();
            if (varNames == null) {
                return null;
            }
            dataSetIdentifierList.clear();
            for (int i = 0; i < varNames.length; i++) {
                VariableType varType = null;
                try {
                    varType = VariableType.getVariableTypeFromInteger(varTypeInts[i]);
                } catch (IllegalArgumentException e) {
                    if (LG.isWarnEnabled()) {
                        LG.warn("Exception typing " + varNames[i] + " has unsupported type " + varTypeInts[i] + ": " + e.getMessage());
                    }
                    varType = SimulationData.getVariableTypeFromLength(mesh, dataSet.getDataLength(varNames[i]));
                }
                Domain domain = Variable.getDomainFromCombinedIdentifier(varNames[i]);
                String varName = Variable.getNameFromCombinedIdentifier(varNames[i]);
                dataSetIdentifierList.addElement(new DataSetIdentifier(varName, varType, domain));
            }
            refreshDataProcessingOutputInfo(outputContext);
            if (dataProcessingOutputInfo != null) {
                for (int i = 0; i < dataProcessingOutputInfo.getVariableNames().length; i++) {
                    if (dataProcessingOutputInfo.getPostProcessDataType(dataProcessingOutputInfo.getVariableNames()[i]).equals(DataProcessingOutputInfo.PostProcessDataType.image)) {
                        dataSetIdentifierList.addElement(new DataSetIdentifier(dataProcessingOutputInfo.getVariableNames()[i], VariableType.POSTPROCESSING, null));
                    }
                }
            }
        }
        // always read functions file since functions might change
        getFunctionDataIdentifiers(outputContext);
    }
    if ((isRulesData || getIsODEData()) && dataSetIdentifierList.size() == 0) {
        ODEDataBlock odeDataBlock = getODEDataBlock();
        if (odeDataBlock == null) {
            throw new DataAccessException("Results are not availabe yet. Please try again later.");
        }
        ODESimData odeSimData = odeDataBlock.getODESimData();
        int colCount = odeSimData.getColumnDescriptionsCount();
        // assume index=0 is time "t"
        int DATA_OFFSET = 1;
        dataSetIdentifierList.clear();
        for (int i = 0; i < (colCount - DATA_OFFSET); i++) {
            String varName = odeSimData.getColumnDescriptions(i + DATA_OFFSET).getDisplayName();
            // TODO domain
            Domain domain = null;
            dataSetIdentifierList.addElement(new DataSetIdentifier(varName, VariableType.NONSPATIAL, domain));
        }
    }
    if (bIsComsol && dataSetIdentifierList.size() == 0) {
        ComsolSimFiles comsolSimFiles = getComsolSimFiles();
        if (comsolSimFiles.simTaskXMLFile != null) {
            try {
                String xmlString = FileUtils.readFileToString(comsolSimFiles.simTaskXMLFile);
                SimulationTask simTask = XmlHelper.XMLToSimTask(xmlString);
                Enumeration<Variable> variablesEnum = simTask.getSimulation().getMathDescription().getVariables();
                while (variablesEnum.hasMoreElements()) {
                    Variable var = variablesEnum.nextElement();
                    if (var instanceof VolVariable) {
                        dataSetIdentifierList.addElement(new DataSetIdentifier(var.getName(), VariableType.VOLUME, var.getDomain()));
                    } else if (var instanceof MemVariable) {
                        dataSetIdentifierList.addElement(new DataSetIdentifier(var.getName(), VariableType.MEMBRANE, var.getDomain()));
                    } else if (var instanceof Function) {
                        VariableType varType = VariableType.UNKNOWN;
                        if (var.getDomain() != null && var.getDomain().getName() != null) {
                            SubDomain subDomain = simTask.getSimulation().getMathDescription().getSubDomain(var.getDomain().getName());
                            if (subDomain instanceof CompartmentSubDomain) {
                                varType = VariableType.VOLUME;
                            } else if (subDomain instanceof MembraneSubDomain) {
                                varType = VariableType.MEMBRANE;
                            } else if (subDomain instanceof FilamentSubDomain) {
                                throw new RuntimeException("filament subdomains not supported");
                            } else if (subDomain instanceof PointSubDomain) {
                                varType = VariableType.POINT_VARIABLE;
                            }
                        }
                        dataSetIdentifierList.addElement(new DataSetIdentifier(var.getName(), varType, var.getDomain()));
                    } else if (var instanceof Constant) {
                        System.out.println("ignoring Constant " + var.getName());
                    } else if (var instanceof InsideVariable) {
                        System.out.println("ignoring InsideVariable " + var.getName());
                    } else if (var instanceof OutsideVariable) {
                        System.out.println("ignoring OutsideVariable " + var.getName());
                    } else {
                        throw new RuntimeException("unexpected variable " + var.getName() + " of type " + var.getClass().getName());
                    }
                }
            } catch (XmlParseException | ExpressionException e) {
                e.printStackTrace();
                throw new RuntimeException("failed to read sim task file, msg: " + e.getMessage(), e);
            }
        }
    }
    DataIdentifier[] dis = new DataIdentifier[dataSetIdentifierList.size()];
    for (int i = 0; i < dataSetIdentifierList.size(); i++) {
        DataSetIdentifier dsi = (DataSetIdentifier) dataSetIdentifierList.elementAt(i);
        String displayName = dsi.getName();
        if (dsi.isFunction()) {
            AnnotatedFunction f = null;
            for (int j = 0; j < annotatedFunctionList.size(); j++) {
                AnnotatedFunction function = (AnnotatedFunction) annotatedFunctionList.elementAt(j);
                if (function.getName().equals(dsi.getName())) {
                    f = function;
                    break;
                }
            }
            if (f != null) {
                displayName = f.getDisplayName();
            }
        }
        dis[i] = new DataIdentifier(dsi.getName(), dsi.getVariableType(), dsi.getDomain(), dsi.isFunction(), displayName);
    }
    return dis;
}
Also used : MembraneSubDomain(cbit.vcell.math.MembraneSubDomain) SimulationTask(cbit.vcell.messaging.server.SimulationTask) InsideVariable(cbit.vcell.math.InsideVariable) VolVariable(cbit.vcell.math.VolVariable) ReservedVariable(cbit.vcell.math.ReservedVariable) MemVariable(cbit.vcell.math.MemVariable) OutsideVariable(cbit.vcell.math.OutsideVariable) Variable(cbit.vcell.math.Variable) VCDataIdentifier(org.vcell.util.document.VCDataIdentifier) VCSimulationDataIdentifier(cbit.vcell.solver.VCSimulationDataIdentifier) ExternalDataIdentifier(org.vcell.util.document.ExternalDataIdentifier) Constant(cbit.vcell.math.Constant) FileNotFoundException(java.io.FileNotFoundException) PointSubDomain(cbit.vcell.math.PointSubDomain) ODESimData(cbit.vcell.solver.ode.ODESimData) InsideVariable(cbit.vcell.math.InsideVariable) ExpressionException(cbit.vcell.parser.ExpressionException) CompartmentSubDomain(cbit.vcell.math.CompartmentSubDomain) FilamentSubDomain(cbit.vcell.math.FilamentSubDomain) SubDomain(cbit.vcell.math.SubDomain) PointSubDomain(cbit.vcell.math.PointSubDomain) MembraneSubDomain(cbit.vcell.math.MembraneSubDomain) Function(cbit.vcell.math.Function) AnnotatedFunction(cbit.vcell.solver.AnnotatedFunction) MemVariable(cbit.vcell.math.MemVariable) DataAccessException(org.vcell.util.DataAccessException) AnnotatedFunction(cbit.vcell.solver.AnnotatedFunction) ComsolSimFiles(org.vcell.vis.io.ComsolSimFiles) VariableType(cbit.vcell.math.VariableType) VolVariable(cbit.vcell.math.VolVariable) XmlParseException(cbit.vcell.xml.XmlParseException) FilamentSubDomain(cbit.vcell.math.FilamentSubDomain) MathException(cbit.vcell.math.MathException) CompartmentSubDomain(cbit.vcell.math.CompartmentSubDomain) OutsideVariable(cbit.vcell.math.OutsideVariable) CompartmentSubDomain(cbit.vcell.math.CompartmentSubDomain) FilamentSubDomain(cbit.vcell.math.FilamentSubDomain) SubDomain(cbit.vcell.math.SubDomain) PointSubDomain(cbit.vcell.math.PointSubDomain) Domain(cbit.vcell.math.Variable.Domain) MembraneSubDomain(cbit.vcell.math.MembraneSubDomain) ZipFile(org.apache.commons.compress.archivers.zip.ZipFile) File(java.io.File)

Example 32 with Variable

use of cbit.vcell.math.Variable in project vcell by virtualcell.

the class NFsimXMLWriter method writeNFsimXML.

public static Element writeNFsimXML(SimulationTask origSimTask, long randomSeed, NFsimSimulationOptions nfsimSimulationOptions, boolean bUseLocationMarks) throws SolverException {
    try {
        System.out.println("VCML ORIGINAL .... START\n" + origSimTask.getSimulation().getMathDescription().getVCML_database() + "\nVCML ORIGINAL .... END\n====================\n");
    } catch (MathException e1) {
        // TODO Auto-generated catch block
        e1.printStackTrace();
    }
    SimulationTask clonedSimTask = null;
    try {
        clonedSimTask = (SimulationTask) BeanUtils.cloneSerializable(origSimTask);
    } catch (Exception eee) {
        throw new SolverException("failed to clone mathDescription while preparing NFSim input: " + eee.getMessage(), eee);
    }
    MathDescription clonedMathDesc = clonedSimTask.getSimulation().getMathDescription();
    if (bUseLocationMarks) {
        try {
            // 
            // get list of Compartment Names (stored in locations).
            // 
            ArrayList<String> locations = new ArrayList<String>();
            Enumeration<Variable> varEnum = clonedMathDesc.getVariables();
            ArrayList<VolumeParticleSpeciesPattern> volumeParticleSpeciesPatterns = new ArrayList<VolumeParticleSpeciesPattern>();
            while (varEnum.hasMoreElements()) {
                Variable var = varEnum.nextElement();
                if (var instanceof VolumeParticleSpeciesPattern) {
                    VolumeParticleSpeciesPattern speciesPattern = (VolumeParticleSpeciesPattern) var;
                    if (!locations.contains(speciesPattern.getLocationName())) {
                        locations.add(speciesPattern.getLocationName());
                    }
                    volumeParticleSpeciesPatterns.add(speciesPattern);
                }
            }
            // 
            for (ParticleMolecularType particleMolecularType : clonedMathDesc.getParticleMolecularTypes()) {
                String pmcLocationName = RbmUtils.SiteStruct;
                String pmcLocationId = particleMolecularType.getName() + "_" + RbmUtils.SiteStruct;
                ParticleMolecularComponent locationComponent = new ParticleMolecularComponent(pmcLocationId, pmcLocationName);
                for (String location : locations) {
                    locationComponent.addComponentStateDefinition(new ParticleComponentStateDefinition(location));
                }
                particleMolecularType.insertMolecularComponent(0, locationComponent);
                String pmcMarkName = RbmUtils.SiteProduct;
                String pmcMarkId = particleMolecularType.getName() + "_" + RbmUtils.SiteProduct;
                ParticleMolecularComponent markComponent = new ParticleMolecularComponent(pmcMarkId, pmcMarkName);
                markComponent.addComponentStateDefinition(new ParticleComponentStateDefinition("0"));
                markComponent.addComponentStateDefinition(new ParticleComponentStateDefinition("1"));
                particleMolecularType.insertMolecularComponent(1, markComponent);
            }
            // 
            for (VolumeParticleSpeciesPattern speciesPattern : volumeParticleSpeciesPatterns) {
                for (ParticleMolecularTypePattern molTypePattern : speciesPattern.getParticleMolecularTypePatterns()) {
                    // 
                    // add location component to pattern ... state=<location>
                    // 
                    {
                        final ParticleMolecularComponent locationComponentDefinition = molTypePattern.getMolecularType().getComponentList().get(0);
                        ParticleMolecularComponentPattern locationPattern = new ParticleMolecularComponentPattern(locationComponentDefinition);
                        ParticleComponentStateDefinition locationStateDefinition = null;
                        for (ParticleComponentStateDefinition stateDef : locationComponentDefinition.getComponentStateDefinitions()) {
                            if (stateDef.getName().equals(speciesPattern.getLocationName())) {
                                locationStateDefinition = stateDef;
                            }
                        }
                        ParticleComponentStatePattern locationStatePattern = new ParticleComponentStatePattern(locationStateDefinition);
                        locationPattern.setComponentStatePattern(locationStatePattern);
                        locationPattern.setBondType(ParticleBondType.None);
                        locationPattern.setBondId(-1);
                        molTypePattern.insertMolecularComponentPattern(0, locationPattern);
                    }
                    // 
                    // add mark component to pattern ... state="0" (for observables and reactants ... later we will clone and use "1" for products).
                    {
                        final ParticleMolecularComponent markComponentDefinition = molTypePattern.getMolecularType().getComponentList().get(1);
                        ParticleMolecularComponentPattern markPattern = new ParticleMolecularComponentPattern(markComponentDefinition);
                        final int clearStateIndex = 0;
                        final int setStateIndex = 1;
                        ParticleComponentStateDefinition markStateClearedDefinition = markComponentDefinition.getComponentStateDefinitions().get(clearStateIndex);
                        ParticleComponentStatePattern markStatePattern = new ParticleComponentStatePattern(markStateClearedDefinition);
                        markPattern.setComponentStatePattern(markStatePattern);
                        markPattern.setBondType(ParticleBondType.None);
                        markPattern.setBondId(-1);
                        molTypePattern.insertMolecularComponentPattern(1, markPattern);
                    }
                }
            }
            // 
            // when processing ParticleJumpProcesses, we add a new "product" species pattern (by cloning the original speciesPattern)
            // and setting the mark site to "1", change name to name+"_PRODUCT", and add to math model if it doesn't already exist.
            // 
            // cloned the "standard" reactant/observable speciesPattern, set the mark for all molecules, and add to mathDesc.
            // 
            CompartmentSubDomain subDomain = (CompartmentSubDomain) clonedMathDesc.getSubDomains().nextElement();
            for (ParticleJumpProcess particleJumpProcess : subDomain.getParticleJumpProcesses()) {
                for (Action action : particleJumpProcess.getActions()) {
                    if (action.getOperation().equals(Action.ACTION_CREATE)) {
                        VolumeParticleSpeciesPattern volumeParticleSpeciesPattern = (VolumeParticleSpeciesPattern) action.getVar();
                        String newSpeciesPatternName = volumeParticleSpeciesPattern.getName() + "_" + particleJumpProcess.getName();
                        VolumeParticleSpeciesPattern productPattern = new VolumeParticleSpeciesPattern(volumeParticleSpeciesPattern, newSpeciesPatternName);
                        // VolumeParticleSpeciesPattern productPattern = (VolumeParticleSpeciesPattern) BeanUtils.cloneSerializable(volumeParticleSpeciesPattern);
                        for (ParticleMolecularTypePattern productMolTypePattern : productPattern.getParticleMolecularTypePatterns()) {
                            ParticleComponentStateDefinition markSet = productMolTypePattern.getMolecularType().getComponentList().get(1).getComponentStateDefinitions().get(1);
                            productMolTypePattern.getMolecularComponentPatternList().get(1).setComponentStatePattern(new ParticleComponentStatePattern(markSet));
                        }
                        System.out.println(productPattern.getName());
                        if (clonedMathDesc.getVariable(productPattern.getName()) == null) {
                            clonedMathDesc.addVariable(productPattern);
                        }
                        action.setVar(productPattern);
                    }
                }
            }
            try {
                System.out.println("===============================\n ----------- VCML HACKED .... START\n" + clonedMathDesc.getVCML_database() + "\nVCML HACKED .... END\n====================\n");
            } catch (MathException e1) {
                // TODO Auto-generated catch block
                e1.printStackTrace();
            }
        } catch (Exception e) {
            throw new SolverException("failed to apply location mark transformation: " + e.getMessage(), e);
        }
    }
    Element sbmlElement = new Element("sbml");
    Element modelElement = new Element("model");
    modelElement.setAttribute("id", "nameless");
    SimulationSymbolTable simulationSymbolTable = new SimulationSymbolTable(clonedSimTask.getSimulation(), clonedSimTask.getSimulationJob().getJobIndex());
    Element listOfParametersElement = getListOfParameters(clonedMathDesc, simulationSymbolTable);
    Element listOfMoleculeTypesElement = getListOfMoleculeTypes(clonedMathDesc);
    Element listOfSpeciesElement = getListOfSpecies(clonedMathDesc, simulationSymbolTable);
    CompartmentSubDomain compartmentSubDomain = (CompartmentSubDomain) clonedMathDesc.getSubDomains().nextElement();
    Element listOfReactionRules = new Element("ListOfReactionRules");
    for (int reactionRuleIndex = 0; reactionRuleIndex < compartmentSubDomain.getParticleJumpProcesses().size(); reactionRuleIndex++) {
        ParticleJumpProcess particleJumpProcess = compartmentSubDomain.getParticleJumpProcesses().get(reactionRuleIndex);
        MathRuleFactory mathRuleFactory = new MathRuleFactory();
        MathRuleEntry rule = mathRuleFactory.createRuleEntry(particleJumpProcess, reactionRuleIndex);
        RuleAnalysisReport report = RuleAnalysis.analyze(rule, true);
        // remember, we have to add RateLaw
        Element reactionRuleElement = RuleAnalysis.getNFSimXML(rule, report);
        // ArrayList<MolecularTypeOfReactionParticipant> currentReactantElementsOfReaction = new ArrayList<MolecularTypeOfReactionParticipant>();
        // ArrayList<ComponentOfMolecularTypeOfReactionParticipant> currentComponentOfReactantElementsOfReaction = new ArrayList<ComponentOfMolecularTypeOfReactionParticipant>();
        // ArrayList<MolecularTypeOfReactionParticipant> currentProductElementsOfReaction = new ArrayList<MolecularTypeOfReactionParticipant>();
        // ArrayList<ComponentOfMolecularTypeOfReactionParticipant> currentComponentOfProductElementsOfReaction = new ArrayList<ComponentOfMolecularTypeOfReactionParticipant>();
        // currentMappingOfReactionParticipants.clear();
        // reactionProductBondSites.clear();
        // reactionReactantBondSites.clear();
        // 
        // Element reactionRuleElement = new Element("ReactionRule");
        // String reactionRuleID = "RR" + (reactionRuleIndex + 1);
        // reactionRuleElement.setAttribute("id",reactionRuleID);
        // reactionRuleElement.setAttribute("name",particleJumpProcess.getName());
        // reactionRuleElement.setAttribute("symmetry_factor","1");
        // reactionRule.resolveBonds();
        // 
        // ArrayList<VolumeParticleSpeciesPattern> selectedPatterns = new ArrayList<VolumeParticleSpeciesPattern>();
        // for (ParticleVariable particleVariable : particleJumpProcess.getParticleVariables()){
        // if (!(particleVariable instanceof VolumeParticleSpeciesPattern)){
        // throw new SolverException("expecting only "+VolumeParticleSpeciesPattern.class.getSimpleName()+"s for "+ParticleJumpProcess.class.getSimpleName()+" "+particleJumpProcess.getName());
        // }
        // selectedPatterns.add((VolumeParticleSpeciesPattern) particleVariable);
        // }
        // ArrayList<VolumeParticleSpeciesPattern> createdPatterns = new ArrayList<VolumeParticleSpeciesPattern>();
        // HashSet<VolumeParticleSpeciesPattern> destroyedPatterns = new HashSet<VolumeParticleSpeciesPattern>();
        // for (Action action : particleJumpProcess.getActions()){
        // if (!(action.getVar() instanceof VolumeParticleSpeciesPattern)){
        // throw new SolverException("expecting only "+VolumeParticleSpeciesPattern.class.getSimpleName()+"s for "+ParticleJumpProcess.class.getSimpleName()+" "+particleJumpProcess.getName());
        // }
        // if (action.getOperation().equals(Action.ACTION_CREATE)){
        // createdPatterns.add((VolumeParticleSpeciesPattern) action.getVar());
        // }else if (action.getOperation().equals(Action.ACTION_DESTROY)){
        // destroyedPatterns.add((VolumeParticleSpeciesPattern) action.getVar());
        // }else{
        // throw new RuntimeException("unexpected action operation "+action.getOperation()+" for jump process "+particleJumpProcess.getName());
        // }
        // }
        // 
        // Element listOfReactantPatternsElement = new Element("ListOfReactantPatterns");
        // for(int reactantPatternIndex=0; reactantPatternIndex < selectedPatterns.size(); reactantPatternIndex++) {
        // VolumeParticleSpeciesPattern reactantSpeciesPattern = selectedPatterns.get(reactantPatternIndex);
        // String reactantPatternID = "RP" + (reactantPatternIndex + 1);
        // patternReactantBondSites.clear();
        // Element reactantPatternElement = getReactionParticipantPattern1(reactionRuleID, reactantPatternID, reactantSpeciesPattern,
        // currentReactantElementsOfReaction, currentComponentOfReactantElementsOfReaction, "ReactantPattern");
        // listOfReactantPatternsElement.addContent(reactantPatternElement);
        // reactionReactantBondSites.addAll(patternReactantBondSites);
        // }
        // reactionRuleElement.addContent(listOfReactantPatternsElement);
        // 
        // Element listOfProductPatternsElement = new Element("ListOfProductPatterns");
        // ArrayList<VolumeParticleSpeciesPattern> productSpeciesPatterns = new ArrayList<VolumeParticleSpeciesPattern>(selectedPatterns);
        // productSpeciesPatterns.removeAll(destroyedPatterns);
        // productSpeciesPatterns.addAll(createdPatterns);
        // // for products, add all "created" species from Actions and all "particles" that are selected but not destroyed
        // for(int productPatternIndex=0; productPatternIndex < productSpeciesPatterns.size(); productPatternIndex++) {
        // VolumeParticleSpeciesPattern productSpeciesPattern = productSpeciesPatterns.get(productPatternIndex);
        // String productPatternID = "PP" + (productPatternIndex + 1);
        // patternProductBondSites.clear();
        // Element productPatternElement = getReactionParticipantPattern1(reactionRuleID, productPatternID, productSpeciesPattern,
        // currentProductElementsOfReaction, currentComponentOfProductElementsOfReaction, "ProductPattern");
        // listOfProductPatternsElement.addContent(productPatternElement);
        // reactionProductBondSites.addAll(patternProductBondSites);
        // }
        // reactionRuleElement.addContent(listOfProductPatternsElement);
        // <RateLaw id="RR1_RateLaw" type="Ele" totalrate="0">
        // <ListOfRateConstants>
        // <RateConstant value="kon"/>
        // </ListOfRateConstants>
        // </RateLaw>
        Element rateLawElement = new Element("RateLaw");
        rateLawElement.setAttribute("id", RuleAnalysis.getID(rule));
        String rateConstantValue = null;
        JumpProcessRateDefinition particleProbabilityRate = particleJumpProcess.getParticleRateDefinition();
        if (particleProbabilityRate.getExpressions().length > 0) {
            JumpProcessRateDefinition particleRateDefinition = particleJumpProcess.getParticleRateDefinition();
            Expression expression = null;
            if (particleRateDefinition instanceof MacroscopicRateConstant) {
                expression = ((MacroscopicRateConstant) particleProbabilityRate).getExpression();
            } else {
                throw new SolverException("ParticleRateDefinition type " + particleRateDefinition.getClass().getSimpleName() + " not supported");
            }
            rateConstantValue = expression.infixBng();
            // all rates constants are being flattened and given reserved names
            Expression substitutedValExpr = null;
            try {
                substitutedValExpr = simulationSymbolTable.substituteFunctions(expression);
            } catch (MathException | ExpressionException e) {
                e.printStackTrace(System.out);
                throw new SolverException("ParticleJumpProcess " + particleJumpProcess.getName() + " substitution failed : exp = \"" + expression.infix() + "\": " + e.getMessage());
            }
            Double value = null;
            try {
                value = substitutedValExpr.evaluateConstant();
                Element parameterElement = new Element("Parameter");
                String id = "K_reserved_" + reactionRuleIndex;
                parameterElement.setAttribute("id", id);
                if (value != null) {
                    parameterElement.setAttribute("type", "Constant");
                    parameterElement.setAttribute("value", value.toString());
                    parameterElement.addContent(new Comment(rateConstantValue));
                    rateConstantValue = id;
                    listOfParametersElement.addContent(parameterElement);
                }
            } catch (ExpressionException e) {
                System.out.println("ParticleJumpProcess " + particleJumpProcess.getName() + " = " + substitutedValExpr.infix() + " does not have a constant value");
            }
        }
        if (isFunction(rateConstantValue, clonedMathDesc, simulationSymbolTable)) {
            rateLawElement.setAttribute("type", "Function");
            rateLawElement.setAttribute("totalrate", "0");
            rateLawElement.setAttribute("name", rateConstantValue);
        } else {
            rateLawElement.setAttribute("type", "Ele");
            rateLawElement.setAttribute("totalrate", "0");
            Element listOfRateConstantsElement = new Element("ListOfRateConstants");
            Element rateConstantElement = new Element("RateConstant");
            // System.out.println(" --- " + particleJumpProcess.getParticleRateDefinition().getExpressions());
            if (particleProbabilityRate.getExpressions().length > 0) {
                rateConstantElement.setAttribute("value", rateConstantValue);
            }
            listOfRateConstantsElement.addContent(rateConstantElement);
            rateLawElement.addContent(listOfRateConstantsElement);
        }
        reactionRuleElement.addContent(rateLawElement);
        // //  <Map>
        // //    <MapItem sourceID="RR1_RP1_M1" targetID="RR1_PP1_M1"/>
        // //    <MapItem sourceID="RR1_RP1_M1_C1" targetID="RR1_PP1_M1_C1"/>
        // //    <MapItem sourceID="RR1_RP1_M1_C2" targetID="RR1_PP1_M1_C2"/>
        // //    <MapItem sourceID="RR1_RP2_M1" targetID="RR1_PP1_M2"/>
        // //    <MapItem sourceID="RR1_RP2_M1_C1" targetID="RR1_PP1_M2_C1"/>
        // //  </Map>
        // Element mapElement = new Element("Map");
        // System.out.println("----------------------------------------------------------------------");
        // for(MolecularTypeOfReactionParticipant p : currentReactantElementsOfReaction) {
        // System.out.println(p.moleculeName + ", " + p.elementID);
        // }
        // for(ComponentOfMolecularTypeOfReactionParticipant c : currentComponentOfReactantElementsOfReaction) {
        // System.out.println(c.moleculeName + ", " + c.componentName + ", " + c.elementID);
        // }
        // System.out.println("----------------------------------------------------------------------");
        // for(MolecularTypeOfReactionParticipant p : currentProductElementsOfReaction) {
        // System.out.println(p.moleculeName + ", " + p.elementID);
        // }
        // for(ComponentOfMolecularTypeOfReactionParticipant c : currentComponentOfProductElementsOfReaction) {
        // System.out.println(c.moleculeName + ", " + c.componentName + ", " + c.elementID);
        // }
        // System.out.println("----------------------------------------------------------------------");
        // 
        // List<MolecularTypeOfReactionParticipant> cloneOfReactants = new ArrayList<MolecularTypeOfReactionParticipant>(currentReactantElementsOfReaction);
        // List<MolecularTypeOfReactionParticipant> cloneOfProducts = new ArrayList<MolecularTypeOfReactionParticipant>(currentProductElementsOfReaction);
        // for(Iterator<MolecularTypeOfReactionParticipant> itReactant = cloneOfReactants.iterator(); itReactant.hasNext();) {	// participants
        // MolecularTypeOfReactionParticipant reactant = itReactant.next();
        // boolean foundProduct = false;
        // for(Iterator<MolecularTypeOfReactionParticipant> itProduct = cloneOfProducts.iterator(); itProduct.hasNext();) {
        // MolecularTypeOfReactionParticipant product = itProduct.next();
        // if(reactant.find(product)) {
        // MappingOfReactionParticipants m = new MappingOfReactionParticipants(reactant.elementID, product.elementID, "");
        // currentMappingOfReactionParticipants.add(m );
        // itProduct.remove();
        // foundProduct = true;
        // break;		// we exit inner loop if we find a match for current reactant
        // }
        // }
        // if(foundProduct == false) {
        // System.out.println("Did not found a match for reactant " + reactant.moleculeName + ", " + reactant.elementID);
        // }
        // itReactant.remove();		// found or not, we remove the reactant
        // }
        // if(!currentProductElementsOfReaction.isEmpty()) {
        // for(MolecularTypeOfReactionParticipant p : currentProductElementsOfReaction) {
        // System.out.println("Did not found a match for product " + p.moleculeName + ", " + p.elementID);
        // }
        // }
        // for(Iterator<ComponentOfMolecularTypeOfReactionParticipant> itReactant = currentComponentOfReactantElementsOfReaction.iterator(); itReactant.hasNext();) {	// components
        // ComponentOfMolecularTypeOfReactionParticipant reactant = itReactant.next();
        // boolean foundProduct = false;
        // for(Iterator<ComponentOfMolecularTypeOfReactionParticipant> itProduct = currentComponentOfProductElementsOfReaction.iterator(); itProduct.hasNext();) {
        // ComponentOfMolecularTypeOfReactionParticipant product = itProduct.next();
        // String state = "";
        // if(reactant.find(product)) {
        // if(!reactant.state.equals(product.state)) {
        // state = product.state;
        // }
        // MappingOfReactionParticipants m = new MappingOfReactionParticipants(reactant.elementID, product.elementID, state);
        // currentMappingOfReactionParticipants.add(m );
        // itProduct.remove();
        // foundProduct = true;
        // break;		// we exit inner loop if we find a match for current reactant
        // }
        // }
        // if(foundProduct == false) {
        // System.out.println("Did not found a match for reactant " + reactant.moleculeName + ", " + reactant.elementID);
        // }
        // itReactant.remove();		// found or not, we remove the reactant
        // }
        // if(!currentComponentOfProductElementsOfReaction.isEmpty()) {
        // for(ComponentOfMolecularTypeOfReactionParticipant p : currentComponentOfProductElementsOfReaction) {
        // System.out.println("Did not found a match for product " + p.moleculeName + ", " + p.elementID);
        // }
        // }
        // for(Iterator<MappingOfReactionParticipants> it = currentMappingOfReactionParticipants.iterator(); it.hasNext();) {
        // MappingOfReactionParticipants m = it.next();
        // Element mapItemElement = new Element("MapItem");
        // mapItemElement.setAttribute("sourceID", m.reactantElementID);
        // mapItemElement.setAttribute("targetID", m.productElementID);
        // mapElement.addContent(mapItemElement);
        // }
        // reactionRuleElement.addContent(mapElement);
        // 
        // //  <ListOfOperations>
        // //      <AddBond site1="RR1_RP1_M1_C1" site2="RR1_RP2_M1_C1"/>
        // //		<StateChange site="RR0_RP0_M0_C2" finalState="Y"/>
        // //  </ListOfOperations>
        // Element listOfOperationsElement = new Element("ListOfOperations");
        // 
        // // AddBond elements
        // // add any bond in the product which is not present in the reactant
        // Iterator<BondSites> it = patternProductBondSites.iterator();
        // while (it.hasNext()) {
        // BondSites bs = it.next();
        // String reactantS1 = MappingOfReactionParticipants.findMatchingReactant(bs.component1, currentMappingOfReactionParticipants);
        // String reactantS2 = MappingOfReactionParticipants.findMatchingReactant(bs.component2, currentMappingOfReactionParticipants);
        // // we check if the bonds in the product existed already in the reactant, in which case they were not "added" in this reaction
        // BondSites candidate = new BondSites(reactantS1, reactantS2);
        // boolean preExistent = false;
        // for(BondSites bsReactant : reactionReactantBondSites) {
        // if(bsReactant.equals(candidate)) {
        // preExistent = true;
        // break;
        // }
        // }
        // if(preExistent == true) {
        // continue;		// we don't add preexisting bonds
        // }
        // Element addBondElement = new Element("AddBond");
        // addBondElement.setAttribute("site1", reactantS1);
        // addBondElement.setAttribute("site2", reactantS2);
        // listOfOperationsElement.addContent(addBondElement);
        // }
        // // StateChange elements
        // for(Iterator<MappingOfReactionParticipants> it1 = currentMappingOfReactionParticipants.iterator(); it1.hasNext();) {
        // MappingOfReactionParticipants m = it1.next();
        // if(!m.componentFinalState.equals("")) {		// state has changed if it's different from ""
        // Element stateChangeElement = new Element("StateChange");
        // stateChangeElement.setAttribute("site", m.reactantElementID);
        // stateChangeElement.setAttribute("finalState", m.componentFinalState);
        // listOfOperationsElement.addContent(stateChangeElement);
        // }
        // }
        // // eliminate all the common entries (molecule types) in reactants and products
        // // what's left in reactants was deleted, what's left in products was added
        // List<MolecularTypeOfReactionParticipant> commonParticipants = new ArrayList<MolecularTypeOfReactionParticipant>();
        // for(Iterator<MolecularTypeOfReactionParticipant> itReactant = currentReactantElementsOfReaction.iterator(); itReactant.hasNext();) {	// participants
        // MolecularTypeOfReactionParticipant reactant = itReactant.next();
        // for(Iterator<MolecularTypeOfReactionParticipant> itProduct = currentProductElementsOfReaction.iterator(); itProduct.hasNext();) {
        // MolecularTypeOfReactionParticipant product = itProduct.next();
        // if(reactant.find(product)) {
        // // commonParticipants contains the reactant molecules with a equivalent molecule in the product (meaning they are not in the "Deleted" category)
        // commonParticipants.add(reactant);
        // itReactant.remove();
        // itProduct.remove();
        // break;		// we exit inner loop if we find a match for current reactant
        // }
        // }
        // }
        // // DeleteBond element
        // // there is no need to mention deletion of bond if the particleSpeciesPattern
        // // or the MolecularType involved in the bond are deleted as well
        // // We only keep those "Deleted" bonds which belong to the molecules (of the reactant) present in commonParticipants
        // // Both components (sites) of the bond need to have their molecules in commonParticipants
        // boolean foundMoleculeForComponent1 = false;
        // boolean foundMoleculeForComponent2 = false;
        // HashSet<BondSites> cloneOfReactantBondSites = new HashSet<BondSites>(patternReactantBondSites);
        // Iterator<BondSites> itbs = cloneOfReactantBondSites.iterator();
        // while (itbs.hasNext()) {
        // BondSites bs = itbs.next();
        // String bondComponent1MoleculeId = BondSites.extractMoleculeId(bs.component1);
        // String bondComponent2MoleculeId = BondSites.extractMoleculeId(bs.component2);
        // for(MolecularTypeOfReactionParticipant commonReactionMoleculeule : commonParticipants) {
        // String commonReactantPatternId = commonReactionMoleculeule.elementID;
        // if(bondComponent1MoleculeId.equals(commonReactantPatternId)) {
        // foundMoleculeForComponent1 = true;
        // }
        // if(bondComponent2MoleculeId.equals(commonReactantPatternId)) {
        // foundMoleculeForComponent2 = true;
        // }
        // }
        // if(!foundMoleculeForComponent1 || !foundMoleculeForComponent2) {
        // // at least one of bond's molecule is not in common, hence we don't need to report the deletion of this bond
        // itbs.remove();
        // }
        // }
        // // the clone has now all the deleted bonds whose molecules have not been deleted
        // itbs = cloneOfReactantBondSites.iterator();
        // while (itbs.hasNext()) {
        // BondSites bs = itbs.next();
        // Element addBondElement = new Element("DeleteBond");
        // addBondElement.setAttribute("site1", bs.component1);
        // addBondElement.setAttribute("site2", bs.component2);
        // listOfOperationsElement.addContent(addBondElement);
        // }
        // // Add MolecularType element
        // for(MolecularTypeOfReactionParticipant molecule : currentProductElementsOfReaction) {
        // System.out.println("created molecule: " + molecule.elementID + "' " + molecule.moleculeName);
        // Element addMolecularTypePatternElement = new Element("Add");
        // addMolecularTypePatternElement.setAttribute("id", molecule.elementID);
        // listOfOperationsElement.addContent(addMolecularTypePatternElement);
        // }
        // // Delete MolecularType element
        // // if the reactant pattern of the molecule being deleted still exists as part of the common, then we only delete the molecule
        // // if the reactant pattern of the molecule being deleted is not as part of the common, then it's gone completely and we delete the reactant pattern
        // ArrayList<String> patternsToDelete = new ArrayList<String>();
        // for(MolecularTypeOfReactionParticipant molecule : currentReactantElementsOfReaction) {
        // String reactantPatternId = molecule.extractReactantPatternId();
        // boolean found = false;
        // for(MolecularTypeOfReactionParticipant common : commonParticipants) {
        // String commonId = common.extractReactantPatternId();
        // if(reactantPatternId.equals(commonId)) {
        // found = true;
        // break;		// some other molecule of this pattern still there, we don't delete the pattern
        // }
        // }
        // if(found == true) {		// some other molecule of this pattern still there, we don't delete the pattern
        // System.out.println("deleted molecule: " + molecule.elementID + "' " + molecule.moleculeName);
        // Element addMolecularTypePatternElement = new Element("Delete");
        // addMolecularTypePatternElement.setAttribute("id", molecule.elementID);
        // addMolecularTypePatternElement.setAttribute("DeleteMolecules", "0");
        // listOfOperationsElement.addContent(addMolecularTypePatternElement);
        // } else {				// no molecule of this pattern left, we delete the pattern
        // if(patternsToDelete.contains(reactantPatternId)) {
        // // nothing to do, we're already deleting this pattern
        // break;
        // } else {
        // patternsToDelete.add(reactantPatternId);
        // System.out.println("deleted pattern: " + reactantPatternId);
        // Element addParticleSpeciesPatternElement = new Element("Delete");
        // addParticleSpeciesPatternElement.setAttribute("id", reactantPatternId);
        // addParticleSpeciesPatternElement.setAttribute("DeleteMolecules", "0");
        // listOfOperationsElement.addContent(addParticleSpeciesPatternElement);
        // }
        // }
        // }
        // reactionRuleElement.addContent(listOfOperationsElement);
        listOfReactionRules.addContent(reactionRuleElement);
    }
    Element listOfObservablesElement = getListOfObservables(clonedMathDesc);
    Element listOfFunctionsElement = getListOfFunctions(clonedMathDesc, simulationSymbolTable);
    modelElement.addContent(listOfParametersElement);
    modelElement.addContent(listOfMoleculeTypesElement);
    modelElement.addContent(listOfSpeciesElement);
    modelElement.addContent(listOfReactionRules);
    modelElement.addContent(listOfObservablesElement);
    modelElement.addContent(listOfFunctionsElement);
    sbmlElement.addContent(modelElement);
    // //		return e1;
    return sbmlElement;
}
Also used : Action(cbit.vcell.math.Action) SimulationTask(cbit.vcell.messaging.server.SimulationTask) Variable(cbit.vcell.math.Variable) MathDescription(cbit.vcell.math.MathDescription) Element(org.jdom.Element) ArrayList(java.util.ArrayList) ParticleMolecularComponent(cbit.vcell.math.ParticleMolecularComponent) ExpressionException(cbit.vcell.parser.ExpressionException) ParticleComponentStateDefinition(cbit.vcell.math.ParticleComponentStateDefinition) ParticleMolecularComponentPattern(cbit.vcell.math.ParticleMolecularComponentPattern) MacroscopicRateConstant(cbit.vcell.math.MacroscopicRateConstant) MathRuleFactory(cbit.vcell.math.MathRuleFactory) ParticleMolecularType(cbit.vcell.math.ParticleMolecularType) RuleAnalysisReport(org.vcell.model.rbm.RuleAnalysisReport) Comment(org.jdom.Comment) JumpProcessRateDefinition(cbit.vcell.math.JumpProcessRateDefinition) ParticleComponentStatePattern(cbit.vcell.math.ParticleComponentStatePattern) ParticleJumpProcess(cbit.vcell.math.ParticleJumpProcess) VolumeParticleSpeciesPattern(cbit.vcell.math.VolumeParticleSpeciesPattern) SimulationSymbolTable(cbit.vcell.solver.SimulationSymbolTable) SolverException(cbit.vcell.solver.SolverException) ExpressionException(cbit.vcell.parser.ExpressionException) MathException(cbit.vcell.math.MathException) ParticleMolecularTypePattern(cbit.vcell.math.ParticleMolecularTypePattern) MathRuleEntry(cbit.vcell.math.MathRuleFactory.MathRuleEntry) Expression(cbit.vcell.parser.Expression) MathException(cbit.vcell.math.MathException) CompartmentSubDomain(cbit.vcell.math.CompartmentSubDomain) SolverException(cbit.vcell.solver.SolverException)

Example 33 with Variable

use of cbit.vcell.math.Variable in project vcell by virtualcell.

the class NFsimXMLWriter method getListOfObservables.

private static Element getListOfObservables(MathDescription mathDesc) throws SolverException {
    Element listOfObservablesElement = new Element("ListOfObservables");
    int observableIndex = 0;
    Enumeration<Variable> enum1 = mathDesc.getVariables();
    while (enum1.hasMoreElements()) {
        Variable var = enum1.nextElement();
        if (var instanceof ParticleObservable) {
            ParticleObservable particleObservable = (ParticleObservable) var;
            Element observableElement = new Element("Observable");
            String observableId = "O" + (observableIndex + 1);
            observableElement.setAttribute("id", observableId);
            observableElement.setAttribute("name", particleObservable.getName());
            observableElement.setAttribute("type", particleObservable.getType().getText());
            Element listOfPatterns = getParticleSpeciesPatternList(observableId, particleObservable);
            observableElement.addContent(listOfPatterns);
            listOfObservablesElement.addContent(observableElement);
            observableIndex++;
        }
    }
    Format format = Format.getPrettyFormat();
    XMLOutputter outp = new XMLOutputter(format);
    String sOurs = outp.outputString(listOfObservablesElement);
    return listOfObservablesElement;
}
Also used : XMLOutputter(org.jdom.output.XMLOutputter) Variable(cbit.vcell.math.Variable) Format(org.jdom.output.Format) Element(org.jdom.Element) ParticleObservable(cbit.vcell.math.ParticleObservable)

Example 34 with Variable

use of cbit.vcell.math.Variable in project vcell by virtualcell.

the class SmoldynFileWriter method init.

private void init() throws SolverException {
    simulation = simTask.getSimulation();
    mathDesc = simulation.getMathDescription();
    simulationSymbolTable = simTask.getSimulationJob().getSimulationSymbolTable();
    particleVariableList = new ArrayList<ParticleVariable>();
    Variable[] variables = simulationSymbolTable.getVariables();
    for (Variable variable : variables) {
        if (variable instanceof ParticleVariable) {
            if (variable.getDomain() == null) {
                throw new SolverException("Particle Variables are required to be defined in a subdomain using syntax Subdomain::Variable.");
            }
            particleVariableList.add((ParticleVariable) variable);
        }
    }
    // write geometry
    Geometry geometry = mathDesc.getGeometry();
    dimension = geometry.getDimension();
    try {
        // clone and resample geometry
        resampledGeometry = (Geometry) BeanUtils.cloneSerializable(geometry);
        GeometrySurfaceDescription geoSurfaceDesc = resampledGeometry.getGeometrySurfaceDescription();
        ISize newSize = simulation.getMeshSpecification().getSamplingSize();
        geoSurfaceDesc.setVolumeSampleSize(newSize);
        geoSurfaceDesc.updateAll();
        bHasNoSurface = geoSurfaceDesc.getSurfaceClasses() == null || geoSurfaceDesc.getSurfaceClasses().length == 0;
    } catch (Exception e) {
        e.printStackTrace();
        throw new SolverException(e.getMessage());
    }
    if (!bGraphicOpenGL) {
        writeMeshFile();
    }
    colors = ColorUtil.generateAutoColor(particleVariableList.size() + resampledGeometry.getGeometrySurfaceDescription().getSurfaceClasses().length, bg, new Integer(5));
}
Also used : Geometry(cbit.vcell.geometry.Geometry) ReservedVariable(cbit.vcell.math.ReservedVariable) VolumeParticleVariable(cbit.vcell.math.VolumeParticleVariable) MembraneParticleVariable(cbit.vcell.math.MembraneParticleVariable) ParticleVariable(cbit.vcell.math.ParticleVariable) Variable(cbit.vcell.math.Variable) GeometrySurfaceDescription(cbit.vcell.geometry.surface.GeometrySurfaceDescription) VolumeParticleVariable(cbit.vcell.math.VolumeParticleVariable) MembraneParticleVariable(cbit.vcell.math.MembraneParticleVariable) ParticleVariable(cbit.vcell.math.ParticleVariable) ISize(org.vcell.util.ISize) SolverException(cbit.vcell.solver.SolverException) ProgrammingException(org.vcell.util.ProgrammingException) GeometryException(cbit.vcell.geometry.GeometryException) IOException(java.io.IOException) DataAccessException(org.vcell.util.DataAccessException) PropertyVetoException(java.beans.PropertyVetoException) DivideByZeroException(cbit.vcell.parser.DivideByZeroException) ImageException(cbit.image.ImageException) ExpressionBindingException(cbit.vcell.parser.ExpressionBindingException) SolverException(cbit.vcell.solver.SolverException) ExpressionException(cbit.vcell.parser.ExpressionException) MathException(cbit.vcell.math.MathException)

Example 35 with Variable

use of cbit.vcell.math.Variable in project vcell by virtualcell.

the class SmoldynFileWriter method writeReactions.

private void writeReactions() throws ExpressionException, MathException {
    printWriter.println("# reactions");
    Enumeration<SubDomain> subdomains = mathDesc.getSubDomains();
    while (subdomains.hasMoreElements()) {
        SubDomain subdomain = subdomains.nextElement();
        for (ParticleJumpProcess pjp : subdomain.getParticleJumpProcesses()) {
            ArrayList<Variable> reactants = new ArrayList<Variable>();
            ArrayList<Variable> products = new ArrayList<Variable>();
            for (Action a : pjp.getActions().toArray(new Action[pjp.getActions().size()])) {
                if (a.getOperation().equals(Action.ACTION_CREATE)) {
                    products.add(a.getVar());
                } else if (a.getOperation().equals(Action.ACTION_DESTROY)) {
                    reactants.add(a.getVar());
                }
            }
            Expression rateDefinition = null;
            JumpProcessRateDefinition jprd = pjp.getParticleRateDefinition();
            if (jprd instanceof MacroscopicRateConstant) {
                rateDefinition = subsituteFlatten(((MacroscopicRateConstant) jprd).getExpression());
            } else if (jprd instanceof InteractionRadius) {
                rateDefinition = subsituteFlatten(((InteractionRadius) jprd).getExpression());
            } else {
                new RuntimeException("The jump process rate definition is not supported");
            }
            if (rateDefinition.isZero()) {
                continue;
            }
            if (mathDesc.isSpatialHybrid()) {
                String[] symbols = rateDefinition.getSymbols();
                if (symbols != null) {
                    if (subdomain instanceof MembraneSubDomain) {
                        rateDefinition = new Expression(FiniteVolumeFileWriter.replaceVolumeVariable(getSimulationTask(), (MembraneSubDomain) subdomain, rateDefinition));
                    }
                }
            } else {
                try {
                    rateDefinition.evaluateConstant();
                } catch (ExpressionException ex) {
                    throw new ExpressionException("reaction rate for jump process " + pjp.getName() + " is not a constant. Constants are required for all reaction rates.");
                }
            }
            // Smoldyn takes maximum 2nd order reaction.
            if (reactants.size() > 2) {
                throw new MathException("VCell spatial stochastic models support up to 2nd order reactions. \n" + "The reaction:" + pjp.getName() + " has more than 2 reactants.");
            }
            if (products.size() > 2) {
                throw new MathException("VCell spatial stochastic models support up to 2nd order reactions. \n" + "The reaction:" + pjp.getName() + " has more than 2 products.");
            }
            String rateDefinitionStr = simulation.getMathDescription().isSpatialHybrid() ? rateDefinition.infix() + ";" : rateDefinition.evaluateConstant() + "";
            if (subdomain instanceof CompartmentSubDomain) {
                // 0th order reaction, product limited to one and we'll let the reaction know where it happens
                if (reactants.size() == 0 && products.size() == 1) {
                    printWriter.print(SmoldynVCellMapper.SmoldynKeyword.reaction_cmpt + " " + subdomain.getName() + " " + pjp.getName() + " ");
                } else {
                    printWriter.print(SmoldynVCellMapper.SmoldynKeyword.reaction + " " + /* + subdomain.getName() + " "*/
                    pjp.getName() + " ");
                }
                writeReactionCommand(reactants, products, subdomain, rateDefinitionStr);
            } else if (subdomain instanceof MembraneSubDomain) {
                // 0th order reaction, product limited to one and it can be on mem or in vol
                if (reactants.size() == 0 && products.size() == 1) {
                    printWriter.print(SmoldynVCellMapper.SmoldynKeyword.reaction_surface + " " + subdomain.getName() + " " + pjp.getName() + " ");
                    writeReactionCommand(reactants, products, subdomain, rateDefinitionStr);
                } else // consuming of a species to nothing, limited to one reactant
                if (reactants.size() == 1 && products.size() == 0) {
                    if (// consuming a mem species in mem reaction
                    getMembraneVariableCount(reactants) == 1) {
                        printWriter.print(SmoldynVCellMapper.SmoldynKeyword.reaction_surface + " " + subdomain.getName() + " " + pjp.getName() + " ");
                        writeReactionCommand(reactants, products, subdomain, rateDefinitionStr);
                    } else // it equals to adsorption, species A from vol adsorbed to mem as again species A, and then we kill the speceis A on mem.
                    if (getVolumeVariableCount(reactants) == 1) {
                        writeRateTransitionCommand(reactants, products, subdomain, rateDefinitionStr);
                        String speciesName = reactants.get(0).getName();
                        String killMolCmd = "cmd " + SmoldynVCellMapper.SmoldynKeyword.E + " " + SmoldynVCellMapper.SmoldynKeyword.killmol + " " + speciesName + "(" + SmoldynKeyword.up + ")";
                        killMolCommands.add(killMolCmd);
                    }
                } else // Use rate command for any membrane reactions with 1 reactant and 1 product
                if ((reactants.size() == 1) && (products.size() == 1)) {
                    // Membrane reaction (1 react to 1 product).
                    if (getMembraneVariableCount(products) == 1 && getMembraneVariableCount(reactants) == 1) {
                        printWriter.print(SmoldynVCellMapper.SmoldynKeyword.reaction_surface + " " + subdomain.getName() + " " + pjp.getName() + " ");
                        writeReactionCommand(reactants, products, subdomain, rateDefinitionStr);
                    } else // Other single molecular reactions
                    {
                        writeRateTransitionCommand(reactants, products, subdomain, rateDefinitionStr);
                    }
                } else // membrane reactions which are not one to one, or 0th order, or consuming species
                {
                    if (// membrane reaction has one membrane bound reactant
                    (getMembraneVariableCount(reactants) == 1)) {
                        printWriter.print(SmoldynVCellMapper.SmoldynKeyword.reaction_surface + " " + subdomain.getName() + " " + pjp.getName() + " ");
                        writeReactionCommand(reactants, products, subdomain, rateDefinitionStr);
                    } else if (// bimolecular membrane reaction
                    getMembraneVariableCount(reactants) == 2) {
                        if (jprd instanceof InteractionRadius) {
                            printWriter.print(SmoldynVCellMapper.SmoldynKeyword.reaction_surface + " " + subdomain.getName() + " " + pjp.getName() + " ");
                            writeReactionByInteractionRadius(reactants, products, subdomain, rateDefinitionStr, pjp.getName());
                        } else {
                            // throw new MathException("Error with reaction: " + pjp.getName() + ".\nVCell Spatial stochastic modeling requires macroscopic or microscopic kinetics for bimolecular membrane reactions.");
                            printWriter.print(SmoldynVCellMapper.SmoldynKeyword.reaction_surface + " " + subdomain.getName() + " " + pjp.getName() + " ");
                            writeReactionCommand(reactants, products, subdomain, rateDefinitionStr);
                        }
                    } else if (getMembraneVariableCount(reactants) == 0) {
                        throw new MathException("Error with reaction: " + pjp.getName() + ".\nIn VCell spatial stochastic modeling, the membrane reaction requires at least one membrane bound reactant.");
                    }
                }
            }
        }
    }
    printWriter.println();
}
Also used : Action(cbit.vcell.math.Action) JumpProcessRateDefinition(cbit.vcell.math.JumpProcessRateDefinition) MembraneSubDomain(cbit.vcell.math.MembraneSubDomain) ReservedVariable(cbit.vcell.math.ReservedVariable) VolumeParticleVariable(cbit.vcell.math.VolumeParticleVariable) MembraneParticleVariable(cbit.vcell.math.MembraneParticleVariable) ParticleVariable(cbit.vcell.math.ParticleVariable) Variable(cbit.vcell.math.Variable) InteractionRadius(cbit.vcell.math.InteractionRadius) ParticleJumpProcess(cbit.vcell.math.ParticleJumpProcess) ArrayList(java.util.ArrayList) ExpressionException(cbit.vcell.parser.ExpressionException) CompartmentSubDomain(cbit.vcell.math.CompartmentSubDomain) SubDomain(cbit.vcell.math.SubDomain) MembraneSubDomain(cbit.vcell.math.MembraneSubDomain) Expression(cbit.vcell.parser.Expression) MathException(cbit.vcell.math.MathException) CompartmentSubDomain(cbit.vcell.math.CompartmentSubDomain) MacroscopicRateConstant(cbit.vcell.math.MacroscopicRateConstant)

Aggregations

Variable (cbit.vcell.math.Variable)110 Expression (cbit.vcell.parser.Expression)71 VolVariable (cbit.vcell.math.VolVariable)64 MemVariable (cbit.vcell.math.MemVariable)48 ReservedVariable (cbit.vcell.math.ReservedVariable)43 MathException (cbit.vcell.math.MathException)38 MembraneRegionVariable (cbit.vcell.math.MembraneRegionVariable)38 VolumeRegionVariable (cbit.vcell.math.VolumeRegionVariable)36 ExpressionException (cbit.vcell.parser.ExpressionException)36 FilamentVariable (cbit.vcell.math.FilamentVariable)35 InsideVariable (cbit.vcell.math.InsideVariable)34 OutsideVariable (cbit.vcell.math.OutsideVariable)34 Function (cbit.vcell.math.Function)33 MathDescription (cbit.vcell.math.MathDescription)33 Constant (cbit.vcell.math.Constant)32 FilamentRegionVariable (cbit.vcell.math.FilamentRegionVariable)29 VolumeParticleVariable (cbit.vcell.math.VolumeParticleVariable)25 MembraneParticleVariable (cbit.vcell.math.MembraneParticleVariable)24 ParticleVariable (cbit.vcell.math.ParticleVariable)24 Vector (java.util.Vector)24