Search in sources :

Example 36 with Function

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

the class ApplicationConstraintsGenerator method steadyStateFromApplication.

/**
 * Insert the method's description here.
 * Creation date: (6/26/01 8:25:55 AM)
 * @return cbit.vcell.constraints.ConstraintContainerImpl
 */
public static ConstraintContainerImpl steadyStateFromApplication(SimulationContext simContext, double tolerance) {
    try {
        ConstraintContainerImpl ccImpl = new ConstraintContainerImpl();
        // ====================
        // add physical limits
        // ====================
        // 
        // no negative concentrations
        // 
        cbit.vcell.model.Model model = simContext.getModel();
        cbit.vcell.model.SpeciesContext[] speciesContexts = model.getSpeciesContexts();
        for (int i = 0; i < speciesContexts.length; i++) {
            ccImpl.addSimpleBound(new SimpleBounds(speciesContexts[i].getName(), new RealInterval(0, Double.POSITIVE_INFINITY), AbstractConstraint.PHYSICAL_LIMIT, "non-negative concentration"));
        }
        for (int i = 0; i < speciesContexts.length; i++) {
            SpeciesContextSpecParameter initParam = (simContext.getReactionContext().getSpeciesContextSpec(speciesContexts[i])).getInitialConditionParameter();
            if (initParam != null) {
                double initialValue = initParam.getExpression().evaluateConstant();
                double lowInitialValue = Math.min(initialValue / tolerance, initialValue * tolerance);
                double highInitialValue = Math.max(initialValue / tolerance, initialValue * tolerance);
                ccImpl.addSimpleBound(new SimpleBounds(speciesContexts[i].getName(), new RealInterval(lowInitialValue, highInitialValue), AbstractConstraint.MODELING_ASSUMPTION, "close to specified \"initialCondition\""));
            }
        }
        // =========================
        // add modeling assumptions
        // =========================
        // 
        // mass action forward and reverse rates should be non-negative
        // 
        cbit.vcell.model.ReactionStep[] reactionSteps = model.getReactionSteps();
        for (int i = 0; i < reactionSteps.length; i++) {
            Kinetics kinetics = reactionSteps[i].getKinetics();
            if (kinetics instanceof MassActionKinetics) {
                Expression forwardRateConstraintExp = new Expression(((MassActionKinetics) kinetics).getForwardRateParameter().getExpression().infix() + ">=0");
                forwardRateConstraintExp = getSteadyStateExpression(forwardRateConstraintExp);
                if (!forwardRateConstraintExp.compareEqual(new Expression(1.0))) {
                    ccImpl.addGeneralConstraint(new GeneralConstraint(forwardRateConstraintExp, AbstractConstraint.MODELING_ASSUMPTION, "non-negative forward rate"));
                }
                Expression reverseRateConstraintExp = new Expression(((MassActionKinetics) kinetics).getReverseRateParameter().getExpression().infix() + ">=0");
                reverseRateConstraintExp = getSteadyStateExpression(reverseRateConstraintExp);
                if (!reverseRateConstraintExp.compareEqual(new Expression(1.0))) {
                    ccImpl.addGeneralConstraint(new GeneralConstraint(reverseRateConstraintExp, AbstractConstraint.MODELING_ASSUMPTION, "non-negative reverse rate"));
                }
            }
            KineticsParameter authoritativeParameter = kinetics.getAuthoritativeParameter();
            Expression kineticRateConstraintExp = new Expression(authoritativeParameter.getName() + "==" + authoritativeParameter.getExpression().infix());
            kineticRateConstraintExp = getSteadyStateExpression(kineticRateConstraintExp);
            if (!kineticRateConstraintExp.compareEqual(new Expression(1.0))) {
                ccImpl.addGeneralConstraint(new GeneralConstraint(kineticRateConstraintExp, AbstractConstraint.MODELING_ASSUMPTION, "definition"));
            }
        }
        // 
        try {
            simContext.setMathDescription(simContext.createNewMathMapping().getMathDescription());
        } catch (Throwable e) {
            e.printStackTrace(System.out);
            throw new RuntimeException("cannot create mathDescription");
        }
        MathDescription mathDesc = simContext.getMathDescription();
        if (mathDesc.getGeometry().getDimension() > 0) {
            throw new RuntimeException("spatial simulations not yet supported");
        }
        CompartmentSubDomain subDomain = (CompartmentSubDomain) mathDesc.getSubDomains().nextElement();
        java.util.Enumeration<Equation> enumEquations = subDomain.getEquations();
        while (enumEquations.hasMoreElements()) {
            Equation equation = (Equation) enumEquations.nextElement();
            Expression rateConstraintExp = new Expression(equation.getRateExpression().infix() + "==0");
            rateConstraintExp = getSteadyStateExpression(rateConstraintExp);
            if (!rateConstraintExp.compareEqual(new Expression(1.0))) {
                // not a trivial constraint (always true)
                ccImpl.addGeneralConstraint(new GeneralConstraint(rateConstraintExp, AbstractConstraint.PHYSICAL_LIMIT, "definition of steady state"));
            }
        }
        // 
        for (int i = 0; i < reactionSteps.length; i++) {
            Kinetics kinetics = reactionSteps[i].getKinetics();
            Kinetics.KineticsParameter[] parameters = kinetics.getKineticsParameters();
            for (int j = 0; j < parameters.length; j++) {
                Expression exp = parameters[j].getExpression();
                if (exp.getSymbols() == null || exp.getSymbols().length == 0) {
                    // 
                    try {
                        double constantValue = exp.evaluateConstant();
                        double lowValue = Math.min(constantValue / tolerance, constantValue * tolerance);
                        double highValue = Math.max(constantValue / tolerance, constantValue * tolerance);
                        RealInterval interval = new RealInterval(lowValue, highValue);
                        ccImpl.addSimpleBound(new SimpleBounds(parameters[j].getName(), interval, AbstractConstraint.MODELING_ASSUMPTION, "parameter close to model default"));
                    } catch (cbit.vcell.parser.ExpressionException e) {
                        System.out.println("error evaluating parameter " + parameters[j].getName() + " in reaction step " + reactionSteps[i].getName());
                    }
                } else {
                    Expression parameterDefinitionExp = new Expression(parameters[j].getName() + "==" + parameters[j].getExpression().infix());
                    ccImpl.addGeneralConstraint(new GeneralConstraint(getSteadyStateExpression(parameterDefinitionExp), AbstractConstraint.MODELING_ASSUMPTION, "parameter definition"));
                }
            }
        }
        ccImpl.addSimpleBound(new SimpleBounds(model.getFARADAY_CONSTANT().getName(), new RealInterval(model.getFARADAY_CONSTANT().getExpression().evaluateConstant()), AbstractConstraint.PHYSICAL_LIMIT, "Faraday's constant"));
        ccImpl.addSimpleBound(new SimpleBounds(model.getTEMPERATURE().getName(), new RealInterval(300), AbstractConstraint.PHYSICAL_LIMIT, "Absolute Temperature Kelvin"));
        ccImpl.addSimpleBound(new SimpleBounds(model.getGAS_CONSTANT().getName(), new RealInterval(model.getGAS_CONSTANT().getExpression().evaluateConstant()), AbstractConstraint.PHYSICAL_LIMIT, "ideal gas constant"));
        ccImpl.addSimpleBound(new SimpleBounds(model.getKMILLIVOLTS().getName(), new RealInterval(model.getKMILLIVOLTS().getExpression().evaluateConstant()), AbstractConstraint.PHYSICAL_LIMIT, "ideal gas constant"));
        // 
        // add K_fluxs
        // 
        java.util.Enumeration<Variable> enumVars = mathDesc.getVariables();
        while (enumVars.hasMoreElements()) {
            Variable var = (Variable) enumVars.nextElement();
            if (var.getName().startsWith("Kflux_") && var instanceof Function) {
                Expression kfluxExp = new Expression(((Function) var).getExpression());
                kfluxExp.bindExpression(mathDesc);
                kfluxExp = MathUtilities.substituteFunctions(kfluxExp, mathDesc);
                kfluxExp = kfluxExp.flatten();
                ccImpl.addSimpleBound(new SimpleBounds(var.getName(), new RealInterval(kfluxExp.evaluateConstant()), AbstractConstraint.MODELING_ASSUMPTION, "flux conversion factor"));
            }
        }
        return ccImpl;
    } catch (cbit.vcell.parser.ExpressionException e) {
        e.printStackTrace(System.out);
        return null;
    } catch (java.beans.PropertyVetoException e) {
        e.printStackTrace(System.out);
        return null;
    }
}
Also used : Variable(cbit.vcell.math.Variable) SimpleBounds(cbit.vcell.constraints.SimpleBounds) MathDescription(cbit.vcell.math.MathDescription) GeneralConstraint(cbit.vcell.constraints.GeneralConstraint) RealInterval(net.sourceforge.interval.ia_math.RealInterval) Function(cbit.vcell.math.Function) KineticsParameter(cbit.vcell.model.Kinetics.KineticsParameter) ConstraintContainerImpl(cbit.vcell.constraints.ConstraintContainerImpl) SpeciesContextSpecParameter(cbit.vcell.mapping.SpeciesContextSpec.SpeciesContextSpecParameter) Equation(cbit.vcell.math.Equation) AbstractConstraint(cbit.vcell.constraints.AbstractConstraint) GeneralConstraint(cbit.vcell.constraints.GeneralConstraint) Expression(cbit.vcell.parser.Expression) CompartmentSubDomain(cbit.vcell.math.CompartmentSubDomain) MassActionKinetics(cbit.vcell.model.MassActionKinetics) MassActionKinetics(cbit.vcell.model.MassActionKinetics) Kinetics(cbit.vcell.model.Kinetics)

