.. _tutorials_uq: .. _guided_pytuq_integration: .. _pytuq_quickstart: UQ with PyTUQ =========== .. admonition:: **Time to Complete**: 30-60 minutes :class: note **What you will learn**: - Install PyTUQ - Run an AMReX+PyTUQ to perform sensitivity analysis Overview -------- AMReX simulations deliver high-fidelity, accurate results for complex physics problems, but the effect on simulation results due to uncertain inputs can require hundreds or thousands of simulations, making comprehensive analysis computationally prohibitive. This tutorial demonstrates how to improve efficiency without sacrificing accuracy by using polynomial chaos expansions to build fast surrogate models that identify which input parameters truly matter. PyTUQ (Python interface to the UQ Toolkit) provides specialized tools for surrogate construction and global sensitivity analysis, enabling rapid parameter space exploration and dimensionality reduction for scientific applications. We demonstrate how to integrate PyTUQ with your AMReX application through a practical workflow; the AMReX-based heat equation tutorial is equipped to perform sensitivity analysis. In these examples we model the heat equation .. math:: \frac{\partial\phi}{\partial t} = D\nabla^2 \phi, with initial condition .. math:: \phi_0 = 1 + A e^{-r^2 / (2V)}, where ``r`` is the distance from the center of the domain, and with uncertain parameters ``diffusion_coeff`` (:math:`D`), ``init_amplitude`` (:math:`A`), and ``init_variance`` (:math:`V`). The outputs of interest are the maximum temperature, mean temperature, standard deviation of temperature, and the temperature at a specified point. .. toctree:: :maxdepth: 1 HeatEquation_UQ_MathematicalDetails Located in ``amrex-tutorials/ExampleCodes/UQ/HeatEquation``, this example illustrates a complete forward UQ workflow from parameter sampling randomized input parameters to perform sensitivity analysis. By understanding this example, you will have a basis for understanding how to adapt this workflow to your own AMReX application. More specifically, you can directly compare/diff ``amrex-tutorials/ExampleCodes/UQ/HeatEquation/main.cpp`` against the original heat equation tutorial ``amrex-tutorials/GuidedTutorials/HeatEquation_Simple/main.cpp`` to see exactly what source code changes are made to the AMReX application in this example. Installation ------------ We now describe the installation and workflow process on a local workstation. First, install pytuq using this script (based on information provided in `pytuq/README.md `_): .. code-block:: bash :caption: Pytuq installation script #!/bin/bash # 1. Clone repositories git clone https://github.com/sandialabs/pytuq.git # 2. Create a conda environment with python (optional, you can add to an existing env) conda create --name pytuq conda activate pytuq conda install python=3.11 # 3. Install PyTUQ and requirements cd pytuq pip install -r requirements.txt pip install . conda install dill # 4. Verify installation conda list | grep pytuq # Should show pytuq 1.0.0 Examples -------- C++ AMReX + PyTUQ (BASH driven) ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ **Prerequisites**: - AMReX and pytuq cloned at the same directory level as amrex-tutorials - pytuq installation described above - GNU Parallel for task management: ``sudo apt-get install parallel`` (mpiexec option also exists in script) .. code-block:: bash :caption: Build AMReX appplication cd /path/to/amrex-tutorials/ExampleCodes/UQ/HeatEquation/ make -j4 .. code-block:: bash :caption: Run sensitivity analysis with bash script ./workflow_uqpc.x Understanding GNU Parallel Workflow Pattern ------------------------------------------- The ``workflow_uqpc.x`` bash script relies on the user augmenting their codes to write outputs of interest to ASCII text files. In this case, the ``main.cpp`` was modified from the ``amrex-tutorials/GuidedTutorials/HeatEquation_Simple/main.cpp`` in the following ways: First, support for parsing ``diffusion_coeff``, ``init_amplitude``, and ``init_variance`` from the input file and command line were added. Second, support for writing outputs of interest to ASCII text files is added. The ``datalog_filename`` input parameter is generated by the bash script to give each simulation output a unique identifier. The ``datalog_int`` parameter gives the user the option to write the outputs of interest at a given time step interval, but in this example, the outputs of interest after the final step are those that matter, and are extracted by the bash script to create a master output file containing a separate set of simulation outputs of interest in each row. The bash script calls PyTUQ scripts that generate an input parameter files for training and testing (``ptrain.txt`` and ``ptest.txt``) based on polynomial chaos settings, and then uses GNU Parallel to run multiple simulations efficiently and collect outputs into results files (``ytrain.txt`` and ``ytest.txt``) that PyTUQ can use for surrogate model fitting. Understanding the Output ------------------------ | `pnames.txt`, `outnames.txt` | -names of input and output parameters | `param_margpc.txt` | -Each row contains the mean and standard deviation for an uncertain input parameter | `qtrain.txt`, `qtest.txt` | -each row is a separate set of normal random variables used to generate uncertain inputs | `ptrain.txt`, `ptest.txt`, `pall.txt` | -each row is a separate set of input parameters for each simulation | `stdoutlog_train##.txt`, `stdoutlog_test##.txt` | -All standard output from AMReX simulations for training and testing data. The numbers refer to separate simulations. | `datalog_train##.txt`, `datalog_test##.txt` | -Specifically-chosen output from AMReX simulations for training and testing data. In this example it reports the outputs (max, mean, standard deviation, and specific cell temperature) at user-specified intervals. | `ytrain.txt`, `ytest.txt`, `yall.txt` | -agglomeration of outputs of interest from all simulations | `results.pk` | -Python pickle file encapsulating the results | `labels.txt` | -list of labels of simulation types (training or testing) for diagnostic/plot generation purposes | `xx__.png` | -scatter plot of 2 inputs for training and testing | `pcoord_1.png`, `pcoord_1_lab_Testing.png`, `pcoord_1_lab_Training.png` | -graphical representation of how inputs correlate to outputs for each individual simulation | `yx_.png`, `yx__log.png` | -scatter plots of output as a function of each input | `yxx_#.png` | -scatter plots of output a function of multiple inputs | `pdf_tri_inputs.png`, `pdf_tri_output.png` | -PDFs of inputs and outputs | `ensemble.png` | -graphical display of all output values | `idm_#_training.png`, `idm_#_testing.png` | -graphical display of output values | `dm_#.png` | -”data vs model” parity plots for each output; compares the predicted values from a surrogate model or approximation against the true or actual values from the full computational model | `fit_s#_training.png`, `fit_s#_testing.png` | -Shows model vs PC approximation for a single simulation. | `pdf_output_#.png`, `pdf_joyplot.png` | -PDF of output variables | `allsens_main.txt`, `sens_main.png` | -raw data and plot for parameter sensitivities | `jsens_#.png`, `Jsens_ave.png` | -joint sensitivities of output with respect to inputs | `sensmat_main.png` | -sensitivity matrix | `pcslices_o#.png` | -polynomial chaos fits | `pccont_o#_d#_d#.png` | -polynomial chaos fits of output variables with respect to two input variables Additional Resources -------------------- **PyTUQ Resources:** - `PyTUQ Documentation `_ - `PyTUQ Examples directory `_ - eebaill, ksargsyan, & Bert Debusschere. (2025). sandialabs/pytuq: v1.0.0z (v1.0.0z). Zenodo. https://doi.org/10.5281/zenodo.17110054 **AMReX Resources:** - `AMReX Documentation `_ **Uncertainty Quantification Theory:** - Ghanem, Roger, David Higdon, and Houman Owhadi, eds. *Handbook of Uncertainty Quantification*. Vol. 6. New York: Springer, 2017. (For workflow, plotting, and analysis specifics)