Search in sources :

Example 41 with Simulation

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

the class FiniteVolumeFileWriter method writeSimulationParamters.

/**
 * @param timeFunction
 * @param startTime
 * @param endTime
 * @param rootFinder
 * @param uniqueRootTimes
 * @param bPrintIterations
 * @throws ExpressionException
 *
 * for testing within scrapbook, see below:
 *
 * try {
 *	edu.northwestern.at.utils.math.rootfinders.MonadicFunctionRootFinder rootFinder =
 *//			new edu.northwestern.at.utils.math.rootfinders.Brent();
 *			new edu.northwestern.at.utils.math.rootfinders.Bisection();
 *//			new edu.northwestern.at.utils.math.rootfinders.NewtonRaphson();
 *//			new edu.northwestern.at.utils.math.rootfinders.Secant();
 *	cbit.vcell.parser.SimpleSymbolTable simpleSymbolTable = new cbit.vcell.parser.SimpleSymbolTable(new String[] { "t" });
 *
 *	cbit.vcell.parser.Expression exp = new cbit.vcell.parser.Expression("t-0.56");
 *
 *	exp.bindExpression(simpleSymbolTable);
 *	java.util.TreeSet<Double> rootTimes = new java.util.TreeSet<Double>();
 *	double startTime = 0.0;
 *	double endTime = 100.0;
 *	System.out.print("exp = '"+ exp.infix() + "'");
 *	long currentTimeMS = System.currentTimeMillis();
 *	cbit.vcell.solvers.FiniteVolumeFileWriter.findAllRoots(exp,startTime,endTime,rootFinder,rootTimes,false);
 *	long finalTimeMS = System.currentTimeMillis();
 *	for (double root : rootTimes){
 *		System.out.println("root = "+root);
 *	}
 *	System.out.println("elapsedTime of computation = "+(finalTimeMS-currentTimeMS)+" ms, found " + rootTimes.size() + " roots (not unique)");
 *
 *}catch (Exception e){
 *	e.printStackTrace(System.out);
 *}
 */
/*public static void findAllRoots(Expression timeFunction, double startTime, double endTime, MonadicFunctionRootFinder rootFinder, TreeSet<Double> uniqueRootTimes, boolean bPrintIterations) throws ExpressionException{
	TreeSet<Double> allRootTimes = new TreeSet<Double>();
	final Expression function_exp = new Expression(timeFunction);
	MonadicFunction valueFunction = new MonadicFunction() {
		double[] values = new double[1];
		public double f(double t) {
			values[0] = t;
			try {
				return function_exp.evaluateVector(values);
			} catch (ExpressionException e) {
				e.printStackTrace();
				throw new RuntimeException("expression exception "+e.getMessage());
			}
		}
	};
	
	final Expression derivative_exp = new Expression(timeFunction.differentiate(ReservedVariable.TIME.getName()));
	MonadicFunction derivativeFunction = new MonadicFunction() {
		double[] values = new double[1];
		public double f(double t) {
			values[0] = t;
			try {
				return derivative_exp.evaluateVector(values);
			} catch (ExpressionException e) {
				e.printStackTrace();
				throw new RuntimeException("expression exception "+e.getMessage());
			}
		}
	};
	
	RootFinderConvergenceTest convergenceTest = new StandardRootFinderConvergenceTest();
	RootFinderIterationInformation iterationInformation = null;
	if (bPrintIterations){
		iterationInformation = new RootFinderIterationInformation() {				
			public void iterationInformation(double x, double fx, double dfx, int currentIteration) {
				System.out.println(currentIteration+") x="+x+", fx="+fx+", dfx="+dfx);
			}
		};
	}
	int NUM_BRACKETS = 1000;
	double simulationTime = endTime - startTime;
	double tolerance = simulationTime/1e10;
	int maxIter = 1000;
	
	for (int i=0;i<NUM_BRACKETS-1;i++){
		double bracketMin = startTime + simulationTime*i/NUM_BRACKETS;
		double bracketMax = startTime + simulationTime*(i+1)/NUM_BRACKETS;
	
		double root = rootFinder.findRoot(bracketMin, bracketMax, tolerance, maxIter, valueFunction, derivativeFunction, convergenceTest, iterationInformation);
		if (root>startTime && root<endTime && valueFunction.f(root)<=tolerance){
			allRootTimes.add(root);
		}
	}
	double uniqueTolerance = tolerance * 100;
	double lastUniqueRoot = Double.NEGATIVE_INFINITY;
	for (double root : allRootTimes){
		if (root-lastUniqueRoot > uniqueTolerance){
			uniqueRootTimes.add(root);
		}
		lastUniqueRoot = root;
	}

}  ---------------------------------JIM's CODE COMMENTTED FOR FUTURE DEVELOPMENT*/
/**
 *# Simulation Parameters
 *SIMULATION_PARAM_BEGIN
 *SOLVER SUNDIALS_PDE_SOLVER 1.0E-7 1.0E-9
 *DISCONTINUITY_TIMES 2 1.0E-4 3.0000000000000003E-4
 *BASE_FILE_NAME c:/Vcell/users/fgao/SimID_31746636_0_
 *ENDING_TIME 4.0E-4
 *KEEP_EVERY ONE_STEP 3
 *KEEP_AT_MOST 1000
 *SIMULATION_PARAM_END
 * @throws MathException
 * @throws ExpressionException
 */
