Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Uncertainty fixes and user experience improvement #1593

Merged
merged 32 commits into from
Jun 7, 2019
Merged

Uncertainty fixes and user experience improvement #1593

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
32 commits
Select commit Hold shift + click to select a range
1b88161
Automatically format rmgpy.tools.uncertainty
mliu49 Apr 11, 2019
ad168ba
Automatically format rmgpy.tools.muq
mliu49 Apr 11, 2019
ea79bcd
Adjust file paths for local uncertainty analysis
mliu49 Jul 10, 2017
e1bbf01
Automatically retrieve internal indices
mliu49 Jun 20, 2018
a27e59a
Return uncertainty results
mliu49 Jun 2, 2018
abc551c
Retrieve uncertainties for uncorrelated global analysis
mliu49 Jan 29, 2019
cfb7607
Perform PCE fitting in ln(x) space
mliu49 Mar 6, 2019
d2c213d
Bugfix in rmgpy.tools.canteraModel
mliu49 Apr 11, 2019
ebdcabb
Add reactionSystemIndex argument to Uncertainty.localAnalysis
mliu49 Apr 11, 2019
11c7cdb
Add function to process and prettify local uncertainty results
mliu49 Apr 11, 2019
3b39af9
Add input block for uncertainty settings
mliu49 Apr 11, 2019
6852f51
Enable automatic uncertainty analysis following RMG job
mliu49 Apr 11, 2019
467f436
Replace labels with species in sensConditions
mliu49 Apr 11, 2019
b7d5fe9
Bugfix: Cast parsed indices to int when loading sensitivity results
mliu49 Apr 11, 2019
cfd07f0
For local uncertainty, save RMG index for reactions
mliu49 Apr 11, 2019
879d7db
Refactor output printing for global uncertainty analysis
mliu49 Jan 29, 2019
da82473
Save extra species created for correlated uncertainty analysis
mliu49 Apr 16, 2019
e3cf9ba
Split run_uncertainty_analysis off of run_model_analysis
mliu49 Apr 19, 2019
1870e87
Move plot_sensitivity to rmgpy.tools.plot
mliu49 Apr 19, 2019
2e55f78
Enable simulate.py to run uncertainty analysis
mliu49 Apr 19, 2019
b8662d4
Make file naming in initializeLog more generic
mliu49 Apr 19, 2019
37dbd53
Fix formatting issue in writing Chemkin files
mliu49 Apr 19, 2019
035e481
Update findParameterSourcesAndAssignUncertainties.ipynb
mliu49 Apr 19, 2019
8494c36
Bugfix: Disable species constraints for training reactions
mliu49 Apr 22, 2019
92bd96b
Check if global analysis parameters are already in list
mliu49 Apr 29, 2019
9502978
Update localUncertainty.ipynb and example models
mliu49 Apr 30, 2019
75a9ff6
Make ln(mole fraction) an option for global uncertainty analysis
mliu49 May 1, 2019
65f3db0
Update globalUncertainty.ipynb
mliu49 May 1, 2019
3d36c5c
Update Cantera simulation IPython notebooks and data files
mliu49 May 1, 2019
3a06b6a
Inform user if local uncertainty is being enabled
mliu49 May 10, 2019
a3ab31d
Add documentation on using uncertainty options in input file
mliu49 May 23, 2019
7e15900
Add note about uncertainty to simulate.py documentation
mliu49 May 23, 2019
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
64 changes: 64 additions & 0 deletions documentation/source/users/rmg/input.rst
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -771,6 +771,70 @@ to turn off pressure dependence for all molecules larger than the given number
of atoms (16 in the above example).


.. _uncertaintyanalysis:

Uncertainty Analysis
====================

It is possible to request automatic uncertainty analysis following the convergence of an RMG simulation by including
an uncertainty options block in the input file::

uncertainty(
localAnalysis=False,
globalAnalysis=False,
uncorrelated=True,
correlated=True,
localNumber=10,
globalNumber=5,
terminationTime=None,
pceRunTime=1800,
logx=True
)

