Search in sources :

Example 6 with OutputTimeSpec

use of cbit.vcell.solver.OutputTimeSpec in project vcell by virtualcell.

the class SundialsSolver method writeLogFile.

private void writeLogFile() throws SolverException {
    String logFile = getBaseName() + LOGFILE_EXTENSION;
    String ideDataFileName = new File(getBaseName() + IDA_DATA_EXTENSION).getName();
    PrintWriter pw = null;
    try {
        pw = new PrintWriter(logFile);
        pw.println(IDA_DATA_IDENTIFIER);
        pw.println(IDA_DATA_FORMAT_ID);
        pw.println(ideDataFileName);
        OutputTimeSpec outputTimeSpec = simTask.getSimulation().getSolverTaskDescription().getOutputTimeSpec();
        if (outputTimeSpec.isDefault()) {
            pw.println(KEEP_MOST + " " + ((DefaultOutputTimeSpec) outputTimeSpec).getKeepAtMost());
        }
        pw.close();
    } catch (FileNotFoundException e) {
        e.printStackTrace();
        throw new SolverException(e.getMessage());
    } finally {
        if (pw != null) {
            pw.close();
        }
    }
}
Also used : DefaultOutputTimeSpec(cbit.vcell.solver.DefaultOutputTimeSpec) OutputTimeSpec(cbit.vcell.solver.OutputTimeSpec) FileNotFoundException(java.io.FileNotFoundException) SolverException(cbit.vcell.solver.SolverException) File(java.io.File) DefaultOutputTimeSpec(cbit.vcell.solver.DefaultOutputTimeSpec) PrintWriter(java.io.PrintWriter)

Example 7 with OutputTimeSpec

use of cbit.vcell.solver.OutputTimeSpec in project vcell by virtualcell.

the class MovingBoundaryFileWriter method gettimeReport.

private Element gettimeReport() {
    Element e = new Element(MBTags.timeReport);
    Element e1 = null;
    OutputTimeSpec outTimeSpec = simTask.getSimulation().getSolverTaskDescription().getOutputTimeSpec();
    if (outTimeSpec instanceof UniformOutputTimeSpec) {
        UniformOutputTimeSpec uniformOutputTimeSpec = (UniformOutputTimeSpec) outTimeSpec;
        double interval = uniformOutputTimeSpec.getOutputTimeStep();
        // hard coded.
        double startTime = 0.0;
        e1 = new Element(MBTags.startTime);
        e1.setText("0");
        e.addContent(e1);
        e1 = new Element(MBTags.interval);
        e1.setText(interval + "");
        e.addContent(e1);
    } else if (outTimeSpec instanceof DefaultOutputTimeSpec) {
        DefaultOutputTimeSpec defaultOutputTimeSpec = (DefaultOutputTimeSpec) outTimeSpec;
        Integer step = defaultOutputTimeSpec.getKeepEvery();
        // hard code
        double startTime = 0.0;
        e1 = new Element(MBTags.startTime);
        e1.setText("0");
        e.addContent(e1);
        e1 = new Element(MBTags.step);
        e1.setText(step.toString());
        e.addContent(e1);
    }
    return e;
}
Also used : DefaultOutputTimeSpec(cbit.vcell.solver.DefaultOutputTimeSpec) UniformOutputTimeSpec(cbit.vcell.solver.UniformOutputTimeSpec) OutputTimeSpec(cbit.vcell.solver.OutputTimeSpec) UniformOutputTimeSpec(cbit.vcell.solver.UniformOutputTimeSpec) Element(org.jdom.Element) DefaultOutputTimeSpec(cbit.vcell.solver.DefaultOutputTimeSpec)

Example 8 with OutputTimeSpec

use of cbit.vcell.solver.OutputTimeSpec in project vcell by virtualcell.

the class XmlHelper method sedmlToBioModel.

