Search in sources :

Example 21 with UniformOutputTimeSpec

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;
}
Also used : UniformOutputTimeSpec(cbit.vcell.solver.UniformOutputTimeSpec) Element(org.jdom.Element)

Example 22 with UniformOutputTimeSpec

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();
        }
    }
}
Also used : QueryHashtable(cbit.sql.QueryHashtable) User(org.vcell.util.document.User) KeyValue(org.vcell.util.document.KeyValue) MathMappingCallback(cbit.vcell.mapping.SimulationContext.MathMappingCallback) UniformOutputTimeSpec(cbit.vcell.solver.UniformOutputTimeSpec) MathDescription(cbit.vcell.math.MathDescription) SimulationContext(cbit.vcell.mapping.SimulationContext) SpeciesContextSpec(cbit.vcell.mapping.SpeciesContextSpec) ServerDocumentManager(cbit.vcell.modeldb.ServerDocumentManager) PropertyVetoException(java.beans.PropertyVetoException) SQLException(java.sql.SQLException) XmlParseException(cbit.vcell.xml.XmlParseException) ImageException(cbit.image.ImageException) IOException(java.io.IOException) DataAccessException(org.vcell.util.DataAccessException) ExpressionException(cbit.vcell.parser.ExpressionException) MappingException(cbit.vcell.mapping.MappingException) GeometryException(cbit.vcell.geometry.GeometryException) UniformOutputTimeSpec(cbit.vcell.solver.UniformOutputTimeSpec) OutputTimeSpec(cbit.vcell.solver.OutputTimeSpec) StringWriter(java.io.StringWriter) BioModel(cbit.vcell.biomodel.BioModel) MathMapping(cbit.vcell.mapping.MathMapping) XMLSource(cbit.vcell.xml.XMLSource) File(java.io.File) PrintWriter(java.io.PrintWriter)

Example 23 with UniformOutputTimeSpec

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);
    }
}
Also used : MathMappingCallbackTaskAdapter(cbit.vcell.mapping.MathMappingCallbackTaskAdapter) UniformOutputTimeSpec(cbit.vcell.solver.UniformOutputTimeSpec) ArrayList(java.util.ArrayList) SimulationContext(cbit.vcell.mapping.SimulationContext) SpeciesContextSpec(cbit.vcell.mapping.SpeciesContextSpec) UniformOutputTimeSpec(cbit.vcell.solver.UniformOutputTimeSpec) OutputTimeSpec(cbit.vcell.solver.OutputTimeSpec) TempSimulation(cbit.vcell.solver.TempSimulation) Simulation(cbit.vcell.solver.Simulation) BioModel(cbit.vcell.biomodel.BioModel) NetworkGenerationRequirements(cbit.vcell.mapping.SimulationContext.NetworkGenerationRequirements) XMLSource(cbit.vcell.xml.XMLSource) AnnotatedFunction(cbit.vcell.solver.AnnotatedFunction)

Example 24 with UniformOutputTimeSpec

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);
    }
}
Also used : KeyValue(org.vcell.util.document.KeyValue) MathDescription(cbit.vcell.math.MathDescription) StructureMapping(cbit.vcell.mapping.StructureMapping) TimeStep(cbit.vcell.solver.TimeStep) ErrorTolerance(cbit.vcell.solver.ErrorTolerance) Structure(cbit.vcell.model.Structure) SimulationJob(cbit.vcell.solver.SimulationJob) UniformOutputTimeSpec(cbit.vcell.solver.UniformOutputTimeSpec) VCLogger(cbit.util.xml.VCLogger) IOException(java.io.IOException) SimulationContext(cbit.vcell.mapping.SimulationContext) IOException(java.io.IOException) Simulation(cbit.vcell.solver.Simulation) BioModel(cbit.vcell.biomodel.BioModel) File(java.io.File) XMLSource(cbit.vcell.xml.XMLSource) VCLogger(cbit.util.xml.VCLogger)

Example 25 with UniformOutputTimeSpec

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);
    }
}
Also used : KeyValue(org.vcell.util.document.KeyValue) SimulationTask(cbit.vcell.messaging.server.SimulationTask) Variable(cbit.vcell.math.Variable) MathDescription(cbit.vcell.math.MathDescription) ReservedSymbol(cbit.vcell.model.Model.ReservedSymbol) CVodeFileWriter(cbit.vcell.solver.ode.CVodeFileWriter) TimeBounds(cbit.vcell.solver.TimeBounds) TimeStep(cbit.vcell.solver.TimeStep) SimulationVersion(org.vcell.util.document.SimulationVersion) ErrorTolerance(cbit.vcell.solver.ErrorTolerance) ODESolverResultSet(cbit.vcell.solver.ode.ODESolverResultSet) Executable(org.vcell.util.exe.Executable) SimulationJob(cbit.vcell.solver.SimulationJob) PrintWriter(java.io.PrintWriter) UniformOutputTimeSpec(cbit.vcell.solver.UniformOutputTimeSpec) VCLogger(cbit.util.xml.VCLogger) IOException(java.io.IOException) SimulationContext(cbit.vcell.mapping.SimulationContext) ExecutableException(org.vcell.util.exe.ExecutableException) XMLStreamException(javax.xml.stream.XMLStreamException) XmlParseException(cbit.vcell.xml.XmlParseException) SolverException(cbit.vcell.solver.SolverException) SbmlException(org.vcell.sbml.SbmlException) IOException(java.io.IOException) SBMLImportException(org.vcell.sbml.vcell.SBMLImportException) Simulation(cbit.vcell.solver.Simulation) BioModel(cbit.vcell.biomodel.BioModel) MathMapping(cbit.vcell.mapping.MathMapping) BioModel(cbit.vcell.biomodel.BioModel) Model(cbit.vcell.model.Model) SolverException(cbit.vcell.solver.SolverException) File(java.io.File) XMLSource(cbit.vcell.xml.XMLSource)

Aggregations

UniformOutputTimeSpec (cbit.vcell.solver.UniformOutputTimeSpec)34 Simulation (cbit.vcell.solver.Simulation)18 OutputTimeSpec (cbit.vcell.solver.OutputTimeSpec)17 TimeBounds (cbit.vcell.solver.TimeBounds)15 DefaultOutputTimeSpec (cbit.vcell.solver.DefaultOutputTimeSpec)14 BioModel (cbit.vcell.biomodel.BioModel)13 SimulationContext (cbit.vcell.mapping.SimulationContext)13 SpeciesContextSpec (cbit.vcell.mapping.SpeciesContextSpec)10 MathDescription (cbit.vcell.math.MathDescription)10 KeyValue (org.vcell.util.document.KeyValue)10 Expression (cbit.vcell.parser.Expression)9 TimeStep (cbit.vcell.solver.TimeStep)9 IOException (java.io.IOException)9 ImageException (cbit.image.ImageException)8 Model (cbit.vcell.model.Model)8 SimulationVersion (org.vcell.util.document.SimulationVersion)8 Geometry (cbit.vcell.geometry.Geometry)7 MathMapping (cbit.vcell.mapping.MathMapping)7 ErrorTolerance (cbit.vcell.solver.ErrorTolerance)7 ExplicitOutputTimeSpec (cbit.vcell.solver.ExplicitOutputTimeSpec)7