Example 37 with Function

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

the class PDEDataViewer method calcAutoAllTimes.

private void calcAutoAllTimes() throws Exception {
    HashSet<String> stateVarNames = null;
    Variable theVariable = null;
    boolean bStateVar = true;
    boolean isFieldData = getPdeDataContext().getVCDataIdentifier() instanceof ExternalDataIdentifier || getPdeDataContext().getVCDataIdentifier() instanceof MergedDataInfo;
    if (isFieldData) {
        // fielddata
        DataIdentifier[] dataids = getPdeDataContext().getDataIdentifiers();
        stateVarNames = new HashSet<>();
        for (int i = 0; i < dataids.length; i++) {
            if (!dataids[i].isFunction()) {
                stateVarNames.add(dataids[i].getName());
            }
        // System.out.println("name:'"+dataids[i].getName()+"' type:"+dataids[i].getVariableType()+" func:"+dataids[i].isFunction());
        }
        bStateVar = !getPdeDataContext().getDataIdentifier().isFunction();
        if (bStateVar) {
            theVariable = new VolVariable(getPdeDataContext().getDataIdentifier().getName(), getPdeDataContext().getDataIdentifier().getDomain());
        } else {
            AnnotatedFunction[] funcs = getPdeDataContext().getFunctions();
            for (int i = 0; i < funcs.length; i++) {
                if (funcs[i].getName().equals(getPdeDataContext().getDataIdentifier().getName())) {
                    theVariable = funcs[i];
                    break;
                }
            }
        }
    } else {
        stateVarNames = getSimulation().getMathDescription().getStateVariableNames();
        theVariable = getSimulation().getMathDescription().getVariable(getPdeDataContext().getVariableName());
        if (theVariable == null) {
            theVariable = ((ClientPDEDataContext) getPdeDataContext()).getDataManager().getOutputContext().getOutputFunction(getPdeDataContext().getVariableName());
        }
        if (theVariable == null) {
            DataProcessingOutputInfo dataProcessingOutputInfo = DataProcessingResultsPanel.getDataProcessingOutputInfo(getPdeDataContext());
            if (dataProcessingOutputInfo != null && Arrays.asList(dataProcessingOutputInfo.getVariableNames()).contains(getPdeDataContext().getVariableName())) {
                // PostProcess Variable
                return;
            }
        }
        bStateVar = stateVarNames.contains(getPdeDataContext().getVariableName());
    }
    if (theVariable == null) {
        throw new Exception("Unexpected Alltimes... selected variable '" + getPdeDataContext().getVariableName() + "' is not stateVariable or OutputFunction");
    }
    if (getPDEDataContextPanel1().getdisplayAdapterService1().getAllTimes()) {
        // min-max over all timepoints (allTimes)
        if (theVariable.isConstant()) {
            getPDEDataContextPanel1().getdisplayAdapterServicePanel1().changeAllTimesButtonText(DisplayAdapterServicePanel.ALL_TIMES__STATE_TEXT);
            double constVal = theVariable.getExpression().evaluateConstant();
            getPDEDataContextPanel1().setFunctionStatisticsRange(new Range(constVal, constVal));
        } else if (bStateVar) {
            getPDEDataContextPanel1().getdisplayAdapterServicePanel1().changeAllTimesButtonText(DisplayAdapterServicePanel.ALL_TIMES__STATE_TEXT);
            ArrayList<VarStatistics> varStatsArr = calcVarStat(getPdeDataContext(), new String[] { theVariable.getName() });
            if (errorAutoAllTimes(varStatsArr != null, (varStatsArr == null ? null : varStatsArr.size() > 0), isFieldData)) {
                // no postprocessinfo
                return;
            }
            FunctionStatistics functionStatistics = new FunctionStatistics(varStatsArr.get(0).minValuesOverTime, varStatsArr.get(0).maxValuesOverTime);
            getPDEDataContextPanel1().setFunctionStatisticsRange(new Range(functionStatistics.getMinOverTime(), functionStatistics.getMaxOverTime()));
        } else if (theVariable instanceof Function) {
            getPDEDataContextPanel1().getdisplayAdapterServicePanel1().changeAllTimesButtonText(DisplayAdapterServicePanel.ALL_TIMES__APPROX_TEXT);
            Function flattened = MathDescription.getFlattenedFunctions(SimulationSymbolTable.createMathSymbolTableFactory(), getSimulation().getMathDescription(), new String[] { theVariable.getName() })[0];
            if (flattened == null) {
                flattened = (Function) theVariable;
            }
            ArrayList<VarStatistics> varStatsArr = calcVarStat(getPdeDataContext(), stateVarNames.toArray(new String[0]));
            if (errorAutoAllTimes(varStatsArr != null, (varStatsArr == null ? null : varStatsArr.size() > 0), isFieldData)) {
                // check for no postprocessinfo
                return;
            }
            if (varStatsArr.size() == stateVarNames.size()) {
                if (getSimulation().getMeshSpecification().getGeometry().getGeometrySurfaceDescription().getRegionImage() == null) {
                    getSimulation().getMeshSpecification().getGeometry().getGeometrySurfaceDescription().updateAll();
                }
                FunctionStatistics functionStatistics = FunctionRangeGenerator.getFunctionStatistics(flattened.getExpression(), varStatsArr.toArray(new VarStatistics[0]), getPdeDataContext().getTimePoints(), getPdeDataContext().getCartesianMesh(), calcInDomainBitSet(), getPdeDataContext().getDataIdentifier().getVariableType());
                getPDEDataContextPanel1().setFunctionStatisticsRange(new Range(functionStatistics.getMinOverTime(), functionStatistics.getMaxOverTime()));
            } else {
                throw new Exception("Unexpectede AllTimes... calculated state var stats size != mathdescr state var size");
            }
        } else {
            throw new Exception("Unexpected AllTimes... not constant, stateVar or function");
        }
    } else {
        // min-max at each timepoint (currTime)
        if (!(theVariable instanceof Function)) {
            getPDEDataContextPanel1().getdisplayAdapterServicePanel1().changeAllTimesButtonText(DisplayAdapterServicePanel.ALL_TIMES__STATE_TEXT);
        } else {
            getPDEDataContextPanel1().getdisplayAdapterServicePanel1().changeAllTimesButtonText(DisplayAdapterServicePanel.ALL_TIMES__APPROX_TEXT);
        }
        getPDEDataContextPanel1().setFunctionStatisticsRange(null);
    }
}
Also used : VolVariable(cbit.vcell.math.VolVariable) ReservedVariable(cbit.vcell.math.ReservedVariable) Variable(cbit.vcell.math.Variable) VCSimulationDataIdentifier(cbit.vcell.solver.VCSimulationDataIdentifier) LocalVCSimulationDataIdentifier(cbit.vcell.client.ClientSimManager.LocalVCSimulationDataIdentifier) ExternalDataIdentifier(org.vcell.util.document.ExternalDataIdentifier) VCDataIdentifier(org.vcell.util.document.VCDataIdentifier) DataIdentifier(cbit.vcell.simdata.DataIdentifier) VolVariable(cbit.vcell.math.VolVariable) ArrayList(java.util.ArrayList) Range(org.vcell.util.Range) Point(java.awt.Point) SinglePoint(cbit.vcell.geometry.SinglePoint) DataAccessException(org.vcell.util.DataAccessException) PropertyVetoException(java.beans.PropertyVetoException) ImageException(cbit.image.ImageException) UserCancelException(org.vcell.util.UserCancelException) MergedDataInfo(cbit.vcell.simdata.MergedDataInfo) AnnotatedFunction(cbit.vcell.solver.AnnotatedFunction) Function(cbit.vcell.math.Function) VarStatistics(cbit.vcell.util.FunctionRangeGenerator.VarStatistics) DataProcessingOutputInfo(cbit.vcell.simdata.DataOperationResults.DataProcessingOutputInfo) ExternalDataIdentifier(org.vcell.util.document.ExternalDataIdentifier) AnnotatedFunction(cbit.vcell.solver.AnnotatedFunction) FunctionStatistics(cbit.vcell.util.FunctionRangeGenerator.FunctionStatistics)