private void writeSimulationParamters() throws ExpressionException, MathException {
    Simulation simulation = simTask.getSimulation();
    SolverTaskDescription solverTaskDesc = simulation.getSolverTaskDescription();
    printWriter.println("# Simulation Parameters");
    printWriter.println(FVInputFileKeyword.SIMULATION_PARAM_BEGIN);
    if (solverTaskDesc.getSolverDescription().equals(SolverDescription.SundialsPDE)) {
        printWriter.print(FVInputFileKeyword.SOLVER + " " + FVInputFileKeyword.SUNDIALS_PDE_SOLVER + " " + solverTaskDesc.getErrorTolerance().getRelativeErrorTolerance() + " " + solverTaskDesc.getErrorTolerance().getAbsoluteErrorTolerance() + " " + solverTaskDesc.getTimeStep().getMaximumTimeStep());
        if (simulation.getMathDescription().hasVelocity()) {
            printWriter.print(" " + solverTaskDesc.getSundialsPdeSolverOptions().getMaxOrderAdvection());
        }
        printWriter.println();
        Vector<Discontinuity> discontinuities = new Vector<Discontinuity>();
        TreeSet<Double> discontinuityTimes = new TreeSet<Double>();
        MathDescription mathDesc = simulation.getMathDescription();
        Enumeration<SubDomain> enum1 = mathDesc.getSubDomains();
        while (enum1.hasMoreElements()) {
            SubDomain sd = enum1.nextElement();
            Enumeration<Equation> enum_equ = sd.getEquations();
            while (enum_equ.hasMoreElements()) {
                Equation equation = enum_equ.nextElement();
                equation.getDiscontinuities(simTask.getSimulationJob().getSimulationSymbolTable(), discontinuities);
            }
        }
        getDiscontinuityTimes(discontinuities, discontinuityTimes);
        if (discontinuityTimes.size() > 0) {
            printWriter.print(FVInputFileKeyword.DISCONTINUITY_TIMES + " " + discontinuityTimes.size());
            for (double d : discontinuityTimes) {
                printWriter.print(" " + d);
            }
            printWriter.println();
        }
    } else if (solverTaskDesc.getSolverDescription().equals(SolverDescription.Chombo)) {
        printWriter.println(FVInputFileKeyword.SOLVER + " " + FVInputFileKeyword.CHOMBO_SEMIIMPLICIT_SOLVER);
    } else if (solverTaskDesc.getSolverDescription().equals(SolverDescription.VCellPetsc)) {
        printWriter.println(FVInputFileKeyword.SOLVER + " " + FVInputFileKeyword.VCELL_PETSC_SOLVER);
    } else {
        printWriter.println(FVInputFileKeyword.SOLVER + " " + FVInputFileKeyword.FV_SOLVER + " " + solverTaskDesc.getErrorTolerance().getRelativeErrorTolerance());
    }
    printWriter.println(FVInputFileKeyword.BASE_FILE_NAME + " " + new File(workingDirectory, simTask.getSimulationJob().getSimulationJobID()).getAbsolutePath());
    if (solverTaskDesc.isParallel() && destinationDirectory != null && !destinationDirectory.equals(workingDirectory)) {
        printWriter.println(FVInputFileKeyword.PRIMARY_DATA_DIR + " " + destinationDirectory.getAbsolutePath());
    }
    printWriter.println(FVInputFileKeyword.ENDING_TIME + " " + solverTaskDesc.getTimeBounds().getEndingTime());
    OutputTimeSpec outputTimeSpec = solverTaskDesc.getOutputTimeSpec();
    if (solverTaskDesc.getSolverDescription().isChomboSolver()) {
        List<TimeInterval> timeIntervalList = solverTaskDesc.getChomboSolverSpec().getTimeIntervalList();
        printWriter.println(FVInputFileKeyword.TIME_INTERVALS + " " + timeIntervalList.size());
        for (TimeInterval ti : timeIntervalList) {
            printWriter.println(ti.getEndingTime() + " " + ti.getTimeStep() + " " + ti.getKeepEvery());
        }
    } else if (solverTaskDesc.getSolverDescription().equals(SolverDescription.SundialsPDE)) {
        if (outputTimeSpec.isDefault()) {
            DefaultOutputTimeSpec defaultOutputTimeSpec = (DefaultOutputTimeSpec) outputTimeSpec;
            printWriter.println(FVInputFileKeyword.KEEP_EVERY + " " + FVInputFileKeyword.ONE_STEP + " " + defaultOutputTimeSpec.getKeepEvery());
            printWriter.println(FVInputFileKeyword.KEEP_AT_MOST + " " + defaultOutputTimeSpec.getKeepAtMost());
        } else {
            printWriter.println(FVInputFileKeyword.TIME_STEP + " " + ((UniformOutputTimeSpec) outputTimeSpec).getOutputTimeStep());
            printWriter.println(FVInputFileKeyword.KEEP_EVERY + " 1");
        }
    } else {
        double defaultTimeStep = solverTaskDesc.getTimeStep().getDefaultTimeStep();
        printWriter.println(FVInputFileKeyword.TIME_STEP + " " + defaultTimeStep);
        int keepEvery = 1;
        if (outputTimeSpec.isDefault()) {
            keepEvery = ((DefaultOutputTimeSpec) outputTimeSpec).getKeepEvery();
        } else if (outputTimeSpec.isUniform()) {
            UniformOutputTimeSpec uots = (UniformOutputTimeSpec) outputTimeSpec;
            double ots = uots.getOutputTimeStep();
            keepEvery = (int) Math.round(ots / defaultTimeStep);
        } else {
            throw new RuntimeException("unexpected OutputTime specification type :" + outputTimeSpec.getClass().getName());
        }
        if (keepEvery <= 0) {
            throw new RuntimeException("Output KeepEvery must be a positive integer. Try to change the output option.");
        }
        printWriter.println(FVInputFileKeyword.KEEP_EVERY + " " + keepEvery);
    }
    ErrorTolerance stopAtSpatiallyUniformErrorTolerance = solverTaskDesc.getStopAtSpatiallyUniformErrorTolerance();
    if (stopAtSpatiallyUniformErrorTolerance != null) {
        printWriter.println(FVInputFileKeyword.CHECK_SPATIALLY_UNIFORM + " " + stopAtSpatiallyUniformErrorTolerance.getAbsoluteErrorTolerance() + " " + stopAtSpatiallyUniformErrorTolerance.getRelativeErrorTolerance());
    }
    printWriter.println(FVInputFileKeyword.SIMULATION_PARAM_END);
    printWriter.println();
}
Also used : Discontinuity(cbit.vcell.parser.Discontinuity) TimeInterval(org.vcell.chombo.TimeInterval) UniformOutputTimeSpec(cbit.vcell.solver.UniformOutputTimeSpec) MathDescription(cbit.vcell.math.MathDescription) MeasureEquation(cbit.vcell.math.MeasureEquation) PdeEquation(cbit.vcell.math.PdeEquation) VolumeRegionEquation(cbit.vcell.math.VolumeRegionEquation) MembraneRegionEquation(cbit.vcell.math.MembraneRegionEquation) Equation(cbit.vcell.math.Equation) CompartmentSubDomain(cbit.vcell.math.CompartmentSubDomain) SubDomain(cbit.vcell.math.SubDomain) MembraneSubDomain(cbit.vcell.math.MembraneSubDomain) DefaultOutputTimeSpec(cbit.vcell.solver.DefaultOutputTimeSpec) OutputTimeSpec(cbit.vcell.solver.OutputTimeSpec) UniformOutputTimeSpec(cbit.vcell.solver.UniformOutputTimeSpec) Simulation(cbit.vcell.solver.Simulation) TreeSet(java.util.TreeSet) ErrorTolerance(cbit.vcell.solver.ErrorTolerance) SolverTaskDescription(cbit.vcell.solver.SolverTaskDescription) Vector(java.util.Vector) File(java.io.File) DefaultOutputTimeSpec(cbit.vcell.solver.DefaultOutputTimeSpec)

