.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples/plot_ishigami.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code .. rst-class:: sphx-glr-example-title .. _sphx_glr_auto_examples_plot_ishigami.py: Example: Ishigami Testfunction ============================== About the model ^^^^^^^^^^^^^^^ This easy tutorial shows the application of pygpc to the Ishigami function, which can be found in the testfunctions section. The model consists of three random variables that are considered as input parameters (*x1*, *x2*, *x3*). The shape parameters of the function are chosen to be :code:`a=7` and :code:`b=0.1`. The model returns an output array with a value *y* for every sampling point. .. GENERATED FROM PYTHON SOURCE LINES 14-82 .. code-block:: default # Windows users have to encapsulate the code into a main function to avoid multiprocessing errors. # def main(): import os import pygpc import numpy as np from collections import OrderedDict import matplotlib # matplotlib.use("Qt5Agg") fn_results = "tmp/example_ishigami" if os.path.exists(fn_results + ".hdf5"): os.remove(fn_results + ".hdf5") if os.path.exists(fn_results + "_val.hdf5"): os.remove(fn_results + "_val.hdf5") if os.path.exists(fn_results + "_mc.hdf5"): os.remove(fn_results + "_mc.hdf5") # define model model = pygpc.testfunctions.Ishigami() # define problem parameters = OrderedDict() parameters["x1"] = pygpc.Beta(pdf_shape=[1, 1], pdf_limits=[-np.pi, np.pi]) parameters["x2"] = pygpc.Beta(pdf_shape=[1, 1], pdf_limits=[-np.pi, np.pi]) parameters["x3"] = pygpc.Beta(pdf_shape=[1, 1], pdf_limits=[-np.pi, np.pi]) parameters["a"] = 7. parameters["b"] = 0.1 parameters_random = OrderedDict() parameters_random["x1"] = pygpc.Beta(pdf_shape=[1, 1], pdf_limits=[-np.pi, np.pi]) parameters_random["x2"] = pygpc.Beta(pdf_shape=[1, 1], pdf_limits=[-np.pi, np.pi]) parameters_random["x3"] = pygpc.Beta(pdf_shape=[1, 1], pdf_limits=[-np.pi, np.pi]) problem = pygpc.Problem(model, parameters) # gPC options options = dict() options["order"] = [15] * problem.dim options["order_max"] = 15 options["order_start"] = 15 options["method"] = 'reg' options["solver"] = "Moore-Penrose" options["interaction_order"] = 2 options["order_max_norm"] = 1.0 options["n_cpu"] = 0 options["eps"] = 0.01 options["fn_results"] = fn_results options["basis_increment_strategy"] = None options["plot_basis"] = False options["n_grid"] = 1300 options["save_session_format"] = ".pkl" options["matrix_ratio"] = 2 options["grid"] = pygpc.Random options["grid_options"] = {"seed": 1} # define algorithm algorithm = pygpc.Static(problem=problem, options=options, grid=None) # Initialize gPC Session session = pygpc.Session(algorithm=algorithm) # run gPC session session, coeffs, results = session.run() .. rst-class:: sphx-glr-script-out .. code-block:: none Performing 1300 simulations! It/Sub-it: 15/2 Performing simulation 0001 from 1300 [ ] 0.1% Total parallel function evaluation: 0.0005440711975097656 sec Determine gPC coefficients using 'Moore-Penrose' solver ... LOOCV 01 from 25 [= ] 4.0% LOOCV 02 from 25 [=== ] 8.0% LOOCV 03 from 25 [==== ] 12.0% LOOCV 04 from 25 [====== ] 16.0% LOOCV 05 from 25 [======== ] 20.0% LOOCV 06 from 25 [========= ] 24.0% LOOCV 07 from 25 [=========== ] 28.0% LOOCV 08 from 25 [============ ] 32.0% LOOCV 09 from 25 [============== ] 36.0% LOOCV 10 from 25 [================ ] 40.0% LOOCV 11 from 25 [================= ] 44.0% LOOCV 12 from 25 [=================== ] 48.0% LOOCV 13 from 25 [==================== ] 52.0% LOOCV 14 from 25 [====================== ] 56.0% LOOCV 15 from 25 [======================== ] 60.0% LOOCV 16 from 25 [========================= ] 64.0% LOOCV 17 from 25 [=========================== ] 68.0% LOOCV 18 from 25 [============================ ] 72.0% LOOCV 19 from 25 [============================== ] 76.0% LOOCV 20 from 25 [================================ ] 80.0% LOOCV 21 from 25 [================================= ] 84.0% LOOCV 22 from 25 [=================================== ] 88.0% LOOCV 23 from 25 [==================================== ] 92.0% LOOCV 24 from 25 [====================================== ] 96.0% LOOCV 25 from 25 [========================================] 100.0% LOOCV computation time: 2.577364444732666 sec -> relative loocv error = 5.062046906871139e-06 .. GENERATED FROM PYTHON SOURCE LINES 83-86 Postprocessing ^^^^^^^^^^^^^^ Postprocess gPC and add results to .hdf5 file .. GENERATED FROM PYTHON SOURCE LINES 86-93 .. code-block:: default pygpc.get_sensitivities_hdf5(fn_gpc=session.fn_results, output_idx=None, calc_sobol=True, calc_global_sens=True, calc_pdf=True, n_samples=int(1e4)) .. rst-class:: sphx-glr-script-out .. code-block:: none > Loading gpc session object: tmp/example_ishigami.pkl > Loading gpc coeffs: tmp/example_ishigami.hdf5 > Adding results to: tmp/example_ishigami.hdf5 .. GENERATED FROM PYTHON SOURCE LINES 94-97 Validation ^^^^^^^^^^ Validate gPC vs original model function .. GENERATED FROM PYTHON SOURCE LINES 97-105 .. code-block:: default pygpc.validate_gpc_plot(session=session, coeffs=coeffs, random_vars=["x1", "x2"], n_grid=[51, 51], output_idx=0, fn_out=session.fn_results + '_val', n_cpu=session.n_cpu) .. image-sg:: /auto_examples/images/sphx_glr_plot_ishigami_001.png :alt: Original model, gPC approximation, Difference (Original vs gPC) :srcset: /auto_examples/images/sphx_glr_plot_ishigami_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 106-107 Validate gPC vs original model function (Monte Carlo) .. GENERATED FROM PYTHON SOURCE LINES 107-114 .. code-block:: default nrmsd = pygpc.validate_gpc_mc(session=session, coeffs=coeffs, n_samples=int(1e4), output_idx=0, n_cpu=session.n_cpu, fn_out=session.fn_results + '_mc') .. image-sg:: /auto_examples/images/sphx_glr_plot_ishigami_002.png :alt: plot ishigami :srcset: /auto_examples/images/sphx_glr_plot_ishigami_002.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 115-117 Sensitivity analysis ^^^^^^^^^^^^^^^^^^^^ .. GENERATED FROM PYTHON SOURCE LINES 117-126 .. code-block:: default sobol, gsens = pygpc.get_sens_summary(fn_results, parameters_random, fn_results + "_sens_summary.txt") pygpc.plot_sens_summary(sobol=sobol, gsens=gsens) # # On Windows subprocesses will import (i.e. execute) the main module at start. # You need to insert an if __name__ == '__main__': guard in the main module to avoid # creating subprocesses recursively. # # if __name__ == '__main__': # main() .. image-sg:: /auto_examples/images/sphx_glr_plot_ishigami_003.png :alt: Normalized Sobol indices :srcset: /auto_examples/images/sphx_glr_plot_ishigami_003.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 19.180 seconds) .. _sphx_glr_download_auto_examples_plot_ishigami.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_ishigami.py ` .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_ishigami.ipynb ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_