Examples with StructureMapping cbit.vcell.mapping.StructureMapping used on opensource projects

Example 16 with StructureMapping

use of cbit.vcell.mapping.StructureMapping in project vcell by virtualcell.

the class MathMapping_4_8 method getMathSymbol0.

/**
 * Substitutes appropriate variables for speciesContext bindings
 *
 * @return cbit.vcell.parser.Expression
 * @param origExp cbit.vcell.parser.Expression
 * @param structureMapping cbit.vcell.mapping.StructureMapping
 */
private String getMathSymbol0(SymbolTableEntry ste, StructureMapping structureMapping) throws MappingException {
    String steName = ste.getName();
    if (ste instanceof Kinetics.KineticsParameter) {
        Integer count = localNameCountHash.get(steName);
        if (count == null) {
            throw new MappingException("KineticsParameter " + steName + " not found in local name count");
        }
        if (count > 1 || steName.equals("J")) {
            return steName + "_" + ste.getNameScope().getName();
        // return getNameScope().getSymbolName(ste);
        } else {
            return steName;
        }
    }
    if (ste instanceof MathMapping_4_8.ProbabilityParameter) {
        // be careful here, to see if we need mangle the reaction name
        MathMapping_4_8.ProbabilityParameter probParm = (MathMapping_4_8.ProbabilityParameter) ste;
        return probParm.getName();
    }
    if (ste instanceof MathMapping_4_8.SpeciesConcentrationParameter) {
        MathMapping_4_8.SpeciesConcentrationParameter concParm = (MathMapping_4_8.SpeciesConcentrationParameter) ste;
        return concParm.getSpeciesContextSpec().getSpeciesContext().getName() + MATH_FUNC_SUFFIX_SPECIES_CONCENTRATION;
    }
    if (ste instanceof MathMapping_4_8.SpeciesCountParameter) {
        MathMapping_4_8.SpeciesCountParameter countParm = (MathMapping_4_8.SpeciesCountParameter) ste;
        return countParm.getSpeciesContextSpec().getSpeciesContext().getName() + MATH_VAR_SUFFIX_SPECIES_COUNT;
    }
    if (ste instanceof MathMapping_4_8.EventAssignmentInitParameter) {
        MathMapping_4_8.EventAssignmentInitParameter eventInitParm = (MathMapping_4_8.EventAssignmentInitParameter) ste;
        return eventInitParm.getName() + MATH_FUNC_SUFFIX_EVENTASSIGN_INIT;
    }
    if (ste instanceof Model.ReservedSymbol) {
        return steName;
    }
    if (ste instanceof Membrane.MembraneVoltage) {
        return steName;
    }
    if (ste instanceof Structure.StructureSize) {
        Structure structure = ((Structure.StructureSize) ste).getStructure();
        StructureMapping.StructureMappingParameter sizeParameter = simContext.getGeometryContext().getStructureMapping(structure).getSizeParameter();
        return getMathSymbol(sizeParameter, structureMapping);
    }
    if (ste instanceof ProxyParameter) {
        ProxyParameter pp = (ProxyParameter) ste;
        return getMathSymbol0(pp.getTarget(), structureMapping);
    }
    // 
    Model model = simContext.getModel();
    if (ste instanceof ModelParameter) {
        ModelParameter mp = (ModelParameter) ste;
        if (simContext.getGeometry().getDimension() == 0) {
            return mp.getName();
        } else {
            if (mp.getExpression().getSymbols() == null) {
                return mp.getName();
            }
            // check if global param variant name exists in globalVarsHash. If so, return it, else, throw exception.
            Hashtable<String, Expression> smVariantsHash = globalParamVariantsHash.get(mp);
            String variantName = mp.getName() + "_" + TokenMangler.fixTokenStrict(structureMapping.getStructure().getName());
            if (smVariantsHash.get(variantName) != null) {
                return variantName;
            } else {
                // global param variant doesn't exist in the hash, so get the substituted expression for global param and
                // gather all symbols (speciesContexts) that do not match with arg 'structureMapping' to display a proper error message.
                Expression expr = null;
                try {
                    expr = substituteGlobalParameters(mp.getExpression());
                } catch (ExpressionException e) {
                    e.printStackTrace(System.out);
                    throw new RuntimeException("Could not substitute expression for global parameter '" + mp.getName() + "' with expression '" + "'" + e.getMessage());
                }
                // find symbols (typically speciesContexts) in 'exp' that do not match with the arg 'structureMapping'
                String[] symbols = expr.getSymbols();
                String msg = "";
                if (symbols != null) {
                    Vector<String> spContextNamesVector = new Vector<String>();
                    for (int j = 0; j < symbols.length; j++) {
                        SpeciesContext sc = model.getSpeciesContext(symbols[j]);
                        if (sc != null) {
                            if (!sc.getStructure().compareEqual(structureMapping.getStructure())) {
                                spContextNamesVector.addElement(sc.getName());
                            }
                        }
                    }
                    for (int i = 0; (spContextNamesVector != null && i < spContextNamesVector.size()); i++) {
                        if (i == 0) {
                            msg += "'" + spContextNamesVector.elementAt(i) + ", ";
                        } else if (i == spContextNamesVector.size() - 1) {
                            msg += spContextNamesVector.elementAt(i) + "'";
                        } else {
                            msg += spContextNamesVector.elementAt(i) + ", ";
                        }
                    }
                }
                throw new RuntimeException("Global parameter '" + mp.getName() + "' is not defined in compartment '" + structureMapping.getStructure().getName() + "', but was referenced in that compartment." + "\n\nExpression '" + mp.getExpression().infix() + "' for global parameter '" + mp.getName() + "' expands to '" + expr.infix() + "' " + "and contains species " + msg + " that is/are not in adjacent compartments.");
            }
        // return (mp.getName()+"_"+TokenMangler.fixTokenStrict(structureMapping.getStructure().getName()));
        }
    }
    if (ste instanceof SpeciesContextSpec.SpeciesContextSpecParameter) {
        SpeciesContextSpec.SpeciesContextSpecParameter scsParm = (SpeciesContextSpec.SpeciesContextSpecParameter) ste;
        if (scsParm.getRole() == SpeciesContextSpec.ROLE_InitialConcentration) {
            return ((SpeciesContextSpec) (scsParm.getNameScope().getScopedSymbolTable())).getSpeciesContext().getName() + MATH_FUNC_SUFFIX_SPECIES_INIT_CONCENTRATION;
        }
        if (scsParm.getRole() == SpeciesContextSpec.ROLE_InitialCount) {
            return ((SpeciesContextSpec) (scsParm.getNameScope().getScopedSymbolTable())).getSpeciesContext().getName() + MATH_FUNC_SUFFIX_SPECIES_INIT_COUNT;
        }
        if (scsParm.getRole() == SpeciesContextSpec.ROLE_DiffusionRate) {
            return ((SpeciesContextSpec) (scsParm.getNameScope().getScopedSymbolTable())).getSpeciesContext().getName() + "_diffusionRate";
        }
        if (scsParm.getRole() == SpeciesContextSpec.ROLE_BoundaryValueXm) {
            return ((SpeciesContextSpec) (scsParm.getNameScope().getScopedSymbolTable())).getSpeciesContext().getName() + "_boundaryXm";
        }
        if (scsParm.getRole() == SpeciesContextSpec.ROLE_BoundaryValueXp) {
            return ((SpeciesContextSpec) (scsParm.getNameScope().getScopedSymbolTable())).getSpeciesContext().getName() + "_boundaryXp";
        }
        if (scsParm.getRole() == SpeciesContextSpec.ROLE_BoundaryValueYm) {
            return ((SpeciesContextSpec) (scsParm.getNameScope().getScopedSymbolTable())).getSpeciesContext().getName() + "_boundaryYm";
        }
        if (scsParm.getRole() == SpeciesContextSpec.ROLE_BoundaryValueYp) {
            return ((SpeciesContextSpec) (scsParm.getNameScope().getScopedSymbolTable())).getSpeciesContext().getName() + "_boundaryYp";
        }
        if (scsParm.getRole() == SpeciesContextSpec.ROLE_BoundaryValueZm) {
            return ((SpeciesContextSpec) (scsParm.getNameScope().getScopedSymbolTable())).getSpeciesContext().getName() + "_boundaryZm";
        }
        if (scsParm.getRole() == SpeciesContextSpec.ROLE_BoundaryValueZp) {
            return ((SpeciesContextSpec) (scsParm.getNameScope().getScopedSymbolTable())).getSpeciesContext().getName() + "_boundaryZp";
        }
        if (scsParm.getRole() == SpeciesContextSpec.ROLE_VelocityX) {
            return ((SpeciesContextSpec) (scsParm.getNameScope().getScopedSymbolTable())).getSpeciesContext().getName() + "_velocityX";
        }
        if (scsParm.getRole() == SpeciesContextSpec.ROLE_VelocityY) {
            return ((SpeciesContextSpec) (scsParm.getNameScope().getScopedSymbolTable())).getSpeciesContext().getName() + "_velocityY";
        }
        if (scsParm.getRole() == SpeciesContextSpec.ROLE_VelocityZ) {
            return ((SpeciesContextSpec) (scsParm.getNameScope().getScopedSymbolTable())).getSpeciesContext().getName() + "_velocityZ";
        }
    }
    if (ste instanceof ElectricalDevice.ElectricalDeviceParameter) {
        ElectricalDevice.ElectricalDeviceParameter edParm = (ElectricalDevice.ElectricalDeviceParameter) ste;
        ElectricalDevice electricalDevice = (ElectricalDevice) edParm.getNameScope().getScopedSymbolTable();
        if (electricalDevice instanceof MembraneElectricalDevice) {
            String nameWithScope = ((MembraneElectricalDevice) electricalDevice).getMembraneMapping().getMembrane().getNameScope().getName();
            if (edParm.getRole() == ElectricalDevice.ROLE_TotalCurrent) {
                return "I_" + nameWithScope;
            }
            if (edParm.getRole() == ElectricalDevice.ROLE_TransmembraneCurrent) {
                return "F_" + nameWithScope;
            }
        // }else if (electricalDevice instanceof CurrentClampElectricalDevice) {
        // if (edParm.getRole()==ElectricalDevice.ROLE_TotalCurrentDensity){
        // return "I_"+((CurrentClampElectricalDevice)electricalDevice).getCurrentClampStimulus().getNameScope().getName();
        // }
        // if (edParm.getRole()==ElectricalDevice.ROLE_TransmembraneCurrentDensity){
        // return "F_"+((CurrentClampElectricalDevice)electricalDevice).getCurrentClampStimulus().getNameScope().getName();
        // }
        // }else if (electricalDevice instanceof VoltageClampElectricalDevice) {
        // if (edParm.getRole()==ElectricalDevice.ROLE_TotalCurrentDensity){
        // return "I_"+((VoltageClampElectricalDevice)electricalDevice).getVoltageClampStimulus().getNameScope().getName();
        // }
        // if (edParm.getRole()==ElectricalDevice.ROLE_TransmembraneCurrentDensity){
        // return "F_"+((VoltageClampElectricalDevice)electricalDevice).getVoltageClampStimulus().getNameScope().getName();
        // }
        }
    }
    if (ste instanceof LocalParameter && ((LocalParameter) ste).getNameScope() instanceof ElectricalStimulus.ElectricalStimulusNameScope) {
        LocalParameter esParm = (LocalParameter) ste;
        String nameWithScope = esParm.getNameScope().getName();
        if (esParm.getRole() == ElectricalStimulusParameterType.TotalCurrent) {
            return "I_" + nameWithScope;
        } else if (esParm.getRole() == ElectricalStimulusParameterType.Voltage) {
            return "V_" + nameWithScope;
        }
    }
    StructureTopology structTopology = model.getStructureTopology();
    if (ste instanceof StructureMapping.StructureMappingParameter) {
        StructureMapping.StructureMappingParameter smParm = (StructureMapping.StructureMappingParameter) ste;
        Structure structure = ((StructureMapping) (smParm.getNameScope().getScopedSymbolTable())).getStructure();
        int role = smParm.getRole();
        if (role == StructureMapping.ROLE_VolumeFraction) {
            return "VolFract_" + (structTopology.getInsideFeature((Membrane) structure)).getNameScope().getName();
        } else {
            String nameWithScope = structure.getNameScope().getName();
            if (role == StructureMapping.ROLE_SurfaceToVolumeRatio) {
                return "SurfToVol_" + nameWithScope;
            } else if (role == StructureMapping.ROLE_InitialVoltage) {
                return smParm.getName();
            } else if (role == StructureMapping.ROLE_SpecificCapacitance) {
                return "C_" + nameWithScope;
            } else if (role == StructureMapping.ROLE_AreaPerUnitArea) {
                return "AreaPerUnitArea_" + nameWithScope;
            } else if (role == StructureMapping.ROLE_AreaPerUnitVolume) {
                return "AreaPerUnitVolume_" + nameWithScope;
            } else if (role == StructureMapping.ROLE_VolumePerUnitArea) {
                return "VolumePerUnitArea_" + nameWithScope;
            } else if (role == StructureMapping.ROLE_VolumePerUnitVolume) {
                return "VolumePerUnitVolume_" + nameWithScope;
            } else if (role == StructureMapping.ROLE_Size) {
                if (simContext.getGeometry().getDimension() == 0) {
                    // if geometry is compartmental, make sure compartment sizes are set if referenced in model.
                    if (smParm.getExpression() == null || smParm.getExpression().isZero()) {
                        throw new MappingException("\nIn non-spatial application '" + getSimulationContext().getName() + "', " + "size of structure '" + structure.getName() + "' must be assigned a " + "positive value if referenced in the model.\n\nPlease go to 'Structure Mapping' tab to check the size.");
                    }
                }
                return "Size_" + nameWithScope;
            }
        }
    }
    // 
    if (ste instanceof SpeciesContext) {
        SpeciesContext sc = (SpeciesContext) ste;
        SpeciesContextMapping scm = getSpeciesContextMapping(sc);
        // 
        if (structureMapping instanceof FeatureMapping) {
            // 
            if (scm.getVariable() != null && !scm.getVariable().getName().equals(steName)) {
                return scm.getVariable().getName();
            }
        // 
        // for reactions within a spatially resolved membrane, may need "_INSIDE" or "_OUTSIDE" for jump condition
        // 
        // if the membrane is distributed, then always use the plain variable.
        // 
        } else if (structureMapping instanceof MembraneMapping) {
            Membrane membrane = ((MembraneMapping) structureMapping).getMembrane();
            // 
            if (sc.getStructure() instanceof Membrane || getResolved(structureMapping) == false) {
                if (scm.getVariable() != null && !(scm.getVariable().getName().equals(steName))) {
                    return scm.getVariable().getName();
                }
            // 
            // if the speciesContext is outside the membrane
            // 
            } else {
                SpeciesContextSpec scs = simContext.getReactionContext().getSpeciesContextSpec(sc);
                if (sc.getStructure() == structTopology.getInsideFeature(membrane) || sc.getStructure() == structTopology.getOutsideFeature(membrane)) {
                    if (getResolved(structureMapping) && !scs.isConstant()) {
                        if (!scs.isDiffusing()) {
                            throw new MappingException("Enable diffusion in Application '" + simContext.getName() + "'. This must be done for any species (e.g '" + sc.getName() + "') in flux reactions.\n\n" + "To save or run simulations, set the diffusion rate to a non-zero " + "value in Initial Conditions or disable those reactions in Specifications->Reactions.");
                        }
                        return scm.getVariable().getName() + (sc.getStructure() == structTopology.getInsideFeature(membrane) ? "_INSIDE" : "_OUTSIDE");
                    } else {
                        return scm.getSpeciesContext().getName();
                    }
                } else {
                    throw new MappingException(sc.getName() + " shouldn't be involved with structure " + structureMapping.getStructure().getName());
                }
            }
        }
    }
    return getNameScope().getSymbolName(ste);
}

