Search in sources :

Example 1 with Parameter

use of cbit.vcell.model.Parameter in project vcell by virtualcell.

the class ParameterMappingPanel method jMenuItemPaste_ActionPerformed.

/**
 * Comment
 */
private void jMenuItemPaste_ActionPerformed(java.awt.event.ActionEvent actionEvent) {
    java.util.Vector<String> pasteDescriptionsV = new java.util.Vector<String>();
    java.util.Vector<Expression> newExpressionsV = new java.util.Vector<Expression>();
    java.util.Vector<ParameterMappingSpec> changedParametersV = new java.util.Vector<ParameterMappingSpec>();
    try {
        if (actionEvent.getSource() == getJMenuItemPaste() || actionEvent.getSource() == getJMenuItemPasteAll()) {
            Object pasteThis = VCellTransferable.getFromClipboard(VCellTransferable.OBJECT_FLAVOR);
            SimulationContext sc = getParameterEstimationTask().getModelOptimizationSpec().getSimulationContext();
            MathSymbolMapping msm = null;
            Exception mathMappingException = null;
            try {
                MathMapping mm = sc.createNewMathMapping();
                msm = mm.getMathSymbolMapping();
            } catch (Exception e) {
                mathMappingException = e;
                e.printStackTrace(System.out);
            }
            // if(msm == null){
            // try{
            // getParameterEstimationTask().refreshMappings();
            // msm = getParameterEstimationTask().getMathSymbolMapping();
            // }catch(Exception e){
            // e.printStackTrace();
            // }
            // }
            int[] rows = null;
            if (actionEvent.getSource() == getJMenuItemPasteAll()) {
                rows = new int[getScrollPaneTable().getRowCount()];
                for (int i = 0; i < rows.length; i += 1) {
                    rows[i] = i;
                }
            } else {
                rows = getScrollPaneTable().getSelectedRows();
            }
            // 
            // Check paste
            // 
            StringBuffer errors = null;
            for (int i = 0; i < rows.length; i += 1) {
                ParameterMappingSpec pms = parameterMappingTableModel.getValueAt(rows[i]);
                try {
                    if (pasteThis instanceof VCellTransferable.ResolvedValuesSelection) {
                        VCellTransferable.ResolvedValuesSelection rvs = (VCellTransferable.ResolvedValuesSelection) pasteThis;
                        for (int j = 0; j < rvs.getPrimarySymbolTableEntries().length; j += 1) {
                            ParameterMappingSpec pasteDestination = null;
                            Parameter clipboardBiologicalParameter = null;
                            if (rvs.getPrimarySymbolTableEntries()[j] instanceof Parameter) {
                                clipboardBiologicalParameter = (Parameter) rvs.getPrimarySymbolTableEntries()[j];
                            } else if (rvs.getAlternateSymbolTableEntries() != null && rvs.getAlternateSymbolTableEntries()[j] instanceof Parameter) {
                                clipboardBiologicalParameter = (Parameter) rvs.getAlternateSymbolTableEntries()[j];
                            }
                            if (clipboardBiologicalParameter == null) {
                                Variable pastedMathVariable = null;
                                if (rvs.getPrimarySymbolTableEntries()[j] instanceof Variable) {
                                    pastedMathVariable = (Variable) rvs.getPrimarySymbolTableEntries()[j];
                                } else if (rvs.getAlternateSymbolTableEntries() != null && rvs.getAlternateSymbolTableEntries()[j] instanceof Variable) {
                                    pastedMathVariable = (Variable) rvs.getAlternateSymbolTableEntries()[j];
                                }
                                if (pastedMathVariable != null) {
                                    if (msm == null) {
                                        throw mathMappingException;
                                    }
                                    Variable localMathVariable = msm.findVariableByName(pastedMathVariable.getName());
                                    if (localMathVariable != null) {
                                        SymbolTableEntry[] localBiologicalSymbolArr = msm.getBiologicalSymbol(localMathVariable);
                                        for (int k = 0; k < localBiologicalSymbolArr.length; k += 1) {
                                            if (localBiologicalSymbolArr[k] == pms.getModelParameter()) {
                                                pasteDestination = pms;
                                                break;
                                            }
                                        }
                                    }
                                }
                            } else {
                                if (pms.getModelParameter().getName().equals(clipboardBiologicalParameter.getName()) && pms.getModelParameter().getClass().equals(clipboardBiologicalParameter.getClass()) && pms.getModelParameter().getNameScope().getName().equals(clipboardBiologicalParameter.getNameScope().getName())) {
                                    pasteDestination = pms;
                                }
                            }
                            if (pasteDestination != null) {
                                changedParametersV.add(pasteDestination);
                                newExpressionsV.add(rvs.getExpressionValues()[j]);
                                pasteDescriptionsV.add(VCellCopyPasteHelper.formatPasteList(pms.getModelParameter().getNameScope().getName(), pms.getModelParameter().getName(), pasteDestination.getCurrent() + "", rvs.getExpressionValues()[j].infix()));
                            }
                        }
                    }
                } catch (Throwable e) {
                    if (errors == null) {
                        errors = new StringBuffer();
                    }
                    errors.append(pms.getModelParameter().getName() + " (" + e.getClass().getName() + ") " + e.getMessage() + "\n");
                }
            }
            if (errors != null) {
                throw new Exception(errors.toString());
            }
        }
    } catch (Throwable e) {
        PopupGenerator.showErrorDialog(this, "Paste failed during pre-check (no changes made).\n" + e.getMessage(), e);
        return;
    }
    // Do paste
    try {
        if (pasteDescriptionsV.size() > 0) {
            String[] pasteDescriptionArr = new String[pasteDescriptionsV.size()];
            pasteDescriptionsV.copyInto(pasteDescriptionArr);
            ParameterMappingSpec[] changedParametersArr = new ParameterMappingSpec[changedParametersV.size()];
            changedParametersV.copyInto(changedParametersArr);
            Expression[] newExpressionsArr = new Expression[newExpressionsV.size()];
            newExpressionsV.copyInto(newExpressionsArr);
            VCellCopyPasteHelper.chooseApplyPaste(this, pasteDescriptionArr, changedParametersArr, newExpressionsArr);
        } else {
            PopupGenerator.showInfoDialog(this, "No paste items match the destination (no changes made).");
        }
    } catch (Throwable e) {
        PopupGenerator.showErrorDialog(this, "Paste Error\n" + e.getMessage(), e);
    }
}
Also used : Variable(cbit.vcell.math.Variable) VCellTransferable(cbit.vcell.desktop.VCellTransferable) SymbolTableEntry(cbit.vcell.parser.SymbolTableEntry) ParameterMappingSpec(cbit.vcell.modelopt.ParameterMappingSpec) SimulationContext(cbit.vcell.mapping.SimulationContext) MathSymbolMapping(cbit.vcell.mapping.MathSymbolMapping) Expression(cbit.vcell.parser.Expression) MathMapping(cbit.vcell.mapping.MathMapping) ModelParameter(cbit.vcell.model.Model.ModelParameter) Parameter(cbit.vcell.model.Parameter) UnresolvedParameter(cbit.vcell.model.Kinetics.UnresolvedParameter)

