use of cbit.vcell.math.Variable in project vcell by virtualcell.
the class NetCDFWriter method getVariableSymbols.
private String[] getVariableSymbols(String[] symbols) {
Simulation simulation = simTask.getSimulation();
Vector<String> vars = new Vector<String>();
if (symbols != null) {
for (int i = 0; i < symbols.length; i++) {
Variable v = simulation.getMathDescription().getVariable(symbols[i]);
if ((v != null) && (v instanceof StochVolVariable)) {
vars.add(symbols[i]);
}
}
return vars.toArray(new String[vars.size()]);
} else
return vars.toArray(new String[0]);
}
use of cbit.vcell.math.Variable in project vcell by virtualcell.
the class MathTestingUtilities method getOutwardNormal.
/**
* 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 Function[] getOutwardNormal(Expression analyticSubVolume, String baseName) throws ExpressionException, MappingException, MathException {
VariableHash varHash = new VariableHash();
Expression[] insideOutsideFunctions = getInsideOutsideFunctions(analyticSubVolume);
StringBuffer normalBufferX = new StringBuffer("0.0");
StringBuffer normalBufferY = new StringBuffer("0.0");
StringBuffer normalBufferZ = new StringBuffer("0.0");
for (int i = 0; i < insideOutsideFunctions.length; i++) {
//
// each one gets a turn being the minimum
//
Function[] functions = getOutwardNormalFromInsideOutsideFunction(insideOutsideFunctions[i], baseName + i);
for (int j = 0; j < functions.length; j++) {
varHash.addVariable(functions[j]);
}
String closestName = baseName + i + "_closest";
StringBuffer closestBuffer = new StringBuffer("1.0");
for (int j = 0; j < insideOutsideFunctions.length; j++) {
if (i != j) {
closestBuffer.append(" && (" + baseName + i + "_distance < " + baseName + j + "_distance)");
}
}
Expression closest = new Expression(closestBuffer.toString());
varHash.addVariable(new Function(closestName, closest, null));
normalBufferX.append(" + (" + baseName + i + "_closest * " + baseName + i + "_Nx)");
normalBufferY.append(" + (" + baseName + i + "_closest * " + baseName + i + "_Ny)");
normalBufferZ.append(" + (" + baseName + i + "_closest * " + baseName + i + "_Nz)");
}
varHash.addVariable(new Function(baseName + "_Nx", new Expression(normalBufferX.toString()), null));
varHash.addVariable(new Function(baseName + "_Ny", new Expression(normalBufferY.toString()), null));
varHash.addVariable(new Function(baseName + "_Nz", new Expression(normalBufferZ.toString()), null));
Variable[] vars = varHash.getAlphabeticallyOrderedVariables();
java.util.Vector<Variable> varList = new java.util.Vector<Variable>(java.util.Arrays.asList(vars));
return (Function[]) BeanUtils.getArray(varList, Function.class);
}
use of cbit.vcell.math.Variable in project vcell by virtualcell.
the class MathTestingUtilities method constructOdesForSensitivity.
//
// This method is used to solve for sensitivity of variables to a given parameter.
// The mathDescription and the sensitivity parameter are passed as arguments.
// New variables and ODEs are constructed according to the rule listed below and are added to the mathDescription.
// The method returns the modified mathDescription.
//
public static MathDescription constructOdesForSensitivity(MathDescription mathDesc, Constant sensParam) throws ExpressionException, MathException, MappingException {
//
// For each ODE :
//
// dX/dt = F(X, P)
//
// (where P is the sensitivity parameter)
//
// we create two other ODEs :
//
// dX_1/dt = F(X_1, P_1) and
//
// dX_2/dt = F(X_2, P_2)
//
// where P_1 = P + epsilon, and
// P_2 = P - epsilon.
//
// We keep the initial conditions for both the new ODEs to be the same,
// i.e., X_1_init = X_2_init.
//
// Then, solving for X_1 & X_2, sensitivity of X wrt P can be computed as :
//
// dX = (X_1 - X_2)
// -- ----------- .
// dP (P_1 - P_2)
//
//
// REMOVE PRINTS AFTER CHECKING !!!
System.out.println(" \n\n------------ Old Math Description -----------------");
System.out.println(mathDesc.getVCML_database());
if (mathDesc.getGeometry().getDimension() > 0) {
throw new RuntimeException("Suppport for Spatial systems not yet implemented.");
}
VariableHash varHash = new VariableHash();
Enumeration<Variable> enumVar = mathDesc.getVariables();
while (enumVar.hasMoreElements()) {
varHash.addVariable(enumVar.nextElement());
}
//
// Get 2 values of senstivity parameter (P + epsilon) & (P - epsilon)
//
Constant epsilon = new Constant("epsilon", new Expression(sensParam.getConstantValue() * 1e-3));
Constant sensParam1 = new Constant(sensParam.getName() + "_1", new Expression(sensParam.getConstantValue() + epsilon.getConstantValue()));
Constant sensParam2 = new Constant(sensParam.getName() + "_2", new Expression(sensParam.getConstantValue() - epsilon.getConstantValue()));
//
// Iterate through each subdomain (only 1 in compartmental case), and each equation in the subdomain
//
Enumeration<SubDomain> subDomainEnum = mathDesc.getSubDomains();
//
// Create a vector of equations to store the 2 equations for each ODE variable in the subdomain.
// Later, add it to the equations list in the subdomain.
//
Vector<Equation> equnsVector = new Vector<Equation>();
Vector<Variable> varsVector = new Vector<Variable>();
Vector<Variable> var1s = new Vector<Variable>();
Vector<Variable> var2s = new Vector<Variable>();
while (subDomainEnum.hasMoreElements()) {
SubDomain subDomain = subDomainEnum.nextElement();
Enumeration<Equation> equationEnum = subDomain.getEquations();
Domain domain = new Domain(subDomain);
while (equationEnum.hasMoreElements()) {
Equation equation = equationEnum.nextElement();
if (equation instanceof OdeEquation) {
OdeEquation odeEquation = (OdeEquation) equation;
// Similar to substituteWithExactSolutions, to bind and substitute functions in the ODE
Expression substitutedRateExp = substituteFunctions(odeEquation.getRateExpression(), mathDesc);
String varName = odeEquation.getVariable().getName();
VolVariable var = new VolVariable(varName, domain);
varsVector.addElement(var);
//
// Create the variable var1, and get the initExpr and rateExpr from the original ODE.
// Substitute the new vars (var1 and param1) in the old initExpr and rateExpr and create a new ODE
//
String varName1 = new String("__" + varName + "_1");
Expression initExpr1 = odeEquation.getInitialExpression();
Expression rateExpr1 = new Expression(substitutedRateExp);
rateExpr1.substituteInPlace(new Expression(varName), new Expression(varName1));
rateExpr1.substituteInPlace(new Expression(sensParam.getName()), new Expression(sensParam1.getName()));
VolVariable var1 = new VolVariable(varName1, domain);
var1s.addElement(var1);
OdeEquation odeEqun1 = new OdeEquation(var1, initExpr1, rateExpr1);
equnsVector.addElement(odeEqun1);
//
// Create the variable var2, and get the initExpr and rateExpr from the original ODE.
// Substitute the new vars (var2 and param2) in the old initExpr and rateExpr and create a new ODE
//
String varName2 = new String("__" + varName + "_2");
Expression initExpr2 = odeEquation.getInitialExpression();
Expression rateExpr2 = new Expression(substitutedRateExp);
rateExpr2.substituteInPlace(new Expression(varName), new Expression(varName2));
rateExpr2.substituteInPlace(new Expression(sensParam.getName()), new Expression(sensParam2.getName()));
VolVariable var2 = new VolVariable(varName2, domain);
var2s.addElement(var2);
OdeEquation odeEqun2 = new OdeEquation(var2, initExpr2, rateExpr2);
equnsVector.addElement(odeEqun2);
//
// Create a function for the sensitivity function expression (X1-X2)/(P1-P2), and save in varHash
//
Expression diffVar = Expression.add(new Expression(var1.getName()), Expression.negate(new Expression(var2.getName())));
Expression diffParam = Expression.add(new Expression(sensParam1.getName()), Expression.negate(new Expression(sensParam2.getName())));
Expression sensitivityExpr = Expression.mult(diffVar, Expression.invert(diffParam));
Function sens_Func = new Function("__sens" + varName + "_wrt_" + sensParam.getName(), sensitivityExpr, domain);
varHash.addVariable(epsilon);
varHash.addVariable(sensParam1);
varHash.addVariable(sensParam2);
varHash.addVariable(var1);
varHash.addVariable(var2);
varHash.addVariable(sens_Func);
} else {
// sensitivity not implemented for PDEs or other equation types.
throw new RuntimeException("SolverTest.constructedExactMath(): equation type " + equation.getClass().getName() + " not yet implemented");
}
}
//
// Need to substitute the new variables in the new ODEs.
// i.e., if Rate Expr for ODE_1 for variable X_1 contains variable Y, variable Z, etc.
// Rate Expr is already substituted with X_1, but it also needs substitute Y with Y_1, Z with Z_1, etc.
// So get the volume variables, from the vectors for vars, var_1s and var_2s
// Substitute the rate expressions for the newly added ODEs in equnsVector.
//
Variable[] vars = (Variable[]) BeanUtils.getArray(varsVector, Variable.class);
Variable[] var_1s = (Variable[]) BeanUtils.getArray(var1s, Variable.class);
Variable[] var_2s = (Variable[]) BeanUtils.getArray(var2s, Variable.class);
Vector<Equation> newEqunsVector = new Vector<Equation>();
for (int i = 0; i < equnsVector.size(); i++) {
Equation equn = equnsVector.elementAt(i);
Expression initEx = equn.getInitialExpression();
Expression rateEx = equn.getRateExpression();
for (int j = 0; j < vars.length; j++) {
if (equn.getVariable().getName().endsWith("_1")) {
rateEx.substituteInPlace(new Expression(vars[j].getName()), new Expression(var_1s[j].getName()));
} else if (equn.getVariable().getName().endsWith("_2")) {
rateEx.substituteInPlace(new Expression(vars[j].getName()), new Expression(var_2s[j].getName()));
}
}
OdeEquation odeEqun = new OdeEquation(equn.getVariable(), initEx, rateEx);
newEqunsVector.addElement(odeEqun);
}
//
for (int i = 0; i < newEqunsVector.size(); i++) {
mathDesc.getSubDomain(subDomain.getName()).addEquation((Equation) newEqunsVector.elementAt(i));
}
//
// FAST SYSTEM
// If the subdomain has a fast system, create a new fast system by substituting the high-low variables/parameters
// in the expressions for the fastInvariants and fastRates and adding them to the fast system.
//
Vector<FastInvariant> invarsVector = new Vector<FastInvariant>();
Vector<FastRate> ratesVector = new Vector<FastRate>();
Enumeration<FastInvariant> fastInvarsEnum = null;
Enumeration<FastRate> fastRatesEnum = null;
// Get the fast invariants and fast rates in the system.
FastSystem fastSystem = subDomain.getFastSystem();
if (fastSystem != null) {
fastInvarsEnum = fastSystem.getFastInvariants();
fastRatesEnum = fastSystem.getFastRates();
//
while (fastInvarsEnum.hasMoreElements()) {
FastInvariant fastInvar = fastInvarsEnum.nextElement();
Expression fastInvarExpr = fastInvar.getFunction();
fastInvarExpr = MathUtilities.substituteFunctions(fastInvarExpr, mathDesc);
Expression fastInvarExpr1 = new Expression(fastInvarExpr);
Expression fastInvarExpr2 = new Expression(fastInvarExpr);
for (int i = 0; i < vars.length; i++) {
fastInvarExpr1.substituteInPlace(new Expression(vars[i].getName()), new Expression(var_1s[i].getName()));
fastInvarExpr2.substituteInPlace(new Expression(vars[i].getName()), new Expression(var_2s[i].getName()));
}
fastInvarExpr1.substituteInPlace(new Expression(sensParam.getName()), new Expression(sensParam1.getName()));
FastInvariant fastInvar1 = new FastInvariant(fastInvarExpr1);
invarsVector.addElement(fastInvar1);
fastInvarExpr2.substituteInPlace(new Expression(sensParam.getName()), new Expression(sensParam2.getName()));
FastInvariant fastInvar2 = new FastInvariant(fastInvarExpr2);
invarsVector.addElement(fastInvar2);
}
// Add the newly created fast invariants to the existing list of fast invariants in the fast system.
for (int i = 0; i < invarsVector.size(); i++) {
FastInvariant inVar = (FastInvariant) invarsVector.elementAt(i);
fastSystem.addFastInvariant(inVar);
}
//
while (fastRatesEnum.hasMoreElements()) {
FastRate fastRate = fastRatesEnum.nextElement();
Expression fastRateExpr = fastRate.getFunction();
fastRateExpr = MathUtilities.substituteFunctions(fastRateExpr, mathDesc);
Expression fastRateExpr1 = new Expression(fastRateExpr);
Expression fastRateExpr2 = new Expression(fastRateExpr);
for (int i = 0; i < vars.length; i++) {
fastRateExpr1.substituteInPlace(new Expression(vars[i].getName()), new Expression(var_1s[i].getName()));
fastRateExpr2.substituteInPlace(new Expression(vars[i].getName()), new Expression(var_2s[i].getName()));
}
fastRateExpr1.substituteInPlace(new Expression(sensParam.getName()), new Expression(sensParam1.getName()));
FastRate fastRate1 = new FastRate(fastRateExpr1);
ratesVector.addElement(fastRate1);
fastRateExpr2.substituteInPlace(new Expression(sensParam.getName()), new Expression(sensParam2.getName()));
FastRate fastRate2 = new FastRate(fastRateExpr2);
ratesVector.addElement(fastRate2);
}
// Add the newly created fast rates to the existing list of fast rates in the fast system.
for (int i = 0; i < ratesVector.size(); i++) {
FastRate rate = (FastRate) ratesVector.elementAt(i);
fastSystem.addFastRate(rate);
}
}
}
// Reset all variables in mathDesc.
mathDesc.setAllVariables(varHash.getAlphabeticallyOrderedVariables());
// REMOVE PRINTS AFTER CHECKING
System.out.println(" \n\n------------ New Math Description -----------------");
System.out.println(mathDesc.getVCML_database());
return mathDesc;
}
use of cbit.vcell.math.Variable in project vcell by virtualcell.
the class ConstructedSolutionTemplate method initialize.
/**
* Insert the method's description here.
* Creation date: (5/12/2003 11:04:13 AM)
* @param mathDesc cbit.vcell.math.MathDescription
*/
private void initialize(cbit.vcell.math.MathDescription mathDesc) {
//
// for all valid combinations of variables/domains ... add a solution template with a default solution.
//
java.util.Vector solutionTemplateList = new java.util.Vector();
Enumeration enumSubDomains = mathDesc.getSubDomains();
while (enumSubDomains.hasMoreElements()) {
SubDomain subDomain = (SubDomain) enumSubDomains.nextElement();
Enumeration enumEquations = subDomain.getEquations();
while (enumEquations.hasMoreElements()) {
Equation equation = (Equation) enumEquations.nextElement();
Variable var = equation.getVariable();
String baseName = " " + var.getName() + "_" + subDomain.getName();
String amplitudeName = baseName + "_A";
String tau1Name = baseName + "_tau1";
String tau2Name = baseName + "_tau2";
if (equation instanceof OdeEquation) {
try {
Expression exp = new Expression(amplitudeName + " * (1.0 + exp(-t/" + tau1Name + ")*sin(2*" + Math.PI + "/" + tau2Name + "*t))");
solutionTemplateList.add(new SolutionTemplate(equation.getVariable().getName(), subDomain.getName(), exp));
} catch (cbit.vcell.parser.ExpressionException e) {
e.printStackTrace(System.out);
throw new RuntimeException(e.getMessage());
}
} else if (equation instanceof PdeEquation) {
try {
Expression exp = new Expression(amplitudeName + " * (1.0 + exp(-t/" + tau1Name + ") + " + tau2Name + "*x)");
solutionTemplateList.add(new SolutionTemplate(equation.getVariable().getName(), subDomain.getName(), exp));
} catch (cbit.vcell.parser.ExpressionException e) {
e.printStackTrace(System.out);
throw new RuntimeException(e.getMessage());
}
}
}
}
this.solutionTemplates = (SolutionTemplate[]) org.vcell.util.BeanUtils.getArray(solutionTemplateList, SolutionTemplate.class);
}
use of cbit.vcell.math.Variable in project vcell by virtualcell.
the class ModelOptimizationMapping method getRemappedReferenceData.
/**
* Gets the constraintData property (cbit.vcell.opt.ConstraintData) value.
* @return The constraintData property value.
* @see #setConstraintData
*/
private ReferenceData getRemappedReferenceData(MathMapping mathMapping) throws MappingException {
if (modelOptimizationSpec.getReferenceData() == null) {
return null;
}
//
// make sure time is mapped
//
ReferenceData refData = modelOptimizationSpec.getReferenceData();
ReferenceDataMappingSpec[] refDataMappingSpecs = modelOptimizationSpec.getReferenceDataMappingSpecs();
RowColumnResultSet rowColResultSet = new RowColumnResultSet();
Vector<SymbolTableEntry> modelObjectList = new Vector<SymbolTableEntry>();
Vector<double[]> dataList = new Vector<double[]>();
//
// find bound columns, (time is always mapped to the first column)
//
int mappedColumnCount = 0;
for (int i = 0; i < refDataMappingSpecs.length; i++) {
SymbolTableEntry modelObject = refDataMappingSpecs[i].getModelObject();
if (modelObject != null) {
int mappedColumnIndex = mappedColumnCount;
if (modelObject instanceof Model.ReservedSymbol && ((ReservedSymbol) modelObject).isTime()) {
mappedColumnIndex = 0;
}
String origRefDataColumnName = refDataMappingSpecs[i].getReferenceDataColumnName();
int origRefDataColumnIndex = refData.findColumn(origRefDataColumnName);
if (origRefDataColumnIndex < 0) {
throw new RuntimeException("reference data column named '" + origRefDataColumnName + "' not found");
}
double[] columnData = refData.getDataByColumn(origRefDataColumnIndex);
if (modelObjectList.contains(modelObject)) {
throw new RuntimeException("multiple reference data columns mapped to same model object '" + modelObject.getName() + "'");
}
modelObjectList.insertElementAt(modelObject, mappedColumnIndex);
dataList.insertElementAt(columnData, mappedColumnIndex);
mappedColumnCount++;
}
}
//
if (modelObjectList.size() == 0) {
throw new RuntimeException("reference data was not associated with model");
}
if (modelObjectList.size() == 1) {
throw new RuntimeException("reference data was not associated with model, must map time and at least one other column");
}
boolean bFoundTimeVar = false;
for (SymbolTableEntry ste : modelObjectList) {
if (ste instanceof Model.ReservedSymbol && ((ReservedSymbol) ste).isTime()) {
bFoundTimeVar = true;
break;
}
}
if (!bFoundTimeVar) {
throw new RuntimeException("must map time column of reference data to model");
}
//
for (int i = 0; i < modelObjectList.size(); i++) {
SymbolTableEntry modelObject = (SymbolTableEntry) modelObjectList.elementAt(i);
try {
// Find by name because MathSybolMapping has different 'objects' than refDataMapping 'objects'
Variable variable = mathMapping.getMathSymbolMapping().findVariableByName(modelObject.getName());
if (variable != null) {
String symbol = variable.getName();
rowColResultSet.addDataColumn(new ODESolverResultSetColumnDescription(symbol));
} else if (modelObject instanceof Model.ReservedSymbol && ((Model.ReservedSymbol) modelObject).isTime()) {
Model.ReservedSymbol time = (Model.ReservedSymbol) modelObject;
String symbol = time.getName();
rowColResultSet.addDataColumn(new ODESolverResultSetColumnDescription(symbol));
}
} catch (MathException | MatrixException | ExpressionException | ModelException e) {
e.printStackTrace();
throw new MappingException(e.getMessage(), e);
}
}
//
// populate data columns (time and rest)
//
double[] weights = new double[rowColResultSet.getColumnDescriptionsCount()];
weights[0] = 1.0;
int numRows = ((double[]) dataList.elementAt(0)).length;
int numColumns = modelObjectList.size();
for (int j = 0; j < numRows; j++) {
double[] row = new double[numColumns];
for (int i = 0; i < numColumns; i++) {
row[i] = ((double[]) dataList.elementAt(i))[j];
if (i > 0) {
weights[i] += row[i] * row[i];
}
}
rowColResultSet.addRow(row);
}
for (int i = 0; i < numColumns; i++) {
if (weights[i] == 0) {
weights[i] = 1;
} else {
weights[i] = 1 / weights[i];
}
}
SimpleReferenceData remappedRefData = new SimpleReferenceData(rowColResultSet, weights);
return remappedRefData;
}
Aggregations