Example 17 with StructureMapping

use of cbit.vcell.mapping.StructureMapping in project vcell by virtualcell.

the class MathMapping_4_8 method refreshKFluxParameters.

/**
 * Insert the method's description here.
 * Creation date: (6/13/2004 8:10:34 AM)
 */
private void refreshKFluxParameters() throws ExpressionException {
    MathMapping_4_8.MathMappingParameter[] newMathMappingParameters = (MathMappingParameter[]) fieldMathMappingParameters.clone();
    // 
    for (int i = 0; i < newMathMappingParameters.length; i++) {
        if (newMathMappingParameters[i].getRole() == PARAMETER_ROLE_KFLUX) {
            newMathMappingParameters = (MathMappingParameter[]) BeanUtils.removeElement(newMathMappingParameters, newMathMappingParameters[i]);
        }
    }
    // 
    // Add new KFlux Parameters
    // 
    StructureMapping[] structureMappings = simContext.getGeometryContext().getStructureMappings();
    StructureTopology structTopology = simContext.getModel().getStructureTopology();
    VCUnitDefinition lengthInverseUnit = simContext.getModel().getX().getUnitDefinition().getInverse();
    for (int i = 0; i < structureMappings.length; i++) {
        if (structureMappings[i] instanceof MembraneMapping && !getResolved(structureMappings[i])) {
            MembraneMapping membraneMapping = (MembraneMapping) structureMappings[i];
            // 
            // add "inside" flux correction
            // 
            Expression insideCorrectionExp = getInsideFluxCorrectionExpression(membraneMapping);
            insideCorrectionExp.bindExpression(this);
            Membrane memMappingMembrane = membraneMapping.getMembrane();
            Feature insideFeature = structTopology.getInsideFeature(memMappingMembrane);
            String membraneNameWithScope = membraneMapping.getNameScope().getName();
            String insideName = "KFlux_" + membraneNameWithScope + "_" + insideFeature.getNameScope().getName();
            KFluxParameter insideKFluxParameter = new KFluxParameter(insideName, insideCorrectionExp, lengthInverseUnit, membraneMapping, insideFeature);
            newMathMappingParameters = (MathMappingParameter[]) BeanUtils.addElement(newMathMappingParameters, insideKFluxParameter);
            // 
            // add "outside" flux correction
            // 
            Expression outsideCorrectionExp = getOutsideFluxCorrectionExpression(simContext, membraneMapping);
            outsideCorrectionExp.bindExpression(this);
            Feature outsideFeature = structTopology.getOutsideFeature(memMappingMembrane);
            String outsideName = "KFlux_" + membraneNameWithScope + "_" + outsideFeature.getNameScope().getName();
            KFluxParameter outsideKFluxParameter = new KFluxParameter(outsideName, outsideCorrectionExp, lengthInverseUnit, membraneMapping, outsideFeature);
            newMathMappingParameters = (MathMappingParameter[]) BeanUtils.addElement(newMathMappingParameters, outsideKFluxParameter);
        }
    }
    try {
        setMathMapppingParameters(newMathMappingParameters);
    } catch (java.beans.PropertyVetoException e) {
        e.printStackTrace(System.out);
        throw new RuntimeException(e.getMessage());
    }
}