Example 42 with Simulation

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

the class FiniteVolumeFileWriter method writeChomboSpec.

private void writeChomboSpec() throws ExpressionException, SolverException, PropertyVetoException, ClassNotFoundException, IOException, GeometryException, ImageException {
    if (!bChomboSolver) {
        return;
    }
    GeometrySpec geometrySpec = resampledGeometry.getGeometrySpec();
    int dimension = geometrySpec.getDimension();
    if (dimension == 1) {
        throw new SolverException(simTask.getSimulation().getSolverTaskDescription().getSolverDescription().getDisplayLabel() + " is only supported for simulations with 2D or 3D geometry.");
    }
    Simulation simulation = getSimulationTask().getSimulation();
    SolverTaskDescription solverTaskDescription = simulation.getSolverTaskDescription();
    ChomboSolverSpec chomboSolverSpec = solverTaskDescription.getChomboSolverSpec();
    printWriter.println(FVInputFileKeyword.CHOMBO_SPEC_BEGIN);
    printWriter.println(FVInputFileKeyword.DIMENSION + " " + geometrySpec.getDimension());
    Extent extent = geometrySpec.getExtent();
    Origin origin = geometrySpec.getOrigin();
    ISize isize = simulation.getMeshSpecification().getSamplingSize();
    switch(geometrySpec.getDimension()) {
        case 2:
            printWriter.println(FVInputFileKeyword.MESH_SIZE + " " + isize.getX() + " " + isize.getY());
            printWriter.println(FVInputFileKeyword.DOMAIN_SIZE + " " + extent.getX() + " " + extent.getY());
            printWriter.println(FVInputFileKeyword.DOMAIN_ORIGIN + " " + origin.getX() + " " + origin.getY());
            break;
        case 3:
            printWriter.println(FVInputFileKeyword.MESH_SIZE + " " + isize.getX() + " " + isize.getY() + " " + isize.getZ());
            printWriter.println(FVInputFileKeyword.DOMAIN_SIZE + " " + extent.getX() + " " + extent.getY() + " " + extent.getZ());
            printWriter.println(FVInputFileKeyword.DOMAIN_ORIGIN + " " + origin.getX() + " " + origin.getY() + " " + origin.getZ());
            break;
    }
    List<CompartmentSubDomain> featureList = new ArrayList<CompartmentSubDomain>();
    Enumeration<SubDomain> enum1 = simulation.getMathDescription().getSubDomains();
    while (enum1.hasMoreElements()) {
        SubDomain sd = enum1.nextElement();
        if (sd instanceof CompartmentSubDomain) {
            featureList.add((CompartmentSubDomain) sd);
        }
    }
    int numFeatures = featureList.size();
    CompartmentSubDomain[] features = featureList.toArray(new CompartmentSubDomain[0]);
    int[] phases = new int[numFeatures];
    Arrays.fill(phases, -1);
    phases[numFeatures - 1] = 0;
    int[] numAssigned = new int[] { 1 };
    assignPhases(features, numFeatures - 1, phases, numAssigned);
    Map<String, Integer> subDomainPhaseMap = new HashMap<String, Integer>();
    for (int i = 0; i < phases.length; ++i) {
        if (phases[i] == -1) {
            throw new SolverException("Failed to assign a phase to CompartmentSubdomain '" + features[i].getName() + "'. It might be caused by too coarsh a mesh.");
        }
        subDomainPhaseMap.put(features[i].getName(), phases[i]);
    }
    SubVolume[] subVolumes = geometrySpec.getSubVolumes();
    if (geometrySpec.hasImage()) {
        Geometry geometry = (Geometry) BeanUtils.cloneSerializable(simulation.getMathDescription().getGeometry());
        Geometry simGeometry = geometry;
        VCImage img = geometry.getGeometrySpec().getImage();
        int factor = Math.max(Math.max(img.getNumX(), img.getNumY()), img.getNumZ()) < 512 ? 2 : 1;
        ISize distanceMapMeshSize = new ISize(img.getNumX() * factor, img.getNumY() * factor, img.getNumZ() * factor);
        Vect3d deltaX = null;
        boolean bCellCentered = false;
        double dx = 0.5;
        double dy = 0.5;
        double dz = 0.5;
        int Nx = distanceMapMeshSize.getX();
        int Ny = distanceMapMeshSize.getY();
        int Nz = distanceMapMeshSize.getZ();
        if (dimension == 2) {
            // pad the 2D image with itself in order to obtain a 3D image used to compute the distance map
            // because the distance map algorithm is 3D only (using distance to triangles)
            byte[] oldPixels = img.getPixels();
            byte[] newPixels = new byte[oldPixels.length * 3];
            System.arraycopy(oldPixels, 0, newPixels, 0, oldPixels.length);
            System.arraycopy(oldPixels, 0, newPixels, oldPixels.length, oldPixels.length);
            System.arraycopy(oldPixels, 0, newPixels, oldPixels.length * 2, oldPixels.length);
            double distX = geometry.getExtent().getX() / img.getNumX();
            double distY = geometry.getExtent().getY() / img.getNumY();
            // we set the distance on the z axis to something that makes sense
            double distZ = Math.max(distX, distY);
            Extent newExtent = new Extent(geometry.getExtent().getX(), geometry.getExtent().getY(), distZ * 3);
            VCImage newImage = new VCImageUncompressed(null, newPixels, newExtent, img.getNumX(), img.getNumY(), 3);
            // copy the pixel classes too
            ArrayList<VCPixelClass> newPixelClasses = new ArrayList<VCPixelClass>();
            for (VCPixelClass origPixelClass : geometry.getGeometrySpec().getImage().getPixelClasses()) {
                SubVolume origSubvolume = geometry.getGeometrySpec().getImageSubVolumeFromPixelValue(origPixelClass.getPixel());
                newPixelClasses.add(new VCPixelClass(null, origSubvolume.getName(), origPixelClass.getPixel()));
            }
            newImage.setPixelClasses(newPixelClasses.toArray(new VCPixelClass[newPixelClasses.size()]));
            simGeometry = new Geometry(geometry, newImage);
            Nz = 3;
        }
        GeometrySpec simGeometrySpec = simGeometry.getGeometrySpec();
        Extent simExtent = simGeometrySpec.getExtent();
        dx = simExtent.getX() / (Nx - 1);
        dy = simExtent.getY() / (Ny - 1);
        dz = simExtent.getZ() / (Nz - 1);
        if (Math.abs(dx - dy) > 0.1 * Math.max(dx, dy)) {
            dx = Math.min(dx, dy);
            dy = dx;
            Nx = (int) (simExtent.getX() / dx + 1);
            Ny = (int) (simExtent.getY() / dx + 1);
            if (dimension == 3) {
                dz = dx;
                Nz = (int) (simExtent.getZ() / dx + 1);
            }
        }
        deltaX = new Vect3d(dx, dy, dz);
        // one more point in each direction
        distanceMapMeshSize = new ISize(Nx + 1, Ny + 1, Nz + 1);
        Extent distanceMapExtent = new Extent(simExtent.getX() + dx, simExtent.getY() + dy, simExtent.getZ() + dz);
        simGeometrySpec.setExtent(distanceMapExtent);
        GeometrySurfaceDescription geoSurfaceDesc = simGeometry.getGeometrySurfaceDescription();
        geoSurfaceDesc.setVolumeSampleSize(distanceMapMeshSize);
        geoSurfaceDesc.updateAll();
        VCImage vcImage = RayCaster.sampleGeometry(simGeometry, distanceMapMeshSize, bCellCentered);
        SubvolumeSignedDistanceMap[] distanceMaps = DistanceMapGenerator.computeDistanceMaps(simGeometry, vcImage, bCellCentered);
        if (dimension == 2) {
            distanceMaps = DistanceMapGenerator.extractMiddleSlice(distanceMaps);
        }
        printWriter.println(FVInputFileKeyword.SUBDOMAINS + " " + simGeometrySpec.getNumSubVolumes() + " " + FVInputFileKeyword.DISTANCE_MAP);
        for (int i = 0; i < subVolumes.length; i++) {
            File distanceMapFile = new File(workingDirectory, getSimulationTask().getSimulationJobID() + "_" + subVolumes[i].getName() + DISTANCE_MAP_FILE_EXTENSION);
            writeDistanceMapFile(deltaX, distanceMaps[i], distanceMapFile);
            int phase = subDomainPhaseMap.get(subVolumes[i].getName());
            printWriter.println(subVolumes[i].getName() + " " + phase + " " + distanceMapFile.getAbsolutePath());
        }
    } else {
        printWriter.println(FVInputFileKeyword.SUBDOMAINS + " " + geometrySpec.getNumSubVolumes());
        Expression[] rvachevExps = convertAnalyticGeometryToRvachevFunction(geometrySpec);
        for (int i = 0; i < subVolumes.length; i++) {
            if (subVolumes[i] instanceof AnalyticSubVolume) {
                String name = subVolumes[i].getName();
                int phase = subDomainPhaseMap.get(name);
                printWriter.println(name + " " + phase + " ");
                printWriter.println(FVInputFileKeyword.IF + " " + rvachevExps[i].infix() + ";");
                printWriter.println(FVInputFileKeyword.USER + " " + ((AnalyticSubVolume) subVolumes[i]).getExpression().infix() + ";");
            }
        }
    }
    printWriter.println(FVInputFileKeyword.MAX_BOX_SIZE + " " + chomboSolverSpec.getMaxBoxSize());
    printWriter.println(FVInputFileKeyword.FILL_RATIO + " " + chomboSolverSpec.getFillRatio());
    printWriter.println(FVInputFileKeyword.RELATIVE_TOLERANCE + " " + simulation.getSolverTaskDescription().getErrorTolerance().getRelativeErrorTolerance());
    printWriter.println(FVInputFileKeyword.SAVE_VCELL_OUTPUT + " " + chomboSolverSpec.isSaveVCellOutput());
    printWriter.println(FVInputFileKeyword.SAVE_CHOMBO_OUTPUT + " " + chomboSolverSpec.isSaveChomboOutput());
    printWriter.println(FVInputFileKeyword.ACTIVATE_FEATURE_UNDER_DEVELOPMENT + " " + chomboSolverSpec.isActivateFeatureUnderDevelopment());
    printWriter.println(FVInputFileKeyword.SMALL_VOLFRAC_THRESHOLD + " " + chomboSolverSpec.getSmallVolfracThreshold());
    printWriter.println(FVInputFileKeyword.BLOCK_FACTOR + " " + chomboSolverSpec.getBlockFactor());
    printWriter.println(FVInputFileKeyword.TAGS_GROW + " " + chomboSolverSpec.getTagsGrow());
    // Refinement
    int numLevels = chomboSolverSpec.getNumRefinementLevels();
    // Refinements #Levels ratio 1, ratio 2, etc
    printWriter.print(FVInputFileKeyword.REFINEMENTS + " " + (numLevels + 1));
    List<Integer> ratios = chomboSolverSpec.getRefineRatioList();
    for (int i : ratios) {
        printWriter.print(" " + i);
    }
    // write last refinement ratio, fake
    printWriter.println(" 2");
    // membrane rois
    List<RefinementRoi> memRios = chomboSolverSpec.getMembraneRefinementRois();
    printWriter.println(FVInputFileKeyword.REFINEMENT_ROIS + " " + RoiType.Membrane + " " + memRios.size());
    for (RefinementRoi roi : memRios) {
        if (roi.getRoiExpression() == null) {
            throw new SolverException("ROI expression cannot be null");
        }
        // level tagsGrow ROIexpression
        printWriter.println(roi.getLevel() + " " + roi.getRoiExpression().infix() + ";");
    }
    List<RefinementRoi> volRios = chomboSolverSpec.getVolumeRefinementRois();
    printWriter.println(FVInputFileKeyword.REFINEMENT_ROIS + " " + RoiType.Volume + " " + volRios.size());
    for (RefinementRoi roi : volRios) {
        if (roi.getRoiExpression() == null) {
            throw new SolverException("ROI expression cannot be null");
        }
        printWriter.println(roi.getLevel() + " " + roi.getRoiExpression().infix() + ";");
    }
    printWriter.println(FVInputFileKeyword.VIEW_LEVEL + " " + chomboSolverSpec.getViewLevel());
    printWriter.println(FVInputFileKeyword.CHOMBO_SPEC_END);
    printWriter.println();
}
Also used : Origin(org.vcell.util.Origin) VCPixelClass(cbit.image.VCPixelClass) GeometrySurfaceDescription(cbit.vcell.geometry.surface.GeometrySurfaceDescription) Extent(org.vcell.util.Extent) HashMap(java.util.HashMap) ISize(org.vcell.util.ISize) ArrayList(java.util.ArrayList) VCImage(cbit.image.VCImage) SubvolumeSignedDistanceMap(cbit.vcell.geometry.surface.SubvolumeSignedDistanceMap) GeometrySpec(cbit.vcell.geometry.GeometrySpec) CompartmentSubDomain(cbit.vcell.math.CompartmentSubDomain) SubDomain(cbit.vcell.math.SubDomain) MembraneSubDomain(cbit.vcell.math.MembraneSubDomain) RefinementRoi(org.vcell.chombo.RefinementRoi) SubVolume(cbit.vcell.geometry.SubVolume) AnalyticSubVolume(cbit.vcell.geometry.AnalyticSubVolume) SolverTaskDescription(cbit.vcell.solver.SolverTaskDescription) VCImageUncompressed(cbit.image.VCImageUncompressed) ChomboSolverSpec(org.vcell.chombo.ChomboSolverSpec) Vect3d(cbit.vcell.render.Vect3d) Geometry(cbit.vcell.geometry.Geometry) Simulation(cbit.vcell.solver.Simulation) Expression(cbit.vcell.parser.Expression) CompartmentSubDomain(cbit.vcell.math.CompartmentSubDomain) SolverException(cbit.vcell.solver.SolverException) File(java.io.File) AnalyticSubVolume(cbit.vcell.geometry.AnalyticSubVolume)

