import java.io.InputStream;\r
import java.io.LineNumberReader;\r
import java.util.ArrayList;\r
+import java.util.ConcurrentModificationException;\r
import java.util.HashMap;\r
import java.util.List;\r
import java.util.regex.Matcher;\r
if(resArray.length == 2 && timesArray.length > 2)\r
timesArray = new double[] {timesArray[0], timesArray[timesArray.length - 1]};\r
ds = new DataSet(variable, timesArray, resArray);\r
- for (DataSet tempds : variables){\r
- if (tempds.name.equals(variable))\r
- // We should never need to go there unless some change in\r
- // the logic of reading variables are made.\r
- // Also helps in seeking memory leaks.\r
- return null;\r
+ try {\r
+ for (DataSet tempds : variables){\r
+ if (tempds.name.equals(variable))\r
+ // We should never need to go there unless some change in\r
+ // the logic of reading variables are made.\r
+ // Also helps in seeking memory leaks.\r
+ System.err.println("Find me and add a comment that I was printed, please.");\r
+ // If it seems that in no case we need to go here, comment this try catch.\r
+ return null;\r
+ }\r
+ } catch (ConcurrentModificationException e) {\r
}\r
variables.add(ds);\r
return ds;\r
"http://www.simantics.org/Sysdyn-1.1" : L0.URI
"http://www.simantics.org/Layer0-1.1" : L0.URI
"http://www.simantics.org/SelectionView-1.2" : L0.URI
- "http://www.simantics.org/Documentation-1.1" : L0.URI // Experimental documentation tool
- "http://www.simantics.org/DocumentWorkbench-1.0" : L0.URI // Experimental documentation tool
+// "http://www.simantics.org/Documentation-1.1" : L0.URI // Experimental documentation tool
+// "http://www.simantics.org/DocumentWorkbench-1.0" : L0.URI // Experimental documentation tool
SYSDYN.SharedFunctionOntology <T L0.Ontology
SYSDYN.SharedModuleOntology <T L0.Ontology
>-- SYSDYN.SensitivityAnalysisExperiment.parameterList --> L0.List <R L0.IsComposedOf : L0.FunctionalRelation
>-- SYSDYN.SensitivityAnalysisExperiment.randomSeed --> L0.Integer <R L0.HasProperty : L0.FunctionalRelation
>-- SYSDYN.SensitivityAnalysisExperiment.resultRefreshRate --> L0.Integer <R L0.HasProperty : L0.FunctionalRelation
+ >-- SYSDYN.SensitivityAnalysisExperiment.numberOfValues --> L0.Integer <R L0.HasProperty : L0.FunctionalRelation
@L0.assert SYSDYN.SensitivityAnalysisExperiment.randomSeed 123
+ @L0.assert SYSDYN.SensitivityAnalysisExperiment.numberOfValues 10
SYSDYN.SensitivityAnalysisExperiment.Parameter <T L0.Entity
>-- SYSDYN.SensitivityAnalysisExperiment.Parameter.propabilityDistribution --> SYSDYN.ProbabilityDistribution <R L0.HasProperty : L0.FunctionalRelation
>-- SYSDYN.SensitivityAnalysisExperiment.Parameter.variable --> L0.String <R L0.HasProperty : L0.FunctionalRelation
>-- SYSDYN.SensitivityAnalysisExperiment.Parameter.indexes --> L0.StringArray <R L0.HasProperty : L0.FunctionalRelation
- >-- SYSDYN.SensitivityAnalysisExperiment.Parameter.numberOfValues --> L0.Integer <R L0.HasProperty : L0.FunctionalRelation
-
-
+
SYSDYN.ProbabilityDistribution <T L0.Entity
SYSDYN.UniformDistribution <T SYSDYN.ProbabilityDistribution
>-- SYSDYN.UniformDistribution.minValue --> L0.Double <R L0.HasProperty : L0.FunctionalRelation
public final Resource SensitivityAnalysisExperiment_ParameterBrowseContext_ParameterLabelRule;\r
public final Resource SensitivityAnalysisExperiment_Parameter_indexes;\r
public final Resource SensitivityAnalysisExperiment_Parameter_indexes_Inverse;\r
- public final Resource SensitivityAnalysisExperiment_Parameter_numberOfValues;\r
- public final Resource SensitivityAnalysisExperiment_Parameter_numberOfValues_Inverse;\r
public final Resource SensitivityAnalysisExperiment_Parameter_propabilityDistribution;\r
public final Resource SensitivityAnalysisExperiment_Parameter_propabilityDistribution_Inverse;\r
public final Resource SensitivityAnalysisExperiment_Parameter_variable;\r
public final Resource SensitivityAnalysisExperiment_Parameter_variable_Inverse;\r
+ public final Resource SensitivityAnalysisExperiment_numberOfValues;\r
+ public final Resource SensitivityAnalysisExperiment_numberOfValues_Inverse;\r
public final Resource SensitivityAnalysisExperiment_parameterList;\r
public final Resource SensitivityAnalysisExperiment_parameterList_Inverse;\r
public final Resource SensitivityAnalysisExperiment_randomSeed;\r
public static final String SensitivityAnalysisExperiment_ParameterBrowseContext_ParameterLabelRule = "http://www.simantics.org/Sysdyn-1.1/SensitivityAnalysisExperiment/ParameterBrowseContext/ParameterLabelRule";\r
public static final String SensitivityAnalysisExperiment_Parameter_indexes = "http://www.simantics.org/Sysdyn-1.1/SensitivityAnalysisExperiment/Parameter/indexes";\r
public static final String SensitivityAnalysisExperiment_Parameter_indexes_Inverse = "http://www.simantics.org/Sysdyn-1.1/SensitivityAnalysisExperiment/Parameter/indexes/Inverse";\r
- public static final String SensitivityAnalysisExperiment_Parameter_numberOfValues = "http://www.simantics.org/Sysdyn-1.1/SensitivityAnalysisExperiment/Parameter/numberOfValues";\r
- public static final String SensitivityAnalysisExperiment_Parameter_numberOfValues_Inverse = "http://www.simantics.org/Sysdyn-1.1/SensitivityAnalysisExperiment/Parameter/numberOfValues/Inverse";\r
public static final String SensitivityAnalysisExperiment_Parameter_propabilityDistribution = "http://www.simantics.org/Sysdyn-1.1/SensitivityAnalysisExperiment/Parameter/propabilityDistribution";\r
public static final String SensitivityAnalysisExperiment_Parameter_propabilityDistribution_Inverse = "http://www.simantics.org/Sysdyn-1.1/SensitivityAnalysisExperiment/Parameter/propabilityDistribution/Inverse";\r
public static final String SensitivityAnalysisExperiment_Parameter_variable = "http://www.simantics.org/Sysdyn-1.1/SensitivityAnalysisExperiment/Parameter/variable";\r
public static final String SensitivityAnalysisExperiment_Parameter_variable_Inverse = "http://www.simantics.org/Sysdyn-1.1/SensitivityAnalysisExperiment/Parameter/variable/Inverse";\r
+ public static final String SensitivityAnalysisExperiment_numberOfValues = "http://www.simantics.org/Sysdyn-1.1/SensitivityAnalysisExperiment/numberOfValues";\r
+ public static final String SensitivityAnalysisExperiment_numberOfValues_Inverse = "http://www.simantics.org/Sysdyn-1.1/SensitivityAnalysisExperiment/numberOfValues/Inverse";\r
public static final String SensitivityAnalysisExperiment_parameterList = "http://www.simantics.org/Sysdyn-1.1/SensitivityAnalysisExperiment/parameterList";\r
public static final String SensitivityAnalysisExperiment_parameterList_Inverse = "http://www.simantics.org/Sysdyn-1.1/SensitivityAnalysisExperiment/parameterList/Inverse";\r
public static final String SensitivityAnalysisExperiment_randomSeed = "http://www.simantics.org/Sysdyn-1.1/SensitivityAnalysisExperiment/randomSeed";\r
SensitivityAnalysisExperiment_ParameterBrowseContext_ParameterLabelRule = getResourceOrNull(graph, URIs.SensitivityAnalysisExperiment_ParameterBrowseContext_ParameterLabelRule);\r
SensitivityAnalysisExperiment_Parameter_indexes = getResourceOrNull(graph, URIs.SensitivityAnalysisExperiment_Parameter_indexes);\r
SensitivityAnalysisExperiment_Parameter_indexes_Inverse = getResourceOrNull(graph, URIs.SensitivityAnalysisExperiment_Parameter_indexes_Inverse);\r
- SensitivityAnalysisExperiment_Parameter_numberOfValues = getResourceOrNull(graph, URIs.SensitivityAnalysisExperiment_Parameter_numberOfValues);\r
- SensitivityAnalysisExperiment_Parameter_numberOfValues_Inverse = getResourceOrNull(graph, URIs.SensitivityAnalysisExperiment_Parameter_numberOfValues_Inverse);\r
SensitivityAnalysisExperiment_Parameter_propabilityDistribution = getResourceOrNull(graph, URIs.SensitivityAnalysisExperiment_Parameter_propabilityDistribution);\r
SensitivityAnalysisExperiment_Parameter_propabilityDistribution_Inverse = getResourceOrNull(graph, URIs.SensitivityAnalysisExperiment_Parameter_propabilityDistribution_Inverse);\r
SensitivityAnalysisExperiment_Parameter_variable = getResourceOrNull(graph, URIs.SensitivityAnalysisExperiment_Parameter_variable);\r
SensitivityAnalysisExperiment_Parameter_variable_Inverse = getResourceOrNull(graph, URIs.SensitivityAnalysisExperiment_Parameter_variable_Inverse);\r
+ SensitivityAnalysisExperiment_numberOfValues = getResourceOrNull(graph, URIs.SensitivityAnalysisExperiment_numberOfValues);\r
+ SensitivityAnalysisExperiment_numberOfValues_Inverse = getResourceOrNull(graph, URIs.SensitivityAnalysisExperiment_numberOfValues_Inverse);\r
SensitivityAnalysisExperiment_parameterList = getResourceOrNull(graph, URIs.SensitivityAnalysisExperiment_parameterList);\r
SensitivityAnalysisExperiment_parameterList_Inverse = getResourceOrNull(graph, URIs.SensitivityAnalysisExperiment_parameterList_Inverse);\r
SensitivityAnalysisExperiment_randomSeed = getResourceOrNull(graph, URIs.SensitivityAnalysisExperiment_randomSeed);\r
Resource parameter = GraphUtils.create2(graph, sr.SensitivityAnalysisExperiment_Parameter,\r
sr.SensitivityAnalysisExperiment_Parameter_propabilityDistribution, distribution,\r
sr.SensitivityAnalysisExperiment_Parameter_variable, "<Write variable name>",\r
- sr.SensitivityAnalysisExperiment_Parameter_numberOfValues, 10,\r
L0.PartOf, experiment);\r
\r
ArrayList<Resource> parameterList = new ArrayList<Resource>();\r
+/*******************************************************************************\r
+ * Copyright (c) 2007, 2011 Association for Decentralized Information Management in\r
+ * Industry THTH ry.\r
+ * All rights reserved. This program and the accompanying materials\r
+ * are made available under the terms of the Eclipse Public License v1.0\r
+ * which accompanies this distribution, and is available at\r
+ * http://www.eclipse.org/legal/epl-v10.html\r
+ *\r
+ * Contributors:\r
+ * VTT Technical Research Centre of Finland - initial API and implementation\r
+ *******************************************************************************/\r
package org.simantics.sysdyn.ui.properties;\r
\r
import java.util.Collection;\r
// Label\r
Composite labelComposite = new Composite(content, SWT.NONE);\r
GridDataFactory.fillDefaults().grab(true, false).span(2, 1).applyTo(labelComposite);\r
- GridLayoutFactory.fillDefaults().numColumns(4).applyTo(labelComposite);\r
+ GridLayoutFactory.fillDefaults().numColumns(6).applyTo(labelComposite);\r
Label label = new Label(labelComposite, SWT.NONE);\r
label.setText("Name");\r
\r
name.setColorProvider(new SysdynBasicColorProvider(new LocalResourceManager(JFaceResources.getResources(), name.getWidget())));\r
GridDataFactory.fillDefaults().grab(true, false).applyTo(name.getWidget());\r
\r
+ label = new Label(labelComposite, SWT.NONE);\r
+ label.setText("Number of runs");\r
+ \r
+ TrackedText n = new TrackedText(labelComposite, support, SWT.BORDER);\r
+ n.setTextFactory(new IntegerPropertyFactory(SysdynResource.URIs.SensitivityAnalysisExperiment_numberOfValues));\r
+ n.addModifyListener(new IntegerPropertyModifier(context, SysdynResource.URIs.SensitivityAnalysisExperiment_numberOfValues));\r
+ n.setInputValidator(new IntegerValidator());\r
+ n.setColorProvider(new SysdynBasicColorProvider(new LocalResourceManager(JFaceResources.getResources(), n.getWidget())));\r
+ GridDataFactory.fillDefaults().hint(80, SWT.DEFAULT).applyTo(n.getWidget());\r
+ \r
label = new Label(labelComposite, SWT.NONE);\r
label.setText("Seed");\r
\r
seed.setColorProvider(new SysdynBasicColorProvider(new LocalResourceManager(JFaceResources.getResources(), seed.getWidget())));\r
GridDataFactory.fillDefaults().hint(80, SWT.DEFAULT).applyTo(seed.getWidget());\r
\r
-\r
// (Ontology-based) GraphExplorer displaying range axis and variables mapped to those axis\r
explorer = new AxisAndVariablesExplorerComposite(ArrayMap.keys(\r
"displaySelectors", "displayFilter", "treeView").values(false, false, true), site, content, support, SWT.FULL_SELECTION | SWT.BORDER | SWT.SINGLE);\r
// variable.setColorProvider(new SysdynBasicColorProvider(new LocalResourceManager(JFaceResources.getResources(), variable.getWidget())));\r
// GridDataFactory.fillDefaults().grab(true, false).applyTo(variable.getWidget());\r
\r
- label = new Label(parameterProperties, SWT.NONE);\r
- label.setText("Number of values:");\r
- GridDataFactory.fillDefaults().align(SWT.END, SWT.CENTER).applyTo(label);\r
-\r
-\r
- TrackedText numValues = new TrackedText(parameterProperties, parameterSupport, SWT.BORDER);\r
- numValues.setTextFactory(new IntegerPropertyFactory(SysdynResource.URIs.SensitivityAnalysisExperiment_Parameter_numberOfValues));\r
- numValues.addModifyListener(new IntegerPropertyModifier(context, SysdynResource.URIs.SensitivityAnalysisExperiment_Parameter_numberOfValues));\r
- numValues.setInputValidator(new IntegerValidator());\r
- numValues.setColorProvider(new SysdynBasicColorProvider(new LocalResourceManager(JFaceResources.getResources(), numValues.getWidget())));\r
- GridDataFactory.fillDefaults().align(SWT.BEGINNING, SWT.CENTER).hint(80, SWT.DEFAULT).applyTo(numValues.getWidget());\r
- \r
label = new Label(parameterProperties, SWT.NONE);\r
label.setText("Distribution:");\r
GridDataFactory.fillDefaults().align(SWT.END, SWT.CENTER).applyTo(label);\r
Resource parameter = GraphUtils.create2(graph, sr.SensitivityAnalysisExperiment_Parameter,\r
sr.SensitivityAnalysisExperiment_Parameter_propabilityDistribution, distribution,\r
sr.SensitivityAnalysisExperiment_Parameter_variable, "<Write variable name>",\r
- sr.SensitivityAnalysisExperiment_Parameter_numberOfValues, 10,\r
L0.PartOf, input);\r
\r
Resource parameterList = graph.getPossibleObject(input, sr.SensitivityAnalysisExperiment_parameterList);\r
*\r
* Contributors:\r
* Semantum Oy - initial API and implementation\r
+ * VTT Technical Research Centre of Finland\r
*******************************************************************************/\r
package org.simantics.sysdyn.ui.properties.widgets.sensitivity;\r
\r
HashMap<String, String> result = new HashMap<String, String>();\r
\r
String variable = graph.getPossibleRelatedValue((Resource)content, SysdynResource.getInstance(graph).SensitivityAnalysisExperiment_Parameter_variable);\r
- Integer n = graph.getPossibleRelatedValue((Resource)content, SysdynResource.getInstance(graph).SensitivityAnalysisExperiment_Parameter_numberOfValues);\r
\r
StringBuilder sb = new StringBuilder();\r
\r
else\r
sb.append("No variable");\r
\r
- sb.append(" (n=");\r
- \r
- if(n != null)\r
- sb.append(n);\r
- else\r
- sb.append("undefined");\r
- \r
- sb.append(")");\r
- \r
result.put(ColumnKeys.SINGLE, sb.toString()) ;\r
\r
return result;\r
*\r
* Contributors:\r
* Semantum Oy - initial API and implementation\r
+ * VTT Technical Research Centre of Finland\r
*******************************************************************************/\r
package org.simantics.sysdyn.adapter.distribution;\r
\r
public interface IDistribution {\r
- \r
- /**\r
- * Return next value of the distribution\r
- * @return\r
- */\r
- public double next();\r
- \r
+\r
/**\r
+ * Map a probability to the random variable.\r
+ * When random numbers are used to generate the Monte Carlo simulation parameters, the probability\r
+ * is itself selected randomly (between (0,1); at the domain endpoints the value is undefined).\r
* \r
- * @return does the distribution have next value;\r
+ * @param randomVariable\r
+ * @return Inverse cumulative distribution function at probability, i.e. the random variable at the\r
+ * point where the cumulative distribution function yields probability.\r
*/\r
- public boolean hasNext();\r
- \r
- /**\r
- * Initialize distribution\r
- */\r
- public void initialize();\r
-\r
+ public double inverseCDF(double probability);\r
}\r
*\r
* Contributors:\r
* Semantum Oy - initial API and implementation\r
+ * VTT Technical Research Centre of Finland\r
*******************************************************************************/\r
package org.simantics.sysdyn.adapter.distribution;\r
\r
-import java.util.Random;\r
-\r
import org.simantics.databoard.Bindings;\r
import org.simantics.db.ReadGraph;\r
import org.simantics.db.Resource;\r
\r
public class Interval implements IDistribution {\r
\r
- private Random random;\r
- \r
private double min = 0;\r
private double max = 10;\r
private double intervalLength = 1;\r
private int numberOfValues = 10;\r
- private int seed = 123;\r
private int index = 0;\r
\r
public Interval(ReadGraph graph, Resource distribution) {\r
Double maxValue = graph.getPossibleRelatedValue(distribution, SR.Interval_maxValue, Bindings.DOUBLE);\r
if(maxValue != null)\r
this.max = maxValue;\r
- \r
- \r
+ \r
Resource parameter = graph.getPossibleObject(distribution, SR.SensitivityAnalysisExperiment_Parameter_propabilityDistribution_Inverse);\r
+ \r
+ Resource experiment = graph.getPossibleObject(parameter, Layer0.getInstance(graph).PartOf);\r
\r
- Integer numberOfValues = graph.getPossibleRelatedValue(parameter, SR.SensitivityAnalysisExperiment_Parameter_numberOfValues, Bindings.INTEGER);\r
+ Integer numberOfValues = graph.getPossibleRelatedValue(experiment, SR.SensitivityAnalysisExperiment_numberOfValues, Bindings.INTEGER);\r
if(numberOfValues != null)\r
this.numberOfValues = numberOfValues;\r
\r
intervalLength = (maxValue - minValue) / (this.numberOfValues - 1);\r
\r
- Resource experiment = graph.getPossibleObject(parameter, Layer0.getInstance(graph).PartOf);\r
-\r
- Integer seed = graph.getPossibleRelatedValue(experiment, SR.SensitivityAnalysisExperiment_randomSeed, Bindings.INTEGER);\r
- if(seed != null)\r
- this.seed = seed;\r
- \r
} catch (DatabaseException e) {\r
\r
}\r
- \r
- initialize();\r
- }\r
- \r
- @Override\r
- public double next() {\r
- double value = min + (intervalLength * index);\r
- index++;\r
- return value;\r
- }\r
-\r
- public Random getRandom() {\r
- return random;\r
}\r
\r
public double getMin() {\r
return numberOfValues;\r
}\r
\r
- public int getSeed() {\r
- return seed;\r
- }\r
- \r
public double getIntervalLength() {\r
return intervalLength;\r
}\r
\r
- @Override\r
- public void initialize() {\r
- random = new Random(this.seed);\r
- index = 0;\r
- }\r
-\r
- @Override\r
- public boolean hasNext() {\r
- return index < numberOfValues;\r
- }\r
+ @Override\r
+ public double inverseCDF(double probability) {\r
+ // This is a bit awkward...\r
+ double value = min + (intervalLength * index);\r
+ index++;\r
+ if (index >= numberOfValues)\r
+ index = 0;\r
+ return value;\r
+ }\r
\r
}\r
*\r
* Contributors:\r
* Semantum Oy - initial API and implementation\r
+ * VTT Technical Research Centre of Finland\r
*******************************************************************************/\r
package org.simantics.sysdyn.adapter.distribution;\r
\r
\r
-import java.util.Random;\r
-\r
import org.simantics.databoard.Bindings;\r
import org.simantics.db.ReadGraph;\r
import org.simantics.db.Resource;\r
import org.simantics.db.exception.DatabaseException;\r
-import org.simantics.layer0.Layer0;\r
import org.simantics.sysdyn.SysdynResource;\r
\r
public class NormalDistribution implements IDistribution {\r
- \r
- private Random random;\r
- \r
+\r
private Double min = null;\r
private Double max = null;\r
+ private double minProbability = 0;\r
+ private double maxProbability = 1;\r
private double mean;\r
private double stdDeviation; \r
- private int numberOfValues = 10;\r
- private int seed = 123;\r
- private int index = 0;\r
- \r
- private static int MAX_TRIES = 100;\r
-\r
\r
public NormalDistribution(ReadGraph graph, Resource distribution) {\r
\r
try {\r
+\r
SysdynResource SR = SysdynResource.getInstance(graph);\r
\r
min = graph.getPossibleRelatedValue(distribution, SR.NormalDistribution_minValue, Bindings.DOUBLE);\r
stdDeviation = graph.getPossibleRelatedValue(distribution, SR.NormalDistribution_stdDeviation, Bindings.DOUBLE);\r
mean = graph.getPossibleRelatedValue(distribution, SR.NormalDistribution_mean, Bindings.DOUBLE);\r
\r
- \r
- Resource parameter = graph.getPossibleObject(distribution, SR.SensitivityAnalysisExperiment_Parameter_propabilityDistribution_Inverse);\r
+ // Determine the max and min probabilities.\r
+ if (min != null)\r
+ minProbability = cdf(min, mean, stdDeviation, false);\r
+ if (max != null)\r
+ maxProbability = cdf(max, mean, stdDeviation, false);\r
\r
- Integer numberOfValues = graph.getPossibleRelatedValue(parameter, SR.SensitivityAnalysisExperiment_Parameter_numberOfValues, Bindings.INTEGER);\r
- if(numberOfValues != null)\r
- this.numberOfValues = numberOfValues;\r
- \r
-\r
- \r
- Resource experiment = graph.getPossibleObject(parameter, Layer0.getInstance(graph).PartOf);\r
-\r
- Integer seed = graph.getPossibleRelatedValue(experiment, SR.SensitivityAnalysisExperiment_randomSeed, Bindings.INTEGER);\r
- if(seed != null)\r
- this.seed = seed;\r
- \r
- \r
- random = new Random(this.seed);\r
- \r
} catch (DatabaseException e) {\r
\r
}\r
}\r
\r
- @Override\r
- public double next() {\r
-\r
- /* Try to find value between min and max, if they are defined.\r
- * If value cannot be found after MAX_TRIES number of tries,\r
- * return the last try. \r
- */\r
- double value = 0;\r
- for(int i = 0; i < MAX_TRIES; i++) {\r
- value = random.nextGaussian() * stdDeviation + mean;\r
- if(accept(value))\r
- break;\r
- }\r
- \r
- index++;\r
- return value;\r
- }\r
- \r
- \r
- private boolean accept(double value) {\r
- if(min != null) {\r
- if(value < min)\r
- return false;\r
- }\r
- \r
- if(max != null) {\r
- if(value > max)\r
- return false;\r
- }\r
- \r
- return true;\r
- }\r
- \r
- public Random getRandom() {\r
- return random;\r
- }\r
-\r
public double getMin() {\r
return min;\r
}\r
public double getMax() {\r
return max;\r
}\r
+ \r
+ @Override\r
+ public double inverseCDF(double probability) {\r
+ // Map probability to [min,max]\r
+ double mappedProbability = minProbability + (maxProbability - minProbability) * probability;\r
+ // Calculate the inverse CDF\r
+ return quantile(mappedProbability, mean, stdDeviation);\r
+ }\r
\r
- public int getNumberOfValues() {\r
- return numberOfValues;\r
- }\r
-\r
- public int getSeed() {\r
- return seed;\r
+ //The following is from BEAST-MCMC (https://code.google.com/p/beast-mcmc/)\r
+ \r
+ /*\r
+ * NormalDistribution.java\r
+ *\r
+ * Copyright (c) 2002-2011 Alexei Drummond, Andrew Rambaut and Marc Suchard\r
+ *\r
+ * This file is part of BEAST.\r
+ * See the NOTICE file distributed with this work for additional\r
+ * information regarding copyright ownership and licensing.\r
+ *\r
+ * BEAST is free software; you can redistribute it and/or modify\r
+ * it under the terms of the GNU Lesser General Public License as\r
+ * published by the Free Software Foundation; either version 2\r
+ * of the License, or (at your option) any later version.\r
+ *\r
+ * BEAST is distributed in the hope that it will be useful,\r
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of\r
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the\r
+ * GNU Lesser General Public License for more details.\r
+ *\r
+ * You should have received a copy of the GNU Lesser General Public\r
+ * License along with BEAST; if not, write to the\r
+ * Free Software Foundation, Inc., 51 Franklin St, Fifth Floor,\r
+ * Boston, MA 02110-1301 USA\r
+ */\r
+ \r
+ /**\r
+ * quantiles (=inverse cumulative density function)\r
+ *\r
+ * @param z argument\r
+ * @param m mean\r
+ * @param sd standard deviation\r
+ * @return icdf at z\r
+ */\r
+ public static double quantile(double z, double m, double sd) {\r
+ return m + Math.sqrt(2.0) * sd * inverseErf(2.0 * z - 1.0);\r
}\r
\r
+ /** A more accurate and faster implementation of the cdf (taken from function pnorm in the R statistical language)\r
+ * This implementation has discrepancies depending on the programming language and system architecture\r
+ * In Java, returned values become zero once z reaches -37.5193 exactly on the machine tested\r
+ * In the other implementation, the returned value 0 at about z = -8\r
+ * In C, this 0 value is reached approximately z = -37.51938\r
+ *\r
+ * Will later need to be optimised for BEAST\r
+ *\r
+ * @param x argument\r
+ * @param mu mean\r
+ * @param sigma standard deviation\r
+ * @param log_p is p logged\r
+ * @return cdf at x\r
+ */\r
+ public static double cdf(double x, double mu, double sigma, boolean log_p) {\r
+ boolean i_tail=false;\r
+ double p, cp = Double.NaN;\r
\r
- @Override\r
- public void initialize() {\r
- random = new Random(this.seed);\r
- index = 0;\r
- }\r
+ if(Double.isNaN(x) || Double.isNaN(mu) || Double.isNaN(sigma)) {\r
+ return Double.NaN;\r
+ }\r
+ if(Double.isInfinite(x) && mu == x) { /* x-mu is NaN */\r
+ return Double.NaN;\r
+ }\r
+ if (sigma <= 0) {\r
+ if(sigma < 0) {\r
+ return Double.NaN;\r
+ }\r
+ return (x < mu) ? 0.0 : 1.0;\r
+ }\r
+ p = (x - mu) / sigma;\r
+ if(Double.isInfinite(p)) {\r
+ return (x < mu) ? 0.0 : 1.0;\r
+ }\r
+ x = p;\r
+ if(Double.isNaN(x)) {\r
+ return Double.NaN;\r
+ }\r
+\r
+ double xden, xnum, temp, del, eps, xsq, y;\r
+ int i;\r
+ boolean lower, upper;\r
+ eps = DBL_EPSILON * 0.5;\r
+ lower = !i_tail;\r
+ upper = i_tail;\r
+\r
+ y = Math.abs(x);\r
+ if (y <= 0.67448975) { /* Normal.quantile(3/4, 1, 0) = 0.67448975 */\r
+ if (y > eps) {\r
+ xsq = x * x;\r
+ xnum = a[4] * xsq;\r
+ xden = xsq;\r
+ for (i = 0; i < 3; i++) {\r
+ xnum = (xnum + a[i]) * xsq;\r
+ xden = (xden + b[i]) * xsq;\r
+ }\r
+ }\r
+ else {\r
+ xnum = xden = 0.0;\r
+ }\r
+ temp = x * (xnum + a[3]) / (xden + b[3]);\r
+ if(lower) {\r
+ p = 0.5 + temp;\r
+ }\r
+ if(upper) {\r
+ cp = 0.5 - temp;\r
+ }\r
+ if(log_p) {\r
+ if(lower) {\r
+ p = Math.log(p);\r
+ }\r
+ if(upper) {\r
+ cp = Math.log(cp);\r
+ }\r
+ }\r
+ }\r
+\r
+\r
+ else if (y <= M_SQRT_32) {\r
+ /* Evaluate pnorm for 0.67448975 = Normal.quantile(3/4, 1, 0) < |x| <= sqrt(32) ~= 5.657 */\r
+\r
+ xnum = c[8] * y;\r
+ xden = y;\r
+ for (i = 0; i < 7; i++) {\r
+ xnum = (xnum + c[i]) * y;\r
+ xden = (xden + d[i]) * y;\r
+ }\r
+ temp = (xnum + c[7]) / (xden + d[7]);\r
+\r
+ //do_del(y);\r
+ //swap_tail;\r
+ //#define do_del(X) \\r
+ xsq = ((int) (y * CUTOFF)) * 1.0 / CUTOFF;\r
+ del = (y - xsq) * (y + xsq);\r
+ if(log_p) {\r
+ p = (-xsq * xsq * 0.5) + (-del * 0.5) + Math.log(temp);\r
+ if((lower && x > 0.0) || (upper && x <= 0.0)) {\r
+ cp = Math.log(1.0-Math.exp(-xsq * xsq * 0.5) * Math.exp(-del * 0.5) * temp);\r
+ }\r
+ }\r
+ else {\r
+ p = Math.exp(-xsq * xsq * 0.5) * Math.exp(-del * 0.5) * temp;\r
+ cp = 1.0 - p;\r
+ }\r
+ //#define swap_tail \\r
+ if (x > 0.0) {\r
+ temp = p;\r
+ if(lower) {\r
+ p = cp;\r
+ }\r
+ cp = temp;\r
+ }\r
+ }\r
+ /* else |x| > sqrt(32) = 5.657 :\r
+ * the next two case differentiations were really for lower=T, log=F\r
+ * Particularly *not* for log_p !\r
+ * Cody had (-37.5193 < x && x < 8.2924) ; R originally had y < 50\r
+ * Note that we do want symmetry(0), lower/upper -> hence use y\r
+ */\r
+ else if(log_p || (lower && -37.5193 < x && x < 8.2924)\r
+ || (upper && -8.2924 < x && x < 37.5193)) {\r
+\r
+ /* Evaluate pnorm for x in (-37.5, -5.657) union (5.657, 37.5) */\r
+ xsq = 1.0 / (x * x);\r
+ xnum = p_[5] * xsq;\r
+ xden = xsq;\r
+ for (i = 0; i < 4; i++) {\r
+ xnum = (xnum + p_[i]) * xsq;\r
+ xden = (xden + q[i]) * xsq;\r
+ }\r
+ temp = xsq * (xnum + p_[4]) / (xden + q[4]);\r
+ temp = (M_1_SQRT_2PI - temp) / y;\r
+\r
+ //do_del(x);\r
+ xsq = ((int) (x * CUTOFF)) * 1.0 / CUTOFF;\r
+ del = (x - xsq) * (x + xsq);\r
+ if(log_p) {\r
+ p = (-xsq * xsq * 0.5) + (-del * 0.5) + Math.log(temp);\r
+ if((lower && x > 0.0) || (upper && x <= 0.0)) {\r
+ cp = Math.log(1.0-Math.exp(-xsq * xsq * 0.5) * Math.exp(-del * 0.5) * temp);\r
+ }\r
+ }\r
+ else {\r
+ p = Math.exp(-xsq * xsq * 0.5) * Math.exp(-del * 0.5) * temp;\r
+ cp = 1.0 - p;\r
+ }\r
+ //swap_tail;\r
+ if (x > 0.0) {\r
+ temp = p;\r
+ if(lower) {\r
+ p = cp;\r
+ }\r
+ cp = temp;\r
+ }\r
+ }\r
+ else { /* no log_p , large x such that probs are 0 or 1 */\r
+ if(x > 0) {\r
+ p = 1.0;\r
+ cp = 0.0;\r
+ }\r
+ else {\r
+ p = 0.0;\r
+ cp = 1.0;\r
+ }\r
+ }\r
+ return p;\r
\r
- @Override\r
- public boolean hasNext() {\r
- return index < numberOfValues;\r
}\r
\r
+ // Private\r
+\r
+ protected double m, sd;\r
+\r
+ private static final double[] a = {\r
+ 2.2352520354606839287,\r
+ 161.02823106855587881,\r
+ 1067.6894854603709582,\r
+ 18154.981253343561249,\r
+ 0.065682337918207449113\r
+ };\r
+ private static final double[] b = {\r
+ 47.20258190468824187,\r
+ 976.09855173777669322,\r
+ 10260.932208618978205,\r
+ 45507.789335026729956\r
+ };\r
+ private static final double[] c = {\r
+ 0.39894151208813466764,\r
+ 8.8831497943883759412,\r
+ 93.506656132177855979,\r
+ 597.27027639480026226,\r
+ 2494.5375852903726711,\r
+ 6848.1904505362823326,\r
+ 11602.651437647350124,\r
+ 9842.7148383839780218,\r
+ 1.0765576773720192317e-8\r
+ };\r
+ private static final double[] d = {\r
+ 22.266688044328115691,\r
+ 235.38790178262499861,\r
+ 1519.377599407554805,\r
+ 6485.558298266760755,\r
+ 18615.571640885098091,\r
+ 34900.952721145977266,\r
+ 38912.003286093271411,\r
+ 19685.429676859990727\r
+ };\r
+ private static final double[] p_ = {\r
+ 0.21589853405795699,\r
+ 0.1274011611602473639,\r
+ 0.022235277870649807,\r
+ 0.001421619193227893466,\r
+ 2.9112874951168792e-5,\r
+ 0.02307344176494017303\r
+ };\r
+ private static final double[] q = {\r
+ 1.28426009614491121,\r
+ 0.468238212480865118,\r
+ 0.0659881378689285515,\r
+ 0.00378239633202758244,\r
+ 7.29751555083966205e-5\r
+ };\r
+\r
+ private static final int CUTOFF = 16; /* Cutoff allowing exact "*" and "/" */\r
+\r
+ private static final double M_SQRT_32 = 5.656854249492380195206754896838; /* The square root of 32 */\r
+ private static final double M_1_SQRT_2PI = 0.398942280401432677939946059934;\r
+ private static final double DBL_EPSILON = 2.2204460492503131e-016;\r
+ \r
+ /*\r
+ * ErrorFunction.java\r
+ *\r
+ * Copyright (C) 2002-2006 Alexei Drummond and Andrew Rambaut\r
+ *\r
+ * This file is part of BEAST.\r
+ * See the NOTICE file distributed with this work for additional\r
+ * information regarding copyright ownership and licensing.\r
+ *\r
+ * BEAST is free software; you can redistribute it and/or modify\r
+ * it under the terms of the GNU Lesser General Public License as\r
+ * published by the Free Software Foundation; either version 2\r
+ * of the License, or (at your option) any later version.\r
+ *\r
+ * BEAST is distributed in the hope that it will be useful,\r
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of\r
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the\r
+ * GNU Lesser General Public License for more details.\r
+ *\r
+ * You should have received a copy of the GNU Lesser General Public\r
+ * License along with BEAST; if not, write to the\r
+ * Free Software Foundation, Inc., 51 Franklin St, Fifth Floor,\r
+ * Boston, MA 02110-1301 USA\r
+ */\r
+ \r
+ /**\r
+ * inverse error function\r
+ *\r
+ * @param z argument\r
+ *\r
+ * @return function value\r
+ */\r
+ public static double inverseErf(double z)\r
+ {\r
+ return pointNormal(0.5*z+0.5)/Math.sqrt(2.0);\r
+ }\r
+ \r
+ // Private\r
+ \r
+ // Returns z so that Prob{x<z}=prob where x ~ N(0,1) and (1e-12) < prob<1-(1e-12)\r
+ private static double pointNormal(double prob)\r
+ {\r
+ // Odeh RE & Evans JO (1974) The percentage points of the normal distribution.\r
+ // Applied Statistics 22: 96-97 (AS70)\r
+ \r
+ // Newer methods:\r
+ // Wichura MJ (1988) Algorithm AS 241: the percentage points of the\r
+ // normal distribution. 37: 477-484.\r
+ // Beasley JD & Springer SG (1977). Algorithm AS 111: the percentage \r
+ // points of the normal distribution. 26: 118-121.\r
+ \r
+ double a0 = -0.322232431088, a1 = -1, a2 = -0.342242088547, a3 = -0.0204231210245;\r
+ double a4 = -0.453642210148e-4, b0 = 0.0993484626060, b1 = 0.588581570495;\r
+ double b2 = 0.531103462366, b3 = 0.103537752850, b4 = 0.0038560700634;\r
+ double y, z = 0, p = prob, p1;\r
+ \r
+ p1 = (p < 0.5 ? p : 1-p);\r
+ if (p1 < 1e-20)\r
+ {\r
+ new IllegalArgumentException("Argument prob out of range");\r
+ }\r
+ \r
+ y = Math.sqrt(Math.log(1/(p1*p1))); \r
+ z = y + ((((y*a4+a3)*y+a2)*y+a1)*y+a0)/((((y*b4+b3)*y+b2)*y+b1)*y+b0);\r
+ return (p < 0.5 ? -z : z);\r
+ }\r
+\r
}\r
*\r
* Contributors:\r
* Semantum Oy - initial API and implementation\r
+ * VTT Technical Research Centre of Finland\r
*******************************************************************************/\r
package org.simantics.sysdyn.adapter.distribution;\r
\r
-import java.util.Random;\r
-\r
import org.simantics.databoard.Bindings;\r
import org.simantics.db.ReadGraph;\r
import org.simantics.db.Resource;\r
import org.simantics.db.exception.DatabaseException;\r
-import org.simantics.layer0.Layer0;\r
import org.simantics.sysdyn.SysdynResource;\r
\r
public class UniformDistribution implements IDistribution {\r
\r
- \r
- private Random random;\r
- \r
private double min;\r
private double max;\r
- private int numberOfValues = 10;\r
- private int seed = 123;\r
- private int index = 0;\r
\r
public UniformDistribution(ReadGraph graph, Resource distribution) {\r
\r
if(maxValue != null)\r
this.max = maxValue;\r
\r
- \r
- Resource parameter = graph.getPossibleObject(distribution, SR.SensitivityAnalysisExperiment_Parameter_propabilityDistribution_Inverse);\r
-\r
- Integer numberOfValues = graph.getPossibleRelatedValue(parameter, SR.SensitivityAnalysisExperiment_Parameter_numberOfValues, Bindings.INTEGER);\r
- if(numberOfValues != null)\r
- this.numberOfValues = numberOfValues;\r
- \r
-\r
- \r
- Resource experiment = graph.getPossibleObject(parameter, Layer0.getInstance(graph).PartOf);\r
-\r
- Integer seed = graph.getPossibleRelatedValue(experiment, SR.SensitivityAnalysisExperiment_randomSeed, Bindings.INTEGER);\r
- if(seed != null)\r
- this.seed = seed;\r
- \r
} catch (DatabaseException e) {\r
\r
}\r
- initialize();\r
- }\r
- \r
- public Random getRandom() {\r
- return random;\r
}\r
\r
public double getMin() {\r
return max;\r
}\r
\r
- public int getNumberOfValues() {\r
- return numberOfValues;\r
- }\r
-\r
- public int getSeed() {\r
- return seed;\r
- }\r
- \r
- \r
-\r
- \r
@Override\r
- public double next() {\r
- index++;\r
- double d = min + ((double)random.nextInt(101) / 100.0) * (max-min); \r
- return d;\r
- }\r
- \r
-\r
- @Override\r
- public void initialize() {\r
- random = new Random(this.seed);\r
- index = 0;\r
- }\r
-\r
- @Override\r
- public boolean hasNext() {\r
- return index < numberOfValues;\r
- }\r
+ public double inverseCDF(double probability) {\r
+ return min + probability * (max - min); \r
+ }\r
\r
}\r
--- /dev/null
+/*******************************************************************************\r
+ * Copyright (c) 2013 Association for Decentralized Information Management in\r
+ * Industry THTH ry.\r
+ * All rights reserved. This program and the accompanying materials\r
+ * are made available under the terms of the Eclipse Public License v1.0\r
+ * which accompanies this distribution, and is available at\r
+ * http://www.eclipse.org/legal/epl-v10.html\r
+ *\r
+ * Contributors:\r
+ * VTT Technical Research Centre of Finland - initial API and implementation\r
+ *******************************************************************************/\r
+package org.simantics.sysdyn.manager;\r
+\r
+import java.util.ArrayList;\r
+import java.util.Random;\r
+\r
+/**\r
+ * Random generator that creates ArrayLists of random Doubles between 0 and 1.\r
+ * \r
+ * @author Tuomas Miettinen\r
+ */\r
+\r
+public class RandomGenerator {\r
+ \r
+ private int seed;\r
+ private int dimension;\r
+ private Random random;\r
+ \r
+ public RandomGenerator(int seed, int dimension) {\r
+ this.dimension = dimension;\r
+ this.seed = seed;\r
+ random = new Random(this.seed);\r
+ }\r
+ \r
+ public ArrayList<Double> next() {\r
+ ArrayList<Double> randoms = new ArrayList<Double>();\r
+ for (int i = 0; i < dimension; ++i)\r
+ randoms.add(random.nextDouble());\r
+ return randoms;\r
+ }\r
+ \r
+ public void initialize() {\r
+ random = new Random(this.seed);\r
+ }\r
+\r
+}\r
*/\r
public class SysdynSensitivityAnalysisExperiment extends SysdynExperiment {\r
\r
- \r
-\r
private ArrayList<MemoryResult> results = null;\r
-\r
+ private int seed = 124;\r
+ private int numberOfRuns = 0;\r
+ \r
public SysdynSensitivityAnalysisExperiment(Resource experiment, Resource model) {\r
super(experiment, model);\r
}\r
\r
- \r
private void findValuesAndRun(List<SensitivityExperimentParameter> parameters, HashMap<String, String> values, HashMap<String, String> experimentParameters) {\r
- findValuesAndRun(parameters, 0, values, experimentParameters);\r
- }\r
- \r
- private void findValuesAndRun(List<SensitivityExperimentParameter> parameters, int index, HashMap<String, String> values, HashMap<String, String> experimentParameters) {\r
- SensitivityExperimentParameter p = parameters.get(index);\r
+ ArrayList<HashMap<String, String>> parameterMaps = new ArrayList<HashMap<String, String>>();\r
+ int parametersSize = parameters.size();\r
+\r
+ // Get seed\r
+ try {\r
+ session.syncRequest(new ReadRequest() {\r
+ @Override\r
+ public void run(ReadGraph graph) throws DatabaseException {\r
+ SysdynResource SR = SysdynResource.getInstance(graph);\r
+ seed = graph.getPossibleRelatedValue(experiment, SR.SensitivityAnalysisExperiment_randomSeed, Bindings.INTEGER);\r
+ numberOfRuns = graph.getPossibleRelatedValue(experiment, SR.SensitivityAnalysisExperiment_numberOfValues, Bindings.INTEGER);\r
+ //refreshRate = numberOfRuns / 20 + 1;\r
+ refreshRate = 1;\r
+ }\r
+ });\r
+ } catch (DatabaseException e) {\r
+ e.printStackTrace();\r
+ }\r
\r
- while(p.getDistribution().hasNext()) {\r
- values.put(p.getFullName(), Double.toString(p.getDistribution().next()));\r
- if((index + 1) < parameters.size()) {\r
- findValuesAndRun(parameters, index + 1, values, experimentParameters);\r
- } else {\r
- // Run with values\r
- runSensitivityRun(values, experimentParameters);\r
- }\r
- }\r
- p.getDistribution().initialize();\r
+ // Generate values for parameters (randomly)\r
+ RandomGenerator randomGenerator = new RandomGenerator(seed, parametersSize);\r
+ \r
+ // Determine the changed paramters and their values for each individual run\r
+ for (int i = 0; i < numberOfRuns; ++i) {\r
+ ArrayList<Double> randoms = randomGenerator.next();// multi-dimensional random\r
+ HashMap<String, String> parameterMap = new HashMap<String, String>(parametersSize);\r
+ for (int j = 0; j < parametersSize; ++j) {\r
+ SensitivityExperimentParameter p = parameters.get(j);\r
+ // Get the random value\r
+ double value = p.getDistribution().inverseCDF(randoms.get(j));\r
+ // Add the new parameter-value-pair \r
+ parameterMap.put(p.getFullName(), Double.toString(value));\r
+ }\r
+ // Add the complete list of one individual run parameter-value-pairs.\r
+ parameterMaps.add(parameterMap);\r
+ }\r
+ \r
+ // Set the parameters for each individual experiment and run them.\r
+ for (HashMap<String, String> parameterSet : parameterMaps) {\r
+ values.putAll(parameterSet);\r
+ runSensitivityRun(values, experimentParameters);\r
+ }\r
}\r
\r
private void runSensitivityRun(HashMap<String, String> values, HashMap<String, String> experimentParameters) {\r
ModelicaManager.printProcessOutput(process, monitor);\r
\r
File resFile = new File(experimentParameters.get(ModelicaManager.RESULT_FILE_NAME));\r
- Thread resultThread = getResultThread(resFile, experimentParameters, monitor, progressMonitor, currentRun % refreshRate == 0);\r
+ boolean updateMonitors = currentRun % refreshRate == 0 || currentRun == numberOfRuns - 1;\r
+ Thread resultThread = getResultThread(resFile, experimentParameters, monitor, progressMonitor, updateMonitors);\r
resultThread.run();\r
\r
process = null;\r
private IModelicaMonitor monitor;\r
private IProgressMonitor progressMonitor;\r
// private HashMap<String, String> experimentParameters;\r
- private int numberOfRuns = 0;\r
private int currentRun = 0;\r
private int refreshRate = 1;\r
private List<SensitivityExperimentParameter> parameters = new ArrayList<SensitivityExperimentParameter>();\r
@Override\r
public Integer perform(ReadGraph graph) throws DatabaseException {\r
SysdynResource SR = SysdynResource.getInstance(graph);\r
- \r
- Resource parameterListResource = graph.getPossibleObject(experiment, SR.SensitivityAnalysisExperiment_parameterList);\r
- List<Resource> parameterResources = ListUtils.toList(graph, parameterListResource);\r
- \r
- Integer numberOfIterations = null;\r
- \r
- for(Resource parameter : parameterResources) {\r
- Integer n = graph.getPossibleRelatedValue(parameter, SR.SensitivityAnalysisExperiment_Parameter_numberOfValues, Bindings.INTEGER);\r
- if(n != null) {\r
- if(numberOfIterations == null)\r
- numberOfIterations = n;\r
- else\r
- numberOfIterations = numberOfIterations * n;\r
- } \r
- }\r
+ Integer numberOfIterations = graph.getPossibleRelatedValue(experiment, SR.SensitivityAnalysisExperiment_numberOfValues, Bindings.INTEGER);\r
\r
return numberOfIterations;\r
}\r
Code Review / simantics / sysdyn.git / commitdiff