Example 18 with StructureMapping

use of cbit.vcell.mapping.StructureMapping in project vcell by virtualcell.

the class StructureSizeSolver method updateAbsoluteStructureSizes.

public static void updateAbsoluteStructureSizes(SimulationContext simContext, Structure struct, double structSize, VCUnitDefinition structSizeUnit) throws Exception {
    StructureMapping[] structMappings = simContext.getGeometryContext().getStructureMappings();
    try {
        StructureTopology structTopology = simContext.getModel().getStructureTopology();
        SolverParameterCollection solverParameters = new SolverParameterCollection();
        ArrayList<String> unknownVars = new ArrayList<String>();
        for (StructureMapping sm : structMappings) {
            System.out.println("Structure " + sm.getStructure().getDisplayName() + " size " + sm.getStructure().getStructureSize().getExpression());
            if (sm.getStructure() instanceof Membrane) {
                MembraneMapping mm = (MembraneMapping) sm;
                StructureMappingParameter svRatioParam = mm.getSurfaceToVolumeParameter();
                StructureMappingParameter volFractParam = mm.getVolumeFractionParameter();
                StructureMappingParameter sizeParam = mm.getSizeParameter();
                solverParameters.add(new SolverParameter(svRatioParam, TokenMangler.mangleToSName("sv_" + mm.getMembrane().getName()), svRatioParam.getExpression().evaluateConstant()));
                solverParameters.add(new SolverParameter(volFractParam, TokenMangler.mangleToSName("vf_" + mm.getMembrane().getName()), volFractParam.getExpression().evaluateConstant()));
            }
            StructureMappingParameter sizeParam = sm.getSizeParameter();
            Double priorKnownValue = null;
            String varName = TokenMangler.mangleToSName("size_" + sm.getStructure().getName());
            if (sizeParam.getExpression() != null) {
                priorKnownValue = sizeParam.getExpression().evaluateConstant();
            } else if (sm.getStructure() == struct) {
                priorKnownValue = structSize;
            } else {
                unknownVars.add(varName);
            }
            solverParameters.add(new SolverParameter(sizeParam, varName, priorKnownValue));
        }
        ArrayList<Expression> expressions = new ArrayList<Expression>();
        for (int i = 0; i < structMappings.length; i++) {
            if (structMappings[i] instanceof MembraneMapping) {
                MembraneMapping membraneMapping = (MembraneMapping) structMappings[i];
                Feature insideFeature = structTopology.getInsideFeature(membraneMapping.getMembrane());
                Feature outsideFeature = structTopology.getOutsideFeature(membraneMapping.getMembrane());
                StructureMappingParameter sizeParameter = membraneMapping.getSizeParameter();
                StructureMappingParameter volFractParameter = membraneMapping.getVolumeFractionParameter();
                StructureMappingParameter surfToVolParameter = membraneMapping.getSurfaceToVolumeParameter();
                // 
                // EC eclosing cyt, which contains er and golgi
                // "(cyt_size+ er_size + golgi_size) * cyt_svRatio - PM_size"  ... implicit equation, exp == 0 is implied
                // 
                Expression sumOfInsideVolumeExp = new Expression(0.0);
                for (int j = 0; j < structMappings.length; j++) {
                    if (structMappings[j] instanceof FeatureMapping && structTopology.enclosedBy(structMappings[j].getStructure(), insideFeature)) {
                        FeatureMapping childFeatureMappingOfInside = ((FeatureMapping) structMappings[j]);
                        sumOfInsideVolumeExp = Expression.add(sumOfInsideVolumeExp, new Expression(solverParameters.getName(childFeatureMappingOfInside.getSizeParameter())));
                    }
                }
                Expression tempExpr = Expression.mult(sumOfInsideVolumeExp, new Expression(solverParameters.getName(surfToVolParameter)));
                tempExpr = Expression.add(tempExpr, new Expression("-" + solverParameters.getName(sizeParameter)));
                expressions.add(tempExpr);
                // 
                // EC eclosing cyt, which contains er and golgi
                // (EC_size + cyt_size + er_size + golgi_size) * cyt_vfRatio - (cyt_size + er_size + golgi_size)  ... implicit equation, exp == 0 is implied
                // 
                Expression sumOfParentVolumeExp = new Expression(0.0);
                for (int j = 0; j < structMappings.length; j++) {
                    if (structMappings[j] instanceof FeatureMapping && structTopology.enclosedBy(structMappings[j].getStructure(), outsideFeature)) {
                        FeatureMapping childFeatureMappingOfParent = ((FeatureMapping) structMappings[j]);
                        sumOfParentVolumeExp = Expression.add(sumOfParentVolumeExp, new Expression(TokenMangler.mangleToSName(solverParameters.getName(childFeatureMappingOfParent.getSizeParameter()))));
                    }
                }
                Expression exp = Expression.mult(sumOfParentVolumeExp, new Expression(solverParameters.getName(volFractParameter)));
                exp = Expression.add(exp, Expression.negate(sumOfInsideVolumeExp));
                expressions.add(exp);
            }
        }
        if (expressions.size() != unknownVars.size()) {
            throw new RuntimeException("number of unknowns is " + unknownVars.size() + ", number of equations is " + expressions.size());
        }
        if (unknownVars.size() == 0 && expressions.size() == 0) {
            StructureMappingParameter sizeParam = simContext.getGeometryContext().getStructureMapping(struct).getSizeParameter();
            sizeParam.setExpression(new Expression(structSize));
            return;
        }
        RationalExp[][] rowColData = new RationalExp[unknownVars.size()][unknownVars.size() + 1];
        for (int row = 0; row < unknownVars.size(); row++) {
            // 
            // verify that there is no "constant" term (without an unknown)
            // 
            // System.out.println("equation("+row+"): "+expressions.get(row).infix());
            Expression constantTerm = new Expression(expressions.get(row));
            for (String var : unknownVars) {
                constantTerm.substituteInPlace(new Expression(var), new Expression(0.0));
            }
            constantTerm = constantTerm.flatten();
            // 
            for (int col = 0; col < unknownVars.size(); col++) {
                Expression equation = new Expression(expressions.get(row));
                String colVariable = unknownVars.get(col);
                Expression deriv = equation.differentiate(colVariable).flatten();
                String[] symbols = deriv.getSymbols();
                if (symbols != null) {
                    for (String symbol : symbols) {
                        if (unknownVars.contains(symbol)) {
                            throw new RuntimeException("equation is not linear in the unknowns");
                        }
                    }
                }
                rowColData[row][col] = RationalExpUtils.getRationalExp(deriv);
            }
            rowColData[row][unknownVars.size()] = RationalExpUtils.getRationalExp(constantTerm).minus();
        }
        RationalExpMatrix rationalExpMatrix = new RationalExpMatrix(rowColData);
        // rationalExpMatrix.show();
        RationalExp[] solutions = rationalExpMatrix.solveLinearExpressions();
        double[] solutionValues = new double[solutions.length];
        for (int i = 0; i < unknownVars.size(); i++) {
            Expression vcSolution = new Expression(solutions[i].infixString());
            String[] symbols = vcSolution.getSymbols();
            if (symbols != null) {
                for (String symbol : symbols) {
                    SolverParameter p = solverParameters.get(symbol);
                    if (p.knownValue == null) {
                        throw new RuntimeException("solution for var " + unknownVars.get(i) + " is a function of unknown var " + p.name);
                    }
                    vcSolution.substituteInPlace(new Expression(symbol), new Expression(p.knownValue.doubleValue()));
                }
            }
            double value = vcSolution.flatten().evaluateConstant();
            // System.out.println(unknownVars.get(i)+" = "+value+" = "+solutions[i].infixString());
            solutionValues[i] = value;
        }
        for (int i = 0; i < unknownVars.size(); i++) {
            SolverParameter p = solverParameters.get(unknownVars.get(i));
            p.parameter.setExpression(new Expression(solutionValues[i]));
        }
        // 
        // set the one known value (if not set already by the gui).
        // 
        StructureMappingParameter sizeParam = simContext.getGeometryContext().getStructureMapping(struct).getSizeParameter();
        sizeParam.setExpression(new Expression(structSize));
        System.out.println("done");
    } catch (ExpressionException e) {
        e.printStackTrace(System.out);
        throw new Exception(e.getMessage());
    }
}