Example 2 with Parameter

use of cbit.vcell.model.Parameter in project vcell by virtualcell.

the class RbmNetworkGenerator method writeParameters.

public static void writeParameters(PrintWriter writer, RbmModelContainer rbmModelContainer, boolean ignoreFunctions) {
    writer.println(BEGIN_PARAMETERS);
    List<Parameter> paramList = rbmModelContainer.getParameterList();
    for (Parameter param : paramList) {
        writer.println(RbmUtils.toBnglString(param, false));
    }
    if (ignoreFunctions) {
        // we cheat and transform all functions into constant parameters
        List<Parameter> functionList = rbmModelContainer.getFunctionList();
        for (Parameter function : functionList) {
            writer.println(RbmUtils.toBnglStringIgnoreExpression(function));
        }
    }
    writer.println(END_PARAMETERS);
    writer.println();
}
Also used : LocalParameter(cbit.vcell.mapping.ParameterContext.LocalParameter) ModelParameter(cbit.vcell.model.Model.ModelParameter) Parameter(cbit.vcell.model.Parameter)

Example 3 with Parameter

use of cbit.vcell.model.Parameter in project vcell by virtualcell.

the class RbmNetworkGenerator method writeFunctions.

public static void writeFunctions(PrintWriter writer, RbmModelContainer rbmModelContainer, boolean ignoreFunctions) {
    if (!ignoreFunctions) {
        writer.println(BEGIN_FUNCTIONS);
        List<Parameter> functionList = rbmModelContainer.getFunctionList();
        if (RbmUtils.haveUnsupportedFunctions(functionList)) {
            writer.println("# Commented out functions correspond to observables that cannot be exported into BNGL yet");
        }
        for (Parameter function : functionList) {
            writer.println(RbmUtils.toBnglString(function, true));
        }
        writer.println(END_FUNCTIONS);
        writer.println();
    }
}
Also used : LocalParameter(cbit.vcell.mapping.ParameterContext.LocalParameter) ModelParameter(cbit.vcell.model.Model.ModelParameter) Parameter(cbit.vcell.model.Parameter)

Example 4 with Parameter

use of cbit.vcell.model.Parameter in project vcell by virtualcell.

the class StochMathMapping method refreshMathDescription.

/**
 * set up a math description based on current simulationContext.
 */