RMG can perform local uncertainty analysis using first-order sensitivity coefficients output by the native RMG solver.
This is enabled by setting ``localAnalysis=True``. Performing local uncertainty analysis requires suitable settings in
the reactor block (see :ref:`reactionsystem`). At minimum, the output species to perform sensitivity analysis on must
be specified, via the ``sensitivity`` argument. RMG will then perform local uncertainty analysis on the same species.
Species and reactions with the largest sensitivity indices will be reported in the log file and output figures.
The number of parameters reported can be adjusted using ``localNumber``.

RMG can also perform global uncertainty analysis, implemented using Cantera [Cantera]_ and the MIT Uncertainty
Quantification (MUQ) [MUQ]_ library. This is enabled by setting ``globalAnalysis=True``. Note that local analysis
is a required prerequisite of running the global analysis (at least for this semi-automatic approach), so
``localAnalysis`` will be enabled regardless of the input file setting. The analysis is performed by allowing the
input parameters with the largest sensitivity indices (as determined from the local uncertainty analysis) to vary
while performing reactor simulations using Cantera. MUQ is used to fit a Polynomial Chaos Expansion (PCE) to the
resulting output surface. The number of input parameters chosen can be adjusted using ``globalNumber``.
Note that this number applies independently to thermo and rate parameters and output species.
For example ``globalNumber=5`` for analysis on a single output species will result in 10 parameters being varied, while
having two output species could result in up to 20 parameters being varied, assuming no overlap in the sensitive input
parameters for each output.

The ``uncorrelated`` and ``correlated`` options refer to two approaches for uncertainty analysis. Uncorrelated means
that all input parameters are considered to be independent, each with their own uncertainty bounds. Thus, the output
uncertainty distribution is determined on the basis that every input parameter could vary within the full range of
its uncertainty bounds. Correlated means that inherent relationships between parameters (such as rate rules for kinetics
or group additivity values for thermochemistry) are accounted for, which reduces the uncertainty space of the input
parameters.

Finally, there are a few miscellaneous options for global uncertainty analysis. The ``terminationTime`` applies for the
reactor simulation. It is only necessary if termination time is not specified in the reactor settings (i.e. only other
termination criteria are used). The ``pceRunTime`` sets the time limit for fitting the PCE to the output surface.
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Define PCE?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I added the abbreviation to the first mention of Polynomial Chaos Expansion earlier in the section (line 806).

Longer run times allow more simulations to be performed, leading to more accurate results. The ``logx`` option toggles
the output parameter space between mole fractions and log mole fractions. Results in mole fraction space are more
physically meaningful, while results in log mole fraction space can be directly compared against local uncertainty
results.

**Important Note:** The current implementation of uncertainty analysis assigns values for input parameter
uncertainties based on the estimation method used by RMG. Actual uncertainties associated with the original data sources
are not used. Thus, the output uncertainties reported by these analyses should be viewed with this in mind.

.. [Cantera] Goodwin, D.G.; Moffat, H.K.; Speth, R.L. Cantera: An object-oriented software toolkit for
chemical kinetics, thermodynamics, and transport processes; http://www.cantera.org
.. [MUQ] Conrad, P.R.; Parno, M.D.; Davis, A.D.; Marzouk, Y.M. MIT Uncertainty Quantification Library (MUQ); http://muq.mit.edu/
.. [Gao2016] Gao, C. W.; Ph.D. Thesis. 2016.


.. _miscellaneousoptions:

Miscellaneous Options
Expand Down
9 changes: 6 additions & 3 deletions documentation/source/users/rmg/modules/simulate.rst
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,8 +1,8 @@
.. _simulate:

***********************************
Simulation and Sensitivity Analysis
***********************************
***********************************************
Simulation and Sensitivity/Uncertainty Analysis
***********************************************

For sensitivity analysis, RMG-Py must be compiled with the DASPK solver, which is done by default but has
some dependency restrictions. (See :ref:`License Restrictions on Dependencies <dependenciesRestrictions>` for more details.)
Expand Down Expand Up @@ -35,3 +35,6 @@ with the file name ``sensitivity_1_SPC_1.csv`` with the first index value indica
the sensitivity analysis is conducted for. Sensitivities to thermo of individual species is also saved as semi normalized sensitivities
dln(C_i)/d(G_j) where the units are given in 1/(kcal mol-1). The sensitivityThreshold is set to some value so that only
sensitivities for dln(C_i)/dln(k_j) > sensitivityThreshold or dlnC_i/d(G_j) > sensitivityThreshold are saved to this file.