Example 19 with StructureMapping

use of cbit.vcell.mapping.StructureMapping in project vcell by virtualcell.

the class StructureSizeSolver method updateRelativeStructureSizes.

/**
 * Insert the method's description here.
 * Creation date: (5/17/2006 10:33:38 AM)
 * @return double[]
 * @param structName java.lang.String
 * @param structSize double
 */
public static void updateRelativeStructureSizes(SimulationContext simContext) throws Exception {
    if (simContext.getGeometry().getDimension() > 0) {
        throw new RuntimeException("not yet supported for spatial applications");
    }
    StructureMapping[] structureMappings = simContext.getGeometryContext().getStructureMappings();
    try {
        // This is rewritten in Feb 2008. Siblings and children are correctly taken into account when calculating the volume fractions.
        StructureTopology structTopology = simContext.getModel().getStructureTopology();
        for (int i = 0; i < structureMappings.length; i++) {
            if (structureMappings[i] instanceof MembraneMapping) {
                // calculate the sum of features' sizes inside this membrane, this is used for calculating both surface volume ratio and volume fraction.
                double sumOfSubFeatures = 0;
                Membrane membrane = ((MembraneMapping) structureMappings[i]).getMembrane();
                Enumeration<Feature> subFeatures = structTopology.getSubFeatures(structTopology.getInsideFeature(membrane));
                while (subFeatures.hasMoreElements()) {
                    Feature feature = subFeatures.nextElement();
                    sumOfSubFeatures = sumOfSubFeatures + simContext.getGeometryContext().getStructureMapping(feature).getSizeParameter().getExpression().evaluateConstant();
                }
                // calculate the sum of features's sizes inside the membrance's parent feature, this is used for calculating the volume fraction.
                double sumOfParentMemSubFeatures = 0;
                Feature parentFeature = structTopology.getOutsideFeature(membrane);
                if (parentFeature != null) {
                    Enumeration<Feature> parentSubFeatures = structTopology.getSubFeatures(parentFeature);
                    while (parentSubFeatures.hasMoreElements()) {
                        Feature feature = parentSubFeatures.nextElement();
                        sumOfParentMemSubFeatures = sumOfParentMemSubFeatures + simContext.getGeometryContext().getStructureMapping(feature).getSizeParameter().getExpression().evaluateConstant();
                    }
                    // set surface volume ratio
                    ((MembraneMapping) structureMappings[i]).getSurfaceToVolumeParameter().setExpression(new Expression(((MembraneMapping) structureMappings[i]).getSizeParameter().getExpression().evaluateConstant() / sumOfSubFeatures));
                    // set volume fraction
                    ((MembraneMapping) structureMappings[i]).getVolumeFractionParameter().setExpression(new Expression(sumOfSubFeatures / sumOfParentMemSubFeatures));
                }
            }
        }
    } catch (NullPointerException e) {
        e.printStackTrace(System.out);
        // DialogUtils.showErrorDialog("structure sizes must all be specified");
        throw new Exception("structure sizes must all be specified");
    } catch (ExpressionException e) {
        e.printStackTrace(System.out);
        throw new Exception(e.getMessage());
    }
}

