use of cbit.vcell.solver.UniformOutputTimeSpec in project vcell by virtualcell.
the class MovingBoundaryFileWriter method getOutputTimeStep.
private Element getOutputTimeStep(SolverTaskDescription std) {
Objects.requireNonNull(std);
if (std.getOutputTimeSpec().isUniform()) {
Element e = new Element(MBTags.outputTimeStep);
double outputTimeStep = ((UniformOutputTimeSpec) std.getOutputTimeSpec()).getOutputTimeStep();
e.setText(outputTimeStep + "");
return e;
}
return null;
}
use of cbit.vcell.solver.UniformOutputTimeSpec in project vcell by virtualcell.
the class StochtestRunService method runOne.
public void runOne() throws IllegalArgumentException, SQLException, DataAccessException, XmlParseException, PropertyVetoException, ExpressionException, MappingException, GeometryException, ImageException, IOException {
StochtestRun stochtestRun = StochtestDbUtils.acceptNextWaitingStochtestRun(conFactory);
String biomodelXML = null;
if (stochtestRun != null) {
String networkGenProbs = null;
try {
User user = new User(PropertyLoader.ADMINISTRATOR_ACCOUNT, new KeyValue(PropertyLoader.ADMINISTRATOR_ID));
ServerDocumentManager serverDocumentManager = new ServerDocumentManager(this.dbServerImpl);
biomodelXML = serverDocumentManager.getBioModelXML(new QueryHashtable(), user, stochtestRun.stochtest.biomodelRef, true);
BioModel bioModel = XmlHelper.XMLToBioModel(new XMLSource(biomodelXML));
bioModel.refreshDependencies();
SimulationContext srcSimContext = null;
for (SimulationContext sc : bioModel.getSimulationContexts()) {
if (sc.getKey().equals(stochtestRun.stochtest.simContextRef)) {
srcSimContext = sc;
}
}
if (srcSimContext == null) {
throw new RuntimeException("cannot find simcontext with key=" + stochtestRun.stochtest.simContextRef);
}
//
for (SpeciesContextSpec scs : srcSimContext.getReactionContext().getSpeciesContextSpecs()) {
scs.setConstant(false);
}
SimulationContext simContext = srcSimContext;
StochtestMathType parentMathType = stochtestRun.parentMathType;
StochtestMathType mathType = stochtestRun.mathType;
if (parentMathType != mathType) {
if (parentMathType == StochtestMathType.nonspatialstochastic && mathType == StochtestMathType.rules) {
simContext = SimulationContext.copySimulationContext(srcSimContext, "generatedRules", false, Application.RULE_BASED_STOCHASTIC);
} else if (parentMathType == StochtestMathType.rules && mathType == StochtestMathType.nonspatialstochastic) {
simContext = SimulationContext.copySimulationContext(srcSimContext, "generatedSSA", false, Application.NETWORK_STOCHASTIC);
} else {
throw new RuntimeException("unexpected copy of simcontext from " + parentMathType + " to " + mathType);
}
bioModel.addSimulationContext(simContext);
}
MathMappingCallback mathMappingCallback = new MathMappingCallback() {
@Override
public void setProgressFraction(float fractionDone) {
}
@Override
public void setMessage(String message) {
}
@Override
public boolean isInterrupted() {
return false;
}
};
MathMapping mathMapping = simContext.createNewMathMapping(mathMappingCallback, NetworkGenerationRequirements.ComputeFullStandardTimeout);
MathDescription mathDesc = mathMapping.getMathDescription(mathMappingCallback);
simContext.setMathDescription(mathDesc);
if (simContext.isInsufficientIterations()) {
networkGenProbs = "insufficientIterations";
} else if (simContext.isInsufficientMaxMolecules()) {
networkGenProbs = "insufficientMaxMolecules";
}
File baseDirectory = StochtestFileUtils.createDirFile(baseDir, stochtestRun);
try {
OutputTimeSpec outputTimeSpec = new UniformOutputTimeSpec(0.5);
double endTime = 10.0;
computeTrials(simContext, stochtestRun, baseDirectory, outputTimeSpec, endTime, numTrials);
StochtestDbUtils.finalizeAcceptedStochtestRun(conFactory, stochtestRun, StochtestRun.StochtestRunStatus.complete, null, networkGenProbs);
} finally {
StochtestFileUtils.clearDir(baseDirectory);
}
} catch (Exception e) {
StochtestDbUtils.finalizeAcceptedStochtestRun(conFactory, stochtestRun, StochtestRun.StochtestRunStatus.failed, e.getMessage(), networkGenProbs);
//
if (biomodelXML != null) {
XmlUtil.writeXMLStringToFile(biomodelXML, new File(baseDir, "stochtestrun_" + stochtestRun.stochtest.key + ".vcml").getPath(), false);
}
//
// write exception trace to .txt file
//
StringWriter stringWriter = new StringWriter();
PrintWriter printWriter = new PrintWriter(stringWriter);
e.printStackTrace(printWriter);
printWriter.flush();
System.out.println(stringWriter.getBuffer().toString());
XmlUtil.writeXMLStringToFile(stringWriter.getBuffer().toString(), new File(baseDir, "stochtestrun_" + stochtestRun.stochtest.key + "_error.txt").getPath(), false);
}
} else {
System.out.println("no jobs waiting");
try {
Thread.sleep(5000);
} catch (InterruptedException e) {
// TODO Auto-generated catch block
e.printStackTrace();
}
}
}
use of cbit.vcell.solver.UniformOutputTimeSpec in project vcell by virtualcell.
the class RuleBasedTest method checkNonspatialStochasticSimContext.
private static void checkNonspatialStochasticSimContext(SimulationContext srcSimContext, File baseDirectory, int numTrials) throws Exception {
if (!srcSimContext.getApplicationType().equals(Application.NETWORK_STOCHASTIC) || srcSimContext.getGeometry().getDimension() != 0) {
throw new RuntimeException("simContext is of type " + srcSimContext.getApplicationType() + " and geometry dimension of " + srcSimContext.getGeometry().getDimension() + ", expecting nonspatial stochastic");
}
BioModel origBioModel = srcSimContext.getBioModel();
BioModel bioModel = XmlHelper.XMLToBioModel(new XMLSource(XmlHelper.bioModelToXML(origBioModel)));
bioModel.refreshDependencies();
// create ODE and RuleBased
SimulationContext newODEApp = SimulationContext.copySimulationContext(srcSimContext, "aUniqueNewODEApp", false, Application.NETWORK_DETERMINISTIC);
SimulationContext newRuleBasedApp = SimulationContext.copySimulationContext(srcSimContext, "aUniqueNewRuleBasedApp", false, Application.RULE_BASED_STOCHASTIC);
newODEApp.setBioModel(bioModel);
newRuleBasedApp.setBioModel(bioModel);
ArrayList<AnnotatedFunction> outputFunctionsList = srcSimContext.getOutputFunctionContext().getOutputFunctionsList();
// OutputContext outputContext = new OutputContext(outputFunctionsList.toArray(new AnnotatedFunction[outputFunctionsList.size()]));
newODEApp.getOutputFunctionContext().setOutputFunctions(outputFunctionsList);
newRuleBasedApp.getOutputFunctionContext().setOutputFunctions(outputFunctionsList);
NetworkGenerationRequirements networkGenerationRequirements = NetworkGenerationRequirements.AllowTruncatedStandardTimeout;
bioModel.addSimulationContext(newODEApp);
newODEApp.refreshMathDescription(new MathMappingCallbackTaskAdapter(null), networkGenerationRequirements);
bioModel.addSimulationContext(newRuleBasedApp);
newRuleBasedApp.refreshMathDescription(new MathMappingCallbackTaskAdapter(null), networkGenerationRequirements);
srcSimContext.refreshMathDescription(new MathMappingCallbackTaskAdapter(null), networkGenerationRequirements);
// Create non-spatialStoch, ODE and RuleBased sims
Simulation nonspatialStochAppNewSim = srcSimContext.addNewSimulation(STOCH_SIM_NAME, /*SimulationOwner.DEFAULT_SIM_NAME_PREFIX*/
new MathMappingCallbackTaskAdapter(null), networkGenerationRequirements);
Simulation newODEAppNewSim = newODEApp.addNewSimulation(ODE_SIM_NAME, new MathMappingCallbackTaskAdapter(null), networkGenerationRequirements);
Simulation newRuleBasedAppNewSim = newRuleBasedApp.addNewSimulation(NFS_SIM_NAME, new MathMappingCallbackTaskAdapter(null), networkGenerationRequirements);
nonspatialStochAppNewSim.setSimulationOwner(srcSimContext);
newODEAppNewSim.setSimulationOwner(newODEApp);
newRuleBasedAppNewSim.setSimulationOwner(newRuleBasedApp);
try {
bioModel.getModel().getSpeciesContexts();
ArrayList<String> varNameList = new ArrayList<String>();
for (SpeciesContextSpec scs : srcSimContext.getReactionContext().getSpeciesContextSpecs()) {
varNameList.add(scs.getSpeciesContext().getName());
}
String[] varNames = varNameList.toArray(new String[0]);
OutputTimeSpec outputTimeSpec = nonspatialStochAppNewSim.getSolverTaskDescription().getOutputTimeSpec();
ArrayList<Double> sampleTimeList = new ArrayList<Double>();
if (outputTimeSpec instanceof UniformOutputTimeSpec) {
double endingTime = nonspatialStochAppNewSim.getSolverTaskDescription().getTimeBounds().getEndingTime();
double dT = ((UniformOutputTimeSpec) outputTimeSpec).getOutputTimeStep();
int currTimeIndex = 0;
while (currTimeIndex * dT <= (endingTime + 1e-8)) {
sampleTimeList.add(currTimeIndex * dT);
currTimeIndex++;
}
}
double[] sampleTimes = new double[sampleTimeList.size()];
for (int i = 0; i < sampleTimes.length; i++) {
sampleTimes[i] = sampleTimeList.get(i);
}
TimeSeriesMultitrialData sampleDataStoch1 = new TimeSeriesMultitrialData("stochastic1", varNames, sampleTimes, numTrials);
TimeSeriesMultitrialData sampleDataStoch2 = new TimeSeriesMultitrialData("stochastic2", varNames, sampleTimes, numTrials);
TimeSeriesMultitrialData sampleDataDeterministic = new TimeSeriesMultitrialData("determinstic", varNames, sampleTimes, 1);
runsolver(nonspatialStochAppNewSim, baseDirectory, numTrials, sampleDataStoch1);
runsolver(newODEAppNewSim, baseDirectory, 1, sampleDataDeterministic);
runsolver(newRuleBasedAppNewSim, baseDirectory, numTrials, sampleDataStoch2);
writeVarDiffData(baseDirectory, sampleDataStoch1, sampleDataStoch2);
writeKolmogorovSmirnovTest(baseDirectory, sampleDataStoch1, sampleDataStoch2);
writeChiSquareTest(baseDirectory, sampleDataStoch1, sampleDataStoch2);
writeData(baseDirectory, sampleDataStoch1);
writeData(baseDirectory, sampleDataStoch2);
writeData(baseDirectory, sampleDataDeterministic);
} finally {
srcSimContext.removeSimulation(nonspatialStochAppNewSim);
newODEApp.removeSimulation(newODEAppNewSim);
newRuleBasedApp.removeSimulation(newRuleBasedAppNewSim);
}
}
use of cbit.vcell.solver.UniformOutputTimeSpec in project vcell by virtualcell.
the class SbmlVcmlConverter method main.
/**
* @param args -import or -export
*/
public static void main(String[] args) {
Logging.init();
ResourceUtil.setNativeLibraryDirectory();
if (args.length < 2 || args.length > 3) {
System.out.println("Usage:\n\t -export path_of_input_file\n\tOR\n\t -import path_of_input_file [-simulate]");
System.exit(1);
}
if (args.length > 2 && args[0].equals("-export")) {
System.out.println("Export cannot have arguments other than input file.");
System.out.println("Usage:\n\t -export path_of_input_file\n\tOR\n\t-import path_of_input_file [-simulate]");
System.exit(1);
}
try {
String pathName = args[1];
// Read xml file (Sbml or Vcml) into stringbuffer, pass the string into VCell's importer/exporter
String xmlString = getXMLString(pathName);
if (args[0].equals("-import")) {
if (args.length > 3 || (args.length == 3 && !args[2].equals("-simulate"))) {
System.out.println(args[2] + " : Unknown arguments for import; please check and re-run import command");
System.out.println("Usage:\n\t -import path_of_input_file [-simulate]");
System.exit(1);
}
boolean bSimulate = false;
if (args.length == 3 && args[2].equals("-simulate")) {
bSimulate = true;
}
// Create a default VCLogger - SBMLImporter needs it
cbit.util.xml.VCLogger logger = new cbit.util.xml.VCLogger() {
@Override
public void sendMessage(Priority p, ErrorType et, String message) {
System.err.println("LOGGER: msgLevel=" + p + ", msgType=" + et + ", " + message);
if (p == VCLogger.Priority.HighPriority) {
throw new RuntimeException("Import failed : " + message);
}
}
public void sendAllMessages() {
}
public boolean hasMessages() {
return false;
}
};
// invoke SBMLImporter, which returns a Biomodel, convert that to VCML using XMLHelper
try {
// import SBML model into VCML, store VCML string in 'fileName.vcml'
// Instantiate an SBMLImporter to get the speciesUnitsHash - to compute the conversion factor from VC->SB species units.
// and import SBML (sbml->bioModel)
org.vcell.sbml.vcell.SBMLImporter sbmlImporter = new org.vcell.sbml.vcell.SBMLImporter(pathName, logger, false);
BioModel bioModel = sbmlImporter.getBioModel();
// Hashtable<String, SBMLImporter.SBVCConcentrationUnits> speciesUnitsHash = sbmlImporter.getSpeciesUnitsHash();
// double timeFactor = sbmlImporter.getSBMLTimeUnitsFactor();
// convert biomodel to vcml and save to file.
String vcmlString = XmlHelper.bioModelToXML(bioModel);
String fileExtensionStr;
if (pathName.endsWith(".xml")) {
fileExtensionStr = ".xml";
} else if (pathName.endsWith(".sbml")) {
fileExtensionStr = ".sbml";
} else {
throw new RuntimeException("Unknown file extension for SBML file name; Exiting SbmlConverter.");
}
String vcmlFileName = pathName.replace(fileExtensionStr, ".vcml");
File vcmlFile = new File(vcmlFileName);
XmlUtil.writeXMLStringToFile(vcmlString, vcmlFile.getAbsolutePath(), true);
// If user doesn't choose to simulate, you are done.
if (!bSimulate) {
return;
}
// Generate math for lone simContext
SimulationContext simContext = (SimulationContext) bioModel.getSimulationContext(0);
MathDescription mathDesc = simContext.createNewMathMapping().getMathDescription();
simContext.setMathDescription(mathDesc);
// Create basic simulation, with IDA solver (set in solve method) and other defaults, and end time 'Te'
org.vcell.util.document.SimulationVersion simVersion = new org.vcell.util.document.SimulationVersion(new KeyValue("12345"), "name", new org.vcell.util.document.User("user", new KeyValue("123")), new org.vcell.util.document.GroupAccessNone(), // versionBranchPointRef
null, // branchID
new java.math.BigDecimal(1.0), new java.util.Date(), org.vcell.util.document.VersionFlag.Archived, "", null);
Simulation sim1 = new Simulation(simVersion, simContext.getMathDescription());
simContext.addSimulation(sim1);
sim1.setName("sim1");
// using a default end time of 10.0 secs, and a forcing default output time step of 0.01.
// If user wants anything different, user can modify the .idaInput file and re-run simulation separately.
// double newEndTime = endTime * timeFactor;
double newEndTime = endTime;
sim1.getSolverTaskDescription().setTimeBounds(new cbit.vcell.solver.TimeBounds(0, newEndTime));
TimeStep timeStep_1 = new TimeStep();
sim1.getSolverTaskDescription().setTimeStep(new TimeStep(timeStep_1.getMinimumTimeStep(), timeStep_1.getDefaultTimeStep(), newEndTime / 10000));
sim1.getSolverTaskDescription().setOutputTimeSpec(new UniformOutputTimeSpec((newEndTime - 0) / numTimeSteps));
sim1.getSolverTaskDescription().setErrorTolerance(new ErrorTolerance(1e-10, 1e-12));
// save the new vcml with generated math and new sim in file named : fileName_IDA_simulation.vcml
String newVcmlString = XmlHelper.bioModelToXML(bioModel);
String newVcmlFileName = vcmlFile.getPath().replace(".vcml", "_IDA_simulation.vcml");
File newVcmlFile = new File(newVcmlFileName);
XmlUtil.writeXMLStringToFile(newVcmlString, newVcmlFile.getAbsolutePath(), true);
// Solve simulation, which also generates and saves the .idainput file and .csv
// SimSpec to get vars to solve and output in .csv
SimSpec simSpec = SimSpec.fromSBML(xmlString);
solveSimulation(new SimulationJob(sim1, 0, null), vcmlFile.getPath(), simSpec);
} catch (Exception e) {
e.printStackTrace(System.err);
// and an annotation saying why it failed.
try {
String fileExtensionStr;
if (pathName.endsWith(".xml")) {
fileExtensionStr = ".xml";
} else if (pathName.endsWith(".sbml")) {
fileExtensionStr = ".sbml";
} else {
throw new RuntimeException("Unknown file extension for SBML file name; Exiting SbmlConverter.");
}
String dummyVcmlFileName = pathName.replace(fileExtensionStr, ".vcml");
File dummyVcmlFile = new File(dummyVcmlFileName);
BioModel dummy_biomodel = new BioModel(null);
String dummyName = dummyVcmlFile.getName().substring(0, dummyVcmlFile.getName().indexOf(".vcml"));
dummy_biomodel.setName(dummyName);
dummy_biomodel.setDescription("SBML Model could not be automatically converted to VCML : " + e.getMessage());
String vcmlString = XmlHelper.bioModelToXML(dummy_biomodel);
XmlUtil.writeXMLStringToFile(vcmlString, dummyVcmlFile.getAbsolutePath(), true);
} catch (Exception e1) {
e.printStackTrace(System.err);
}
}
} else if (args[0].equals("-export")) {
try {
BioModel bioModel = XmlHelper.XMLToBioModel(new XMLSource(xmlString));
for (int i = 0; i < bioModel.getSimulationContexts().length; i++) {
SimulationContext simContext = bioModel.getSimulationContext(i);
if (simContext.getGeometry().getDimension() == 0 && !simContext.isStoch()) {
if (!simContext.getGeometryContext().isAllSizeSpecifiedPositive()) {
Structure structure = simContext.getModel().getStructure(0);
double structureSize = 1.0;
StructureMapping structMapping = simContext.getGeometryContext().getStructureMapping(structure);
StructureSizeSolver.updateAbsoluteStructureSizes(simContext, structure, structureSize, structMapping.getSizeParameter().getUnitDefinition());
}
// Export the application itself, with default overrides
String sbmlString = XmlHelper.exportSBML(bioModel, 2, 4, 0, false, simContext, null);
String filePath = pathName.substring(0, pathName.lastIndexOf("\\") + 1);
String sbmlFileName = TokenMangler.mangleToSName(bioModel.getName() + "_" + i);
File sbmlFile = new File(filePath + sbmlFileName + ".xml");
XmlUtil.writeXMLStringToFile(sbmlString, sbmlFile.getAbsolutePath(), true);
// Now export each simulation in the application that has overrides
Simulation[] simulations = bioModel.getSimulations(simContext);
for (int j = 0; j < simulations.length; j++) {
if (simulations[j].getMathOverrides().hasOverrides()) {
// Check for parameter scan and create simJob to pass into exporter
for (int k = 0; k < simulations[j].getScanCount(); k++) {
SimulationJob simJob = new SimulationJob(simulations[j], k, null);
sbmlString = XmlHelper.exportSBML(bioModel, 2, 4, 0, false, simContext, simJob);
String fileName = TokenMangler.mangleToSName(sbmlFileName + "_" + j + "_" + k);
sbmlFile = new File(filePath + fileName + ".xml");
XmlUtil.writeXMLStringToFile(sbmlString, sbmlFile.getAbsolutePath(), true);
}
}
}
}
}
System.out.println("Successfully translated model : " + pathName);
} catch (Exception e) {
e.printStackTrace(System.err);
}
}
System.exit(0);
} catch (IOException e) {
e.printStackTrace(System.err);
System.exit(1);
}
}
use of cbit.vcell.solver.UniformOutputTimeSpec in project vcell by virtualcell.
the class VCellSBMLSolver method solve.
public File solve(String filePrefix, File outDir, String sbmlFileName, SimSpec testSpec) throws IOException, SolverException, SbmlException {
try {
cbit.util.xml.VCLogger sbmlImportLogger = new LocalLogger();
//
// Instantiate an SBMLImporter to get the speciesUnitsHash - to compute the conversion factor from VC->SB species units.
// and import SBML (sbml->bioModel)
BioModel bioModel = importSBML(sbmlFileName, sbmlImportLogger, false);
// Hashtable<String, SBMLImporter.SBVCConcentrationUnits> speciesUnitsHash = sbmlImporter.getSpeciesUnitsHash();
// double timeFactor = sbmlImporter.getSBMLTimeUnitsFactor();
String vcml_1 = XmlHelper.bioModelToXML(bioModel);
SBMLUtils.writeStringToFile(vcml_1, new File(outDir, filePrefix + ".vcml").getAbsolutePath(), true);
if (bRoundTrip) {
// Round trip the bioModel (bioModel->sbml->bioModel).
// save imported "bioModel" as VCML
// String vcml_1 = XmlHelper.bioModelToXML(bioModel);
// SBMLUtils.writeStringToFile(vcml_1, new File(outDir,filePrefix+".vcml").getAbsolutePath());
// export bioModel as sbml and save
// String vcml_sbml = cbit.vcell.xml.XmlHelper.exportSBML(bioModel, 2, 1, bioModel.getSimulationContexts(0).getName());
// SimulationJob simJob = new SimulationJob(bioModel.getSimulations(bioModel.getSimulationContexts(0))[0], null, 0);
String vcml_sbml = cbit.vcell.xml.XmlHelper.exportSBML(bioModel, 2, 1, 0, false, bioModel.getSimulationContext(0), null);
SBMLUtils.writeStringToFile(vcml_sbml, new File(outDir, filePrefix + ".vcml.sbml").getAbsolutePath(), true);
// re-import bioModel from exported sbml
XMLSource vcml_sbml_Src = new XMLSource(vcml_sbml);
BioModel newBioModel = (BioModel) XmlHelper.importSBML(sbmlImportLogger, vcml_sbml_Src, false);
String vcml_sbml_vcml = XmlHelper.bioModelToXML(newBioModel);
SBMLUtils.writeStringToFile(vcml_sbml_vcml, new File(outDir, filePrefix + ".vcml.sbml.vcml").getAbsolutePath(), true);
// have rest of code use the round-tripped biomodel
bioModel = newBioModel;
}
//
// select only Application, generate math, and create a single Simulation.
//
SimulationContext simContext = bioModel.getSimulationContext(0);
MathMapping mathMapping = simContext.createNewMathMapping();
MathDescription mathDesc = mathMapping.getMathDescription();
String vcml = mathDesc.getVCML();
try (PrintWriter pw = new PrintWriter("vcmlTrace.txt")) {
pw.println(vcml);
}
simContext.setMathDescription(mathDesc);
SimulationVersion simVersion = new SimulationVersion(new KeyValue("100"), "unnamed", null, null, null, null, null, null, null, null);
Simulation sim = new Simulation(simVersion, mathDesc);
sim.setName("unnamed");
// if time factor from SBML is not 1 (i.e., it is not in secs but in minutes or hours), convert endTime to min/hr as : endTime*timeFactor
// double endTime = testSpec.getEndTime()*timeFactor;
double endTime = testSpec.getEndTime();
sim.getSolverTaskDescription().setTimeBounds(new TimeBounds(0, endTime));
TimeStep timeStep = new TimeStep();
sim.getSolverTaskDescription().setTimeStep(new TimeStep(timeStep.getMinimumTimeStep(), timeStep.getDefaultTimeStep(), endTime / 10000));
sim.getSolverTaskDescription().setOutputTimeSpec(new UniformOutputTimeSpec((endTime - 0) / testSpec.getNumTimeSteps()));
sim.getSolverTaskDescription().setErrorTolerance(new ErrorTolerance(1e-10, 1e-12));
// sim.getSolverTaskDescription().setErrorTolerance(new cbit.vcell.solver.ErrorTolerance(1e-10, 1e-12));
// Generate .idaInput string
/* IDAFileWriter idaFileWriter = new IDAFileWriter(sim);
File idaInputFile = new File(filePathName.replace(".vcml", ".idaInput"));
PrintWriter idaPW = new java.io.PrintWriter(idaInputFile);
idaFileWriter.writeInputFile(idaPW);
idaPW.close();
// use the idastandalone solver
File idaOutputFile = new File(filePathName.replace(".vcml", ".ida"));
Executable executable = new Executable("IDAStandalone " + idaInputFile + " " + idaOutputFile);
executable.start();
*/
// Generate .cvodeInput string
File cvodeFile = new File(outDir, filePrefix + SimDataConstants.CVODEINPUT_DATA_EXTENSION);
PrintWriter cvodePW = new java.io.PrintWriter(cvodeFile);
SimulationJob simJob = new SimulationJob(sim, 0, null);
SimulationTask simTask = new SimulationTask(simJob, 0);
CVodeFileWriter cvodeFileWriter = new CVodeFileWriter(cvodePW, simTask);
cvodeFileWriter.write();
cvodePW.close();
// use the cvodeStandalone solver
File cvodeOutputFile = new File(outDir, filePrefix + SimDataConstants.IDA_DATA_EXTENSION);
String executableName = null;
try {
executableName = SolverUtilities.getExes(SolverDescription.CVODE)[0].getAbsolutePath();
} catch (IOException e) {
throw new RuntimeException("failed to get executable for solver " + SolverDescription.CVODE.getDisplayLabel() + ": " + e.getMessage(), e);
}
Executable executable = new Executable(new String[] { executableName, cvodeFile.getAbsolutePath(), cvodeOutputFile.getAbsolutePath() });
executable.start();
// get the result
ODESolverResultSet odeSolverResultSet = getODESolverResultSet(simJob, cvodeOutputFile.getPath());
//
// print header
//
File outputFile = new File(outDir, filePrefix + ".vcell.csv");
java.io.PrintStream outputStream = new java.io.PrintStream(new java.io.BufferedOutputStream(new java.io.FileOutputStream(outputFile)));
outputStream.print("time");
for (int i = 0; i < testSpec.getVarsList().length; i++) {
outputStream.print("," + testSpec.getVarsList()[i]);
}
outputStream.println();
//
// extract data for time and species
//
double[][] data = new double[testSpec.getVarsList().length + 1][];
int column = odeSolverResultSet.findColumn("t");
data[0] = odeSolverResultSet.extractColumn(column);
int origDataLength = data[0].length;
for (int i = 0; i < testSpec.getVarsList().length; i++) {
column = odeSolverResultSet.findColumn(testSpec.getVarsList()[i]);
if (column == -1) {
Variable var = simJob.getSimulationSymbolTable().getVariable(testSpec.getVarsList()[i]);
data[i + 1] = new double[data[0].length];
if (var instanceof cbit.vcell.math.Constant) {
double value = ((cbit.vcell.math.Constant) var).getExpression().evaluateConstant();
for (int j = 0; j < data[i + 1].length; j++) {
data[i + 1][j] = value;
}
} else {
throw new RuntimeException("Did not find " + testSpec.getVarsList()[i] + " in simulation");
}
} else {
data[i + 1] = odeSolverResultSet.extractColumn(column);
}
}
//
// for each time, print row
//
int index = 0;
double[] sampleTimes = new double[testSpec.getNumTimeSteps() + 1];
for (int i = 0; i <= testSpec.getNumTimeSteps(); i++) {
sampleTimes[i] = endTime * i / testSpec.getNumTimeSteps();
}
Model vcModel = bioModel.getModel();
ReservedSymbol kMole = vcModel.getKMOLE();
for (int i = 0; i < sampleTimes.length; i++) {
//
while (true) {
//
if (index == odeSolverResultSet.getRowCount() - 1) {
if (data[0][index] == sampleTimes[i]) {
break;
} else {
throw new RuntimeException("sampleTime does not match at last time point");
}
}
//
if (data[0][index + 1] > sampleTimes[i]) {
break;
}
//
// sampleTime must be later in our data list.
//
index++;
}
// if data[0][index] == sampleTime no need to interpolate
if (data[0][index] == sampleTimes[i]) {
// if timeFactor is not 1.0, time is not in seconds (mins or hrs); if timeFactor is 60, divide sampleTime/60; if it is 3600, divide sampleTime/3600.
// if (timeFactor != 1.0) {
// outputStream.print(data[0][index]/timeFactor);
// } else {
outputStream.print(data[0][index]);
// }
for (int j = 0; j < testSpec.getVarsList().length; j++) {
// SBMLImporter.SBVCConcentrationUnits spConcUnits = speciesUnitsHash.get(testSpec.getVarsList()[j]);
// if (spConcUnits != null) {
// VCUnitDefinition sbunits = spConcUnits.getSBConcentrationUnits();
// VCUnitDefinition vcunits = spConcUnits.getVCConcentrationUnits();
// SBMLUnitParameter unitFactor = SBMLUtils.getConcUnitFactor("spConcParam", vcunits, sbunits, kMole);
// outputStream.print("," + data[j + 1][index] * unitFactor.getExpression().evaluateConstant()); //earlier, hack unitfactor = 0.000001
// earlier, hack unitfactor = 0.000001
outputStream.print("," + data[j + 1][index]);
// }
}
// System.out.println("No interpolation needed!");
outputStream.println();
} else {
// if data[0][index] < sampleTime, must interpolate
double fraction = (sampleTimes[i] - data[0][index]) / (data[0][index + 1] - data[0][index]);
// if timeFactor is not 1.0, time is not in seconds (mins or hrs); if timeFactor is 60, divide sampleTime/60; if it is 3600, divide sampleTime/3600.
// if (timeFactor != 1.0) {
// outputStream.print(sampleTimes[i]/timeFactor);
// } else {
outputStream.print(sampleTimes[i]);
// }
for (int j = 0; j < testSpec.getVarsList().length; j++) {
double interpolatedValue = 0.0;
double[] speciesVals = null;
double[] times = null;
// Currently using 2nd order interpolation
if (index == 0) {
// can only do 1st order interpolation
times = new double[] { data[0][index], data[0][index + 1] };
speciesVals = new double[] { data[j + 1][index], data[j + 1][index + 1] };
interpolatedValue = MathTestingUtilities.taylorInterpolation(sampleTimes[i], times, speciesVals);
} else if (index >= 1 && index <= origDataLength - 3) {
double val_1 = Math.abs(sampleTimes[i] - data[0][index - 1]);
double val_2 = Math.abs(sampleTimes[i] - data[0][index + 2]);
if (val_1 < val_2) {
times = new double[] { data[0][index - 1], data[0][index], data[0][index + 1] };
speciesVals = new double[] { data[j + 1][index - 1], data[j + 1][index], data[j + 1][index + 1] };
} else {
times = new double[] { data[0][index], data[0][index + 1], data[0][index + 2] };
speciesVals = new double[] { data[j + 1][index], data[j + 1][index + 1], data[j + 1][index + 2] };
}
interpolatedValue = MathTestingUtilities.taylorInterpolation(sampleTimes[i], times, speciesVals);
} else {
times = new double[] { data[0][index - 1], data[0][index], data[0][index + 1] };
speciesVals = new double[] { data[j + 1][index - 1], data[j + 1][index], data[j + 1][index + 1] };
interpolatedValue = MathTestingUtilities.taylorInterpolation(sampleTimes[i], times, speciesVals);
}
// // Currently using 1st order interpolation
// times = new double[] { data[0][index], data[0][index+1] };
// speciesVals = new double[] { data[j+1][index], data[j+1][index+1] };
// interpolatedValue = taylorInterpolation(sampleTimes[i], times, speciesVals);
// interpolatedValue = interpolatedValue * unitFactor.getExpression().evaluateConstant(); //earlier, hack unitfactor = 0.000001
// System.out.println("Sample time: " + sampleTimes[i] + ", between time[" + index + "]=" + data[0][index]+" and time["+(index+1)+"]="+(data[0][index+1])+", interpolated = "+interpolatedValue);
outputStream.print("," + interpolatedValue);
}
outputStream.println();
}
}
outputStream.close();
return outputFile;
} catch (RuntimeException e) {
e.printStackTrace(System.out);
// rethrow without losing context
throw e;
} catch (Exception e) {
e.printStackTrace(System.out);
throw new SolverException(e.getMessage(), e);
}
}
Aggregations