Example 43 with Simulation

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

the class FiniteVolumeFileWriter method writeSerialParameterScans.

private void writeSerialParameterScans() throws ExpressionException, DivideByZeroException {
    Simulation simulation = simTask.getSimulation();
    if (!simulation.isSerialParameterScan()) {
        return;
    }
    int scanCount = simulation.getScanCount();
    printWriter.println("# Serial Scan Parameters");
    String[] scanParameters = simulation.getMathOverrides().getOverridenConstantNames();
    printWriter.println("SERIAL_SCAN_PARAMETER_BEGIN " + scanParameters.length);
    for (int i = 0; i < scanParameters.length; i++) {
        String scanParameter = scanParameters[i];
        printWriter.println(scanParameter);
    }
    printWriter.println("SERIAL_SCAN_PARAMETER_END");
    printWriter.println();
    printWriter.println("# Parameter Scan Values");
    printWriter.println("SERIAL_SCAN_PARAMETER_VALUE_BEGIN " + scanCount);
    for (int i = 0; i < scanCount; i++) {
        for (int j = 0; j < scanParameters.length; j++) {
            String scanParameter = scanParameters[j];
            Expression exp = simulation.getMathOverrides().getActualExpression(scanParameter, i);
            double value = exp.evaluateConstant();
            printWriter.print(value + " ");
        }
        printWriter.println();
    }
    printWriter.println("SERIAL_SCAN_PARAMETER_VALUE_END");
    printWriter.println();
}
Also used : Simulation(cbit.vcell.solver.Simulation) Expression(cbit.vcell.parser.Expression)