Example 20 with StructureMapping

use of cbit.vcell.mapping.StructureMapping in project vcell by virtualcell.

the class StochMathMapping_4_8 method refreshSpeciesContextMappings.

/**
 * Insert the method's description here.
 * Creation date: (10/26/2006 11:47:26 AM)
 * @exception cbit.vcell.parser.ExpressionException The exception description.
 * @exception cbit.vcell.mapping.MappingException The exception description.
 * @exception cbit.vcell.math.MathException The exception description.
 */
private void refreshSpeciesContextMappings() throws cbit.vcell.parser.ExpressionException, MappingException, cbit.vcell.math.MathException {
    // 
    // create a SpeciesContextMapping for each speciesContextSpec.
    // 
    // set initialExpression from SpeciesContextSpec.
    // set diffusing5
    // set variable (only if "Constant" or "Function", else leave it as null)-----why commented?
    // 
    // 
    // have to put geometric paras into mathsymbolmapping, since species initial condition needs the volume size symbol.
    // and the parameters later on were added into contants or functions in refreshMathDescription()
    // 
    StructureMapping[] structureMappings = getSimulationContext().getGeometryContext().getStructureMappings();
    for (int i = 0; i < structureMappings.length; i++) {
        StructureMapping sm = structureMappings[i];
        StructureMapping.StructureMappingParameter parm = sm.getParameterFromRole(StructureMapping.ROLE_Size);
        getMathSymbol(parm, sm);
    }
    getSpeciesContextMappingList().removeAllElements();
    SpeciesContextSpec[] speciesContextSpecs = getSimulationContext().getReactionContext().getSpeciesContextSpecs();
    for (int i = 0; i < speciesContextSpecs.length; i++) {
        SpeciesContextSpec scs = speciesContextSpecs[i];
        SpeciesContextMapping scm = new SpeciesContextMapping(scs.getSpeciesContext());
        scm.setPDERequired(getSimulationContext().isPDERequired(scs.getSpeciesContext()));
        scm.setHasEventAssignment(getSimulationContext().hasEventAssignment(scs.getSpeciesContext()));
        scm.setHasHybridReaction(false);
        for (ReactionSpec reactionSpec : getSimulationContext().getReactionContext().getReactionSpecs()) {
            if (!reactionSpec.isExcluded() && reactionSpec.hasHybrid(getSimulationContext(), scs.getSpeciesContext())) {
                scm.setHasHybridReaction(true);
            }
        }
        // scm.setAdvecting(isAdvectionRequired(scs.getSpeciesContext()));
        if (scs.isConstant()) {
            SpeciesContextSpec.SpeciesContextSpecParameter initCountParm = scs.getInitialCountParameter();
            SpeciesContextSpec.SpeciesContextSpecParameter initConcParm = scs.getInitialConcentrationParameter();
            Expression initCondInCount = null;
            // initial condition is concentration
            if (initConcParm != null && initConcParm.getExpression() != null) {
                initCondInCount = getExpressionConcToAmt(new Expression(initConcParm, getNameScope()), speciesContextSpecs[i].getSpeciesContext());
            } else {
                initCondInCount = new Expression(initCountParm, getNameScope());
            }
            // initCondInCount.bindExpression(this);
            initCondInCount = getSubstitutedExpr(initCondInCount, true, true);
            scm.setDependencyExpression(initCondInCount);
        }
        // 
        // test if participant in fast reaction step, request elimination if possible
        // 
        scm.setFastParticipant(false);
        ReactionSpec[] reactionSpecs = getSimulationContext().getReactionContext().getReactionSpecs();
        for (int j = 0; j < reactionSpecs.length; j++) {
            ReactionSpec reactionSpec = reactionSpecs[j];
            if (reactionSpec.isExcluded()) {
                continue;
            }
            ReactionStep rs = reactionSpec.getReactionStep();
            if (rs instanceof SimpleReaction && rs.countNumReactionParticipants(scs.getSpeciesContext()) > 0) {
                if (reactionSpec.isFast()) {
                    scm.setFastParticipant(true);
                }
            }
        }
        getSpeciesContextMappingList().addElement(scm);
    }
}