Example 38 with Function

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

the class XmlReader method getKinetics.

/**
 * This method returns a Kinetics object from a XML Element based on the value of the kinetics type attribute.
 * Creation date: (3/19/2001 4:42:04 PM)
 * @return cbit.vcell.model.Kinetics
 * @param param org.jdom.Element
 */
private Kinetics getKinetics(Element param, ReactionStep reaction, Model model) throws XmlParseException {
    VariableHash varHash = new VariableHash();
    addResevedSymbols(varHash, model);
    String type = param.getAttributeValue(XMLTags.KineticsTypeAttrTag);
    Kinetics newKinetics = null;
    try {
        if (type.equalsIgnoreCase(XMLTags.KineticsTypeGeneralKinetics)) {
            // create a general kinetics
            newKinetics = new GeneralKinetics(reaction);
        } else if (type.equalsIgnoreCase(XMLTags.KineticsTypeGeneralCurrentKinetics)) {
            // Create GeneralCurrentKinetics
            newKinetics = new GeneralCurrentKinetics(reaction);
        } else if (type.equalsIgnoreCase(XMLTags.KineticsTypeMassAction) && reaction instanceof SimpleReaction) {
            // create a Mass Action kinetics
            newKinetics = new MassActionKinetics((SimpleReaction) reaction);
        } else if (type.equalsIgnoreCase(XMLTags.KineticsTypeNernst) && reaction instanceof FluxReaction) {
            // create NernstKinetics
            newKinetics = new NernstKinetics((FluxReaction) reaction);
        } else if (type.equalsIgnoreCase(XMLTags.KineticsTypeGHK) && reaction instanceof FluxReaction) {
            // create GHKKinetics
            newKinetics = new GHKKinetics((FluxReaction) reaction);
        } else if (type.equalsIgnoreCase(XMLTags.KineticsTypeHMM_Irr) && reaction instanceof SimpleReaction) {
            // create HMM_IrrKinetics
            newKinetics = new HMM_IRRKinetics((SimpleReaction) reaction);
        } else if (type.equalsIgnoreCase(XMLTags.KineticsTypeHMM_Rev) && reaction instanceof SimpleReaction) {
            // create HMM_RevKinetics
            newKinetics = new HMM_REVKinetics((SimpleReaction) reaction);
        } else if (type.equalsIgnoreCase(XMLTags.KineticsTypeGeneralTotal_oldname)) {
            // create GeneralTotalKinetics
            newKinetics = new GeneralLumpedKinetics(reaction);
        } else if (type.equalsIgnoreCase(XMLTags.KineticsTypeGeneralLumped)) {
            // create GeneralLumpedKinetics
            newKinetics = new GeneralLumpedKinetics(reaction);
        } else if (type.equalsIgnoreCase(XMLTags.KineticsTypeGeneralCurrentLumped)) {
            // create GeneralCurrentLumpedKinetics
            newKinetics = new GeneralCurrentLumpedKinetics(reaction);
        } else if (type.equalsIgnoreCase(XMLTags.KineticsTypeGeneralPermeability) && reaction instanceof FluxReaction) {
            // create GeneralPermeabilityKinetics
            newKinetics = new GeneralPermeabilityKinetics((FluxReaction) reaction);
        } else if (type.equalsIgnoreCase(XMLTags.KineticsTypeMacroscopic_Irr) && reaction instanceof SimpleReaction) {
            // create Macroscopic_IRRKinetics
            newKinetics = new Macroscopic_IRRKinetics((SimpleReaction) reaction);
        } else if (type.equalsIgnoreCase(XMLTags.KineticsTypeMicroscopic_Irr) && reaction instanceof SimpleReaction) {
            // create Microscopic_IRRKinetics
            newKinetics = new Microscopic_IRRKinetics((SimpleReaction) reaction);
        } else {
            throw new XmlParseException("Unknown kinetics type: " + type);
        }
    } catch (ExpressionException e) {
        e.printStackTrace();
        throw new XmlParseException("Error creating the kinetics for reaction: " + reaction.getName(), e);
    }
    try {
        // transaction begin flag ... yeah, this is a hack
        newKinetics.reading(true);
        // Read all of the parameters
        List<Element> list = param.getChildren(XMLTags.ParameterTag, vcNamespace);
        // add constants that may be used in kinetics.
        // VariableHash varHash = getVariablesHash();
        ArrayList<String> reserved = new ArrayList<String>();
        ReservedSymbol[] reservedSymbols = reaction.getModel().getReservedSymbols();
        for (ReservedSymbol rs : reservedSymbols) {
            reserved.add(rs.getName());
        }
        try {
            if (reaction.getStructure() instanceof Membrane) {
                Membrane membrane = (Membrane) reaction.getStructure();
                varHash.addVariable(new Constant(membrane.getMembraneVoltage().getName(), new Expression(0.0)));
                reserved.add(membrane.getMembraneVoltage().getName());
            }
            // 
            // add Reactants, Products, and Catalysts (ReactionParticipants)
            // 
            ReactionParticipant[] rp = reaction.getReactionParticipants();
            for (int i = 0; i < rp.length; i++) {
                varHash.addVariable(new Constant(rp[i].getName(), new Expression(0.0)));
            }
        } catch (MathException e) {
            e.printStackTrace(System.out);
            throw new XmlParseException("error reordering parameters according to dependencies: ", e);
        }
        // 
        for (Element xmlParam : list) {
            String paramName = unMangle(xmlParam.getAttributeValue(XMLTags.NameAttrTag));
            String role = xmlParam.getAttributeValue(XMLTags.ParamRoleAttrTag);
            String paramExpStr = xmlParam.getText();
            Expression paramExp = unMangleExpression(paramExpStr);
            try {
                if (varHash.getVariable(paramName) == null) {
                    varHash.addVariable(new Function(paramName, paramExp, null));
                } else {
                    if (reserved.contains(paramName)) {
                        varHash.removeVariable(paramName);
                        varHash.addVariable(new Function(paramName, paramExp, null));
                    }
                }
            } catch (MathException e) {
                e.printStackTrace(System.out);
                throw new XmlParseException("error reordering parameters according to dependencies: ", e);
            }
            Kinetics.KineticsParameter tempParam = null;
            if (!role.equals(XMLTags.ParamRoleUserDefinedTag)) {
                tempParam = newKinetics.getKineticsParameterFromRole(Kinetics.getParamRoleFromDefaultDesc(role));
            } else {
                continue;
            }
            // hack for bringing in General Total kinetics without breaking.
            if (tempParam == null && newKinetics instanceof GeneralLumpedKinetics) {
                if (role.equals(Kinetics.GTK_AssumedCompartmentSize_oldname) || role.equals(Kinetics.GTK_ReactionRate_oldname) || role.equals(Kinetics.GTK_CurrentDensity_oldname)) {
                    continue;
                } else if (role.equals(VCMODL.TotalRate_oldname)) {
                    tempParam = newKinetics.getKineticsParameterFromRole(Kinetics.ROLE_LumpedReactionRate);
                }
            }
            // hack from bringing in chargeValence parameters without breaking
            if (tempParam == null && Kinetics.getParamRoleFromDefaultDesc(role) == Kinetics.ROLE_ChargeValence) {
                tempParam = newKinetics.getChargeValenceParameter();
            }
            if (tempParam == null) {
                throw new XmlParseException("parameter with role '" + role + "' not found in kinetics type '" + type + "'");
            }
            // 
            if (!tempParam.getName().equals(paramName)) {
                Kinetics.KineticsParameter multNameParam = newKinetics.getKineticsParameter(paramName);
                int n = 0;
                while (multNameParam != null) {
                    String tempName = paramName + "_" + n++;
                    newKinetics.renameParameter(paramName, tempName);
                    multNameParam = newKinetics.getKineticsParameter(tempName);
                }
                newKinetics.renameParameter(tempParam.getName(), paramName);
            }
        }
        // 
        // create unresolved parameters for all unresolved symbols
        // 
        String unresolvedSymbol = varHash.getFirstUnresolvedSymbol();
        while (unresolvedSymbol != null) {
            try {
                // will turn into an UnresolvedParameter.
                varHash.addVariable(new Function(unresolvedSymbol, new Expression(0.0), null));
            } catch (MathException e) {
                e.printStackTrace(System.out);
                throw new XmlParseException(e);
            }
            newKinetics.addUnresolvedParameter(unresolvedSymbol);
            unresolvedSymbol = varHash.getFirstUnresolvedSymbol();
        }
        Variable[] sortedVariables = varHash.getTopologicallyReorderedVariables();
        ModelUnitSystem modelUnitSystem = reaction.getModel().getUnitSystem();
        for (int i = sortedVariables.length - 1; i >= 0; i--) {
            if (sortedVariables[i] instanceof Function) {
                Function paramFunction = (Function) sortedVariables[i];
                Element xmlParam = null;
                for (int j = 0; j < list.size(); j++) {
                    Element tempParam = (Element) list.get(j);
                    if (paramFunction.getName().equals(unMangle(tempParam.getAttributeValue(XMLTags.NameAttrTag)))) {
                        xmlParam = tempParam;
                        break;
                    }
                }
                if (xmlParam == null) {
                    // must have been an unresolved parameter
                    continue;
                }
                String symbol = xmlParam.getAttributeValue(XMLTags.VCUnitDefinitionAttrTag);
                VCUnitDefinition unit = null;
                if (symbol != null) {
                    unit = modelUnitSystem.getInstance(symbol);
                }
                Kinetics.KineticsParameter tempParam = newKinetics.getKineticsParameter(paramFunction.getName());
                if (tempParam == null) {
                    newKinetics.addUserDefinedKineticsParameter(paramFunction.getName(), paramFunction.getExpression(), unit);
                } else {
                    newKinetics.setParameterValue(tempParam, paramFunction.getExpression());
                    tempParam.setUnitDefinition(unit);
                }
            }
        }
    } catch (java.beans.PropertyVetoException e) {
        e.printStackTrace(System.out);
        throw new XmlParseException("Exception while setting parameters for Reaction : " + reaction.getName(), e);
    } catch (ExpressionException e) {
        e.printStackTrace(System.out);
        throw new XmlParseException("Exception while settings parameters for Reaction : " + reaction.getName(), e);
    } finally {
        newKinetics.reading(false);
    }
    return newKinetics;
}
Also used : FilamentVariable(cbit.vcell.math.FilamentVariable) OutsideVariable(cbit.vcell.math.OutsideVariable) StochVolVariable(cbit.vcell.math.StochVolVariable) RandomVariable(cbit.vcell.math.RandomVariable) VolumeRandomVariable(cbit.vcell.math.VolumeRandomVariable) VolumeParticleVariable(cbit.vcell.math.VolumeParticleVariable) VolumeRegionVariable(cbit.vcell.math.VolumeRegionVariable) InsideVariable(cbit.vcell.math.InsideVariable) VolVariable(cbit.vcell.math.VolVariable) MembraneRegionVariable(cbit.vcell.math.MembraneRegionVariable) PointVariable(cbit.vcell.math.PointVariable) MembraneRandomVariable(cbit.vcell.math.MembraneRandomVariable) MembraneParticleVariable(cbit.vcell.math.MembraneParticleVariable) ParticleVariable(cbit.vcell.math.ParticleVariable) MemVariable(cbit.vcell.math.MemVariable) FilamentRegionVariable(cbit.vcell.math.FilamentRegionVariable) Variable(cbit.vcell.math.Variable) KineticsParameter(cbit.vcell.model.Kinetics.KineticsParameter) VariableHash(cbit.vcell.math.VariableHash) HMM_IRRKinetics(cbit.vcell.model.HMM_IRRKinetics) MacroscopicRateConstant(cbit.vcell.math.MacroscopicRateConstant) Constant(cbit.vcell.math.Constant) Element(org.jdom.Element) ReservedSymbol(cbit.vcell.model.Model.ReservedSymbol) NernstKinetics(cbit.vcell.model.NernstKinetics) ArrayList(java.util.ArrayList) FluxReaction(cbit.vcell.model.FluxReaction) GeneralKinetics(cbit.vcell.model.GeneralKinetics) GeneralLumpedKinetics(cbit.vcell.model.GeneralLumpedKinetics) ExpressionException(cbit.vcell.parser.ExpressionException) PropertyVetoException(java.beans.PropertyVetoException) GeneralCurrentKinetics(cbit.vcell.model.GeneralCurrentKinetics) AnnotatedFunction(cbit.vcell.solver.AnnotatedFunction) Function(cbit.vcell.math.Function) GeneralCurrentLumpedKinetics(cbit.vcell.model.GeneralCurrentLumpedKinetics) Membrane(cbit.vcell.model.Membrane) Macroscopic_IRRKinetics(cbit.vcell.model.Macroscopic_IRRKinetics) GeneralPermeabilityKinetics(cbit.vcell.model.GeneralPermeabilityKinetics) ModelUnitSystem(cbit.vcell.model.ModelUnitSystem) SimpleReaction(cbit.vcell.model.SimpleReaction) GHKKinetics(cbit.vcell.model.GHKKinetics) HMM_REVKinetics(cbit.vcell.model.HMM_REVKinetics) VCUnitDefinition(cbit.vcell.units.VCUnitDefinition) Expression(cbit.vcell.parser.Expression) MathException(cbit.vcell.math.MathException) MassActionKinetics(cbit.vcell.model.MassActionKinetics) Macroscopic_IRRKinetics(cbit.vcell.model.Macroscopic_IRRKinetics) Kinetics(cbit.vcell.model.Kinetics) NernstKinetics(cbit.vcell.model.NernstKinetics) GeneralKinetics(cbit.vcell.model.GeneralKinetics) GeneralLumpedKinetics(cbit.vcell.model.GeneralLumpedKinetics) GeneralPermeabilityKinetics(cbit.vcell.model.GeneralPermeabilityKinetics) GHKKinetics(cbit.vcell.model.GHKKinetics) MassActionKinetics(cbit.vcell.model.MassActionKinetics) Microscopic_IRRKinetics(cbit.vcell.model.Microscopic_IRRKinetics) GeneralCurrentKinetics(cbit.vcell.model.GeneralCurrentKinetics) HMM_REVKinetics(cbit.vcell.model.HMM_REVKinetics) HMM_IRRKinetics(cbit.vcell.model.HMM_IRRKinetics) GeneralCurrentLumpedKinetics(cbit.vcell.model.GeneralCurrentLumpedKinetics) Microscopic_IRRKinetics(cbit.vcell.model.Microscopic_IRRKinetics) ReactionParticipant(cbit.vcell.model.ReactionParticipant)