public static VCDocument sedmlToBioModel(VCLogger transLogger, ExternalDocInfo externalDocInfo, SedML sedml, AbstractTask selectedTask) throws Exception {
    if (sedml.getModels().isEmpty()) {
        return null;
    }
    VCDocument doc = null;
    try {
        // extract the path only from the sedml file
        String fullPath = FileUtils.getFullPath(externalDocInfo.getFile().getAbsolutePath());
        // Namespace namespace = sedml.getNamespace();
        // iterate through all the elements and show them at the console
        List<org.jlibsedml.Model> mmm = sedml.getModels();
        for (Model mm : mmm) {
            System.out.println(mm.toString());
        }
        List<org.jlibsedml.Simulation> sss = sedml.getSimulations();
        for (org.jlibsedml.Simulation ss : sss) {
            System.out.println(ss.toString());
        }
        List<AbstractTask> ttt = sedml.getTasks();
        for (AbstractTask tt : ttt) {
            System.out.println(tt.toString());
        }
        List<DataGenerator> ddd = sedml.getDataGenerators();
        for (DataGenerator dd : ddd) {
            System.out.println(dd.toString());
        }
        List<Output> ooo = sedml.getOutputs();
        for (Output oo : ooo) {
            System.out.println(oo.toString());
        }
        KisaoTerm sedmlKisao = null;
        // this will become the vCell simulation
        org.jlibsedml.Simulation sedmlSimulation = null;
        // the "original" model referred to by the task
        org.jlibsedml.Model sedmlOriginalModel = null;
        String sedmlOriginalModelName = null;
        if (selectedTask == null) {
            // no task, just pick the Model and find its sbml file
            sedmlOriginalModelName = SEDMLUtil.getName(mmm.get(0));
        } else {
            if (selectedTask instanceof Task) {
                sedmlOriginalModel = sedml.getModelWithId(selectedTask.getModelReference());
                sedmlSimulation = sedml.getSimulation(selectedTask.getSimulationReference());
            } else if (selectedTask instanceof RepeatedTask) {
                RepeatedTask rt = (RepeatedTask) selectedTask;
                assert (rt.getSubTasks().size() == 1);
                // first (and only) subtask
                SubTask st = rt.getSubTasks().entrySet().iterator().next().getValue();
                String taskId = st.getTaskId();
                AbstractTask t = sedml.getTaskWithId(taskId);
                // get model and simulation from subtask
                sedmlOriginalModel = sedml.getModelWithId(t.getModelReference());
                sedmlSimulation = sedml.getSimulation(t.getSimulationReference());
            } else {
                throw new RuntimeException("Unexpected task " + selectedTask);
            }
            sedmlOriginalModelName = sedmlOriginalModel.getId();
            sedmlKisao = KisaoOntology.getInstance().getTermById(sedmlSimulation.getAlgorithm().getKisaoID());
        }
        // UniformTimeCourse [initialTime=0.0, numberOfPoints=1000, outputEndTime=1.0, outputStartTime=0.0,
        // Algorithm [kisaoID=KISAO:0000019], getId()=SimSlow]
        // identify the vCell solvers that would match best the sedml solver kisao id
        List<SolverDescription> solverDescriptions = new ArrayList<>();
        for (SolverDescription sd : SolverDescription.values()) {
            KisaoTerm solverKisaoTerm = KisaoOntology.getInstance().getTermById(sd.getKisao());
            if (solverKisaoTerm == null) {
                break;
            }
            boolean isExactlySame = solverKisaoTerm.equals(sedmlKisao);
            if (isExactlySame && !solverKisaoTerm.isObsolete()) {
                // we make a list with all the solvers that match the kisao
                solverDescriptions.add(sd);
            }
        }
        // from the list of vcell solvers that match the sedml kisao we select the ones that have a matching time step
        SolverDescription solverDescription = null;
        for (SolverDescription sd : solverDescriptions) {
            if (true) {
                solverDescription = sd;
                break;
            }
        }
        // find out everything else we need about the application we're going to use,
        // some of the info will be needed when we parse the sbml file
        boolean bSpatial = false;
        Application appType = Application.NETWORK_DETERMINISTIC;
        Set<SolverDescription.SolverFeature> sfList = solverDescription.getSupportedFeatures();
        for (SolverDescription.SolverFeature sf : sfList) {
            switch(sf) {
                case Feature_Rulebased:
                    appType = Application.RULE_BASED_STOCHASTIC;
                    break;
                case Feature_Stochastic:
                    appType = Application.NETWORK_STOCHASTIC;
                    break;
                case Feature_Deterministic:
                    appType = Application.NETWORK_DETERMINISTIC;
                    break;
                case Feature_Spatial:
                    bSpatial = true;
                    break;
                default:
                    break;
            }
        }
        // -------------------------------------------------------------------------------------------
        // extract bioModel name from sedx (or sedml) file
        String bioModelName = FileUtils.getBaseName(externalDocInfo.getFile().getAbsolutePath());
        // if we have repeated task, we ignore them, we just use the normal resolvers for archive and changes
        // once the application and simulation are built, we iterate through the repeated tasks and
        // add math overrides to the simulation for each repeated task
        ArchiveComponents ac = null;
        if (externalDocInfo.getFile().getPath().toLowerCase().endsWith("sedx") || externalDocInfo.getFile().getPath().toLowerCase().endsWith("omex")) {
            ac = Libsedml.readSEDMLArchive(new FileInputStream(externalDocInfo.getFile().getPath()));
        }
        ModelResolver resolver = new ModelResolver(sedml);
        if (ac != null) {
            resolver.add(new ArchiveModelResolver(ac));
        }
        resolver.add(new FileModelResolver());
        resolver.add(new RelativeFileModelResolver(fullPath));
        String newMdl = resolver.getModelString(sedmlOriginalModel);
        // sbmlSource with all the changes applied
        XMLSource sbmlSource = new XMLSource(newMdl);
        doc = XmlHelper.importSBML(transLogger, sbmlSource, bSpatial);
        BioModel bioModel = (BioModel) doc;
        bioModel.setName(bioModelName);
        // we already have an application loaded from the sbml file, with initial conditions and stuff
        // which may be not be suitable because the sedml kisao may need a different app type
        // so we do a "copy as" to the right type and then delete the original we loaded from the sbml file
        // the new application we're making from the old one
        SimulationContext newSimulationContext = null;
        if (bioModel.getSimulationContexts().length == 1) {
            SimulationContext oldSimulationContext = bioModel.getSimulationContext(0);
            newSimulationContext = SimulationContext.copySimulationContext(oldSimulationContext, sedmlOriginalModelName, bSpatial, appType);
            bioModel.removeSimulationContext(oldSimulationContext);
            bioModel.addSimulationContext(newSimulationContext);
        } else {
            // length == 0
            newSimulationContext = bioModel.addNewSimulationContext(sedmlOriginalModelName, appType);
        }
        // making the new vCell simulation based on the sedml simulation
        newSimulationContext.refreshDependencies();
        MathMappingCallback callback = new MathMappingCallbackTaskAdapter(null);
        newSimulationContext.refreshMathDescription(callback, NetworkGenerationRequirements.ComputeFullStandardTimeout);
        Simulation newSimulation = new Simulation(newSimulationContext.getMathDescription());
        newSimulation.setName(SEDMLUtil.getName(sedmlSimulation));
        // TODO: make sure that everything has proper names
        // we check the repeated tasks, if any, and add to the list of math overrides
        // if(selectedTask instanceof RepeatedTask) {
        // for(Change change : ((RepeatedTask) selectedTask).getChanges()) {
        // if(!(change instanceof SetValue)) {
        // throw new RuntimeException("Only 'SetValue' changes are supported for repeated tasks.");
        // }
        // SetValue setValue = (SetValue)change;
        // // TODO: extract target from XPath
        // // ......
        // //
        // String target = "s0";	// for now we just use a hardcoded thing
        // ConstantArraySpec cas;
        // Range range = ((RepeatedTask) selectedTask).getRange(setValue.getRangeReference());
        // if(range instanceof UniformRange) {
        // cas = ConstantArraySpec.createIntervalSpec(target, ((UniformRange) range).getStart(), ((UniformRange) range).getEnd(),
        // range.getNumElements(), ((UniformRange) range).getType() == UniformRange.UniformType.LOG ? true : false);
        // } else if(range instanceof VectorRange) {
        // //    				List<String> constants = new ArrayList<> ();
        // //    				for(int i=0; i<range.getNumElements(); i++) {
        // //    					constants.add(new Constant(i+"", new Expression(range.getElementAt(i))));
        // //    				}
        // //    				cas = ConstantArraySpec.createListSpec(target, constants);
        // 
        // } else {
        // throw new RuntimeException("Only 'Uniform Range' and 'Vector Range' are supported at this time.");
        // }
        // 
        // }
        // }
        // we identify the type of sedml simulation (uniform time course, etc)
        // and set the vCell simulation parameters accordingly
        SolverTaskDescription simTaskDesc = newSimulation.getSolverTaskDescription();
        TimeBounds timeBounds = new TimeBounds();
        TimeStep timeStep = new TimeStep();
        double outputTimeStep = 0.1;
        if (sedmlSimulation instanceof UniformTimeCourse) {
            // we translate initial time to zero, we provide output for the duration of the simulation
            // because we can't select just an interval the way the SEDML simulation can
            double initialTime = ((UniformTimeCourse) sedmlSimulation).getInitialTime();
            double outputStartTime = ((UniformTimeCourse) sedmlSimulation).getOutputStartTime();
            double outputEndTime = ((UniformTimeCourse) sedmlSimulation).getOutputEndTime();
            double outputNumberOfPoints = ((UniformTimeCourse) sedmlSimulation).getNumberOfPoints();
            outputTimeStep = (outputEndTime - outputStartTime) / outputNumberOfPoints;
            timeBounds = new TimeBounds(0, outputEndTime - initialTime);
        } else if (sedmlSimulation instanceof OneStep) {
        // for anything other than UniformTimeCourse we just ignore
        } else if (sedmlSimulation instanceof SteadyState) {
        } else {
        }
        OutputTimeSpec outputTimeSpec = new UniformOutputTimeSpec(outputTimeStep);
        simTaskDesc.setTimeBounds(timeBounds);
        simTaskDesc.setTimeStep(timeStep);
        simTaskDesc.setOutputTimeSpec(outputTimeSpec);
        newSimulation.setSolverTaskDescription(simTaskDesc);
        bioModel.addSimulation(newSimulation);
        newSimulation.refreshDependencies();
    } catch (Exception e) {
        e.printStackTrace();
        throw new RuntimeException("Unable to initialize bioModel for the given selection.");
    }
    return doc;
}
Also used : Task(org.jlibsedml.Task) RepeatedTask(org.jlibsedml.RepeatedTask) SimulationTask(cbit.vcell.messaging.server.SimulationTask) AbstractTask(org.jlibsedml.AbstractTask) SubTask(org.jlibsedml.SubTask) ArrayList(java.util.ArrayList) ArchiveModelResolver(org.jlibsedml.execution.ArchiveModelResolver) FileModelResolver(org.jlibsedml.execution.FileModelResolver) RelativeFileModelResolver(org.vcell.sedml.RelativeFileModelResolver) ModelResolver(org.jlibsedml.execution.ModelResolver) SteadyState(org.jlibsedml.SteadyState) OutputTimeSpec(cbit.vcell.solver.OutputTimeSpec) UniformOutputTimeSpec(cbit.vcell.solver.UniformOutputTimeSpec) SolverTaskDescription(cbit.vcell.solver.SolverTaskDescription) SubTask(org.jlibsedml.SubTask) VCDocument(org.vcell.util.document.VCDocument) MathMappingCallback(cbit.vcell.mapping.SimulationContext.MathMappingCallback) UniformOutputTimeSpec(cbit.vcell.solver.UniformOutputTimeSpec) FileInputStream(java.io.FileInputStream) KisaoTerm(org.jlibsedml.modelsupport.KisaoTerm) ArchiveModelResolver(org.jlibsedml.execution.ArchiveModelResolver) DataGenerator(org.jlibsedml.DataGenerator) UniformTimeCourse(org.jlibsedml.UniformTimeCourse) Application(cbit.vcell.mapping.SimulationContext.Application) AbstractTask(org.jlibsedml.AbstractTask) SolverDescription(cbit.vcell.solver.SolverDescription) RelativeFileModelResolver(org.vcell.sedml.RelativeFileModelResolver) OneStep(org.jlibsedml.OneStep) TimeBounds(cbit.vcell.solver.TimeBounds) TimeStep(cbit.vcell.solver.TimeStep) ArchiveComponents(org.jlibsedml.ArchiveComponents) RepeatedTask(org.jlibsedml.RepeatedTask) Output(org.jlibsedml.Output) MathMappingCallbackTaskAdapter(cbit.vcell.mapping.MathMappingCallbackTaskAdapter) FileModelResolver(org.jlibsedml.execution.FileModelResolver) RelativeFileModelResolver(org.vcell.sedml.RelativeFileModelResolver) SimulationContext(cbit.vcell.mapping.SimulationContext) XMLStreamException(javax.xml.stream.XMLStreamException) SbmlException(org.vcell.sbml.SbmlException) SBMLException(org.sbml.jsbml.SBMLException) IOException(java.io.IOException) ExpressionException(cbit.vcell.parser.ExpressionException) Simulation(cbit.vcell.solver.Simulation) Model(org.jlibsedml.Model) BioModel(cbit.vcell.biomodel.BioModel) MathModel(cbit.vcell.mathmodel.MathModel) BioModel(cbit.vcell.biomodel.BioModel) Model(org.jlibsedml.Model)