@Override
protected void refreshMathDescription() throws MappingException, MatrixException, MathException, ExpressionException, ModelException {
    // use local variable instead of using getter all the time.
    SimulationContext simContext = getSimulationContext();
    GeometryClass geometryClass = simContext.getGeometry().getGeometrySpec().getSubVolumes()[0];
    Domain domain = new Domain(geometryClass);
    // local structure mapping list
    StructureMapping[] structureMappings = simContext.getGeometryContext().getStructureMappings();
    // We have to check if all the reactions are able to tranform to stochastic jump processes before generating the math.
    String stochChkMsg = simContext.getModel().isValidForStochApp();
    if (!(stochChkMsg.equals(""))) {
        throw new ModelException("Problem updating math description: " + simContext.getName() + "\n" + stochChkMsg);
    }
    simContext.checkValidity();
    // 
    if (simContext.getGeometry().getDimension() > 0) {
        throw new MappingException("nonspatial stochastic math mapping requires 0-dimensional geometry");
    }
    // 
    for (int i = 0; i < structureMappings.length; i++) {
        if (structureMappings[i] instanceof MembraneMapping) {
            if (((MembraneMapping) structureMappings[i]).getCalculateVoltage()) {
                throw new MappingException("electric potential not yet supported for particle models");
            }
        }
    }
    // 
    // fail if any events
    // 
    BioEvent[] bioEvents = simContext.getBioEvents();
    if (bioEvents != null && bioEvents.length > 0) {
        throw new MappingException("events not yet supported for particle-based models");
    }
    // 
    // verify that all structures are mapped to subvolumes and all subvolumes are mapped to a structure
    // 
    Structure[] structures = simContext.getGeometryContext().getModel().getStructures();
    for (int i = 0; i < structures.length; i++) {
        StructureMapping sm = simContext.getGeometryContext().getStructureMapping(structures[i]);
        if (sm == null || (sm instanceof FeatureMapping && ((FeatureMapping) sm).getGeometryClass() == null)) {
            throw new MappingException("model structure '" + structures[i].getName() + "' not mapped to a geometry subVolume");
        }
        if (sm != null && (sm instanceof MembraneMapping) && ((MembraneMapping) sm).getVolumeFractionParameter() != null) {
            Expression volFractExp = ((MembraneMapping) sm).getVolumeFractionParameter().getExpression();
            try {
                if (volFractExp != null) {
                    double volFract = volFractExp.evaluateConstant();
                    if (volFract >= 1.0) {
                        throw new MappingException("model structure '" + (getSimulationContext().getModel().getStructureTopology().getInsideFeature(((MembraneMapping) sm).getMembrane()).getName() + "' has volume fraction >= 1.0"));
                    }
                }
            } catch (ExpressionException e) {
                e.printStackTrace(System.out);
            }
        }
    }
    SubVolume[] subVolumes = simContext.getGeometryContext().getGeometry().getGeometrySpec().getSubVolumes();
    for (int i = 0; i < subVolumes.length; i++) {
        Structure[] mappedStructures = simContext.getGeometryContext().getStructuresFromGeometryClass(subVolumes[i]);
        if (mappedStructures == null || mappedStructures.length == 0) {
            throw new MappingException("geometry subVolume '" + subVolumes[i].getName() + "' not mapped from a model structure");
        }
    }
    // 
    // gather only those reactionSteps that are not "excluded"
    // 
    ReactionSpec[] reactionSpecs = simContext.getReactionContext().getReactionSpecs();
    Vector<ReactionStep> rsList = new Vector<ReactionStep>();
    for (int i = 0; i < reactionSpecs.length; i++) {
        if (!reactionSpecs[i].isExcluded()) {
            rsList.add(reactionSpecs[i].getReactionStep());
        }
    }
    // 
    for (ReactionStep reactionStep : rsList) {
        Kinetics.UnresolvedParameter[] unresolvedParameters = reactionStep.getKinetics().getUnresolvedParameters();
        if (unresolvedParameters != null && unresolvedParameters.length > 0) {
            StringBuffer buffer = new StringBuffer();
            for (int j = 0; j < unresolvedParameters.length; j++) {
                if (j > 0) {
                    buffer.append(", ");
                }
                buffer.append(unresolvedParameters[j].getName());
            }
            throw new MappingException("In Application '" + simContext.getName() + "', " + reactionStep.getDisplayType() + " '" + reactionStep.getName() + "' contains unresolved identifier(s): " + buffer);
        }
    }
    // 
    // create new MathDescription (based on simContext's previous MathDescription if possible)
    // 
    MathDescription oldMathDesc = simContext.getMathDescription();
    mathDesc = null;
    if (oldMathDesc != null) {
        if (oldMathDesc.getVersion() != null) {
            mathDesc = new MathDescription(oldMathDesc.getVersion());
        } else {
            mathDesc = new MathDescription(oldMathDesc.getName());
        }
    } else {
        mathDesc = new MathDescription(simContext.getName() + "_generated");
    }
    // 
    // temporarily place all variables in a hashtable (before binding) and discarding duplicates
    // 
    VariableHash varHash = new VariableHash();
    // 
    // conversion factors
    // 
    Model model = simContext.getModel();
    varHash.addVariable(new Constant(getMathSymbol(model.getKMOLE(), null), getIdentifierSubstitutions(model.getKMOLE().getExpression(), model.getKMOLE().getUnitDefinition(), null)));
    varHash.addVariable(new Constant(getMathSymbol(model.getN_PMOLE(), null), getIdentifierSubstitutions(model.getN_PMOLE().getExpression(), model.getN_PMOLE().getUnitDefinition(), null)));
    varHash.addVariable(new Constant(getMathSymbol(model.getPI_CONSTANT(), null), getIdentifierSubstitutions(model.getPI_CONSTANT().getExpression(), model.getPI_CONSTANT().getUnitDefinition(), null)));
    varHash.addVariable(new Constant(getMathSymbol(model.getFARADAY_CONSTANT(), null), getIdentifierSubstitutions(model.getFARADAY_CONSTANT().getExpression(), model.getFARADAY_CONSTANT().getUnitDefinition(), null)));
    varHash.addVariable(new Constant(getMathSymbol(model.getFARADAY_CONSTANT_NMOLE(), null), getIdentifierSubstitutions(model.getFARADAY_CONSTANT_NMOLE().getExpression(), model.getFARADAY_CONSTANT_NMOLE().getUnitDefinition(), null)));
    varHash.addVariable(new Constant(getMathSymbol(model.getGAS_CONSTANT(), null), getIdentifierSubstitutions(model.getGAS_CONSTANT().getExpression(), model.getGAS_CONSTANT().getUnitDefinition(), null)));
    varHash.addVariable(new Constant(getMathSymbol(model.getTEMPERATURE(), null), getIdentifierSubstitutions(new Expression(simContext.getTemperatureKelvin()), model.getTEMPERATURE().getUnitDefinition(), null)));
    Enumeration<SpeciesContextMapping> enum1 = getSpeciesContextMappings();
    while (enum1.hasMoreElements()) {
        SpeciesContextMapping scm = enum1.nextElement();
        if (scm.getVariable() instanceof StochVolVariable) {
            varHash.addVariable(scm.getVariable());
        }
    }
    // deals with model parameters
    ModelParameter[] modelParameters = simContext.getModel().getModelParameters();
    for (int j = 0; j < modelParameters.length; j++) {
        Expression expr = getSubstitutedExpr(modelParameters[j].getExpression(), true, false);
        expr = getIdentifierSubstitutions(expr, modelParameters[j].getUnitDefinition(), geometryClass);
        varHash.addVariable(newFunctionOrConstant(getMathSymbol(modelParameters[j], geometryClass), expr, geometryClass));
    }
    // added July 2009, ElectricalStimulusParameter electric mapping tab
    ElectricalStimulus[] elecStimulus = simContext.getElectricalStimuli();
    if (elecStimulus.length > 0) {
        throw new MappingException("Modles with electrophysiology are not supported for stochastic applications.");
    }
    for (int j = 0; j < structureMappings.length; j++) {
        if (structureMappings[j] instanceof MembraneMapping) {
            MembraneMapping memMapping = (MembraneMapping) structureMappings[j];
            Parameter initialVoltageParm = memMapping.getInitialVoltageParameter();
            try {
                Expression exp = initialVoltageParm.getExpression();
                exp.evaluateConstant();
                varHash.addVariable(newFunctionOrConstant(getMathSymbol(memMapping.getMembrane().getMembraneVoltage(), memMapping.getGeometryClass()), getIdentifierSubstitutions(memMapping.getInitialVoltageParameter().getExpression(), memMapping.getInitialVoltageParameter().getUnitDefinition(), memMapping.getGeometryClass()), memMapping.getGeometryClass()));
            } catch (ExpressionException e) {
                e.printStackTrace(System.out);
                throw new MappingException("Membrane initial voltage: " + initialVoltageParm.getName() + " cannot be evaluated as constant.");
            }
        }
    }
    // 
    for (ReactionStep rs : rsList) {
        if (rs.getKinetics() instanceof LumpedKinetics) {
            throw new RuntimeException("Lumped Kinetics not yet supported for Stochastic Math Generation");
        }
        Kinetics.KineticsParameter[] parameters = rs.getKinetics().getKineticsParameters();
        for (KineticsParameter parameter : parameters) {
            // 
            if ((parameter.getRole() == Kinetics.ROLE_CurrentDensity) && (parameter.getExpression() == null || parameter.getExpression().isZero())) {
                continue;
            }
            // 
            // don't add rate, we'll do it later when creating the jump processes
            // 
            // if (parameter.getRole() == Kinetics.ROLE_ReactionRate) {
            // continue;
            // }
            // 
            // don't add mass action reverse parameter if irreversible
            // 
            // if (!rs.isReversible() && parameters[i].getRole() == Kinetics.ROLE_KReverse){
            // continue;
            // }
            Expression expr = getSubstitutedExpr(parameter.getExpression(), true, false);
            varHash.addVariable(newFunctionOrConstant(getMathSymbol(parameter, geometryClass), getIdentifierSubstitutions(expr, parameter.getUnitDefinition(), geometryClass), geometryClass));
        }
    }
    // the parameter "Size" is already put into mathsymbolmapping in refreshSpeciesContextMapping()
    for (int i = 0; i < structureMappings.length; i++) {
        StructureMapping sm = structureMappings[i];
        StructureMapping.StructureMappingParameter parm = sm.getParameterFromRole(StructureMapping.ROLE_Size);
        if (parm.getExpression() != null) {
            try {
                double value = parm.getExpression().evaluateConstant();
                varHash.addVariable(new Constant(getMathSymbol(parm, sm.getGeometryClass()), new Expression(value)));
            } catch (ExpressionException e) {
                // varHash.addVariable(new Function(getMathSymbol0(parm,sm),getIdentifierSubstitutions(parm.getExpression(),parm.getUnitDefinition(),sm)));
                e.printStackTrace(System.out);
                throw new MappingException("Size of structure:" + sm.getNameScope().getName() + " cannot be evaluated as constant.");
            }
        }
    }
    SpeciesContextSpec[] speciesContextSpecs = getSimulationContext().getReactionContext().getSpeciesContextSpecs();
    addInitialConditions(domain, speciesContextSpecs, varHash);
    // 
    // constant species (either function or constant)
    // 
    enum1 = getSpeciesContextMappings();
    while (enum1.hasMoreElements()) {
        SpeciesContextMapping scm = (SpeciesContextMapping) enum1.nextElement();
        if (scm.getVariable() instanceof Constant) {
            varHash.addVariable(scm.getVariable());
        }
    }
    // 
    if (simContext.getGeometryContext().getGeometry() != null) {
        try {
            mathDesc.setGeometry(simContext.getGeometryContext().getGeometry());
        } catch (java.beans.PropertyVetoException e) {
            e.printStackTrace(System.out);
            throw new MappingException("failure setting geometry " + e.getMessage());
        }
    } else {
        throw new MappingException("Geometry must be defined in Application " + simContext.getName());
    }
    // 
    // create subDomains
    // 
    SubVolume subVolume = simContext.getGeometry().getGeometrySpec().getSubVolumes()[0];
    SubDomain subDomain = new CompartmentSubDomain(subVolume.getName(), 0);
    mathDesc.addSubDomain(subDomain);
    // 
    // functions: species which is not a variable, but has dependency expression
    // 
    enum1 = getSpeciesContextMappings();
    while (enum1.hasMoreElements()) {
        SpeciesContextMapping scm = (SpeciesContextMapping) enum1.nextElement();
        if (scm.getVariable() == null && scm.getDependencyExpression() != null) {
            StructureMapping sm = simContext.getGeometryContext().getStructureMapping(scm.getSpeciesContext().getStructure());
            Expression exp = scm.getDependencyExpression();
            exp.bindExpression(this);
            SpeciesCountParameter spCountParam = getSpeciesCountParameter(scm.getSpeciesContext());
            varHash.addVariable(new Function(getMathSymbol(spCountParam, sm.getGeometryClass()), getIdentifierSubstitutions(exp, spCountParam.getUnitDefinition(), sm.getGeometryClass()), domain));
        }
    }
    addJumpProcesses(varHash, geometryClass, subDomain);
    // 
    for (int i = 0; i < fieldMathMappingParameters.length; i++) {
        if (fieldMathMappingParameters[i] instanceof UnitFactorParameter) {
            varHash.addVariable(newFunctionOrConstant(getMathSymbol(fieldMathMappingParameters[i], geometryClass), getIdentifierSubstitutions(fieldMathMappingParameters[i].getExpression(), fieldMathMappingParameters[i].getUnitDefinition(), geometryClass), fieldMathMappingParameters[i].getGeometryClass()));
        }
    }
    // 
    // set Variables to MathDescription all at once with the order resolved by "VariableHash"
    // 
    mathDesc.setAllVariables(varHash.getAlphabeticallyOrderedVariables());
    // 
    // set up variable initial conditions in subDomain
    // 
    SpeciesContextSpec[] scSpecs = simContext.getReactionContext().getSpeciesContextSpecs();
    for (int i = 0; i < speciesContextSpecs.length; i++) {
        // get stochastic variable by name
        SpeciesCountParameter spCountParam = getSpeciesCountParameter(speciesContextSpecs[i].getSpeciesContext());
        StructureMapping sm = simContext.getGeometryContext().getStructureMapping(speciesContextSpecs[i].getSpeciesContext().getStructure());
        String varName = getMathSymbol(spCountParam, sm.getGeometryClass());
        StochVolVariable var = (StochVolVariable) mathDesc.getVariable(varName);
        // stochastic use initial number of particles
        SpeciesContextSpec.SpeciesContextSpecParameter initParm = scSpecs[i].getInitialCountParameter();
        // stochastic variables initial expression.
        if (initParm != null) {
            VarIniCondition varIni = null;
            if (!scSpecs[i].isConstant() && getSimulationContext().isRandomizeInitCondition()) {
                varIni = new VarIniPoissonExpectedCount(var, new Expression(getMathSymbol(initParm, sm.getGeometryClass())));
            } else {
                varIni = new VarIniCount(var, new Expression(getMathSymbol(initParm, sm.getGeometryClass())));
            }
            subDomain.addVarIniCondition(varIni);
        }
    }
    // 
    for (int i = 0; i < fieldMathMappingParameters.length; i++) {
        if (fieldMathMappingParameters[i] instanceof UnitFactorParameter) {
            Variable variable = newFunctionOrConstant(getMathSymbol(fieldMathMappingParameters[i], geometryClass), getIdentifierSubstitutions(fieldMathMappingParameters[i].getExpression(), fieldMathMappingParameters[i].getUnitDefinition(), geometryClass), fieldMathMappingParameters[i].getGeometryClass());
            if (mathDesc.getVariable(variable.getName()) == null) {
                mathDesc.addVariable(variable);
            }
        }
        if (fieldMathMappingParameters[i] instanceof ObservableCountParameter) {
            Variable variable = newFunctionOrConstant(getMathSymbol(fieldMathMappingParameters[i], geometryClass), getIdentifierSubstitutions(fieldMathMappingParameters[i].getExpression(), fieldMathMappingParameters[i].getUnitDefinition(), geometryClass), fieldMathMappingParameters[i].getGeometryClass());
            if (mathDesc.getVariable(variable.getName()) == null) {
                mathDesc.addVariable(variable);
            }
        }
    }
    if (!mathDesc.isValid()) {
        System.out.println(mathDesc.getVCML_database());
        throw new MappingException("generated an invalid mathDescription: " + mathDesc.getWarning());
    }
}
Also used : VarIniCondition(cbit.vcell.math.VarIniCondition) LumpedKinetics(cbit.vcell.model.LumpedKinetics) MathDescription(cbit.vcell.math.MathDescription) ExpressionException(cbit.vcell.parser.ExpressionException) PropertyVetoException(java.beans.PropertyVetoException) CompartmentSubDomain(cbit.vcell.math.CompartmentSubDomain) SubDomain(cbit.vcell.math.SubDomain) KineticsParameter(cbit.vcell.model.Kinetics.KineticsParameter) SubVolume(cbit.vcell.geometry.SubVolume) Vector(java.util.Vector) ModelException(cbit.vcell.model.ModelException) ModelParameter(cbit.vcell.model.Model.ModelParameter) CompartmentSubDomain(cbit.vcell.math.CompartmentSubDomain) ReactionStep(cbit.vcell.model.ReactionStep) CompartmentSubDomain(cbit.vcell.math.CompartmentSubDomain) SubDomain(cbit.vcell.math.SubDomain) Domain(cbit.vcell.math.Variable.Domain) GeometryClass(cbit.vcell.geometry.GeometryClass) StochVolVariable(cbit.vcell.math.StochVolVariable) Variable(cbit.vcell.math.Variable) VariableHash(cbit.vcell.math.VariableHash) Constant(cbit.vcell.math.Constant) VarIniPoissonExpectedCount(cbit.vcell.math.VarIniPoissonExpectedCount) Function(cbit.vcell.math.Function) Structure(cbit.vcell.model.Structure) StochVolVariable(cbit.vcell.math.StochVolVariable) VarIniCount(cbit.vcell.math.VarIniCount) Expression(cbit.vcell.parser.Expression) Model(cbit.vcell.model.Model) Parameter(cbit.vcell.model.Parameter) KineticsParameter(cbit.vcell.model.Kinetics.KineticsParameter) ModelParameter(cbit.vcell.model.Model.ModelParameter)