Example 39 with Function

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

the class CopasiOptimizationSolver method getOdeSolverResultSet.

private static ODESolverResultSet getOdeSolverResultSet(RowColumnResultSet rcResultSet, SimulationSymbolTable simSymbolTable, String[] parameterNames, double[] parameterValues) {
    // 
    // get simulation results - copy from RowColumnResultSet into OdeSolverResultSet
    // 
    ODESolverResultSet odeSolverResultSet = new ODESolverResultSet();
    for (int i = 0; i < rcResultSet.getDataColumnCount(); i++) {
        odeSolverResultSet.addDataColumn(new ODESolverResultSetColumnDescription(rcResultSet.getColumnDescriptions(i).getName()));
    }
    for (int i = 0; i < rcResultSet.getRowCount(); i++) {
        odeSolverResultSet.addRow(rcResultSet.getRow(i));
    }
    // 
    // add appropriate Function columns to result set
    // 
    Function[] functions = simSymbolTable.getFunctions();
    for (int i = 0; i < functions.length; i++) {
        if (SimulationSymbolTable.isFunctionSaved(functions[i])) {
            Expression exp1 = new Expression(functions[i].getExpression());
            try {
                exp1 = simSymbolTable.substituteFunctions(exp1).flatten();
                // 
                for (int j = 0; parameterNames != null && j < parameterNames.length; j++) {
                    exp1.substituteInPlace(new Expression(parameterNames[j]), new Expression(parameterValues[j]));
                }
            } catch (MathException e) {
                e.printStackTrace(System.out);
                throw new RuntimeException("Substitute function failed on function " + functions[i].getName() + " " + e.getMessage());
            } catch (ExpressionException e) {
                e.printStackTrace(System.out);
                throw new RuntimeException("Substitute function failed on function " + functions[i].getName() + " " + e.getMessage());
            }
            try {
                FunctionColumnDescription cd = new FunctionColumnDescription(exp1.flatten(), functions[i].getName(), null, functions[i].getName(), false);
                odeSolverResultSet.addFunctionColumn(cd);
            } catch (ExpressionException e) {
                e.printStackTrace(System.out);
            }
        }
    }
    return odeSolverResultSet;
}
Also used : Function(cbit.vcell.math.Function) Expression(cbit.vcell.parser.Expression) MathException(cbit.vcell.math.MathException) ODESolverResultSet(cbit.vcell.solver.ode.ODESolverResultSet) ODESolverResultSetColumnDescription(cbit.vcell.math.ODESolverResultSetColumnDescription) FunctionColumnDescription(cbit.vcell.math.FunctionColumnDescription) ExpressionException(cbit.vcell.parser.ExpressionException)