Example 44 with Simulation

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

the class FiniteVolumeFileWriter method writeCompartments.

/**
 *COMPARTMENT_BEGIN nucleus
 *
 *BOUNDARY_CONDITIONS value value value value
 *
 *EQUATION_BEGIN rfB
 *INITIAL _VCell_FieldData_0;
 *RATE ( - (50.0 * rfB * ((x > -5.0) && (x < 5.0) && (y > -5.0) && (y < 5.0))) + (0.02 * ( - rB - rfB + (20.0 * ((x > -3.0) && (x < 3.0) && (y > -5.0) && (y < 5.0))) + _VCell_FieldData_0) * rf) - (0.1 * rfB));
 *EQUATION_END
 *
 *EQUATION_BEGIN r
 *INITIAL 5.0;
 *RATE ( - ((0.02 * r * ( - rB - rfB + (20.0 * ((x > -3.0) && (x < 3.0) && (y > -5.0) && (y < 5.0))) + _VCell_FieldData_0)) - (0.1 * rB)) + (50.0 * rf * ((x > -5.0) && (x < 5.0) && (y > -5.0) && (y < 5.0))));
 *DIFFUSION 10.0;
 *BOUNDARY_XM 5.0;
 *BOUNDARY_XP 5.0;
 *BOUNDARY_YM 5.0;
 *BOUNDARY_YP 5.0;
 *EQUATION_END
 *
 *COMPARTMENT_END
 * @throws ExpressionException
 * @throws MathException
 */