Example 5 with Parameter

use of cbit.vcell.model.Parameter in project vcell by virtualcell.

the class StochMathMapping method addJumpProcesses.

private void addJumpProcesses(VariableHash varHash, GeometryClass geometryClass, SubDomain subDomain) throws ExpressionException, ModelException, MappingException, MathException {
    // set up jump processes
    // get all the reactions from simulation context
    // ReactionSpec[] reactionSpecs = simContext.getReactionContext().getReactionSpecs();---need to take a look here!
    ModelUnitSystem modelUnitSystem = getSimulationContext().getModel().getUnitSystem();
    ReactionSpec[] reactionSpecs = getSimulationContext().getReactionContext().getReactionSpecs();
    for (ReactionSpec reactionSpec : reactionSpecs) {
        if (reactionSpec.isExcluded()) {
            continue;
        }
        // get the reaction
        ReactionStep reactionStep = reactionSpec.getReactionStep();
        Kinetics kinetics = reactionStep.getKinetics();
        // probability parameter from modelUnitSystem
        VCUnitDefinition probabilityParamUnit = modelUnitSystem.getStochasticSubstanceUnit().divideBy(modelUnitSystem.getTimeUnit());
        // Different ways to deal with simple reactions and flux reactions
        if (// simple reactions
        reactionStep instanceof SimpleReaction) {
            // check the reaction rate law to see if we need to decompose a reaction(reversible) into two jump processes.
            // rate constants are important in calculating the probability rate.
            // for Mass Action, we use KForward and KReverse,
            // for General Kinetics we parse reaction rate J to see if it is in Mass Action form.
            Expression forwardRate = null;
            Expression reverseRate = null;
            if (kinetics.getKineticsDescription().equals(KineticsDescription.MassAction) || kinetics.getKineticsDescription().equals(KineticsDescription.General)) {
                Expression rateExp = new Expression(kinetics.getKineticsParameterFromRole(Kinetics.ROLE_ReactionRate), reactionStep.getNameScope());
                Parameter forwardRateParameter = null;
                Parameter reverseRateParameter = null;
                if (kinetics.getKineticsDescription().equals(KineticsDescription.MassAction)) {
                    forwardRateParameter = kinetics.getKineticsParameterFromRole(Kinetics.ROLE_KForward);
                    reverseRateParameter = kinetics.getKineticsParameterFromRole(Kinetics.ROLE_KReverse);
                }
                MassActionSolver.MassActionFunction maFunc = MassActionSolver.solveMassAction(forwardRateParameter, reverseRateParameter, rateExp, reactionStep);
                if (maFunc.getForwardRate() == null && maFunc.getReverseRate() == null) {
                    throw new MappingException("Cannot generate stochastic math mapping for the reaction:" + reactionStep.getName() + "\nLooking for the rate function according to the form of k1*Reactant1^Stoir1*Reactant2^Stoir2...-k2*Product1^Stoip1*Product2^Stoip2.");
                } else {
                    if (maFunc.getForwardRate() != null) {
                        forwardRate = maFunc.getForwardRate();
                    }
                    if (maFunc.getReverseRate() != null) {
                        reverseRate = maFunc.getReverseRate();
                    }
                }
            } else // if it's macro/microscopic kinetics, we'll have them set up as reactions with only forward rate.
            if (kinetics.getKineticsDescription().equals(KineticsDescription.Macroscopic_irreversible) || kinetics.getKineticsDescription().equals(KineticsDescription.Microscopic_irreversible)) {
                Expression Kon = getIdentifierSubstitutions(new Expression(reactionStep.getKinetics().getKineticsParameterFromRole(Kinetics.ROLE_KOn), getNameScope()), reactionStep.getKinetics().getKineticsParameterFromRole(Kinetics.ROLE_Binding_Radius).getUnitDefinition(), geometryClass);
                if (Kon != null) {
                    Expression KonCopy = new Expression(Kon);
                    try {
                        MassActionSolver.substituteParameters(KonCopy, true).evaluateConstant();
                        forwardRate = new Expression(Kon);
                    } catch (ExpressionException e) {
                        throw new MathException(VCellErrorMessages.getMassActionSolverMessage(reactionStep.getName(), "Problem with Kon parameter in " + reactionStep.getName() + ":  '" + KonCopy.infix() + "', " + e.getMessage()));
                    }
                } else {
                    throw new MathException(VCellErrorMessages.getMassActionSolverMessage(reactionStep.getName(), "Kon parameter of " + reactionStep.getName() + " is null."));
                }
            }
            boolean isForwardRatePresent = false;
            boolean isReverseRatePresent = false;
            if (forwardRate != null) {
                isForwardRatePresent = true;
            }
            if (reverseRate != null) {
                isReverseRatePresent = true;
            }
            // we process it as forward reaction
            if ((isForwardRatePresent)) /*|| ((forwardRate == null) && (reverseRate == null))*/
            {
                // get jump process name
                String jpName = TokenMangler.mangleToSName(reactionStep.getName());
                // get probability
                Expression exp = null;
                // reactions are of mass action form
                exp = getProbabilityRate(reactionStep, forwardRate, true);
                ProbabilityParameter probParm = null;
                try {
                    probParm = addProbabilityParameter(PARAMETER_PROBABILITYRATE_PREFIX + jpName, exp, PARAMETER_ROLE_P, probabilityParamUnit, reactionStep);
                } catch (PropertyVetoException pve) {
                    pve.printStackTrace();
                    throw new MappingException(pve.getMessage());
                }
                // add probability to function or constant
                varHash.addVariable(newFunctionOrConstant(getMathSymbol(probParm, geometryClass), getIdentifierSubstitutions(exp, probabilityParamUnit, geometryClass), geometryClass));
                JumpProcess jp = new JumpProcess(jpName, new Expression(getMathSymbol(probParm, geometryClass)));
                // actions
                ReactionParticipant[] reacPart = reactionStep.getReactionParticipants();
                for (int j = 0; j < reacPart.length; j++) {
                    Action action = null;
                    SpeciesCountParameter spCountParam = getSpeciesCountParameter(reacPart[j].getSpeciesContext());
                    if (reacPart[j] instanceof Reactant) {
                        // check if the reactant is a constant. If the species is a constant, there will be no action taken on this species
                        if (// not a constant
                        !simContext.getReactionContext().getSpeciesContextSpec(reacPart[j].getSpeciesContext()).isConstant()) {
                            int stoi = ((Reactant) reacPart[j]).getStoichiometry();
                            action = Action.createIncrementAction(varHash.getVariable(getMathSymbol(spCountParam, geometryClass)), new Expression(-stoi));
                            jp.addAction(action);
                        }
                    } else if (reacPart[j] instanceof Product) {
                        // check if the product is a constant. If the product is a constant, there will be no action taken on this species
                        if (// not a constant
                        !simContext.getReactionContext().getSpeciesContextSpec(reacPart[j].getSpeciesContext()).isConstant()) {
                            int stoi = ((Product) reacPart[j]).getStoichiometry();
                            action = Action.createIncrementAction(varHash.getVariable(getMathSymbol(spCountParam, geometryClass)), new Expression(stoi));
                            jp.addAction(action);
                        }
                    }
                }
                // add jump process to compartment subDomain
                subDomain.addJumpProcess(jp);
            }
            if (// one more jump process for a reversible reaction
            isReverseRatePresent) {
                // get jump process name
                String jpName = TokenMangler.mangleToSName(reactionStep.getName()) + PARAMETER_PROBABILITY_RATE_REVERSE_SUFFIX;
                Expression exp = null;
                // reactions are mass actions
                exp = getProbabilityRate(reactionStep, reverseRate, false);
                ProbabilityParameter probRevParm = null;
                try {
                    probRevParm = addProbabilityParameter(PARAMETER_PROBABILITYRATE_PREFIX + jpName, exp, PARAMETER_ROLE_P_reverse, probabilityParamUnit, reactionStep);
                } catch (PropertyVetoException pve) {
                    pve.printStackTrace();
                    throw new MappingException(pve.getMessage());
                }
                // add probability to function or constant
                varHash.addVariable(newFunctionOrConstant(getMathSymbol(probRevParm, geometryClass), getIdentifierSubstitutions(exp, probabilityParamUnit, geometryClass), geometryClass));
                JumpProcess jp = new JumpProcess(jpName, new Expression(getMathSymbol(probRevParm, geometryClass)));
                // actions
                ReactionParticipant[] reacPart = reactionStep.getReactionParticipants();
                for (int j = 0; j < reacPart.length; j++) {
                    Action action = null;
                    SpeciesCountParameter spCountParam = getSpeciesCountParameter(reacPart[j].getSpeciesContext());
                    if (reacPart[j] instanceof Reactant) {
                        // check if the reactant is a constant. If the species is a constant, there will be no action taken on this species
                        if (// not a constant
                        !simContext.getReactionContext().getSpeciesContextSpec(reacPart[j].getSpeciesContext()).isConstant()) {
                            int stoi = ((Reactant) reacPart[j]).getStoichiometry();
                            action = Action.createIncrementAction(varHash.getVariable(getMathSymbol(spCountParam, geometryClass)), new Expression(stoi));
                            jp.addAction(action);
                        }
                    } else if (reacPart[j] instanceof Product) {
                        // check if the product is a constant. If the product is a constant, there will be no action taken on this species
                        if (// not a constant
                        !simContext.getReactionContext().getSpeciesContextSpec(reacPart[j].getSpeciesContext()).isConstant()) {
                            int stoi = ((Product) reacPart[j]).getStoichiometry();
                            action = Action.createIncrementAction(varHash.getVariable(getMathSymbol(spCountParam, geometryClass)), new Expression(-stoi));
                            jp.addAction(action);
                        }
                    }
                }
                // add jump process to compartment subDomain
                subDomain.addJumpProcess(jp);
            }
        // end of if(isForwardRateNonZero), if(isReverseRateNonRate)
        } else if (// flux reactions
        reactionStep instanceof FluxReaction) {
            // we could set jump processes for general flux rate in forms of p1*Sout + p2*Sin
            if (kinetics.getKineticsDescription().equals(KineticsDescription.General) || kinetics.getKineticsDescription().equals(KineticsDescription.GeneralPermeability)) {
                Expression fluxRate = new Expression(kinetics.getKineticsParameterFromRole(Kinetics.ROLE_ReactionRate), reactionStep.getNameScope());
                // we have to pass the math description para to flux solver, coz somehow math description in simulation context is not updated.
                // forward and reverse rate parameters may be null
                Parameter forwardRateParameter = null;
                Parameter reverseRateParameter = null;
                if (kinetics.getKineticsDescription().equals(KineticsDescription.GeneralPermeability)) {
                    forwardRateParameter = kinetics.getKineticsParameterFromRole(Kinetics.ROLE_Permeability);
                    reverseRateParameter = kinetics.getKineticsParameterFromRole(Kinetics.ROLE_Permeability);
                }
                MassActionSolver.MassActionFunction fluxFunc = MassActionSolver.solveMassAction(forwardRateParameter, reverseRateParameter, fluxRate, (FluxReaction) reactionStep);
                // create jump process for forward flux if it exists.
                Expression rsStructureSize = new Expression(reactionStep.getStructure().getStructureSize(), getNameScope());
                VCUnitDefinition probRateUnit = modelUnitSystem.getStochasticSubstanceUnit().divideBy(modelUnitSystem.getAreaUnit()).divideBy(modelUnitSystem.getTimeUnit());
                Expression rsRateUnitFactor = getUnitFactor(probRateUnit.divideBy(modelUnitSystem.getFluxReactionUnit()));
                if (fluxFunc.getForwardRate() != null && !fluxFunc.getForwardRate().isZero()) {
                    Expression rate = fluxFunc.getForwardRate();
                    // get species expression (depend on structure, if mem: Species/mem_Size, if vol: species*KMOLE/vol_size)
                    if (fluxFunc.getReactants().size() != 1) {
                        throw new MappingException("Flux " + reactionStep.getName() + " should have only one reactant.");
                    }
                    SpeciesContext scReactant = fluxFunc.getReactants().get(0).getSpeciesContext();
                    Expression scConcExpr = new Expression(getSpeciesConcentrationParameter(scReactant), getNameScope());
                    Expression probExp = Expression.mult(rate, rsRateUnitFactor, rsStructureSize, scConcExpr);
                    // jump process name
                    // +"_reverse";
                    String jpName = TokenMangler.mangleToSName(reactionStep.getName());
                    ProbabilityParameter probParm = null;
                    try {
                        probParm = addProbabilityParameter(PARAMETER_PROBABILITYRATE_PREFIX + jpName, probExp, PARAMETER_ROLE_P, probabilityParamUnit, reactionStep);
                    } catch (PropertyVetoException pve) {
                        pve.printStackTrace();
                        throw new MappingException(pve.getMessage());
                    }
                    // add probability to function or constant
                    String ms = getMathSymbol(probParm, geometryClass);
                    Expression is = getIdentifierSubstitutions(probExp, probabilityParamUnit, geometryClass);
                    Variable nfoc = newFunctionOrConstant(ms, is, geometryClass);
                    varHash.addVariable(nfoc);
                    JumpProcess jp = new JumpProcess(jpName, new Expression(getMathSymbol(probParm, geometryClass)));
                    // actions
                    Action action = null;
                    SpeciesContext sc = fluxFunc.getReactants().get(0).getSpeciesContext();
                    if (!simContext.getReactionContext().getSpeciesContextSpec(sc).isConstant()) {
                        SpeciesCountParameter spCountParam = getSpeciesCountParameter(sc);
                        action = Action.createIncrementAction(varHash.getVariable(getMathSymbol(spCountParam, geometryClass)), new Expression(-1));
                        jp.addAction(action);
                    }
                    sc = fluxFunc.getProducts().get(0).getSpeciesContext();
                    if (!simContext.getReactionContext().getSpeciesContextSpec(sc).isConstant()) {
                        SpeciesCountParameter spCountParam = getSpeciesCountParameter(sc);
                        action = Action.createIncrementAction(varHash.getVariable(getMathSymbol(spCountParam, geometryClass)), new Expression(1));
                        jp.addAction(action);
                    }
                    subDomain.addJumpProcess(jp);
                }
                // create jump process for reverse flux if it exists.
                if (fluxFunc.getReverseRate() != null && !fluxFunc.getReverseRate().isZero()) {
                    // jump process name
                    String jpName = TokenMangler.mangleToSName(reactionStep.getName()) + PARAMETER_PROBABILITY_RATE_REVERSE_SUFFIX;
                    Expression rate = fluxFunc.getReverseRate();
                    // get species expression (depend on structure, if mem: Species/mem_Size, if vol: species*KMOLE/vol_size)
                    if (fluxFunc.getProducts().size() != 1) {
                        throw new MappingException("Flux " + reactionStep.getName() + " should have only one product.");
                    }
                    SpeciesContext scProduct = fluxFunc.getProducts().get(0).getSpeciesContext();
                    Expression scConcExpr = new Expression(getSpeciesConcentrationParameter(scProduct), getNameScope());
                    Expression probExp = Expression.mult(rate, rsRateUnitFactor, rsStructureSize, scConcExpr);
                    ProbabilityParameter probRevParm = null;
                    try {
                        probRevParm = addProbabilityParameter(PARAMETER_PROBABILITYRATE_PREFIX + jpName, probExp, PARAMETER_ROLE_P_reverse, probabilityParamUnit, reactionStep);
                    } catch (PropertyVetoException pve) {
                        pve.printStackTrace();
                        throw new MappingException(pve.getMessage());
                    }
                    // add probability to function or constant
                    varHash.addVariable(newFunctionOrConstant(getMathSymbol(probRevParm, geometryClass), getIdentifierSubstitutions(probExp, probabilityParamUnit, geometryClass), geometryClass));
                    JumpProcess jp = new JumpProcess(jpName, new Expression(getMathSymbol(probRevParm, geometryClass)));
                    // actions
                    Action action = null;
                    SpeciesContext sc = fluxFunc.getReactants().get(0).getSpeciesContext();
                    if (!simContext.getReactionContext().getSpeciesContextSpec(sc).isConstant()) {
                        SpeciesCountParameter spCountParam = getSpeciesCountParameter(sc);
                        action = Action.createIncrementAction(varHash.getVariable(getMathSymbol(spCountParam, geometryClass)), new Expression(1));
                        jp.addAction(action);
                    }
                    sc = fluxFunc.getProducts().get(0).getSpeciesContext();
                    if (!simContext.getReactionContext().getSpeciesContextSpec(sc).isConstant()) {
                        SpeciesCountParameter spCountParam = getSpeciesCountParameter(sc);
                        action = Action.createIncrementAction(varHash.getVariable(getMathSymbol(spCountParam, geometryClass)), new Expression(-1));
                        jp.addAction(action);
                    }
                    subDomain.addJumpProcess(jp);
                }
            }
        }
    // end of if (simplereaction)...else if(fluxreaction)
    }
// end of reaction step loop
}
Also used : Action(cbit.vcell.math.Action) StochVolVariable(cbit.vcell.math.StochVolVariable) Variable(cbit.vcell.math.Variable) Product(cbit.vcell.model.Product) FluxReaction(cbit.vcell.model.FluxReaction) SpeciesContext(cbit.vcell.model.SpeciesContext) Reactant(cbit.vcell.model.Reactant) ExpressionException(cbit.vcell.parser.ExpressionException) JumpProcess(cbit.vcell.math.JumpProcess) ModelUnitSystem(cbit.vcell.model.ModelUnitSystem) SimpleReaction(cbit.vcell.model.SimpleReaction) PropertyVetoException(java.beans.PropertyVetoException) VCUnitDefinition(cbit.vcell.units.VCUnitDefinition) Expression(cbit.vcell.parser.Expression) MathException(cbit.vcell.math.MathException) ReactionStep(cbit.vcell.model.ReactionStep) Parameter(cbit.vcell.model.Parameter) KineticsParameter(cbit.vcell.model.Kinetics.KineticsParameter) ModelParameter(cbit.vcell.model.Model.ModelParameter) Kinetics(cbit.vcell.model.Kinetics) LumpedKinetics(cbit.vcell.model.LumpedKinetics) MassActionSolver(cbit.vcell.model.MassActionSolver) ReactionParticipant(cbit.vcell.model.ReactionParticipant)