Example 9 with OutputTimeSpec

use of cbit.vcell.solver.OutputTimeSpec in project vcell by virtualcell.

the class SimulationWorkspace method getEstimatedNumTimePointsForStoch.

private static long getEstimatedNumTimePointsForStoch(SimulationSymbolTable simSymbolTable) {
    Simulation sim = simSymbolTable.getSimulation();
    SolverTaskDescription solverTaskDescription = sim.getSolverTaskDescription();
    TimeBounds tb = solverTaskDescription.getTimeBounds();
    double startTime = tb.getStartingTime();
    double endTime = tb.getEndingTime();
    OutputTimeSpec tSpec = solverTaskDescription.getOutputTimeSpec();
    // hybrid G_E and G_M are fixed time step methods using uniform output time spec
    if (tSpec.isUniform()) {
        double outputTimeStep = ((UniformOutputTimeSpec) tSpec).getOutputTimeStep();
        return (long) ((endTime - startTime) / outputTimeStep);
    }
    double maxProbability = 0;
    SubDomain subDomain = sim.getMathDescription().getSubDomains().nextElement();
    List<VarIniCondition> varInis = subDomain.getVarIniConditions();
    // get all the probability expressions
    ArrayList<Expression> probList = new ArrayList<Expression>();
    for (JumpProcess jp : subDomain.getJumpProcesses()) {
        probList.add(jp.getProbabilityRate());
    }
    // loop through probability expressions
    for (int i = 0; i < probList.size(); i++) {
        try {
            Expression pExp = new Expression(probList.get(i));
            pExp.bindExpression(simSymbolTable);
            pExp = simSymbolTable.substituteFunctions(pExp);
            pExp = pExp.flatten();
            String[] symbols = pExp.getSymbols();
            // substitute stoch vars with it's initial condition expressions
            if (symbols != null) {
                for (int j = 0; symbols != null && j < symbols.length; j++) {
                    for (int k = 0; k < varInis.size(); k++) {
                        if (symbols[j].equals(varInis.get(k).getVar().getName())) {
                            pExp.substituteInPlace(new Expression(symbols[j]), new Expression(varInis.get(k).getIniVal()));
                            break;
                        }
                    }
                }
            }
            pExp = simSymbolTable.substituteFunctions(pExp);
            pExp = pExp.flatten();
            double val = pExp.evaluateConstant();
            if (maxProbability < val) {
                maxProbability = val;
            }
        } catch (ExpressionBindingException e) {
            System.out.println("Cannot estimate the total time points for stochastic simulation!! Due to the reason below...");
            e.printStackTrace();
        } catch (ExpressionException ex) {
            System.out.println("Cannot estimate the total time points for stochastic simulation!! Due to the reason below...");
            ex.printStackTrace();
        } catch (MathException e) {
            System.out.println("Cannot estimate the total time points for stochastic simulation!! Due to the reason below...");
            e.printStackTrace();
        }
    }
    int keepEvery = 1;
    if (tSpec.isDefault()) {
        keepEvery = ((DefaultOutputTimeSpec) tSpec).getKeepEvery();
    }
    // points = (endt-startt)/(t*keepEvery) = (endt - startt)/(keepEvery*1/prob)
    long estimatedPoints = Math.round((tb.getEndingTime() - tb.getStartingTime()) * maxProbability / keepEvery) + 1;
    return estimatedPoints;
}
Also used : VarIniCondition(cbit.vcell.math.VarIniCondition) UniformOutputTimeSpec(cbit.vcell.solver.UniformOutputTimeSpec) ArrayList(java.util.ArrayList) ExpressionBindingException(cbit.vcell.parser.ExpressionBindingException) ExpressionException(cbit.vcell.parser.ExpressionException) SubDomain(cbit.vcell.math.SubDomain) TimeBounds(cbit.vcell.solver.TimeBounds) DefaultOutputTimeSpec(cbit.vcell.solver.DefaultOutputTimeSpec) UniformOutputTimeSpec(cbit.vcell.solver.UniformOutputTimeSpec) OutputTimeSpec(cbit.vcell.solver.OutputTimeSpec) Simulation(cbit.vcell.solver.Simulation) Expression(cbit.vcell.parser.Expression) MathException(cbit.vcell.math.MathException) JumpProcess(cbit.vcell.math.JumpProcess) SolverTaskDescription(cbit.vcell.solver.SolverTaskDescription)