private void writeCompartments() throws ExpressionException, MathException {
    Simulation simulation = simTask.getSimulation();
    MathDescription mathDesc = simulation.getMathDescription();
    Enumeration<SubDomain> enum1 = mathDesc.getSubDomains();
    while (enum1.hasMoreElements()) {
        SubDomain sd = enum1.nextElement();
        if (sd instanceof CompartmentSubDomain) {
            CompartmentSubDomain csd = (CompartmentSubDomain) sd;
            printWriter.println("COMPARTMENT_BEGIN " + csd.getName());
            printWriter.println();
            writeCompartment_VarContext(csd);
            writeFastSystem(csd);
            printWriter.println("COMPARTMENT_END");
            printWriter.println();
        }
    }
    printWriter.println();
}
Also used : CompartmentSubDomain(cbit.vcell.math.CompartmentSubDomain) SubDomain(cbit.vcell.math.SubDomain) MembraneSubDomain(cbit.vcell.math.MembraneSubDomain) Simulation(cbit.vcell.solver.Simulation) MathDescription(cbit.vcell.math.MathDescription) CompartmentSubDomain(cbit.vcell.math.CompartmentSubDomain)

Example 45 with Simulation

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

the class FiniteVolumeFileWriter method writeMembranes.

/**
 *MEMBRANE_BEGIN subVolume0_subVolume1_membrane subVolume0 subVolume1
 *
 *BOUNDARY_CONDITIONS value value value value
 *
 *EQUATION_BEGIN varMem
 *INITIAL (x > 0.75);
 *REACTION 0.0;
 *DIFFUSION 1.0;
 *BOUNDARY_XM (x > 0.75);
 *BOUNDARY_XP (x > 0.75);
 *BOUNDARY_YM (x > 0.75);
 *BOUNDARY_YP (x > 0.75);
 *EQUATION_END
 *
 *JUMP_CONDITION_BEGIN varVol
 *INFLUX 0.0;
 *OUTFLUX 0.0;
 *JUMP_CONDITION_END
 *
 *MEMBRANE_END
 * @throws ExpressionException
 * @throws MathException
 */
