use of cbit.vcell.parser.ExpressionException in project vcell by virtualcell.
the class RulebasedMathMapping method refreshMathDescription.
/**
* This method was created in VisualAge.
*/
@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 transform 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("rule-based particle math mapping not implemented for spatial geometry - dimension >= 1");
}
//
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 reactionRules that are not "excluded"
//
ArrayList<ReactionRule> rrList = new ArrayList<ReactionRule>();
for (ReactionRuleSpec reactionRuleSpec : simContext.getReactionContext().getReactionRuleSpecs()) {
if (!reactionRuleSpec.isExcluded()) {
rrList.add(reactionRuleSpec.getReactionRule());
}
}
//
for (ReactionRule reactionRule : rrList) {
UnresolvedParameter[] unresolvedParameters = reactionRule.getKineticLaw().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() + "', " + reactionRule.getDisplayType() + " '" + reactionRule.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 (ReactionRule reactionRule : rrList) {
// if (reactionRule.getKineticLaw() instanceof LumpedKinetics){
// throw new RuntimeException("Lumped Kinetics not yet supported for RuleBased Modeling");
// }
LocalParameter[] parameters = reactionRule.getKineticLaw().getLocalParameters();
for (LocalParameter parameter : parameters) {
//
if ((parameter.getRole() == RbmKineticLawParameterType.RuleRate)) {
continue;
}
//
if (!reactionRule.isReversible() && parameter.getRole() == RbmKineticLawParameterType.MassActionReverseRate) {
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);
//
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);
//
// define all molecules and unique species patterns (add molecules to mathDesc and speciesPatterns to varHash).
//
HashMap<SpeciesPattern, VolumeParticleSpeciesPattern> speciesPatternMap = addSpeciesPatterns(domain, rrList);
HashSet<VolumeParticleSpeciesPattern> uniqueParticleSpeciesPatterns = new HashSet<>(speciesPatternMap.values());
for (VolumeParticleSpeciesPattern volumeParticleSpeciesPattern : uniqueParticleSpeciesPatterns) {
varHash.addVariable(volumeParticleSpeciesPattern);
}
//
// define observables (those explicitly declared and those corresponding to seed species.
//
List<ParticleObservable> observables = addObservables(geometryClass, domain, speciesPatternMap);
for (ParticleObservable particleObservable : observables) {
varHash.addVariable(particleObservable);
}
try {
addParticleJumpProcesses(varHash, geometryClass, subDomain, speciesPatternMap);
} catch (PropertyVetoException e1) {
e1.printStackTrace();
throw new MappingException(e1.getMessage(), e1);
}
//
for (int i = 0; i < fieldMathMappingParameters.length; i++) {
if (fieldMathMappingParameters[i] instanceof UnitFactorParameter || fieldMathMappingParameters[i] instanceof ObservableConcentrationParameter) {
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());
//
for (SpeciesContext sc : model.getSpeciesContexts()) {
if (!sc.hasSpeciesPattern()) {
throw new MappingException("species " + sc.getName() + " has no molecular pattern");
}
VolumeParticleSpeciesPattern volumeParticleSpeciesPattern = speciesPatternMap.get(sc.getSpeciesPattern());
ArrayList<ParticleInitialCondition> particleInitialConditions = new ArrayList<ParticleProperties.ParticleInitialCondition>();
// initial conditions from scs
SpeciesContextSpec scs = simContext.getReactionContext().getSpeciesContextSpec(sc);
Parameter initialCountParameter = scs.getInitialCountParameter();
Expression e = getIdentifierSubstitutions(new Expression(initialCountParameter, getNameScope()), initialCountParameter.getUnitDefinition(), geometryClass);
particleInitialConditions.add(new ParticleInitialConditionCount(e, new Expression(0.0), new Expression(0.0), new Expression(0.0)));
ParticleProperties particleProperies = new ParticleProperties(volumeParticleSpeciesPattern, new Expression(0.0), new Expression(0.0), new Expression(0.0), new Expression(0.0), particleInitialConditions);
subDomain.addParticleProperties(particleProperies);
}
//
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 ObservableConcentrationParameter) {
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());
}
}
use of cbit.vcell.parser.ExpressionException in project vcell by virtualcell.
the class RulebasedTransformer method transform.
private void transform(SimulationContext originalSimContext, SimulationContext transformedSimulationContext, ArrayList<ModelEntityMapping> entityMappings, MathMappingCallback mathMappingCallback) throws PropertyVetoException {
Model newModel = transformedSimulationContext.getModel();
Model originalModel = originalSimContext.getModel();
ModelEntityMapping em = null;
// list of rules created from the reactions; we apply the symmetry factor computed by bionetgen only to these
Set<ReactionRule> fromReactions = new HashSet<>();
for (SpeciesContext newSpeciesContext : newModel.getSpeciesContexts()) {
final SpeciesContext originalSpeciesContext = originalModel.getSpeciesContext(newSpeciesContext.getName());
// map new and old species contexts
em = new ModelEntityMapping(originalSpeciesContext, newSpeciesContext);
entityMappings.add(em);
if (newSpeciesContext.hasSpeciesPattern()) {
// it's perfect already and can't be improved
continue;
}
try {
MolecularType newmt = newModel.getRbmModelContainer().createMolecularType();
newModel.getRbmModelContainer().addMolecularType(newmt, false);
MolecularTypePattern newmtp_sc = new MolecularTypePattern(newmt);
SpeciesPattern newsp_sc = new SpeciesPattern();
newsp_sc.addMolecularTypePattern(newmtp_sc);
newSpeciesContext.setSpeciesPattern(newsp_sc);
RbmObservable newo = new RbmObservable(newModel, "O0_" + newmt.getName() + "_tot", newSpeciesContext.getStructure(), RbmObservable.ObservableType.Molecules);
MolecularTypePattern newmtp_ob = new MolecularTypePattern(newmt);
SpeciesPattern newsp_ob = new SpeciesPattern();
newsp_ob.addMolecularTypePattern(newmtp_ob);
newo.addSpeciesPattern(newsp_ob);
newModel.getRbmModelContainer().addObservable(newo);
// map new observable to old species context
em = new ModelEntityMapping(originalSpeciesContext, newo);
entityMappings.add(em);
} catch (ModelException e) {
e.printStackTrace();
throw new RuntimeException("unable to transform species context: " + e.getMessage());
} catch (PropertyVetoException e) {
// TODO Auto-generated catch block
e.printStackTrace();
}
}
ReactionSpec[] reactionSpecs = transformedSimulationContext.getReactionContext().getReactionSpecs();
for (ReactionSpec reactionSpec : reactionSpecs) {
if (reactionSpec.isExcluded()) {
// we create rules only from those reactions which are not excluded
continue;
}
ReactionStep rs = reactionSpec.getReactionStep();
String name = rs.getName();
String mangled = TokenMangler.fixTokenStrict(name);
mangled = newModel.getReactionName(mangled);
Kinetics k = rs.getKinetics();
if (!(k instanceof MassActionKinetics)) {
throw new RuntimeException("Only Mass Action Kinetics supported at this time, reaction \"" + rs.getName() + "\" uses kinetic law type \"" + rs.getKinetics().getName() + "\"");
}
boolean bReversible = rs.isReversible();
ReactionRule rr = new ReactionRule(newModel, mangled, rs.getStructure(), bReversible);
fromReactions.add(rr);
MassActionKinetics massActionKinetics = (MassActionKinetics) k;
List<Reactant> rList = rs.getReactants();
List<Product> pList = rs.getProducts();
// counting the stoichiometry - 2A+B means 3 reactants
int numReactants = 0;
for (Reactant r : rList) {
numReactants += r.getStoichiometry();
if (numReactants > 2) {
String message = "NFSim doesn't support more than 2 reactants within a reaction: " + name;
throw new RuntimeException(message);
}
}
int numProducts = 0;
for (Product p : pList) {
numProducts += p.getStoichiometry();
if (bReversible && numProducts > 2) {
String message = "NFSim doesn't support more than 2 products within a reversible reaction: " + name;
throw new RuntimeException(message);
}
}
RateLawType rateLawType = RateLawType.MassAction;
RbmKineticLaw kineticLaw = new RbmKineticLaw(rr, rateLawType);
try {
String forwardRateName = massActionKinetics.getForwardRateParameter().getName();
Expression forwardRateExp = massActionKinetics.getForwardRateParameter().getExpression();
String reverseRateName = massActionKinetics.getReverseRateParameter().getName();
Expression reverseRateExp = massActionKinetics.getReverseRateParameter().getExpression();
LocalParameter fR = kineticLaw.getLocalParameter(RbmKineticLawParameterType.MassActionForwardRate);
fR.setName(forwardRateName);
LocalParameter rR = kineticLaw.getLocalParameter(RbmKineticLawParameterType.MassActionReverseRate);
rR.setName(reverseRateName);
if (rs.hasReactant()) {
kineticLaw.setParameterValue(fR, forwardRateExp, true);
}
if (rs.hasProduct()) {
kineticLaw.setParameterValue(rR, reverseRateExp, true);
}
//
for (KineticsParameter reaction_p : massActionKinetics.getKineticsParameters()) {
if (reaction_p.getRole() == Kinetics.ROLE_UserDefined) {
LocalParameter rule_p = kineticLaw.getLocalParameter(reaction_p.getName());
if (rule_p == null) {
//
// after lazy parameter creation we didn't find a user-defined rule parameter with this same name.
//
// there must be a global symbol with the same name, that the local reaction parameter has overridden.
//
ParameterContext.LocalProxyParameter rule_proxy_parameter = null;
for (ProxyParameter proxyParameter : kineticLaw.getProxyParameters()) {
if (proxyParameter.getName().equals(reaction_p.getName())) {
rule_proxy_parameter = (LocalProxyParameter) proxyParameter;
}
}
if (rule_proxy_parameter != null) {
// we want to convert to local
boolean bConvertToGlobal = false;
kineticLaw.convertParameterType(rule_proxy_parameter, bConvertToGlobal);
} else {
// could find neither local parameter nor proxy parameter
throw new RuntimeException("user defined parameter " + reaction_p.getName() + " from reaction " + rs.getName() + " didn't map to a reactionRule parameter");
}
} else if (rule_p.getRole() == RbmKineticLawParameterType.UserDefined) {
kineticLaw.setParameterValue(rule_p, reaction_p.getExpression(), true);
rule_p.setUnitDefinition(reaction_p.getUnitDefinition());
} else {
throw new RuntimeException("user defined parameter " + reaction_p.getName() + " from reaction " + rs.getName() + " mapped to a reactionRule parameter with unexpected role " + rule_p.getRole().getDescription());
}
}
}
} catch (ExpressionException e) {
e.printStackTrace();
throw new RuntimeException("Problem attempting to set RbmKineticLaw expression: " + e.getMessage());
}
rr.setKineticLaw(kineticLaw);
KineticsParameter[] kpList = k.getKineticsParameters();
ModelParameter[] mpList = rs.getModel().getModelParameters();
ModelParameter mp = rs.getModel().getModelParameter(kpList[0].getName());
ReactionParticipant[] rpList = rs.getReactionParticipants();
for (ReactionParticipant p : rpList) {
if (p instanceof Reactant) {
int stoichiometry = p.getStoichiometry();
for (int i = 0; i < stoichiometry; i++) {
SpeciesPattern speciesPattern = new SpeciesPattern(rs.getModel(), p.getSpeciesContext().getSpeciesPattern());
ReactantPattern reactantPattern = new ReactantPattern(speciesPattern, p.getStructure());
rr.addReactant(reactantPattern);
}
} else if (p instanceof Product) {
int stoichiometry = p.getStoichiometry();
for (int i = 0; i < stoichiometry; i++) {
SpeciesPattern speciesPattern = new SpeciesPattern(rs.getModel(), p.getSpeciesContext().getSpeciesPattern());
ProductPattern productPattern = new ProductPattern(speciesPattern, p.getStructure());
rr.addProduct(productPattern);
}
}
}
// commented code below is probably obsolete, we verify (above) in the reaction the number of participants,
// no need to do it again in the corresponding rule
// if(rr.getReactantPatterns().size() > 2) {
// String message = "NFSim doesn't support more than 2 reactants within a reaction: " + name;
// throw new RuntimeException(message);
// }
// if(rr.getProductPatterns().size() > 2) {
// String message = "NFSim doesn't support more than 2 products within a reaction: " + name;
// throw new RuntimeException(message);
// }
newModel.removeReactionStep(rs);
newModel.getRbmModelContainer().addReactionRule(rr);
}
for (ReactionRuleSpec rrs : transformedSimulationContext.getReactionContext().getReactionRuleSpecs()) {
if (rrs == null) {
continue;
}
ReactionRule rr = rrs.getReactionRule();
if (rrs.isExcluded()) {
// delete those rules which are disabled (excluded) in the Specifications / Reaction table
newModel.getRbmModelContainer().removeReactionRule(rr);
continue;
}
}
// now that we generated the rules we can delete the reaction steps they're coming from
for (ReactionStep rs : newModel.getReactionSteps()) {
newModel.removeReactionStep(rs);
}
try {
// we invoke bngl just for the purpose of generating the xml file, which we'll then use to extract the symmetry factor
generateNetwork(transformedSimulationContext, fromReactions, mathMappingCallback);
} catch (ClassNotFoundException | IOException e) {
// TODO Auto-generated catch block
e.printStackTrace();
}
System.out.println("Finished RuleBased Transformer.");
}
use of cbit.vcell.parser.ExpressionException in project vcell by virtualcell.
the class RulebasedTransformer method parseBngOutput.
private void parseBngOutput(SimulationContext simContext, Set<ReactionRule> fromReactions, BNGOutput bngOutput) {
Model model = simContext.getModel();
Document bngNFSimXMLDocument = bngOutput.getNFSimXMLDocument();
bngRootElement = bngNFSimXMLDocument.getRootElement();
Element modelElement = bngRootElement.getChild("model", Namespace.getNamespace("http://www.sbml.org/sbml/level3"));
Element listOfReactionRulesElement = modelElement.getChild("ListOfReactionRules", Namespace.getNamespace("http://www.sbml.org/sbml/level3"));
List<Element> reactionRuleChildren = new ArrayList<Element>();
reactionRuleChildren = listOfReactionRulesElement.getChildren("ReactionRule", Namespace.getNamespace("http://www.sbml.org/sbml/level3"));
for (Element reactionRuleElement : reactionRuleChildren) {
ReactionRuleAnalysisReport rar = new ReactionRuleAnalysisReport();
boolean bForward = true;
String rule_id_str = reactionRuleElement.getAttributeValue("id");
ruleElementMap.put(rule_id_str, reactionRuleElement);
String rule_name_str = reactionRuleElement.getAttributeValue("name");
String symmetry_factor_str = reactionRuleElement.getAttributeValue("symmetry_factor");
Double symmetry_factor_double = Double.parseDouble(symmetry_factor_str);
ReactionRule rr = null;
if (rule_name_str.startsWith("_reverse_")) {
bForward = false;
rule_name_str = rule_name_str.substring("_reverse_".length());
rr = model.getRbmModelContainer().getReactionRule(rule_name_str);
RbmKineticLaw rkl = rr.getKineticLaw();
try {
if (symmetry_factor_double != 1.0 && fromReactions.contains(rr)) {
Expression expression = rkl.getLocalParameterValue(RbmKineticLawParameterType.MassActionReverseRate);
expression = Expression.div(expression, new Expression(symmetry_factor_double));
rkl.setLocalParameterValue(RbmKineticLawParameterType.MassActionReverseRate, expression);
}
} catch (ExpressionException | PropertyVetoException exc) {
exc.printStackTrace();
throw new RuntimeException("Unexpected transform exception: " + exc.getMessage());
}
rulesReverseMap.put(rr, rar);
} else {
rr = model.getRbmModelContainer().getReactionRule(rule_name_str);
RbmKineticLaw rkl = rr.getKineticLaw();
try {
if (symmetry_factor_double != 1.0 && fromReactions.contains(rr)) {
Expression expression = rkl.getLocalParameterValue(RbmKineticLawParameterType.MassActionForwardRate);
expression = Expression.div(expression, new Expression(symmetry_factor_double));
rkl.setLocalParameterValue(RbmKineticLawParameterType.MassActionForwardRate, expression);
}
} catch (ExpressionException | PropertyVetoException exc) {
exc.printStackTrace();
throw new RuntimeException("Unexpected transform exception: " + exc.getMessage());
}
rulesForwardMap.put(rr, rar);
}
rar.symmetryFactor = symmetry_factor_double;
System.out.println("rule id=" + rule_id_str + ", name=" + rule_name_str + ", symmetry factor=" + symmetry_factor_str);
keyMap.put(rule_id_str, rr);
Element listOfReactantPatternsElement = reactionRuleElement.getChild("ListOfReactantPatterns", Namespace.getNamespace("http://www.sbml.org/sbml/level3"));
List<Element> reactantPatternChildren = new ArrayList<Element>();
reactantPatternChildren = listOfReactantPatternsElement.getChildren("ReactantPattern", Namespace.getNamespace("http://www.sbml.org/sbml/level3"));
for (int i = 0; i < reactantPatternChildren.size(); i++) {
Element reactantPatternElement = reactantPatternChildren.get(i);
String pattern_id_str = reactantPatternElement.getAttributeValue("id");
ReactionRuleParticipant p = bForward ? rr.getReactantPattern(i) : rr.getProductPattern(i);
SpeciesPattern sp = p.getSpeciesPattern();
System.out.println(" reactant id=" + pattern_id_str + ", name=" + sp.toString());
keyMap.put(pattern_id_str, sp);
// list of molecules
extractMolecules(sp, model, reactantPatternElement);
// list of bonds (not implemented)
extractBonds(reactantPatternElement);
}
Element listOfProductPatternsElement = reactionRuleElement.getChild("ListOfProductPatterns", Namespace.getNamespace("http://www.sbml.org/sbml/level3"));
List<Element> productPatternChildren = new ArrayList<Element>();
productPatternChildren = listOfProductPatternsElement.getChildren("ProductPattern", Namespace.getNamespace("http://www.sbml.org/sbml/level3"));
for (int i = 0; i < productPatternChildren.size(); i++) {
Element productPatternElement = productPatternChildren.get(i);
String pattern_id_str = productPatternElement.getAttributeValue("id");
ReactionRuleParticipant p = bForward ? rr.getProductPattern(i) : rr.getReactantPattern(i);
SpeciesPattern sp = p.getSpeciesPattern();
System.out.println(" product id=" + pattern_id_str + ", name=" + sp.toString());
keyMap.put(pattern_id_str, sp);
extractMolecules(sp, model, productPatternElement);
extractBonds(productPatternElement);
}
// extract the Map for this rule
Element listOfMapElement = reactionRuleElement.getChild("Map", Namespace.getNamespace("http://www.sbml.org/sbml/level3"));
List<Element> mapChildren = new ArrayList<Element>();
mapChildren = listOfMapElement.getChildren("MapItem", Namespace.getNamespace("http://www.sbml.org/sbml/level3"));
for (Element mapElement : mapChildren) {
String target_id_str = mapElement.getAttributeValue("targetID");
String source_id_str = mapElement.getAttributeValue("sourceID");
System.out.println("Map: target=" + target_id_str + " source=" + source_id_str);
RbmObject target_object = keyMap.get(target_id_str);
RbmObject source_object = keyMap.get(source_id_str);
if (target_object == null)
System.out.println("!!! Missing map target " + target_id_str);
if (source_object == null)
System.out.println("!!! Missing map source " + source_id_str);
if (source_object != null) {
// target_object may be null
System.out.println(" target=" + target_object + " source=" + source_object);
Pair<RbmObject, RbmObject> mapEntry = new Pair<RbmObject, RbmObject>(target_object, source_object);
rar.objmappingList.add(mapEntry);
}
Pair<String, String> idmapEntry = new Pair<String, String>(target_id_str, source_id_str);
rar.idmappingList.add(idmapEntry);
}
// ListOfOperations
Element listOfOperationsElement = reactionRuleElement.getChild("ListOfOperations", Namespace.getNamespace("http://www.sbml.org/sbml/level3"));
List<Element> operationsChildren = new ArrayList<Element>();
operationsChildren = listOfOperationsElement.getChildren("StateChange", Namespace.getNamespace("http://www.sbml.org/sbml/level3"));
System.out.println("ListOfOperations");
for (Element operationsElement : operationsChildren) {
String finalState_str = operationsElement.getAttributeValue("finalState");
String site_str = operationsElement.getAttributeValue("site");
RbmObject site_object = keyMap.get(site_str);
if (site_object == null)
System.out.println("!!! Missing map object " + site_str);
if (site_object != null) {
System.out.println(" finalState=" + finalState_str + " site=" + site_object);
StateChangeOperation sco = new StateChangeOperation(finalState_str, site_str, site_object);
rar.operationsList.add(sco);
}
}
operationsChildren = listOfOperationsElement.getChildren("AddBond", Namespace.getNamespace("http://www.sbml.org/sbml/level3"));
for (Element operationsElement : operationsChildren) {
String site1_str = operationsElement.getAttributeValue("site1");
String site2_str = operationsElement.getAttributeValue("site2");
RbmObject site1_object = keyMap.get(site1_str);
RbmObject site2_object = keyMap.get(site2_str);
if (site1_object == null)
System.out.println("!!! Missing map object " + site1_str);
if (site2_object == null)
System.out.println("!!! Missing map object " + site2_str);
if (site1_object != null && site2_object != null) {
System.out.println(" site1=" + site1_object + " site2=" + site2_object);
AddBondOperation abo = new AddBondOperation(site1_str, site2_str, site1_object, site2_object);
rar.operationsList.add(abo);
}
}
operationsChildren = listOfOperationsElement.getChildren("DeleteBond", Namespace.getNamespace("http://www.sbml.org/sbml/level3"));
for (Element operationsElement : operationsChildren) {
String site1_str = operationsElement.getAttributeValue("site1");
String site2_str = operationsElement.getAttributeValue("site2");
RbmObject site1_object = keyMap.get(site1_str);
RbmObject site2_object = keyMap.get(site2_str);
if (site1_object == null)
System.out.println("!!! Missing map object " + site1_str);
if (site2_object == null)
System.out.println("!!! Missing map object " + site2_str);
if (site1_object != null && site2_object != null) {
System.out.println(" site1=" + site1_object + " site2=" + site2_object);
DeleteBondOperation dbo = new DeleteBondOperation(site1_str, site2_str, site1_object, site2_object);
rar.operationsList.add(dbo);
}
}
operationsChildren = listOfOperationsElement.getChildren("Add", Namespace.getNamespace("http://www.sbml.org/sbml/level3"));
for (Element operationsElement : operationsChildren) {
String id_str = operationsElement.getAttributeValue("id");
RbmObject id_object = keyMap.get(id_str);
if (id_object == null)
System.out.println("!!! Missing map object " + id_str);
if (id_object != null) {
System.out.println(" id=" + id_str);
AddOperation ao = new AddOperation(id_str, id_object);
rar.operationsList.add(ao);
}
}
operationsChildren = listOfOperationsElement.getChildren("Delete", Namespace.getNamespace("http://www.sbml.org/sbml/level3"));
for (Element operationsElement : operationsChildren) {
String id_str = operationsElement.getAttributeValue("id");
String delete_molecules_str = operationsElement.getAttributeValue("DeleteMolecules");
RbmObject id_object = keyMap.get(id_str);
if (id_object == null)
System.out.println("!!! Missing map object " + id_str);
int delete_molecules_int = 0;
if (delete_molecules_str != null) {
delete_molecules_int = Integer.parseInt(delete_molecules_str);
}
if (id_object != null) {
System.out.println(" id=" + id_str + ", DeleteMolecules=" + delete_molecules_str);
DeleteOperation dop = new DeleteOperation(id_str, id_object, delete_molecules_int);
rar.operationsList.add(dop);
}
}
}
System.out.println("done parsing xml file");
}
use of cbit.vcell.parser.ExpressionException in project vcell by virtualcell.
the class MathTestingUtilities method getInsideOutsideFunctions.
/**
* Insert the method's description here.
* Creation date: (1/23/2003 10:30:23 PM)
* @return cbit.vcell.parser.Expression
* @param analyticSubDomainExp cbit.vcell.parser.Expression
*/
public static Expression[] getInsideOutsideFunctions(Expression analyticSubDomainExp) throws ExpressionException, MappingException {
Expression analyticExp = new Expression(analyticSubDomainExp);
analyticExp.bindExpression(null);
analyticExp = analyticExp.flatten();
java.util.Stack<Expression> unparsedExpressionStack = new java.util.Stack<Expression>();
java.util.Vector<Expression> expList = new java.util.Vector<Expression>();
unparsedExpressionStack.push(analyticExp);
while (unparsedExpressionStack.size() > 0) {
Expression exp = unparsedExpressionStack.pop();
if (exp.isRelational()) {
ExpressionTerm expTerm = exp.extractTopLevelTerm();
if (expTerm.getOperator().equals("<") || expTerm.getOperator().equals("<=")) {
expList.add(new Expression(expTerm.getOperands()[0].infix() + "-" + expTerm.getOperands()[1].infix()));
} else if (expTerm.getOperator().equals(">") || expTerm.getOperator().equals(">=")) {
expList.add(new Expression(expTerm.getOperands()[1].infix() + "-" + expTerm.getOperands()[0].infix()));
} else {
throw new ExpressionException("relational expression '" + exp + "' is not an inequality");
}
} else if (exp.isLogical()) {
ExpressionTerm expTerm = exp.extractTopLevelTerm();
for (int i = 0; i < expTerm.getOperands().length; i++) {
unparsedExpressionStack.push(expTerm.getOperands()[i]);
}
} else {
throw new ExpressionException("expression '" + exp + "' is neither relational nor logical, bad analytic geometry");
}
}
return (Expression[]) BeanUtils.getArray(expList, Expression.class);
}
use of cbit.vcell.parser.ExpressionException in project vcell by virtualcell.
the class FVSolverStandalone method writeVCGAndResampleFieldData.
protected void writeVCGAndResampleFieldData() throws SolverException {
fireSolverStarting(SimulationMessage.MESSAGE_SOLVEREVENT_STARTING_PROC_GEOM);
try {
// write subdomains file
SubdomainInfo.write(new File(getSaveDirectory(), baseName + SimDataConstants.SUBDOMAINS_FILE_SUFFIX), simTask.getSimulation().getMathDescription());
PrintWriter pw = new PrintWriter(new FileWriter(new File(getSaveDirectory(), baseName + SimDataConstants.VCG_FILE_EXTENSION)));
GeometryFileWriter.write(pw, getResampledGeometry());
pw.close();
FieldDataIdentifierSpec[] argFieldDataIDSpecs = simTask.getSimulationJob().getFieldDataIdentifierSpecs();
if (argFieldDataIDSpecs != null && argFieldDataIDSpecs.length > 0) {
fireSolverStarting(SimulationMessage.MESSAGE_SOLVEREVENT_STARTING_RESAMPLE_FD);
FieldFunctionArguments psfFieldFunc = null;
Variable var = simTask.getSimulationJob().getSimulationSymbolTable().getVariable(Simulation.PSF_FUNCTION_NAME);
if (var != null) {
FieldFunctionArguments[] ffas = FieldUtilities.getFieldFunctionArguments(var.getExpression());
if (ffas == null || ffas.length == 0) {
throw new DataAccessException("Point Spread Function " + Simulation.PSF_FUNCTION_NAME + " can only be a single field function.");
} else {
Expression newexp;
try {
newexp = new Expression(ffas[0].infix());
if (!var.getExpression().compareEqual(newexp)) {
throw new DataAccessException("Point Spread Function " + Simulation.PSF_FUNCTION_NAME + " can only be a single field function.");
}
psfFieldFunc = ffas[0];
} catch (ExpressionException e) {
e.printStackTrace();
throw new DataAccessException(e.getMessage());
}
}
}
boolean[] bResample = new boolean[argFieldDataIDSpecs.length];
Arrays.fill(bResample, true);
for (int i = 0; i < argFieldDataIDSpecs.length; i++) {
argFieldDataIDSpecs[i].getFieldFuncArgs().getTime().bindExpression(simTask.getSimulationJob().getSimulationSymbolTable());
if (argFieldDataIDSpecs[i].getFieldFuncArgs().equals(psfFieldFunc)) {
bResample[i] = false;
}
}
int numMembraneElements = getResampledGeometry().getGeometrySurfaceDescription().getSurfaceCollection().getTotalPolygonCount();
CartesianMesh simpleMesh = CartesianMesh.createSimpleCartesianMesh(getResampledGeometry().getOrigin(), getResampledGeometry().getExtent(), simTask.getSimulation().getMeshSpecification().getSamplingSize(), getResampledGeometry().getGeometrySurfaceDescription().getRegionImage());
String secondarySimDataDir = PropertyLoader.getProperty(PropertyLoader.secondarySimDataDirInternalProperty, null);
DataSetControllerImpl dsci = new DataSetControllerImpl(null, getSaveDirectory().getParentFile(), secondarySimDataDir == null ? null : new File(secondarySimDataDir));
dsci.writeFieldFunctionData(null, argFieldDataIDSpecs, bResample, simpleMesh, simResampleInfoProvider, numMembraneElements, HESM_OVERWRITE_AND_CONTINUE);
}
} catch (Exception e) {
throw new SolverException(e.getMessage());
}
}
Aggregations