Example 10 with OutputTimeSpec

use of cbit.vcell.solver.OutputTimeSpec in project vcell by virtualcell.

the class Generate2DExpModelOpAbstract method generateModel.

public final GeneratedModelResults generateModel(double deltaX, double bleachRadius, double cellRadius, double bleachDuration, double bleachRate, double postbleachDelay, double postbleachDuration, double psfSigma, double outputTimeStep, double primaryDiffusionRate, double primaryFraction, double bleachMonitorRate, double secondaryDiffusionRate, double secondaryFraction, String extracellularName, String cytosolName, Context context) throws PropertyVetoException, ExpressionException, GeometryException, ImageException, ModelException, MappingException, MathException, MatrixException {
    double domainSize = 2.2 * cellRadius;
    Extent extent = new Extent(domainSize, domainSize, 1.0);
    Origin origin = new Origin(-extent.getX() / 2.0, -extent.getY() / 2.0, -extent.getZ() / 2.0);
    String EXTRACELLULAR_NAME = extracellularName;
    String CYTOSOL_NAME = cytosolName;
    AnalyticSubVolume cytosolSubVolume = new AnalyticSubVolume(CYTOSOL_NAME, new Expression("pow(x,2)+pow(y,2)<pow(" + cellRadius + ",2)"));
    AnalyticSubVolume extracellularSubVolume = new AnalyticSubVolume(EXTRACELLULAR_NAME, new Expression(1.0));
    Geometry geometry = new Geometry("geometry", 2);
    geometry.getGeometrySpec().setExtent(extent);
    geometry.getGeometrySpec().setOrigin(origin);
    geometry.getGeometrySpec().addSubVolume(extracellularSubVolume);
    geometry.getGeometrySpec().addSubVolume(cytosolSubVolume, true);
    geometry.getGeometrySurfaceDescription().updateAll();
    BioModel bioModel = new BioModel(null);
    bioModel.setName("unnamed");
    Model model = new Model("model");
    bioModel.setModel(model);
    model.addFeature(EXTRACELLULAR_NAME);
    Feature extracellular = (Feature) model.getStructure(EXTRACELLULAR_NAME);
    model.addFeature(CYTOSOL_NAME);
    Feature cytosol = (Feature) model.getStructure(CYTOSOL_NAME);
    SpeciesContext immobileSC = model.createSpeciesContext(cytosol);
    SpeciesContext primarySC = model.createSpeciesContext(cytosol);
    SpeciesContext secondarySC = model.createSpeciesContext(cytosol);
    // 
    // common bleaching rate for all species
    // 
    double bleachStart = 10 * outputTimeStep - bleachDuration - postbleachDelay;
    double bleachEnd = bleachStart + bleachDuration;
    Expression bleachRateExp = createBleachExpression(bleachRadius, bleachRate, bleachMonitorRate, bleachStart, bleachEnd);
    {
        SimpleReaction immobileBWM = model.createSimpleReaction(cytosol);
        GeneralKinetics immobileBWMKinetics = new GeneralKinetics(immobileBWM);
        immobileBWM.setKinetics(immobileBWMKinetics);
        immobileBWM.addReactant(immobileSC, 1);
        immobileBWMKinetics.getReactionRateParameter().setExpression(Expression.mult(bleachRateExp, new Expression(immobileSC.getName())));
    }
    {
        SimpleReaction primaryBWM = model.createSimpleReaction(cytosol);
        GeneralKinetics primaryBWMKinetics = new GeneralKinetics(primaryBWM);
        primaryBWM.setKinetics(primaryBWMKinetics);
        primaryBWM.addReactant(primarySC, 1);
        primaryBWMKinetics.getReactionRateParameter().setExpression(Expression.mult(bleachRateExp, new Expression(primarySC.getName())));
    }
    {
        SimpleReaction secondaryBWM = model.createSimpleReaction(cytosol);
        GeneralKinetics secondaryBWMKinetics = new GeneralKinetics(secondaryBWM);
        secondaryBWM.setKinetics(secondaryBWMKinetics);
        secondaryBWM.addReactant(secondarySC, 1);
        secondaryBWMKinetics.getReactionRateParameter().setExpression(Expression.mult(bleachRateExp, new Expression(secondarySC.getName())));
    }
    // create simulation context
    SimulationContext simContext = bioModel.addNewSimulationContext("simContext", SimulationContext.Application.NETWORK_DETERMINISTIC);
    simContext.setGeometry(geometry);
    FeatureMapping cytosolFeatureMapping = (FeatureMapping) simContext.getGeometryContext().getStructureMapping(cytosol);
    FeatureMapping extracellularFeatureMapping = (FeatureMapping) simContext.getGeometryContext().getStructureMapping(extracellular);
    SubVolume cytSubVolume = geometry.getGeometrySpec().getSubVolume(CYTOSOL_NAME);
    SubVolume exSubVolume = geometry.getGeometrySpec().getSubVolume(EXTRACELLULAR_NAME);
    // unused? SurfaceClass pmSurfaceClass = geometry.getGeometrySurfaceDescription().getSurfaceClass(exSubVolume, cytSubVolume);
    cytosolFeatureMapping.setGeometryClass(cytSubVolume);
    extracellularFeatureMapping.setGeometryClass(exSubVolume);
    cytosolFeatureMapping.getUnitSizeParameter().setExpression(new Expression(1.0));
    extracellularFeatureMapping.getUnitSizeParameter().setExpression(new Expression(1.0));
    double fixedFraction = 1.0 - primaryFraction - secondaryFraction;
    SpeciesContextSpec immobileSCS = simContext.getReactionContext().getSpeciesContextSpec(immobileSC);
    immobileSCS.getInitialConditionParameter().setExpression(new Expression(fixedFraction));
    immobileSCS.getDiffusionParameter().setExpression(new Expression(0.0));
    SpeciesContextSpec primarySCS = simContext.getReactionContext().getSpeciesContextSpec(primarySC);
    primarySCS.getInitialConditionParameter().setExpression(new Expression(primaryFraction));
    primarySCS.getDiffusionParameter().setExpression(new Expression(primaryDiffusionRate));
    SpeciesContextSpec secondarySCS = simContext.getReactionContext().getSpeciesContextSpec(secondarySC);
    secondarySCS.getInitialConditionParameter().setExpression(new Expression(secondaryFraction));
    secondarySCS.getDiffusionParameter().setExpression(new Expression(secondaryDiffusionRate));
    simContext.getMicroscopeMeasurement().addFluorescentSpecies(immobileSC);
    simContext.getMicroscopeMeasurement().addFluorescentSpecies(primarySC);
    simContext.getMicroscopeMeasurement().addFluorescentSpecies(secondarySC);
    simContext.getMicroscopeMeasurement().setConvolutionKernel(new GaussianConvolutionKernel(new Expression(psfSigma), new Expression(psfSigma)));
    MathMapping mathMapping = simContext.createNewMathMapping();
    MathDescription mathDesc = mathMapping.getMathDescription();
    simContext.setMathDescription(mathDesc);
    User owner = context.getDefaultOwner();
    int meshSize = (int) (domainSize / deltaX);
    if (meshSize % 2 == 0) {
        // want an odd-sized mesh in x and y ... so centered at the origin.
        meshSize = meshSize + 1;
    }
    TimeBounds timeBounds = new TimeBounds(0.0, postbleachDuration);
    // 
    // simulation to use for data generation (output time steps as recorded by the microscope)
    // 
    double bleachBlackoutBegin = bleachStart - postbleachDelay;
    double bleachBlackoutEnd = bleachEnd + postbleachDelay;
    // ArrayList<Double> times = new ArrayList<Double>();
    // double time = 0;
    // while (time<=timeBounds.getEndingTime()){
    // if (time<=bleachBlackoutBegin || time>bleachBlackoutEnd){
    // // postbleachDelay is the time it takes to switch the filters.
    // times.add(time);
    // }
    // time += outputTimeStep.getData();
    // }
    // double[] timeArray = new double[times.size()];
    // for (int i=0;i<timeArray.length;i++){
    // timeArray[i] = times.get(i);
    // }
    // OutputTimeSpec fakeDataSimOutputTimeSpec = new ExplicitOutputTimeSpec(timeArray);
    OutputTimeSpec fakeDataSimOutputTimeSpec = new UniformOutputTimeSpec(outputTimeStep);
    KeyValue fakeDataSimKey = context.createNewKeyValue();
    SimulationVersion fakeDataSimVersion = new SimulationVersion(fakeDataSimKey, "fakeDataSim", owner, new GroupAccessNone(), new KeyValue("0"), new BigDecimal(0), new Date(), VersionFlag.Current, "", null);
    Simulation fakeDataSim = new Simulation(fakeDataSimVersion, mathDesc);
    simContext.addSimulation(fakeDataSim);
    fakeDataSim.getSolverTaskDescription().setTimeBounds(timeBounds);
    fakeDataSim.getMeshSpecification().setSamplingSize(new ISize(meshSize, meshSize, 1));
    fakeDataSim.getSolverTaskDescription().setSolverDescription(SolverDescription.SundialsPDE);
    fakeDataSim.getSolverTaskDescription().setOutputTimeSpec(fakeDataSimOutputTimeSpec);
    // 
    // simulation to use for viewing the protocol (output time steps to understand the physics)
    // 
    KeyValue fullExperimentSimKey = context.createNewKeyValue();
    SimulationVersion fullExperimentSimVersion = new SimulationVersion(fullExperimentSimKey, "fullExperiment", owner, new GroupAccessNone(), new KeyValue("0"), new BigDecimal(0), new Date(), VersionFlag.Current, "", null);
    Simulation fullExperimentSim = new Simulation(fullExperimentSimVersion, mathDesc);
    simContext.addSimulation(fullExperimentSim);
    OutputTimeSpec fullExperimentOutputTimeSpec = new UniformOutputTimeSpec(outputTimeStep / 10.0);
    fullExperimentSim.getSolverTaskDescription().setTimeBounds(timeBounds);
    fullExperimentSim.getMeshSpecification().setSamplingSize(new ISize(meshSize, meshSize, 1));
    fullExperimentSim.getSolverTaskDescription().setSolverDescription(SolverDescription.SundialsPDE);
    fullExperimentSim.getSolverTaskDescription().setOutputTimeSpec(fullExperimentOutputTimeSpec);
    GeneratedModelResults results = new GeneratedModelResults();
    results.bioModel_2D = bioModel;
    results.simulation_2D = fakeDataSim;
    results.bleachBlackoutBeginTime = bleachBlackoutBegin;
    results.bleachBlackoutEndTime = bleachBlackoutEnd;
    return results;
}
Also used : Origin(org.vcell.util.Origin) User(org.vcell.util.document.User) KeyValue(org.vcell.util.document.KeyValue) Extent(org.vcell.util.Extent) MathDescription(cbit.vcell.math.MathDescription) ISize(org.vcell.util.ISize) SpeciesContext(cbit.vcell.model.SpeciesContext) GeneralKinetics(cbit.vcell.model.GeneralKinetics) SpeciesContextSpec(cbit.vcell.mapping.SpeciesContextSpec) Feature(cbit.vcell.model.Feature) TimeBounds(cbit.vcell.solver.TimeBounds) GroupAccessNone(org.vcell.util.document.GroupAccessNone) UniformOutputTimeSpec(cbit.vcell.solver.UniformOutputTimeSpec) OutputTimeSpec(cbit.vcell.solver.OutputTimeSpec) SimulationVersion(org.vcell.util.document.SimulationVersion) FeatureMapping(cbit.vcell.mapping.FeatureMapping) SubVolume(cbit.vcell.geometry.SubVolume) AnalyticSubVolume(cbit.vcell.geometry.AnalyticSubVolume) SimpleReaction(cbit.vcell.model.SimpleReaction) UniformOutputTimeSpec(cbit.vcell.solver.UniformOutputTimeSpec) SimulationContext(cbit.vcell.mapping.SimulationContext) GaussianConvolutionKernel(cbit.vcell.mapping.MicroscopeMeasurement.GaussianConvolutionKernel) BigDecimal(java.math.BigDecimal) Date(java.util.Date) Geometry(cbit.vcell.geometry.Geometry) Simulation(cbit.vcell.solver.Simulation) Expression(cbit.vcell.parser.Expression) BioModel(cbit.vcell.biomodel.BioModel) BioModel(cbit.vcell.biomodel.BioModel) Model(cbit.vcell.model.Model) MathMapping(cbit.vcell.mapping.MathMapping) AnalyticSubVolume(cbit.vcell.geometry.AnalyticSubVolume)