private void writeMembranes() throws ExpressionException, MathException {
    Simulation simulation = simTask.getSimulation();
    MathDescription mathDesc = simulation.getMathDescription();
    Enumeration<SubDomain> enum1 = mathDesc.getSubDomains();
    while (enum1.hasMoreElements()) {
        SubDomain sd = enum1.nextElement();
        if (sd instanceof MembraneSubDomain) {
            MembraneSubDomain msd = (MembraneSubDomain) sd;
            printWriter.println("MEMBRANE_BEGIN " + msd.getName() + " " + msd.getInsideCompartment().getName() + " " + msd.getOutsideCompartment().getName());
            printWriter.println();
            writeMembrane_VarContext(msd);
            writeMembrane_jumpConditions(msd);
            writeFastSystem(msd);
            printWriter.println("MEMBRANE_END");
            printWriter.println();
        }
    }
}
Also used : CompartmentSubDomain(cbit.vcell.math.CompartmentSubDomain) SubDomain(cbit.vcell.math.SubDomain) MembraneSubDomain(cbit.vcell.math.MembraneSubDomain) MembraneSubDomain(cbit.vcell.math.MembraneSubDomain) Simulation(cbit.vcell.solver.Simulation) MathDescription(cbit.vcell.math.MathDescription)

Aggregations

Simulation (cbit.vcell.solver.Simulation)195 SimulationContext (cbit.vcell.mapping.SimulationContext)57 BioModel (cbit.vcell.biomodel.BioModel)53 MathDescription (cbit.vcell.math.MathDescription)48 KeyValue (org.vcell.util.document.KeyValue)33 Geometry (cbit.vcell.geometry.Geometry)29 MathModel (cbit.vcell.mathmodel.MathModel)27 Expression (cbit.vcell.parser.Expression)26 DataAccessException (org.vcell.util.DataAccessException)26 File (java.io.File)25 ExpressionException (cbit.vcell.parser.ExpressionException)24 IOException (java.io.IOException)24 SimulationJob (cbit.vcell.solver.SimulationJob)23 ArrayList (java.util.ArrayList)23 PropertyVetoException (java.beans.PropertyVetoException)20 UniformOutputTimeSpec (cbit.vcell.solver.UniformOutputTimeSpec)18 XMLSource (cbit.vcell.xml.XMLSource)18 SimulationTask (cbit.vcell.messaging.server.SimulationTask)17 TimeBounds (cbit.vcell.solver.TimeBounds)16 BigString (org.vcell.util.BigString)16