Uncertainty analysis can also be requested via input file options. For more details, see :ref:`uncertaintyanalysis`.
The results of the analysis will be printed in the ``simulate.log`` file which is generated.
170 changes: 52 additions & 118 deletions ipython/canteraSensitivityComparison.ipynb
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -11,38 +11,31 @@
"Cantera version >= 2.3.0 is required"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Let's check what version of Cantera you have installed"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"import cantera\n",
"print cantera.__version__"
"print cantera.__version__ # Check Cantera version"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"from rmgpy.chemkin import *\n",
"from rmgpy.tools.canteraModel import *\n",
"from rmgpy.tools.plot import parseCSVData\n",
"import shutil\n",
"\n",
"from IPython.display import display, Image\n",
"\n",
"from rmgpy.chemkin import loadChemkinFile\n",
"from rmgpy.species import Species\n",
"from IPython.display import display, Image"
"from rmgpy.tools.canteraModel import Cantera, getRMGSpeciesFromUserSpecies\n",
"from rmgpy.tools.plot import SimulationPlot, ReactionSensitivityPlot, parseCSVData\n",
"from rmgpy.tools.simulate import run_simulation"
]
},
{
Expand All @@ -55,14 +48,12 @@
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": true
},
"metadata": {},
"outputs": [],
"source": [
"speciesList, reactionList = loadChemkinFile('data/minimal_model/chem_annotated.inp',\n",
" 'data/minimal_model/species_dictionary.txt',\n",
" 'data/minimal_model/tran.dat')"
"speciesList, reactionList = loadChemkinFile('./data/ethane_model/chem_annotated.inp',\n",
" './data/ethane_model/species_dictionary.txt',\n",
" './data/ethane_model/tran.dat')"
]
},
{
Expand All @@ -75,9 +66,7 @@
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"# Find the species: ethane and methane\n",
Expand All @@ -91,17 +80,15 @@
"#reactorTypeList = ['IdealGasReactor']\n",
"reactorTypeList = ['IdealGasConstPressureTemperatureReactor']\n",
"molFracList=[{ethane: 1}]\n",
"Tlist = ([1300],'K')#,1500,2000],'K')\n",
"Plist = ([1],'atm')\n",
"Tlist = ([1300], 'K')\n",
"Plist = ([1], 'bar')\n",
"reactionTimeList = ([0.5], 'ms')"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"# Create cantera object, loading in the species and reactions\n",
Expand All @@ -115,23 +102,22 @@
"#job.loadChemkinModel('data/minimal_model/chem_annotated.inp',transportFile='data/minimal_model/tran.dat')\n",
"\n",
"# Generate the conditions based on the settings we declared earlier\n",
"job.generateConditions(reactorTypeList, reactionTimeList, molFracList, Tlist, Plist)\n",
"# Simulate and plot\n",
"alldata = job.simulate()\n",
"job.plot(alldata)"
"job.generateConditions(reactorTypeList, reactionTimeList, molFracList, Tlist, Plist)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"# Simulate and plot\n",
"alldata = job.simulate()\n",
"job.plot(alldata)\n",
"\n",
"# Show the plots in the ipython notebook\n",
"for i, condition in enumerate(job.conditions):\n",
" print 'Cantera Simulation: Condition {0} Mole Fractions'.format(i+1)\n",
" print 'Cantera Simulation: Condition {0} Species Mole Fractions'.format(i+1)\n",
" display(Image(filename=\"temp/{0}_mole_fractions.png\".format(i+1)))\n",
" \n",
" print 'Cantera Simulation: Condition {0} Ethane Reaction Sensitivity'.format(i+1)\n",
Expand All @@ -142,111 +128,59 @@
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": true
"scrolled": true
},
"outputs": [],
"source": [
"# Let's compare against the same simulation in RMG\n",
"# Create an input file\n",
"\n",
"input = '''\n",
"database(\n",
" thermoLibraries = ['primaryThermoLibrary'],\n",
" reactionLibraries = [],\n",
" seedMechanisms = [],\n",
" kineticsDepositories = 'default',\n",
" kineticsFamilies = 'default',\n",
" kineticsEstimator = 'rate rules',\n",
")\n",
"\n",
"species(\n",
" label = \"ethane\",\n",
" structure = SMILES(\"CC\"))\n",
"\n",
"species(\n",
" label = \"methane\",\n",
" structure = SMILES(\"C\"))\n",
"# Copy example input file to temp folder\n",
"shutil.copy('./data/ethane_model/input.py', './temp')\n",
"\n",
"simpleReactor(\n",
" temperature = (1300,\"K\"), \n",
" pressure = (1,\"atm\"),\n",
" initialMoleFractions={\n",
" \"ethane\": 1\n",
" },\n",
" terminationTime = (0.5,\"ms\"),\n",
" sensitivity=['ethane','methane']\n",
")\n",
"\n",
"model(\n",
" toleranceMoveToCore = 0.04,\n",
")\n",
"\n",
"options(\n",
" saveSimulationProfiles = True,\n",
")\n",
"\n",
"'''\n",
"f = open('temp/temp_input.py', 'wb')\n",
"f.write(input)\n",
"f.close()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"from rmgpy.tools.sensitivity import runSensitivity\n",
"runSensitivity('temp/temp_input.py', 'data/minimal_model/chem_annotated.inp', 'data/minimal_model/species_dictionary.txt')"
"# We can run the same simulation using RMG's native solver\n",
"run_simulation(\n",
" './temp/input.py',\n",
" './data/ethane_model/chem_annotated.inp',\n",
" './data/ethane_model/species_dictionary.txt',\n",
")"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"print 'RMG Native Simulation: Species Mole Fractions'\n",
"display(Image(filename=\"temp/solver/simulation_1_19.png\"))\n",
"display(Image(filename=\"./temp/solver/simulation_1_27.png\"))\n",
"\n",
"print 'RMG Native Simulation: Ethane Reaction Sensitivity'\n",
"display(Image(filename=\"temp/solver/sensitivity_1_SPC_1_reactions.png\"))\n"
"display(Image(filename=\"./temp/solver/sensitivity_1_SPC_1_reactions.png\"))"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"# Let's also compare against the same simulation and sensitivity analysis that was conducted in CHEMKIN\n",
"# and saved as a .csv file\n",
"time, dataList = parseCSVData('data/minimal_model/chemkin_mole_fractions.csv')\n",
"SimulationPlot(xVar=time,yVar=dataList).plot('temp/chemkin_mole_fractions.png')\n",
"print 'CHEMKIN Simulation: Mole Fractions'\n",
"display(Image(filename=\"temp/chemkin_mole_fractions.png\"))"
"time, dataList = parseCSVData('./data/ethane_model/chemkin_mole_fractions.csv')\n",
"SimulationPlot(xVar=time, yVar=dataList, numSpecies=10).plot('./temp/chemkin_mole_fractions.png')\n",
"print 'CHEMKIN Simulation: Species Mole Fractions'\n",
"display(Image(filename=\"./temp/chemkin_mole_fractions.png\"))\n",
"\n",
"time, dataList = parseCSVData('./data/ethane_model/chemkin_sensitivity_ethane.csv')\n",
"ReactionSensitivityPlot(xVar=time, yVar=dataList, numReactions=10).barplot('./temp/chemkin_sensitivity_ethane.png')\n",
"print 'CHEMKIN Simulation: Ethane Reaction Sensitivity'\n",
"display(Image(filename=\"./temp/chemkin_sensitivity_ethane.png\"))"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"time, dataList = parseCSVData('data/minimal_model/chemkin_sensitivity_ethane.csv')\n",
"ReactionSensitivityPlot(xVar=time,yVar=dataList).barplot('temp/chemkin_ethane_sensitivity.png')\n",
"print 'CHEMKIN Simulation: Ethane Reaction Sensitivity'\n",
"display(Image(filename=\"temp/chemkin_ethane_sensitivity.png\"))"
]
"source": []
}
],
"metadata": {
Expand All @@ -265,9 +199,9 @@
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython2",
"version": "2.7.11"
"version": "2.7.15"
}
},
"nbformat": 4,
"nbformat_minor": 0
"nbformat_minor": 1
}
Loading