Aggregations

Parameter (cbit.vcell.model.Parameter)54 Expression (cbit.vcell.parser.Expression)33 ModelParameter (cbit.vcell.model.Model.ModelParameter)29 LocalParameter (cbit.vcell.mapping.ParameterContext.LocalParameter)19 KineticsParameter (cbit.vcell.model.Kinetics.KineticsParameter)19 SpeciesContext (cbit.vcell.model.SpeciesContext)19 ExpressionException (cbit.vcell.parser.ExpressionException)18 Model (cbit.vcell.model.Model)15 VCUnitDefinition (cbit.vcell.units.VCUnitDefinition)15 PropertyVetoException (java.beans.PropertyVetoException)15 ReactionStep (cbit.vcell.model.ReactionStep)14 SpeciesContextSpecParameter (cbit.vcell.mapping.SpeciesContextSpec.SpeciesContextSpecParameter)13 ModelUnitSystem (cbit.vcell.model.ModelUnitSystem)12 SpeciesContextSpec (cbit.vcell.mapping.SpeciesContextSpec)11 Structure (cbit.vcell.model.Structure)11 ArrayList (java.util.ArrayList)10 StructureMapping (cbit.vcell.mapping.StructureMapping)9 MathDescription (cbit.vcell.math.MathDescription)9 ProxyParameter (cbit.vcell.model.ProxyParameter)9 SymbolTableEntry (cbit.vcell.parser.SymbolTableEntry)9