Aggregations

OutputTimeSpec (cbit.vcell.solver.OutputTimeSpec)21 UniformOutputTimeSpec (cbit.vcell.solver.UniformOutputTimeSpec)18 DefaultOutputTimeSpec (cbit.vcell.solver.DefaultOutputTimeSpec)13 Simulation (cbit.vcell.solver.Simulation)9 ExplicitOutputTimeSpec (cbit.vcell.solver.ExplicitOutputTimeSpec)7 SolverTaskDescription (cbit.vcell.solver.SolverTaskDescription)7 TimeBounds (cbit.vcell.solver.TimeBounds)6 BioModel (cbit.vcell.biomodel.BioModel)5 SimulationContext (cbit.vcell.mapping.SimulationContext)5 PropertyVetoException (java.beans.PropertyVetoException)5 ArrayList (java.util.ArrayList)5 SpeciesContextSpec (cbit.vcell.mapping.SpeciesContextSpec)4 ExpressionException (cbit.vcell.parser.ExpressionException)4 SolverDescription (cbit.vcell.solver.SolverDescription)4 File (java.io.File)4 IOException (java.io.IOException)4 MathMappingCallbackTaskAdapter (cbit.vcell.mapping.MathMappingCallbackTaskAdapter)3 MathDescription (cbit.vcell.math.MathDescription)3 SubDomain (cbit.vcell.math.SubDomain)3 Expression (cbit.vcell.parser.Expression)3