use of cbit.vcell.model.Kinetics.KineticsParameter in project vcell by virtualcell.
the class ParticleMathMapping method refreshMathDescription.
/**
* This method was created in VisualAge.
*/
private void refreshMathDescription() throws MappingException, MatrixException, MathException, ExpressionException, ModelException {
getSimulationContext().checkValidity();
if (getSimulationContext().getGeometry().getDimension() == 0) {
throw new MappingException("particle math mapping requires spatial geometry - dimension >= 1");
}
StructureMapping[] structureMappings = getSimulationContext().getGeometryContext().getStructureMappings();
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 = getSimulationContext().getBioEvents();
if (bioEvents != null && bioEvents.length > 0) {
throw new MappingException("events not yet supported for particle-based models");
}
//
// gather only those reactionSteps that are not "excluded"
//
ReactionSpec[] reactionSpecs = getSimulationContext().getReactionContext().getReactionSpecs();
Vector<ReactionStep> rsList = new Vector<ReactionStep>();
for (int i = 0; i < reactionSpecs.length; i++) {
if (reactionSpecs[i].isExcluded() == false) {
if (reactionSpecs[i].isFast()) {
throw new MappingException("fast reactions not supported for particle models");
}
rsList.add(reactionSpecs[i].getReactionStep());
}
}
ReactionStep[] reactionSteps = new ReactionStep[rsList.size()];
rsList.copyInto(reactionSteps);
//
for (int i = 0; i < reactionSteps.length; i++) {
Kinetics.UnresolvedParameter[] unresolvedParameters = reactionSteps[i].getKinetics().getUnresolvedParameters();
if (unresolvedParameters != null && unresolvedParameters.length > 0) {
StringBuffer buffer = new StringBuffer();
for (int j = 0; j < unresolvedParameters.length; j++) {
if (j > 0) {
buffer.append(", ");
}
buffer.append(unresolvedParameters[j].getName());
}
throw new MappingException(reactionSteps[i].getDisplayType() + " '" + reactionSteps[i].getName() + "' contains unresolved identifier(s): " + buffer);
}
}
//
// temporarily place all variables in a hashtable (before binding) and discarding duplicates (check for equality)
//
VariableHash varHash = new VariableHash();
// //
// // verify that all structures are mapped to geometry classes and all geometry classes are mapped to a structure
// //
// Structure structures[] = getSimulationContext().getGeometryContext().getModel().getStructures();
// for (int i = 0; i < structures.length; i++){
// StructureMapping sm = getSimulationContext().getGeometryContext().getStructureMapping(structures[i]);
// if (sm==null || (sm.getGeometryClass() == null)){
// throw new MappingException("model structure '"+structures[i].getName()+"' not mapped to a geometry subdomain");
// }
// if (sm.getUnitSizeParameter()!=null){
// Expression unitSizeExp = sm.getUnitSizeParameter().getExpression();
// if(unitSizeExp != null)
// {
// try {
// double unitSize = unitSizeExp.evaluateConstant();
// if (unitSize != 1.0){
// throw new MappingException("model structure '"+sm.getStructure().getName()+"' unit size = "+unitSize+" != 1.0 ... partial volume or surface mapping not yet supported for particles");
// }
// }catch (ExpressionException e){
// e.printStackTrace(System.out);
// throw new MappingException("couldn't evaluate unit size for model structure '"+sm.getStructure().getName()+"' : "+e.getMessage());
// }
// }
// }
// }
// {
// GeometryClass[] geometryClass = getSimulationContext().getGeometryContext().getGeometry().getGeometryClasses();
// for (int i = 0; i < geometryClass.length; i++){
// Structure[] mappedStructures = getSimulationContext().getGeometryContext().getStructuresFromGeometryClass(geometryClass[i]);
// if (mappedStructures==null || mappedStructures.length==0){
// throw new MappingException("geometryClass '"+geometryClass[i].getName()+"' not mapped from a model structure");
// }
// }
// }
// deals with model parameters
Model model = getSimulationContext().getModel();
ModelUnitSystem modelUnitSystem = model.getUnitSystem();
ModelParameter[] modelParameters = model.getModelParameters();
// populate in globalParameterVariants hashtable
for (int j = 0; j < modelParameters.length; j++) {
Expression modelParamExpr = modelParameters[j].getExpression();
GeometryClass geometryClass = getDefaultGeometryClass(modelParamExpr);
modelParamExpr = getIdentifierSubstitutions(modelParamExpr, modelParameters[j].getUnitDefinition(), geometryClass);
varHash.addVariable(newFunctionOrConstant(getMathSymbol(modelParameters[j], geometryClass), modelParamExpr, geometryClass));
}
//
// create new MathDescription (based on simContext's previous MathDescription if possible)
//
MathDescription oldMathDesc = getSimulationContext().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(getSimulationContext().getName() + "_generated");
}
//
// volume particle variables
//
Enumeration<SpeciesContextMapping> enum1 = getSpeciesContextMappings();
while (enum1.hasMoreElements()) {
SpeciesContextMapping scm = enum1.nextElement();
if (scm.getVariable() instanceof ParticleVariable) {
if (!(mathDesc.getVariable(scm.getVariable().getName()) instanceof ParticleVariable)) {
varHash.addVariable(scm.getVariable());
}
}
}
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(getSimulationContext().getTemperatureKelvin()), model.getTEMPERATURE().getUnitDefinition(), null)));
//
for (int j = 0; j < structureMappings.length; j++) {
if (structureMappings[j] instanceof MembraneMapping) {
MembraneMapping membraneMapping = (MembraneMapping) structureMappings[j];
GeometryClass geometryClass = membraneMapping.getGeometryClass();
//
// don't calculate voltage, still may need it though
//
Parameter initialVoltageParm = membraneMapping.getInitialVoltageParameter();
Variable voltageFunction = newFunctionOrConstant(getMathSymbol(membraneMapping.getMembrane().getMembraneVoltage(), geometryClass), getIdentifierSubstitutions(initialVoltageParm.getExpression(), initialVoltageParm.getUnitDefinition(), geometryClass), geometryClass);
varHash.addVariable(voltageFunction);
varHash.addVariable(newFunctionOrConstant(getMathSymbol(membraneMapping.getMembrane().getMembraneVoltage(), membraneMapping.getGeometryClass()), getIdentifierSubstitutions(membraneMapping.getInitialVoltageParameter().getExpression(), membraneMapping.getInitialVoltageParameter().getUnitDefinition(), membraneMapping.getGeometryClass()), membraneMapping.getGeometryClass()));
}
}
//
for (int j = 0; j < reactionSteps.length; j++) {
ReactionStep rs = reactionSteps[j];
if (getSimulationContext().getReactionContext().getReactionSpec(rs).isExcluded()) {
continue;
}
Kinetics.KineticsParameter[] parameters = rs.getKinetics().getKineticsParameters();
GeometryClass geometryClass = null;
if (rs.getStructure() != null) {
geometryClass = getSimulationContext().getGeometryContext().getStructureMapping(rs.getStructure()).getGeometryClass();
}
if (parameters != null) {
for (int i = 0; i < parameters.length; i++) {
// Reaction rate, currentDensity, LumpedCurrent and null parameters are not going to displayed in the particle math description.
if (((parameters[i].getRole() == Kinetics.ROLE_CurrentDensity) || (parameters[i].getRole() == Kinetics.ROLE_LumpedCurrent) || (parameters[i].getRole() == Kinetics.ROLE_ReactionRate)) || (parameters[i].getExpression() == null)) {
continue;
}
varHash.addVariable(newFunctionOrConstant(getMathSymbol(parameters[i], geometryClass), getIdentifierSubstitutions(parameters[i].getExpression(), parameters[i].getUnitDefinition(), geometryClass), geometryClass));
}
}
}
//
// initial constants (either function or constant)
//
SpeciesContextSpec[] speciesContextSpecs = getSimulationContext().getReactionContext().getSpeciesContextSpecs();
for (int i = 0; i < speciesContextSpecs.length; i++) {
SpeciesContextSpecParameter initParm = null;
Expression initExpr = null;
if (getSimulationContext().isUsingConcentration()) {
initParm = speciesContextSpecs[i].getParameterFromRole(SpeciesContextSpec.ROLE_InitialConcentration);
initExpr = new Expression(initParm.getExpression());
// if (speciesContextSpecs[i].getSpeciesContext().getStructure() instanceof Feature) {
// initExpr = Expression.div(initExpr, new Expression(model.getKMOLE, getNameScope())).flatten();
// }
} else {
initParm = speciesContextSpecs[i].getParameterFromRole(SpeciesContextSpec.ROLE_InitialCount);
initExpr = new Expression(initParm.getExpression());
}
if (initExpr != null) {
StructureMapping sm = getSimulationContext().getGeometryContext().getStructureMapping(speciesContextSpecs[i].getSpeciesContext().getStructure());
String[] symbols = initExpr.getSymbols();
// Check if 'initExpr' has other speciesContexts in its expression, need to replace it with 'spContext_init'
for (int j = 0; symbols != null && j < symbols.length; j++) {
// if symbol is a speciesContext, replacing it with a reference to initial condition for that speciesContext.
SpeciesContext spC = null;
SymbolTableEntry ste = initExpr.getSymbolBinding(symbols[j]);
if (ste instanceof SpeciesContextSpecProxyParameter) {
SpeciesContextSpecProxyParameter spspp = (SpeciesContextSpecProxyParameter) ste;
if (spspp.getTarget() instanceof SpeciesContext) {
spC = (SpeciesContext) spspp.getTarget();
SpeciesContextSpec spcspec = getSimulationContext().getReactionContext().getSpeciesContextSpec(spC);
SpeciesContextSpecParameter spCInitParm = spcspec.getParameterFromRole(SpeciesContextSpec.ROLE_InitialConcentration);
// if initConc param expression is null, try initCount
if (spCInitParm.getExpression() == null) {
spCInitParm = spcspec.getParameterFromRole(SpeciesContextSpec.ROLE_InitialCount);
}
// need to get init condn expression, but can't get it from getMathSymbol() (mapping between bio and math), hence get it as below.
Expression scsInitExpr = new Expression(spCInitParm, getNameScope());
// scsInitExpr.bindExpression(this);
initExpr.substituteInPlace(new Expression(spC.getName()), scsInitExpr);
}
}
}
// now create the appropriate function for the current speciesContextSpec.
varHash.addVariable(newFunctionOrConstant(getMathSymbol(initParm, sm.getGeometryClass()), getIdentifierSubstitutions(initExpr, initParm.getUnitDefinition(), sm.getGeometryClass()), sm.getGeometryClass()));
}
}
//
for (int i = 0; i < speciesContextSpecs.length; i++) {
SpeciesContextSpec.SpeciesContextSpecParameter diffParm = speciesContextSpecs[i].getParameterFromRole(SpeciesContextSpec.ROLE_DiffusionRate);
if (diffParm != null) {
StructureMapping sm = getSimulationContext().getGeometryContext().getStructureMapping(speciesContextSpecs[i].getSpeciesContext().getStructure());
varHash.addVariable(newFunctionOrConstant(getMathSymbol(diffParm, sm.getGeometryClass()), getIdentifierSubstitutions(diffParm.getExpression(), diffParm.getUnitDefinition(), sm.getGeometryClass()), sm.getGeometryClass()));
}
}
//
for (int i = 0; i < speciesContextSpecs.length; i++) {
SpeciesContextSpec.SpeciesContextSpecParameter bc_xm = speciesContextSpecs[i].getParameterFromRole(SpeciesContextSpec.ROLE_BoundaryValueXm);
StructureMapping sm = getSimulationContext().getGeometryContext().getStructureMapping(speciesContextSpecs[i].getSpeciesContext().getStructure());
if (bc_xm != null && (bc_xm.getExpression() != null)) {
varHash.addVariable(newFunctionOrConstant(getMathSymbol(bc_xm, sm.getGeometryClass()), getIdentifierSubstitutions(bc_xm.getExpression(), bc_xm.getUnitDefinition(), sm.getGeometryClass()), sm.getGeometryClass()));
}
SpeciesContextSpec.SpeciesContextSpecParameter bc_xp = speciesContextSpecs[i].getParameterFromRole(SpeciesContextSpec.ROLE_BoundaryValueXp);
if (bc_xp != null && (bc_xp.getExpression() != null)) {
varHash.addVariable(newFunctionOrConstant(getMathSymbol(bc_xp, sm.getGeometryClass()), getIdentifierSubstitutions(bc_xp.getExpression(), bc_xp.getUnitDefinition(), sm.getGeometryClass()), sm.getGeometryClass()));
}
SpeciesContextSpec.SpeciesContextSpecParameter bc_ym = speciesContextSpecs[i].getParameterFromRole(SpeciesContextSpec.ROLE_BoundaryValueYm);
if (bc_ym != null && (bc_ym.getExpression() != null)) {
varHash.addVariable(newFunctionOrConstant(getMathSymbol(bc_ym, sm.getGeometryClass()), getIdentifierSubstitutions(bc_ym.getExpression(), bc_ym.getUnitDefinition(), sm.getGeometryClass()), sm.getGeometryClass()));
}
SpeciesContextSpec.SpeciesContextSpecParameter bc_yp = speciesContextSpecs[i].getParameterFromRole(SpeciesContextSpec.ROLE_BoundaryValueYp);
if (bc_yp != null && (bc_yp.getExpression() != null)) {
varHash.addVariable(newFunctionOrConstant(getMathSymbol(bc_yp, sm.getGeometryClass()), getIdentifierSubstitutions(bc_yp.getExpression(), bc_yp.getUnitDefinition(), sm.getGeometryClass()), sm.getGeometryClass()));
}
SpeciesContextSpec.SpeciesContextSpecParameter bc_zm = speciesContextSpecs[i].getParameterFromRole(SpeciesContextSpec.ROLE_BoundaryValueZm);
if (bc_zm != null && (bc_zm.getExpression() != null)) {
varHash.addVariable(newFunctionOrConstant(getMathSymbol(bc_zm, sm.getGeometryClass()), getIdentifierSubstitutions(bc_zm.getExpression(), bc_zm.getUnitDefinition(), sm.getGeometryClass()), sm.getGeometryClass()));
}
SpeciesContextSpec.SpeciesContextSpecParameter bc_zp = speciesContextSpecs[i].getParameterFromRole(SpeciesContextSpec.ROLE_BoundaryValueZp);
if (bc_zp != null && (bc_zp.getExpression() != null)) {
varHash.addVariable(newFunctionOrConstant(getMathSymbol(bc_zp, sm.getGeometryClass()), getIdentifierSubstitutions(bc_zp.getExpression(), bc_zp.getUnitDefinition(), sm.getGeometryClass()), sm.getGeometryClass()));
}
}
//
for (int i = 0; i < speciesContextSpecs.length; i++) {
SpeciesContextSpec.SpeciesContextSpecParameter advection_velX = speciesContextSpecs[i].getParameterFromRole(SpeciesContextSpec.ROLE_VelocityX);
StructureMapping sm = getSimulationContext().getGeometryContext().getStructureMapping(speciesContextSpecs[i].getSpeciesContext().getStructure());
GeometryClass geometryClass = sm.getGeometryClass();
if (advection_velX != null && (advection_velX.getExpression() != null)) {
varHash.addVariable(newFunctionOrConstant(getMathSymbol(advection_velX, geometryClass), getIdentifierSubstitutions(advection_velX.getExpression(), advection_velX.getUnitDefinition(), geometryClass), geometryClass));
}
SpeciesContextSpec.SpeciesContextSpecParameter advection_velY = speciesContextSpecs[i].getParameterFromRole(SpeciesContextSpec.ROLE_VelocityY);
if (advection_velY != null && (advection_velY.getExpression() != null)) {
varHash.addVariable(newFunctionOrConstant(getMathSymbol(advection_velY, geometryClass), getIdentifierSubstitutions(advection_velY.getExpression(), advection_velY.getUnitDefinition(), geometryClass), geometryClass));
}
SpeciesContextSpec.SpeciesContextSpecParameter advection_velZ = speciesContextSpecs[i].getParameterFromRole(SpeciesContextSpec.ROLE_VelocityZ);
if (advection_velZ != null && (advection_velZ.getExpression() != null)) {
varHash.addVariable(newFunctionOrConstant(getMathSymbol(advection_velZ, geometryClass), getIdentifierSubstitutions(advection_velZ.getExpression(), advection_velZ.getUnitDefinition(), geometryClass), geometryClass));
}
}
//
// constant species (either function or constant)
//
enum1 = getSpeciesContextMappings();
while (enum1.hasMoreElements()) {
SpeciesContextMapping scm = (SpeciesContextMapping) enum1.nextElement();
if (scm.getVariable() instanceof Constant) {
varHash.addVariable(scm.getVariable());
}
}
//
// conversion factors
//
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.getKMILLIVOLTS(), null), getIdentifierSubstitutions(model.getKMILLIVOLTS().getExpression(), model.getKMILLIVOLTS().getUnitDefinition(), null)));
varHash.addVariable(new Constant(getMathSymbol(model.getK_GHK(), null), getIdentifierSubstitutions(model.getK_GHK().getExpression(), model.getK_GHK().getUnitDefinition(), null)));
//
for (int i = 0; i < structureMappings.length; i++) {
StructureMapping sm = structureMappings[i];
if (getSimulationContext().getGeometry().getDimension() == 0) {
StructureMappingParameter sizeParm = sm.getSizeParameter();
if (sizeParm != null && sizeParm.getExpression() != null) {
varHash.addVariable(newFunctionOrConstant(getMathSymbol(sizeParm, sm.getGeometryClass()), getIdentifierSubstitutions(sizeParm.getExpression(), sizeParm.getUnitDefinition(), sm.getGeometryClass()), sm.getGeometryClass()));
} else {
if (sm instanceof MembraneMapping) {
MembraneMapping mm = (MembraneMapping) sm;
StructureMappingParameter volFrac = mm.getVolumeFractionParameter();
if (volFrac != null && volFrac.getExpression() != null) {
varHash.addVariable(newFunctionOrConstant(getMathSymbol(volFrac, sm.getGeometryClass()), getIdentifierSubstitutions(volFrac.getExpression(), volFrac.getUnitDefinition(), sm.getGeometryClass()), sm.getGeometryClass()));
}
StructureMappingParameter surfToVol = mm.getSurfaceToVolumeParameter();
if (surfToVol != null && surfToVol.getExpression() != null) {
varHash.addVariable(newFunctionOrConstant(getMathSymbol(surfToVol, sm.getGeometryClass()), getIdentifierSubstitutions(surfToVol.getExpression(), surfToVol.getUnitDefinition(), sm.getGeometryClass()), sm.getGeometryClass()));
}
}
}
} else {
Parameter parm = sm.getParameterFromRole(StructureMapping.ROLE_AreaPerUnitArea);
if (parm != null && parm.getExpression() != null && sm.getGeometryClass() instanceof SurfaceClass) {
varHash.addVariable(newFunctionOrConstant(getMathSymbol(parm, sm.getGeometryClass()), getIdentifierSubstitutions(parm.getExpression(), parm.getUnitDefinition(), sm.getGeometryClass()), sm.getGeometryClass()));
}
parm = sm.getParameterFromRole(StructureMapping.ROLE_AreaPerUnitVolume);
if (parm != null && parm.getExpression() != null && sm.getGeometryClass() instanceof SubVolume) {
varHash.addVariable(newFunctionOrConstant(getMathSymbol(parm, sm.getGeometryClass()), getIdentifierSubstitutions(parm.getExpression(), parm.getUnitDefinition(), sm.getGeometryClass()), sm.getGeometryClass()));
}
parm = sm.getParameterFromRole(StructureMapping.ROLE_VolumePerUnitArea);
if (parm != null && parm.getExpression() != null && sm.getGeometryClass() instanceof SurfaceClass) {
varHash.addVariable(newFunctionOrConstant(getMathSymbol(parm, sm.getGeometryClass()), getIdentifierSubstitutions(parm.getExpression(), parm.getUnitDefinition(), sm.getGeometryClass()), sm.getGeometryClass()));
}
parm = sm.getParameterFromRole(StructureMapping.ROLE_VolumePerUnitVolume);
if (parm != null && parm.getExpression() != null && sm.getGeometryClass() instanceof SubVolume) {
varHash.addVariable(newFunctionOrConstant(getMathSymbol(parm, sm.getGeometryClass()), getIdentifierSubstitutions(parm.getExpression(), parm.getUnitDefinition(), sm.getGeometryClass()), sm.getGeometryClass()));
}
}
}
//
// functions
//
enum1 = getSpeciesContextMappings();
while (enum1.hasMoreElements()) {
SpeciesContextMapping scm = (SpeciesContextMapping) enum1.nextElement();
if (scm.getVariable() == null && scm.getDependencyExpression() != null) {
StructureMapping sm = getSimulationContext().getGeometryContext().getStructureMapping(scm.getSpeciesContext().getStructure());
Variable dependentVariable = newFunctionOrConstant(getMathSymbol(scm.getSpeciesContext(), sm.getGeometryClass()), getIdentifierSubstitutions(scm.getDependencyExpression(), scm.getSpeciesContext().getUnitDefinition(), sm.getGeometryClass()), sm.getGeometryClass());
dependentVariable.setDomain(new Domain(sm.getGeometryClass()));
varHash.addVariable(dependentVariable);
}
}
//
for (int i = 0; i < fieldMathMappingParameters.length; i++) {
if (fieldMathMappingParameters[i] instanceof UnitFactorParameter) {
GeometryClass geometryClass = fieldMathMappingParameters[i].getGeometryClass();
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());
//
if (getSimulationContext().getGeometryContext().getGeometry() != null) {
try {
mathDesc.setGeometry(getSimulationContext().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");
}
//
// create subdomains (volume and surfaces)
//
GeometryClass[] geometryClasses = getSimulationContext().getGeometryContext().getGeometry().getGeometryClasses();
for (int k = 0; k < geometryClasses.length; k++) {
if (geometryClasses[k] instanceof SubVolume) {
SubVolume subVolume = (SubVolume) geometryClasses[k];
//
// get priority of subDomain
//
// now does not have to match spatial feature, *BUT* needs to be unique
int priority = k;
//
// create subDomain
//
CompartmentSubDomain subDomain = new CompartmentSubDomain(subVolume.getName(), priority);
mathDesc.addSubDomain(subDomain);
//
// assign boundary condition types
//
StructureMapping[] mappedSMs = getSimulationContext().getGeometryContext().getStructureMappings(subVolume);
FeatureMapping mappedFM = null;
for (int i = 0; i < mappedSMs.length; i++) {
if (mappedSMs[i] instanceof FeatureMapping) {
if (mappedFM != null) {
lg.warn("WARNING:::: MathMapping.refreshMathDescription() ... assigning boundary condition types not unique");
}
mappedFM = (FeatureMapping) mappedSMs[i];
}
}
if (mappedFM != null) {
subDomain.setBoundaryConditionXm(mappedFM.getBoundaryConditionTypeXm());
subDomain.setBoundaryConditionXp(mappedFM.getBoundaryConditionTypeXp());
if (getSimulationContext().getGeometry().getDimension() > 1) {
subDomain.setBoundaryConditionYm(mappedFM.getBoundaryConditionTypeYm());
subDomain.setBoundaryConditionYp(mappedFM.getBoundaryConditionTypeYp());
}
if (getSimulationContext().getGeometry().getDimension() > 2) {
subDomain.setBoundaryConditionZm(mappedFM.getBoundaryConditionTypeZm());
subDomain.setBoundaryConditionZp(mappedFM.getBoundaryConditionTypeZp());
}
}
} else if (geometryClasses[k] instanceof SurfaceClass) {
SurfaceClass surfaceClass = (SurfaceClass) geometryClasses[k];
// determine membrane inside and outside subvolume
// this preserves backward compatibility so that membrane subdomain
// inside and outside correspond to structure hierarchy when present
Pair<SubVolume, SubVolume> ret = DiffEquMathMapping.computeBoundaryConditionSource(model, simContext, surfaceClass);
SubVolume innerSubVolume = ret.one;
SubVolume outerSubVolume = ret.two;
//
// create subDomain
//
CompartmentSubDomain outerCompartment = mathDesc.getCompartmentSubDomain(outerSubVolume.getName());
CompartmentSubDomain innerCompartment = mathDesc.getCompartmentSubDomain(innerSubVolume.getName());
MembraneSubDomain memSubDomain = new MembraneSubDomain(innerCompartment, outerCompartment, surfaceClass.getName());
mathDesc.addSubDomain(memSubDomain);
}
}
//
// create Particle Contexts for all Particle Variables
//
Enumeration<SpeciesContextMapping> enumSCM = getSpeciesContextMappings();
Expression unitFactor = getUnitFactor(modelUnitSystem.getStochasticSubstanceUnit().divideBy(modelUnitSystem.getVolumeSubstanceUnit()));
while (enumSCM.hasMoreElements()) {
SpeciesContextMapping scm = enumSCM.nextElement();
SpeciesContext sc = scm.getSpeciesContext();
StructureMapping sm = getSimulationContext().getGeometryContext().getStructureMapping(sc.getStructure());
SpeciesContextSpec scs = getSimulationContext().getReactionContext().getSpeciesContextSpec(sc);
if (scm.getVariable() instanceof ParticleVariable && scm.getDependencyExpression() == null) {
ParticleVariable particleVariable = (ParticleVariable) scm.getVariable();
//
// initial distribution of particles
//
ArrayList<ParticleInitialCondition> particleInitialConditions = new ArrayList<ParticleInitialCondition>();
ParticleInitialCondition pic = null;
if (getSimulationContext().isUsingConcentration()) {
Expression initialDistribution = scs.getInitialConcentrationParameter().getExpression() == null ? null : new Expression(getMathSymbol(scs.getInitialConcentrationParameter(), sm.getGeometryClass()));
if (particleVariable instanceof VolumeParticleVariable) {
initialDistribution = Expression.mult(initialDistribution, unitFactor);
}
pic = new ParticleInitialConditionConcentration(initialDistribution);
} else {
Expression initialCount = scs.getInitialCountParameter().getExpression() == null ? null : new Expression(getMathSymbol(scs.getInitialCountParameter(), sm.getGeometryClass()));
if (initialCount == null) {
throw new MappingException("initialCount not defined for speciesContext " + scs.getSpeciesContext().getName());
}
Expression locationX = new Expression("u");
Expression locationY = new Expression("u");
Expression locationZ = new Expression("u");
pic = new ParticleInitialConditionCount(initialCount, locationX, locationY, locationZ);
}
particleInitialConditions.add(pic);
//
// diffusion
//
Expression diffusion = new Expression(getMathSymbol(scs.getDiffusionParameter(), sm.getGeometryClass()));
Expression driftXExp = null;
if (scs.getVelocityXParameter().getExpression() != null) {
driftXExp = new Expression(getMathSymbol(scs.getVelocityXParameter(), sm.getGeometryClass()));
} else {
SpatialQuantity[] velX_quantities = scs.getVelocityQuantities(QuantityComponent.X);
if (velX_quantities.length > 0) {
int numRegions = simContext.getGeometry().getGeometrySurfaceDescription().getGeometricRegions(sm.getGeometryClass()).length;
if (velX_quantities.length == 1 && numRegions == 1) {
driftXExp = new Expression(getMathSymbol(velX_quantities[0], sm.getGeometryClass()));
} else {
throw new MappingException("multiple advection velocities enabled set for multiple volume domains ");
}
}
}
Expression driftYExp = null;
if (scs.getVelocityYParameter().getExpression() != null) {
driftYExp = new Expression(getMathSymbol(scs.getVelocityYParameter(), sm.getGeometryClass()));
} else {
SpatialQuantity[] velY_quantities = scs.getVelocityQuantities(QuantityComponent.Y);
if (velY_quantities.length > 0) {
int numRegions = simContext.getGeometry().getGeometrySurfaceDescription().getGeometricRegions(sm.getGeometryClass()).length;
if (velY_quantities.length == 1 && numRegions == 1) {
driftYExp = new Expression(getMathSymbol(velY_quantities[0], sm.getGeometryClass()));
} else {
throw new MappingException("multiple advection velocities enabled set for multiple volume domains ");
}
}
}
Expression driftZExp = null;
if (scs.getVelocityZParameter().getExpression() != null) {
driftZExp = new Expression(getMathSymbol(scs.getVelocityZParameter(), sm.getGeometryClass()));
} else {
SpatialQuantity[] velZ_quantities = scs.getVelocityQuantities(QuantityComponent.Z);
if (velZ_quantities.length > 0) {
int numRegions = simContext.getGeometry().getGeometrySurfaceDescription().getGeometricRegions(sm.getGeometryClass()).length;
if (velZ_quantities.length == 1 && numRegions == 1) {
driftZExp = new Expression(getMathSymbol(velZ_quantities[0], sm.getGeometryClass()));
} else {
throw new MappingException("multiple advection velocities enabled set for multiple volume domains ");
}
}
}
ParticleProperties particleProperties = new ParticleProperties(particleVariable, diffusion, driftXExp, driftYExp, driftZExp, particleInitialConditions);
GeometryClass myGC = sm.getGeometryClass();
if (myGC == null) {
throw new MappingException("Application '" + getSimulationContext().getName() + "'\nGeometry->StructureMapping->(" + sm.getStructure().getTypeName() + ")'" + sm.getStructure().getName() + "' must be mapped to geometry domain.\n(see 'Problems' tab)");
}
SubDomain subDomain = mathDesc.getSubDomain(myGC.getName());
subDomain.addParticleProperties(particleProperties);
}
}
for (ReactionStep reactionStep : reactionSteps) {
Kinetics kinetics = reactionStep.getKinetics();
StructureMapping sm = getSimulationContext().getGeometryContext().getStructureMapping(reactionStep.getStructure());
GeometryClass reactionStepGeometryClass = sm.getGeometryClass();
SubDomain subdomain = mathDesc.getSubDomain(reactionStepGeometryClass.getName());
KineticsParameter reactionRateParameter = null;
if (kinetics instanceof LumpedKinetics) {
reactionRateParameter = ((LumpedKinetics) kinetics).getLumpedReactionRateParameter();
} else {
reactionRateParameter = ((DistributedKinetics) kinetics).getReactionRateParameter();
}
// macroscopic_irreversible/Microscopic_irreversible for bimolecular membrane reactions. They will NOT go through MassAction solver.
if (kinetics.getKineticsDescription().equals(KineticsDescription.Macroscopic_irreversible) || kinetics.getKineticsDescription().equals(KineticsDescription.Microscopic_irreversible)) {
Expression radiusExp = getIdentifierSubstitutions(reactionStep.getKinetics().getKineticsParameterFromRole(Kinetics.ROLE_Binding_Radius).getExpression(), modelUnitSystem.getBindingRadiusUnit(), reactionStepGeometryClass);
if (radiusExp != null) {
Expression expCopy = new Expression(radiusExp);
try {
MassActionSolver.substituteParameters(expCopy, true).evaluateConstant();
} catch (ExpressionException e) {
throw new MathException(VCellErrorMessages.getMassActionSolverMessage(reactionStep.getName(), "Problem in binding radius of " + reactionStep.getName() + ": '" + radiusExp.infix() + "', " + e.getMessage()));
}
} else {
throw new MathException(VCellErrorMessages.getMassActionSolverMessage(reactionStep.getName(), "Binding radius of " + reactionStep.getName() + " is null."));
}
List<ParticleVariable> reactantParticles = new ArrayList<ParticleVariable>();
List<ParticleVariable> productParticles = new ArrayList<ParticleVariable>();
List<Action> forwardActions = new ArrayList<Action>();
for (ReactionParticipant rp : reactionStep.getReactionParticipants()) {
SpeciesContext sc = rp.getSpeciesContext();
SpeciesContextSpec scs = getSimulationContext().getReactionContext().getSpeciesContextSpec(sc);
GeometryClass scGeometryClass = getSimulationContext().getGeometryContext().getStructureMapping(sc.getStructure()).getGeometryClass();
String varName = getMathSymbol(sc, scGeometryClass);
Variable var = mathDesc.getVariable(varName);
if (var instanceof ParticleVariable) {
ParticleVariable particle = (ParticleVariable) var;
if (rp instanceof Reactant) {
reactantParticles.add(particle);
if (!scs.isConstant() && !scs.isForceContinuous()) {
for (int i = 0; i < Math.abs(rp.getStoichiometry()); i++) {
if (radiusExp != null) {
forwardActions.add(Action.createDestroyAction(particle));
}
}
}
} else if (rp instanceof Product) {
productParticles.add(particle);
if (!scs.isConstant() && !scs.isForceContinuous()) {
for (int i = 0; i < Math.abs(rp.getStoichiometry()); i++) {
if (radiusExp != null) {
forwardActions.add(Action.createCreateAction(particle));
}
}
}
}
} else {
throw new MappingException("particle variable '" + varName + "' not found");
}
}
JumpProcessRateDefinition bindingRadius = new InteractionRadius(radiusExp);
// get jump process name
String jpName = TokenMangler.mangleToSName(reactionStep.getName());
// only for NFSim/Rules for now.
ProcessSymmetryFactor processSymmetryFactor = null;
if (forwardActions.size() > 0) {
ParticleJumpProcess forwardProcess = new ParticleJumpProcess(jpName, reactantParticles, bindingRadius, forwardActions, processSymmetryFactor);
subdomain.addParticleJumpProcess(forwardProcess);
}
} else // other type of reactions
{
/* check the reaction rate law to see if we need to decompose a reaction(reversible) into two jump processes.
rate constants are important in calculating the probability rate.
for Mass Action, we use KForward and KReverse,
for General Kinetics we parse reaction rate J to see if it is in Mass Action form.
*/
Expression forwardRate = null;
Expression reverseRate = null;
// Using the MassActionFunction to write out the math description
MassActionSolver.MassActionFunction maFunc = null;
if (kinetics.getKineticsDescription().equals(KineticsDescription.MassAction) || kinetics.getKineticsDescription().equals(KineticsDescription.General) || kinetics.getKineticsDescription().equals(KineticsDescription.GeneralPermeability)) {
Expression rateExp = kinetics.getKineticsParameterFromRole(Kinetics.ROLE_ReactionRate).getExpression();
Parameter forwardRateParameter = null;
Parameter reverseRateParameter = null;
if (kinetics.getKineticsDescription().equals(KineticsDescription.MassAction)) {
forwardRateParameter = kinetics.getKineticsParameterFromRole(Kinetics.ROLE_KForward);
reverseRateParameter = kinetics.getKineticsParameterFromRole(Kinetics.ROLE_KReverse);
} else if (kinetics.getKineticsDescription().equals(KineticsDescription.GeneralPermeability)) {
forwardRateParameter = kinetics.getKineticsParameterFromRole(Kinetics.ROLE_Permeability);
reverseRateParameter = kinetics.getKineticsParameterFromRole(Kinetics.ROLE_Permeability);
}
maFunc = MassActionSolver.solveMassAction(forwardRateParameter, reverseRateParameter, rateExp, reactionStep);
if (maFunc.getForwardRate() == null && maFunc.getReverseRate() == null) {
throw new MappingException("Cannot generate stochastic math mapping for the reaction:" + reactionStep.getName() + "\nLooking for the rate function according to the form of k1*Reactant1^Stoir1*Reactant2^Stoir2...-k2*Product1^Stoip1*Product2^Stoip2.");
} else {
if (maFunc.getForwardRate() != null) {
forwardRate = maFunc.getForwardRate();
}
if (maFunc.getReverseRate() != null) {
reverseRate = maFunc.getReverseRate();
}
}
}
if (maFunc != null) {
// if the reaction has forward rate (Mass action,HMMs), or don't have either forward or reverse rate (some other rate laws--like general)
// we process it as forward reaction
List<ParticleVariable> reactantParticles = new ArrayList<ParticleVariable>();
List<ParticleVariable> productParticles = new ArrayList<ParticleVariable>();
List<Action> forwardActions = new ArrayList<Action>();
List<Action> reverseActions = new ArrayList<Action>();
List<ReactionParticipant> reactants = maFunc.getReactants();
List<ReactionParticipant> products = maFunc.getProducts();
for (ReactionParticipant rp : reactants) {
SpeciesContext sc = rp.getSpeciesContext();
SpeciesContextSpec scs = getSimulationContext().getReactionContext().getSpeciesContextSpec(sc);
GeometryClass scGeometryClass = getSimulationContext().getGeometryContext().getStructureMapping(sc.getStructure()).getGeometryClass();
String varName = getMathSymbol(sc, scGeometryClass);
Variable var = mathDesc.getVariable(varName);
if (var instanceof ParticleVariable) {
ParticleVariable particle = (ParticleVariable) var;
reactantParticles.add(particle);
if (!scs.isConstant() && !scs.isForceContinuous()) {
for (int i = 0; i < Math.abs(rp.getStoichiometry()); i++) {
if (forwardRate != null) {
forwardActions.add(Action.createDestroyAction(particle));
}
if (reverseRate != null) {
reverseActions.add(Action.createCreateAction(particle));
}
}
}
} else {
throw new MappingException("particle variable '" + varName + "' not found");
}
}
for (ReactionParticipant rp : products) {
SpeciesContext sc = rp.getSpeciesContext();
SpeciesContextSpec scs = getSimulationContext().getReactionContext().getSpeciesContextSpec(sc);
GeometryClass scGeometryClass = getSimulationContext().getGeometryContext().getStructureMapping(sc.getStructure()).getGeometryClass();
String varName = getMathSymbol(sc, scGeometryClass);
Variable var = mathDesc.getVariable(varName);
if (var instanceof ParticleVariable) {
ParticleVariable particle = (ParticleVariable) var;
productParticles.add(particle);
if (!scs.isConstant() && !scs.isForceContinuous()) {
for (int i = 0; i < Math.abs(rp.getStoichiometry()); i++) {
if (forwardRate != null) {
forwardActions.add(Action.createCreateAction(particle));
}
if (reverseRate != null) {
reverseActions.add(Action.createDestroyAction(particle));
}
}
}
} else {
throw new MappingException("particle variable '" + varName + "' not found");
}
}
//
// There are two unit conversions required:
//
// 1) convert entire reaction rate from vcell reaction units to Smoldyn units (molecules/lengthunit^dim/timeunit)
// (where dim is 2 for membrane reactions and 3 for volume reactions)
//
// for forward rates:
// 2) convert each reactant from Smoldyn units (molecules/lengthunit^dim) to VCell units
// (where dim is 2 for membrane reactants and 3 for volume reactants)
//
// or
//
// for reverse rates:
// 2) convert each product from Smoldyn units (molecules/lengthunit^dim) to VCell units
// (where dim is 2 for membrane products and 3 for volume products)
//
RationalNumber reactionLocationDim = new RationalNumber(reactionStep.getStructure().getDimension());
VCUnitDefinition timeUnit = modelUnitSystem.getTimeUnit();
VCUnitDefinition smoldynReactionSizeUnit = modelUnitSystem.getLengthUnit().raiseTo(reactionLocationDim);
VCUnitDefinition smoldynSubstanceUnit = modelUnitSystem.getStochasticSubstanceUnit();
VCUnitDefinition smoldynReactionRateUnit = smoldynSubstanceUnit.divideBy(smoldynReactionSizeUnit).divideBy(timeUnit);
VCUnitDefinition vcellReactionRateUnit = reactionRateParameter.getUnitDefinition();
VCUnitDefinition reactionUnitFactor = smoldynReactionRateUnit.divideBy(vcellReactionRateUnit);
if (forwardRate != null) {
VCUnitDefinition smoldynReactantsUnit = modelUnitSystem.getInstance_DIMENSIONLESS();
// start with factor to translate entire reaction rate.
VCUnitDefinition forwardUnitFactor = reactionUnitFactor;
//
for (ReactionParticipant reactant : maFunc.getReactants()) {
VCUnitDefinition vcellReactantUnit = reactant.getSpeciesContext().getUnitDefinition();
boolean bForceContinuous = simContext.getReactionContext().getSpeciesContextSpec(reactant.getSpeciesContext()).isForceContinuous();
VCUnitDefinition smoldynReactantUnit = null;
if (bForceContinuous) {
// reactant is continuous (vcell units)
smoldynReactantUnit = reactant.getSpeciesContext().getUnitDefinition();
} else {
// reactant is a particle (smoldyn units)
RationalNumber reactantLocationDim = new RationalNumber(reactant.getStructure().getDimension());
VCUnitDefinition smoldynReactantSize = modelUnitSystem.getLengthUnit().raiseTo(reactantLocationDim);
smoldynReactantUnit = smoldynSubstanceUnit.divideBy(smoldynReactantSize);
}
// keep track of units of all reactants
smoldynReactantsUnit = smoldynReactantsUnit.multiplyBy(smoldynReactantUnit);
RationalNumber reactantStoichiometry = new RationalNumber(reactant.getStoichiometry());
VCUnitDefinition reactantUnitFactor = (vcellReactantUnit.divideBy(smoldynReactantUnit)).raiseTo(reactantStoichiometry);
// accumulate unit factors for all reactants
forwardUnitFactor = forwardUnitFactor.multiplyBy(reactantUnitFactor);
}
forwardRate = Expression.mult(forwardRate, getUnitFactor(forwardUnitFactor));
VCUnitDefinition smoldynExpectedForwardRateUnit = smoldynReactionRateUnit.divideBy(smoldynReactantsUnit);
// get probability
Expression exp = getIdentifierSubstitutions(forwardRate, smoldynExpectedForwardRateUnit, reactionStepGeometryClass).flatten();
JumpProcessRateDefinition partRateDef = new MacroscopicRateConstant(exp);
// create particle jump process
String jpName = TokenMangler.mangleToSName(reactionStep.getName());
// only for NFSim/Rules for now.
ProcessSymmetryFactor processSymmetryFactor = null;
if (forwardActions.size() > 0) {
ParticleJumpProcess forwardProcess = new ParticleJumpProcess(jpName, reactantParticles, partRateDef, forwardActions, processSymmetryFactor);
subdomain.addParticleJumpProcess(forwardProcess);
}
}
// end of forward rate not null
if (reverseRate != null) {
VCUnitDefinition smoldynProductsUnit = modelUnitSystem.getInstance_DIMENSIONLESS();
// start with factor to translate entire reaction rate.
VCUnitDefinition reverseUnitFactor = reactionUnitFactor;
//
for (ReactionParticipant product : maFunc.getProducts()) {
VCUnitDefinition vcellProductUnit = product.getSpeciesContext().getUnitDefinition();
boolean bForceContinuous = simContext.getReactionContext().getSpeciesContextSpec(product.getSpeciesContext()).isForceContinuous();
VCUnitDefinition smoldynProductUnit = null;
if (bForceContinuous) {
smoldynProductUnit = product.getSpeciesContext().getUnitDefinition();
} else {
RationalNumber productLocationDim = new RationalNumber(product.getStructure().getDimension());
VCUnitDefinition smoldynProductSize = modelUnitSystem.getLengthUnit().raiseTo(productLocationDim);
smoldynProductUnit = smoldynSubstanceUnit.divideBy(smoldynProductSize);
}
// keep track of units of all products
smoldynProductsUnit = smoldynProductsUnit.multiplyBy(smoldynProductUnit);
RationalNumber productStoichiometry = new RationalNumber(product.getStoichiometry());
VCUnitDefinition productUnitFactor = (vcellProductUnit.divideBy(smoldynProductUnit)).raiseTo(productStoichiometry);
// accumulate unit factors for all products
reverseUnitFactor = reverseUnitFactor.multiplyBy(productUnitFactor);
}
reverseRate = Expression.mult(reverseRate, getUnitFactor(reverseUnitFactor));
VCUnitDefinition smoldynExpectedReverseRateUnit = smoldynReactionRateUnit.divideBy(smoldynProductsUnit);
// get probability
Expression exp = getIdentifierSubstitutions(reverseRate, smoldynExpectedReverseRateUnit, reactionStepGeometryClass).flatten();
JumpProcessRateDefinition partProbRate = new MacroscopicRateConstant(exp);
// get jump process name
String jpName = TokenMangler.mangleToSName(reactionStep.getName() + "_reverse");
// only for NFSim/Rules for now.
ProcessSymmetryFactor processSymmetryFactor = null;
if (reverseActions.size() > 0) {
ParticleJumpProcess reverseProcess = new ParticleJumpProcess(jpName, productParticles, partProbRate, reverseActions, processSymmetryFactor);
subdomain.addParticleJumpProcess(reverseProcess);
}
}
// end of reverse rate not null
}
// end of maFunc not null
}
// end of reaction step for loop
}
//
for (int i = 0; i < fieldMathMappingParameters.length; i++) {
if (fieldMathMappingParameters[i] instanceof UnitFactorParameter) {
GeometryClass geometryClass = fieldMathMappingParameters[i].getGeometryClass();
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()) {
lg.warn(mathDesc.getVCML_database());
throw new MappingException("generated an invalid mathDescription: " + mathDesc.getWarning());
}
if (lg.isDebugEnabled()) {
System.out.println("]]]]]]]]]]]]]]]]]]]]]] VCML string begin ]]]]]]]]]]]]]]]]]]]]]]]]]]]]]]]]]]]]]");
System.out.println(mathDesc.getVCML());
System.out.println("]]]]]]]]]]]]]]]]]]]]]] VCML string end ]]]]]]]]]]]]]]]]]]]]]]]]]]]]]]]]]]]]]");
}
}
use of cbit.vcell.model.Kinetics.KineticsParameter 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.model.Kinetics.KineticsParameter in project vcell by virtualcell.
the class ReactionStep method setKinetics.
/**
* Sets the kinetics property (cbit.vcell.model.Kinetics) value.
* @param kinetics The new value for the property.
* @see #getKinetics
*/
public void setKinetics(Kinetics kinetics) {
Kinetics oldValue = fieldKinetics;
//
if (oldValue != null) {
oldValue.removePropertyChangeListener(this);
oldValue.removePropertyChangeListener(this);
removePropertyChangeListener(oldValue);
removePropertyChangeListener(oldValue);
}
fieldKinetics = kinetics;
if (kinetics != null) {
kinetics.removePropertyChangeListener(this);
kinetics.addPropertyChangeListener(this);
removePropertyChangeListener(kinetics);
addPropertyChangeListener(kinetics);
}
//
try {
KineticsParameter chargeValenceParameter = getKinetics().getChargeValenceParameter();
if (kinetics.getKineticsDescription().isElectrical()) {
if (getPhysicsOptions() == PHYSICS_MOLECULAR_ONLY) {
if (chargeValenceParameter != null && !chargeValenceParameter.getExpression().isZero()) {
setPhysicsOptions(PHYSICS_MOLECULAR_AND_ELECTRICAL);
} else {
setPhysicsOptions(PHYSICS_ELECTRICAL_ONLY);
}
}
} else {
if (getPhysicsOptions() == PHYSICS_ELECTRICAL_ONLY) {
if (chargeValenceParameter != null && !chargeValenceParameter.getExpression().isZero()) {
setPhysicsOptions(PHYSICS_MOLECULAR_AND_ELECTRICAL);
} else {
setPhysicsOptions(PHYSICS_MOLECULAR_ONLY);
}
}
}
if (getKinetics() != null) {
if (chargeValenceParameter != null) {
if (chargeValenceParameter.getExpression().isZero()) {
chargeValenceParameter.setExpression(new Expression(1.0));
}
}
}
} catch (PropertyVetoException e) {
e.printStackTrace(System.out);
}
firePropertyChange(PROPERTY_NAME_KINETICS, oldValue, kinetics);
}
use of cbit.vcell.model.Kinetics.KineticsParameter in project vcell by virtualcell.
the class ReactionCartoonTool method detailsDeleteSpecies.
private static DeleteSpeciesInfo detailsDeleteSpecies(Component requester, SpeciesContext[] speciesContextArr, ReactionStep[] toBeDeletedReactStepArr, ReactionCartoon rxCartoon) throws Exception, UserCancelException {
if (speciesContextArr == null || speciesContextArr.length == 0) {
return null;
}
// Warn user that there may be some BioModel components that reference speciesContext to be removed
// Get ReactionParticipant list
Collection<Shape> rxPartColl = rxCartoon.getShapes();
HashMap<SpeciesContext, HashSet<ReactionParticipant>> rxPartHashMap = new HashMap<SpeciesContext, HashSet<ReactionParticipant>>();
for (Shape objShape : rxPartColl) {
if (objShape instanceof ReactionParticipantShape) {
ReactionParticipant objReactionParticipant = ((ReactionParticipantShape) objShape).getReactionParticipant();
if (Arrays.asList(speciesContextArr).contains(objReactionParticipant.getSpeciesContext())) {
if (!rxPartHashMap.containsKey(objReactionParticipant.getSpeciesContext())) {
rxPartHashMap.put(objReactionParticipant.getSpeciesContext(), new HashSet<ReactionParticipant>());
}
if (!rxPartHashMap.get(objReactionParticipant.getSpeciesContext()).contains(objReactionParticipant)) {
rxPartHashMap.get(objReactionParticipant.getSpeciesContext()).add(objReactionParticipant);
}
}
}
}
int reactionParticipantCount = 0;
for (HashSet<ReactionParticipant> objReactPart : rxPartHashMap.values()) {
reactionParticipantCount += objReactPart.size();
}
// find bioModel and get SymbolTable references to SpeciesContext
BioModelEditor bioModelEditor = (BioModelEditor) BeanUtils.findTypeParentOfComponent(requester, BioModelEditor.class);
if (bioModelEditor == null) {
throw new Exception("Error deleting Speciescontext, Can't find BiomodelEditor");
}
BioModel bioModel = bioModelEditor.getBioModelWindowManager().getVCDocument();
HashMap<SpeciesContext, HashSet<SymbolTableEntry>> referencingSymbolsHashMap = new HashMap<SpeciesContext, HashSet<SymbolTableEntry>>();
for (int i = 0; i < speciesContextArr.length; i++) {
List<SymbolTableEntry> referencingSymbolsList = bioModel.findReferences(speciesContextArr[i]);
if (referencingSymbolsList != null && referencingSymbolsList.size() > 0) {
if (!referencingSymbolsHashMap.containsKey(speciesContextArr[i])) {
referencingSymbolsHashMap.put(speciesContextArr[i], new HashSet<SymbolTableEntry>());
}
referencingSymbolsHashMap.get(speciesContextArr[i]).addAll(referencingSymbolsList);
}
}
int referencingSymbolsCount = 0;
for (HashSet<SymbolTableEntry> objSimTableEntry : referencingSymbolsHashMap.values()) {
referencingSymbolsCount += objSimTableEntry.size();
}
// Warn user about delete
HashMap<SpeciesContext, Boolean> bUnresolvableHashMap = new HashMap<SpeciesContext, Boolean>();
for (int i = 0; i < speciesContextArr.length; i++) {
bUnresolvableHashMap.put(speciesContextArr[i], Boolean.FALSE);
}
String[][] rowData = null;
if (rxPartHashMap.size() == 0 && referencingSymbolsHashMap.size() == 0) {
rowData = new String[speciesContextArr.length][4];
for (int i = 0; i < speciesContextArr.length; i++) {
rowData[i][0] = speciesContextArr[i].getName();
rowData[i][1] = "";
rowData[i][2] = "";
rowData[i][3] = "";
}
Arrays.sort(rowData, new Comparator<String[]>() {
@Override
public int compare(String[] o1, String[] o2) {
return o1[0].compareToIgnoreCase(o2[0]);
}
});
} else {
// find SpeciesContext that had no reference warnings
Vector<SpeciesContext> speciesContextNoReferences = new Vector<SpeciesContext>();
for (int i = 0; i < speciesContextArr.length; i++) {
if (!rxPartHashMap.containsKey(speciesContextArr[i]) && !referencingSymbolsHashMap.containsKey(speciesContextArr[i])) {
speciesContextNoReferences.add(speciesContextArr[i]);
}
}
rowData = new String[reactionParticipantCount + referencingSymbolsCount + speciesContextNoReferences.size()][4];
int count = 0;
for (SpeciesContext objSpeciesContext : speciesContextNoReferences) {
rowData[count][0] = objSpeciesContext.getName();
rowData[count][1] = "";
rowData[count][2] = "";
rowData[count][3] = "";
count++;
}
for (SpeciesContext objSpeciesContext : rxPartHashMap.keySet()) {
Iterator<ReactionParticipant> iterRxPart = rxPartHashMap.get(objSpeciesContext).iterator();
while (iterRxPart.hasNext()) {
rowData[count][0] = objSpeciesContext.getName();
rowData[count][1] = "";
rowData[count][2] = "Reaction Diagram stoichiometry '" + iterRxPart.next().getReactionStep().getName() + "'";
rowData[count][3] = "";
count++;
}
}
for (SpeciesContext objSpeciesContext : referencingSymbolsHashMap.keySet()) {
Iterator<SymbolTableEntry> iterSymbolTableEntry = referencingSymbolsHashMap.get(objSpeciesContext).iterator();
while (iterSymbolTableEntry.hasNext()) {
rowData[count][0] = objSpeciesContext.getName();
rowData[count][1] = "";
SymbolTableEntry objSymbolTableEntry = iterSymbolTableEntry.next();
boolean bKineticsParameter = objSymbolTableEntry instanceof KineticsParameter;
if (bKineticsParameter) {
KineticsParameter kp = (KineticsParameter) objSymbolTableEntry;
boolean isOK = kp.isRegenerated();
for (int i = 0; toBeDeletedReactStepArr != null && i < toBeDeletedReactStepArr.length; i++) {
if (toBeDeletedReactStepArr[i] == kp.getKinetics().getReactionStep()) {
// OK to delete this Speciescontext if were deleting the reaction that contained the reference
isOK = true;
}
}
rowData[count][1] = (isOK ? "" : RXSPECIES_ERROR);
bUnresolvableHashMap.put(objSpeciesContext, bUnresolvableHashMap.get(objSpeciesContext) || !isOK);
}
boolean bReaction = objSymbolTableEntry.getNameScope() instanceof ReactionNameScope;
rowData[count][2] = (bReaction ? "Reaction" : objSymbolTableEntry.getNameScope().getClass().getName()) + "( " + objSymbolTableEntry.getNameScope().getName() + " )";
rowData[count][3] = (bKineticsParameter ? "Parameter" : objSymbolTableEntry.getClass().getName()) + "( " + objSymbolTableEntry.getName() + " )";
count++;
}
}
// for (SymbolTableEntry referencingSTE : referencingSymbols) {
// System.out.println("REFERENCE "+referencingSTE.getClass().getName()+"("+referencingSTE.getName()+") nameScope "+referencingSTE.getNameScope().getClass().getName()+"("+referencingSTE.getNameScope().getName()+")");
// }
Arrays.sort(rowData, new Comparator<String[]>() {
@Override
public int compare(String[] o1, String[] o2) {
return o1[0].compareToIgnoreCase(o2[0]);
}
});
}
if (rowData == null || rowData.length == 0) {
return null;
}
return new DeleteSpeciesInfo(rxPartHashMap, bUnresolvableHashMap, rowData);
}
use of cbit.vcell.model.Kinetics.KineticsParameter in project vcell by virtualcell.
the class BioCartoonTool method pasteReactionSteps0.
/**
* pasteReactionSteps0 : does the actual pasting. First called with a cloned model, to track issues. If user still wants to proceed, the paste
* is performed on the original model.
*
* Insert the method's description here.
* Creation date: (5/10/2003 3:55:25 PM)
* @param pasteToModel cbit.vcell.model.Model
* @param pasteToStructure cbit.vcell.model.Structure
* @param bNew boolean
*/
private static final PasteHelper pasteReactionSteps0(HashMap<String, HashMap<ReactionParticipant, Structure>> rxPartMapStructure, Component parent, IssueContext issueContext, ReactionStep[] copyFromRxSteps, Model pasteToModel, Structure pasteToStructure, boolean bNew, /*boolean bUseDBSpecies,*/
UserResolvedRxElements userResolvedRxElements) throws Exception {
HashMap<BioModelEntityObject, BioModelEntityObject> reactionsAndSpeciesContexts = new HashMap<>();
if (copyFromRxSteps == null || copyFromRxSteps.length == 0 || pasteToModel == null || pasteToStructure == null) {
throw new IllegalArgumentException("CartoonTool.pasteReactionSteps Error " + (copyFromRxSteps == null || copyFromRxSteps.length == 0 ? "reactionStepsArr empty " : "") + (pasteToModel == null ? "model is null " : "") + (pasteToStructure == null ? "struct is null " : ""));
}
if (!pasteToModel.contains(pasteToStructure)) {
throw new IllegalArgumentException("CartoonTool.pasteReactionSteps model " + pasteToModel.getName() + " does not contain structure " + pasteToStructure.getName());
}
// Check PasteToModel has preferred targets if set
if (userResolvedRxElements != null) {
for (int i = 0; i < userResolvedRxElements.toSpeciesArr.length; i++) {
if (userResolvedRxElements.toSpeciesArr[i] != null) {
if (!pasteToModel.contains(userResolvedRxElements.toSpeciesArr[i])) {
throw new RuntimeException("PasteToModel does not contain preferred Species " + userResolvedRxElements.toSpeciesArr[i]);
}
}
if (userResolvedRxElements.toStructureArr[i] != null) {
if (!pasteToModel.contains(userResolvedRxElements.toStructureArr[i])) {
throw new RuntimeException("PasteToModel does not contain preferred Structure " + userResolvedRxElements.toStructureArr[i]);
}
}
}
}
int counter = 0;
Structure currentStruct = pasteToStructure;
String copiedStructName = copyFromRxSteps[counter].getStructure().getName();
StructureTopology structTopology = (copyFromRxSteps[counter].getModel() == null ? pasteToModel.getStructureTopology() : copyFromRxSteps[counter].getModel().getStructureTopology());
IdentityHashMap<Species, Species> speciesHash = new IdentityHashMap<Species, Species>();
IdentityHashMap<SpeciesContext, SpeciesContext> speciesContextHash = new IdentityHashMap<SpeciesContext, SpeciesContext>();
Vector<Issue> issueVector = new Vector<Issue>();
do {
// create a new reaction, instead of cloning the old one; set struc
ReactionStep copyFromReactionStep = copyFromRxSteps[counter];
String newName = copyFromReactionStep.getName();
while (pasteToModel.getReactionStep(newName) != null) {
newName = org.vcell.util.TokenMangler.getNextEnumeratedToken(newName);
}
ReactionStep newReactionStep = null;
if (copyFromReactionStep instanceof SimpleReaction) {
newReactionStep = new SimpleReaction(pasteToModel, currentStruct, newName, copyFromReactionStep.isReversible());
} else if (copyFromReactionStep instanceof FluxReaction && currentStruct instanceof Membrane) {
newReactionStep = new FluxReaction(pasteToModel, (Membrane) currentStruct, null, newName, copyFromReactionStep.isReversible());
}
pasteToModel.addReactionStep(newReactionStep);
reactionsAndSpeciesContexts.put(newReactionStep, copyFromReactionStep);
Structure toRxnStruct = newReactionStep.getStructure();
Structure fromRxnStruct = copyFromReactionStep.getStructure();
if (!fromRxnStruct.getClass().equals(pasteToStructure.getClass())) {
throw new Exception("Cannot copy reaction from " + fromRxnStruct.getTypeName() + " to " + pasteToStructure.getTypeName() + ".");
}
// add appropriate reactionParticipants to newReactionStep.
StructureTopology toStructureTopology = pasteToModel.getStructureTopology();
ReactionParticipant[] copyFromRxParticipantArr = copyFromReactionStep.getReactionParticipants();
if (rxPartMapStructure == null) {
// null during 'issues' trial
rxPartMapStructure = new HashMap<String, HashMap<ReactionParticipant, Structure>>();
}
if (rxPartMapStructure.get(copyFromReactionStep.getName()) == null) {
// Ask user to assign species to compartments for each reaction to be pasted
rxPartMapStructure.put(copyFromReactionStep.getName(), askUserResolveMembraneConnections(parent, pasteToModel.getStructures(), currentStruct, fromRxnStruct, toRxnStruct, copyFromRxParticipantArr, toStructureTopology, structTopology));
}
for (int i = 0; i < copyFromRxParticipantArr.length; i += 1) {
Structure pasteToStruct = currentStruct;
// if(toRxnStruct instanceof Membrane){
pasteToStruct = rxPartMapStructure.get(copyFromReactionStep.getName()).get(copyFromRxParticipantArr[i]);
// if(pasteToStruct == null){
// for(ReactionParticipant myRXPart:rxPartMapStructure.get(copyFromReactionStep.getName()).keySet()){
// if(myRXPart.getSpeciesContext().getName().equals(copyFromRxParticipantArr[i].getSpeciesContext().getName())){
// pasteToStruct = rxPartMapStructure.get(copyFromReactionStep.getName()).get(myRXPart);
// break;
// }
// }
// }
// }
// this adds the speciesContexts and species (if any) to the model)
String rootSC = ReactionCartoonTool.speciesContextRootFinder(copyFromRxParticipantArr[i].getSpeciesContext());
SpeciesContext newSc = null;
SpeciesContext[] matchSC = pasteToModel.getSpeciesContexts();
for (int j = 0; matchSC != null && j < matchSC.length; j++) {
String matchRoot = ReactionCartoonTool.speciesContextRootFinder(matchSC[j]);
if (matchRoot != null && matchRoot.equals(rootSC) && matchSC[j].getStructure().getName().equals(pasteToStruct.getName())) {
newSc = matchSC[j];
reactionsAndSpeciesContexts.put(newSc, matchSC[j]);
break;
}
}
if (newSc == null) {
newSc = pasteSpecies(parent, copyFromRxParticipantArr[i].getSpecies(), rootSC, pasteToModel, pasteToStruct, bNew, /*bUseDBSpecies,*/
speciesHash, UserResolvedRxElements.getPreferredReactionElement(userResolvedRxElements, copyFromRxParticipantArr[i]));
reactionsAndSpeciesContexts.put(newSc, copyFromRxParticipantArr[i].getSpeciesContext());
}
// record the old-new speciesContexts (reactionparticipants) in the IdHashMap, this is useful, esp for 'Paste new', while replacing proxyparams.
SpeciesContext oldSc = copyFromRxParticipantArr[i].getSpeciesContext();
if (speciesContextHash.get(oldSc) == null) {
speciesContextHash.put(oldSc, newSc);
}
if (copyFromRxParticipantArr[i] instanceof Reactant) {
newReactionStep.addReactionParticipant(new Reactant(null, newReactionStep, newSc, copyFromRxParticipantArr[i].getStoichiometry()));
} else if (copyFromRxParticipantArr[i] instanceof Product) {
newReactionStep.addReactionParticipant(new Product(null, newReactionStep, newSc, copyFromRxParticipantArr[i].getStoichiometry()));
} else if (copyFromRxParticipantArr[i] instanceof Catalyst) {
newReactionStep.addCatalyst(newSc);
}
}
// // If 'newReactionStep' is a fluxRxn, set its fluxCarrier
// if (newReactionStep instanceof FluxReaction) {
// if (fluxCarrierSp != null) {
// ((FluxReaction)newReactionStep).setFluxCarrier(fluxCarrierSp, pasteToModel);
// } else {
// throw new RuntimeException("Could not set FluxCarrier species for the flux reaction to be pasted");
// }
// }
// For each kinetic parameter expression for new kinetics, replace the proxyParams from old kinetics with proxyParams in new kinetics
// i.e., if the proxyParams are speciesContexts, replace with corresponding speciesContext in newReactionStep;
// if the proxyParams are structureSizes or MembraneVoltages, replace with corresponding structure quantity in newReactionStep
Kinetics oldKinetics = copyFromReactionStep.getKinetics();
KineticsParameter[] oldKps = oldKinetics.getKineticsParameters();
KineticsProxyParameter[] oldKprps = oldKinetics.getProxyParameters();
Hashtable<String, Expression> paramExprHash = new Hashtable<String, Expression>();
for (int i = 0; oldKps != null && i < oldKps.length; i++) {
Expression newExpression = new Expression(oldKps[i].getExpression());
for (int j = 0; oldKprps != null && j < oldKprps.length; j++) {
// check if kinetic proxy parameter is in kinetic parameter expression
if (newExpression.hasSymbol(oldKprps[j].getName())) {
SymbolTableEntry ste = oldKprps[j].getTarget();
Model pasteFromModel = copyFromReactionStep.getModel();
if (ste instanceof SpeciesContext) {
// if newRxnStruct is a feature/membrane, get matching spContexts from old reaction and replace them in new rate expr.
SpeciesContext oldSC = (SpeciesContext) ste;
SpeciesContext newSC = speciesContextHash.get(oldSC);
if (newSC == null) {
// check if oldSc is present in paste-model; if not, add it.
if (!pasteToModel.equals(pasteFromModel)) {
if (pasteToModel.getSpeciesContext(oldSC.getName()) == null) {
// if paste-model has oldSc struct, paste it there,
Structure newSCStruct = pasteToModel.getStructure(oldSC.getStructure().getName());
if (newSCStruct != null) {
newSC = pasteSpecies(parent, oldSC.getSpecies(), null, pasteToModel, newSCStruct, bNew, /*bUseDBSpecies,*/
speciesHash, UserResolvedRxElements.getPreferredReactionElement(userResolvedRxElements, oldSC));
speciesContextHash.put(oldSC, newSC);
} else {
// oldStruct wasn't found in paste-model, paste it in newRxnStruct and add warning to issues list
newSC = pasteSpecies(parent, oldSC.getSpecies(), null, pasteToModel, toRxnStruct, bNew, /*bUseDBSpecies,*/
speciesHash, UserResolvedRxElements.getPreferredReactionElement(userResolvedRxElements, oldSC));
speciesContextHash.put(oldSC, newSC);
Issue issue = new Issue(oldSC, issueContext, IssueCategory.CopyPaste, "SpeciesContext '" + oldSC.getSpecies().getCommonName() + "' was not found in compartment '" + oldSC.getStructure().getName() + "' in the model; the species was added to the compartment '" + toRxnStruct.getName() + "' where the reaction was pasted.", Issue.SEVERITY_WARNING);
issueVector.add(issue);
}
}
}
// if models are the same and newSc is null, then oldSc is not a rxnParticipant. Leave it as is in the expr.
}
if (newSC != null) {
reactionsAndSpeciesContexts.put(newSC, oldSC);
newExpression.substituteInPlace(new Expression(ste.getName()), new Expression(newSC.getName()));
}
// SpeciesContext sc = null;
// Species newSp = model.getSpecies(oldSc.getSpecies().getCommonName());
// if (oldSc.getStructure() == (oldRxnStruct)) {
// sc = model.getSpeciesContext(newSp, newRxnStruct);
// } else {
// if (newRxnStruct instanceof Membrane) {
// // for a membrane, we need to make sure that inside-outside spContexts used are appropriately replaced.
// if (oldSc.getStructure() == ((Membrane)oldRxnStruct).getOutsideFeature()) {
// // old speciesContext is outside (old) membrane, new spContext should be outside new membrane
// sc = model.getSpeciesContext(newSp, ((Membrane)newRxnStruct).getOutsideFeature());
// } else if (oldSc.getStructure() == ((Membrane)oldRxnStruct).getInsideFeature()) {
// // old speciesContext is inside (old) membrane, new spContext should be inside new membrane
// sc = model.getSpeciesContext(newSp, ((Membrane)newRxnStruct).getInsideFeature());
// }
// }
// }
// if (sc != null) {
// newExpression.substituteInPlace(new Expression(ste.getName()), new Expression(sc.getName()));
// }
} else if (ste instanceof StructureSize) {
Structure str = ((StructureSize) ste).getStructure();
// if the structure size used is same as the structure in which the reaction is present, change the structSize to appropriate new struct
if (str.compareEqual(fromRxnStruct)) {
newExpression.substituteInPlace(new Expression(ste.getName()), new Expression(toRxnStruct.getStructureSize().getName()));
} else {
if (fromRxnStruct instanceof Membrane) {
if (str.equals(structTopology.getOutsideFeature((Membrane) fromRxnStruct))) {
newExpression.substituteInPlace(new Expression(ste.getName()), new Expression(structTopology.getOutsideFeature((Membrane) toRxnStruct).getStructureSize().getName()));
} else if (str.equals(structTopology.getInsideFeature((Membrane) fromRxnStruct))) {
newExpression.substituteInPlace(new Expression(ste.getName()), new Expression(structTopology.getInsideFeature((Membrane) toRxnStruct).getStructureSize().getName()));
}
}
}
} else if (ste instanceof MembraneVoltage) {
Membrane membr = ((MembraneVoltage) ste).getMembrane();
// if the MembraneVoltage used is same as that of the membrane in which the reaction is present, change the MemVoltage
if ((fromRxnStruct instanceof Membrane) && (membr.compareEqual(fromRxnStruct))) {
newExpression.substituteInPlace(new Expression(ste.getName()), new Expression(((Membrane) toRxnStruct).getMembraneVoltage().getName()));
}
} else if (ste instanceof ModelParameter) {
// see if model has this global parameter (if rxn is being pasted into another model, it won't)
if (!pasteToModel.equals(pasteFromModel)) {
ModelParameter oldMp = (ModelParameter) ste;
ModelParameter mp = pasteToModel.getModelParameter(oldMp.getName());
boolean bNonNumeric = false;
String newMpName = oldMp.getName();
if (mp != null) {
// new model has a model parameter with same name - are they the same param?
if (!mp.getExpression().equals(oldMp.getExpression())) {
// no, they are not the same param, so mangle the 'ste' name and add as global in the other model
while (pasteToModel.getModelParameter(newMpName) != null) {
newMpName = TokenMangler.getNextEnumeratedToken(newMpName);
}
// if expression if numeric, add it as such. If not, set it to 0.0 and add it as global
Expression exp = oldMp.getExpression();
if (!exp.flatten().isNumeric()) {
exp = new Expression(0.0);
bNonNumeric = true;
}
ModelParameter newMp = pasteToModel.new ModelParameter(newMpName, exp, Model.ROLE_UserDefined, oldMp.getUnitDefinition());
String annotation = "Copied from model : " + pasteFromModel.getNameScope();
newMp.setModelParameterAnnotation(annotation);
pasteToModel.addModelParameter(newMp);
// if global param name had to be changed, make sure newExpr is updated as well.
if (!newMpName.equals(oldMp.getName())) {
newExpression.substituteInPlace(new Expression(oldMp.getName()), new Expression(newMpName));
}
}
} else {
// no global param with same name was found in other model, so add it to other model.
// if expression if numeric, add it as such. If not, set it to 0.0 and add it as global
Expression exp = oldMp.getExpression();
if (!exp.flatten().isNumeric()) {
exp = new Expression(0.0);
bNonNumeric = true;
}
ModelParameter newMp = pasteToModel.new ModelParameter(newMpName, exp, Model.ROLE_UserDefined, oldMp.getUnitDefinition());
String annotation = "Copied from model : " + pasteFromModel.getNameScope();
newMp.setModelParameterAnnotation(annotation);
pasteToModel.addModelParameter(newMp);
}
// if a non-numeric parameter was encountered in the old model, it was added as a numeric (0.0), warn user of change.
if (bNonNumeric) {
Issue issue = new Issue(oldMp, issueContext, IssueCategory.CopyPaste, "Global parameter '" + oldMp.getName() + "' was non-numeric; it has been added " + "as global parameter '" + newMpName + "' in the new model with value = 0.0. " + "Please update its value, if required, before using it.", Issue.SEVERITY_WARNING);
issueVector.add(issue);
}
}
}
}
// end - if newExpr.hasSymbol(ProxyParam)
}
// now if store <param names, new expression> in hashTable
if (paramExprHash.get(oldKps[i].getName()) == null) {
paramExprHash.put(oldKps[i].getName(), newExpression);
}
}
// end for - oldKps (old kinetic parameters)
// use this new expression to generate 'vcml' for the (new) kinetics (easier way to transfer all kinetic parameters)
String newKineticsStr = oldKinetics.writeTokensWithReplacingProxyParams(paramExprHash);
// convert the kinetics 'vcml' to tokens.
CommentStringTokenizer kineticsTokens = new CommentStringTokenizer(newKineticsStr);
// skip the first token;
kineticsTokens.nextToken();
// second token is the kinetic type; use this to create a dummy kinetics
String kineticType = kineticsTokens.nextToken();
Kinetics newkinetics = KineticsDescription.fromVCMLKineticsName(kineticType).createKinetics(newReactionStep);
// use the remaining tokens to construct the new kinetics
newkinetics.fromTokens(newKineticsStr);
// bind newkinetics to newReactionStep and add it to newReactionStep
newkinetics.bind(newReactionStep);
newReactionStep.setKinetics(newkinetics);
counter += 1;
if (counter == copyFromRxSteps.length) {
break;
}
if (!copiedStructName.equals(fromRxnStruct.getName())) {
if (currentStruct instanceof Feature) {
currentStruct = structTopology.getMembrane((Feature) currentStruct);
} else if (currentStruct instanceof Membrane) {
currentStruct = structTopology.getInsideFeature((Membrane) currentStruct);
}
}
copiedStructName = fromRxnStruct.getName();
} while (true);
return new PasteHelper(issueVector, rxPartMapStructure, reactionsAndSpeciesContexts);
}
Aggregations