Example 40 with Function

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

the class OptXmlWriter method getObjectiveFunctionXML.

public static Element getObjectiveFunctionXML(OptimizationSpec optimizationSpec) {
    Element objectiveFunctionElement = new Element(OptXmlTags.ObjectiveFunction_Tag);
    ObjectiveFunction objectiveFunction = optimizationSpec.getObjectiveFunction();
    if (objectiveFunction instanceof ExplicitObjectiveFunction) {
        ExplicitObjectiveFunction explicitObjectiveFunction = (ExplicitObjectiveFunction) objectiveFunction;
        objectiveFunctionElement.setAttribute(OptXmlTags.ObjectiveFunctionType_Attr, OptXmlTags.ObjectiveFunctionType_Attr_Explicit);
        Element expressionElement = new Element(OptXmlTags.Expression_Tag);
        expressionElement.addContent(explicitObjectiveFunction.getExpression().infix());
        objectiveFunctionElement.addContent(expressionElement);
    } else if (objectiveFunction instanceof ExplicitFitObjectiveFunction) {
        ExplicitFitObjectiveFunction explicitFitObjectiveFunction = (ExplicitFitObjectiveFunction) objectiveFunction;
        objectiveFunctionElement.setAttribute(OptXmlTags.ObjectiveFunctionType_Attr, OptXmlTags.ObjectiveFunctionType_Attr_ExplicitFit);
        // add expression list
        Element expressionListElement = new Element(OptXmlTags.ExpressionList_Tag);
        objectiveFunctionElement.addContent(expressionListElement);
        ExplicitFitObjectiveFunction.ExpressionDataPair[] expDataPairs = explicitFitObjectiveFunction.getExpressionDataPairs();
        for (int i = 0; i < expDataPairs.length; i++) {
            Element expressionElement = new Element(OptXmlTags.Expression_Tag);
            expressionElement.setAttribute(OptXmlTags.ExpRefDataColumnIndex_Attr, expDataPairs[i].getReferenceDataIndex() + "");
            expressionElement.addContent(expDataPairs[i].getFitFunction().infix());
            expressionListElement.addContent(expressionElement);
        }
        if (explicitFitObjectiveFunction.getReferenceData() instanceof SimpleReferenceData) {
            Element dataElement = getDataXML(explicitFitObjectiveFunction.getReferenceData());
            objectiveFunctionElement.addContent(dataElement);
        } else {
            throw new OptimizationException("Optimization XML Writer exports SimpleReferenceData only.");
        }
    } else if (objectiveFunction instanceof PdeObjectiveFunction) {
        PdeObjectiveFunction pdeObjectiveFunction = (PdeObjectiveFunction) objectiveFunction;
        objectiveFunctionElement.setAttribute(OptXmlTags.ObjectiveFunctionType_Attr, OptXmlTags.ObjectiveFunctionType_Attr_PDE);
        SpatialReferenceData refData = pdeObjectiveFunction.getReferenceData();
        Element dataElement = getDataXML(refData);
        objectiveFunctionElement.addContent(dataElement);
        Element modelElement = getModelXML(pdeObjectiveFunction, optimizationSpec.getParameterNames());
        objectiveFunctionElement.addContent(modelElement);
        Simulation tempSim = new Simulation(pdeObjectiveFunction.getMathDescription());
        SimulationSymbolTable simSymbolTable = new SimulationSymbolTable(tempSim, 0);
        // 
        for (int i = 1; i < refData.getNumVariables(); i++) {
            Element modelMappingElement = new Element(OptXmlTags.ModelMapping_Tag);
            modelMappingElement.setAttribute(OptXmlTags.ModelMappingDataColumn_Attr, refData.getVariableNames()[i]);
            modelMappingElement.setAttribute(OptXmlTags.ModelMappingWeight_Attr, String.valueOf(refData.getVariableWeights()[i]));
            Expression mappingExpression = null;
            try {
                Variable var = pdeObjectiveFunction.getMathDescription().getVariable(refData.getVariableNames()[i]);
                if (var instanceof Function) {
                    Expression exp = new Expression(var.getExpression());
                    exp.bindExpression(simSymbolTable);
                    exp = MathUtilities.substituteFunctions(exp, simSymbolTable);
                    mappingExpression = exp.flatten();
                } else {
                    mappingExpression = new Expression(var.getName());
                }
            } catch (ExpressionException e) {
                e.printStackTrace();
                throw new OptimizationException(e.getMessage());
            }
            modelMappingElement.addContent(mappingExpression.infix());
            objectiveFunctionElement.addContent(modelMappingElement);
        }
    }
    return objectiveFunctionElement;
}
Also used : OptimizationException(cbit.vcell.opt.OptimizationException) ReservedVariable(cbit.vcell.math.ReservedVariable) Variable(cbit.vcell.math.Variable) Element(org.jdom.Element) PdeObjectiveFunction(cbit.vcell.opt.PdeObjectiveFunction) ExplicitFitObjectiveFunction(cbit.vcell.opt.ExplicitFitObjectiveFunction) SimulationSymbolTable(cbit.vcell.solver.SimulationSymbolTable) ExplicitFitObjectiveFunction(cbit.vcell.opt.ExplicitFitObjectiveFunction) ExplicitObjectiveFunction(cbit.vcell.opt.ExplicitObjectiveFunction) ObjectiveFunction(cbit.vcell.opt.ObjectiveFunction) PdeObjectiveFunction(cbit.vcell.opt.PdeObjectiveFunction) SimpleReferenceData(cbit.vcell.opt.SimpleReferenceData) Constraint(cbit.vcell.opt.Constraint) ExpressionException(cbit.vcell.parser.ExpressionException) ExplicitFitObjectiveFunction(cbit.vcell.opt.ExplicitFitObjectiveFunction) Function(cbit.vcell.math.Function) ExplicitObjectiveFunction(cbit.vcell.opt.ExplicitObjectiveFunction) ObjectiveFunction(cbit.vcell.opt.ObjectiveFunction) PdeObjectiveFunction(cbit.vcell.opt.PdeObjectiveFunction) Simulation(cbit.vcell.solver.Simulation) Expression(cbit.vcell.parser.Expression) ExplicitObjectiveFunction(cbit.vcell.opt.ExplicitObjectiveFunction) SpatialReferenceData(cbit.vcell.opt.SpatialReferenceData)

Aggregations

Function (cbit.vcell.math.Function)52 Expression (cbit.vcell.parser.Expression)42 Variable (cbit.vcell.math.Variable)32 ExpressionException (cbit.vcell.parser.ExpressionException)26 Constant (cbit.vcell.math.Constant)24 VolVariable (cbit.vcell.math.VolVariable)22 CompartmentSubDomain (cbit.vcell.math.CompartmentSubDomain)18 MathException (cbit.vcell.math.MathException)16 MathDescription (cbit.vcell.math.MathDescription)15 InsideVariable (cbit.vcell.math.InsideVariable)13 MemVariable (cbit.vcell.math.MemVariable)13 OutsideVariable (cbit.vcell.math.OutsideVariable)13 ReservedVariable (cbit.vcell.math.ReservedVariable)13 Domain (cbit.vcell.math.Variable.Domain)13 MembraneRegionVariable (cbit.vcell.math.MembraneRegionVariable)12 SubDomain (cbit.vcell.math.SubDomain)12 VolumeRegionVariable (cbit.vcell.math.VolumeRegionVariable)12 Vector (java.util.Vector)12 Element (org.jdom.Element)11 SubVolume (cbit.vcell.geometry.SubVolume)9