Merge pull request #207 from SURGroup/Development

Development

docs/code/stochastic_processes/karhunen_loeve_2d/README.rst

@@ -0,0 +1,4 @@
+Karhunen Loeve Expansion 2D Examples
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+

docs/code/stochastic_processes/karhunen_loeve_2d/plot_karhunen_loeve_2d.py

@@ -0,0 +1,72 @@
+"""
+
+Karhunen Loeve Expansion 2 Dimesional
+=================================================================
+
+In this example, the KL Expansion is used to generate 2 dimensional stochastic fields from a prescribed Autocorrelation
+Function. This example illustrates how to use the :class:`.KarhunenLoeveExpansion2D` class for a 2 dimensional
+random field and compare the statistics of the generated random field with the expected values.
+
+"""
+
+# %% md
+#
+# Import the necessary libraries. Here we import standard libraries such as numpy and matplotlib, but also need to
+# import the :class:`.KarhunenLoeveExpansionTwoDimension` class from the :py:mod:`stochastic_processes` module of UQpy.
+
+# %%
+from matplotlib import pyplot as plt
+from UQpy.stochastic_process import KarhunenLoeveExpansion2D
+import numpy as np
+
+plt.style.use('seaborn')
+
+# %% md
+#
+# The input parameters necessary for the generation of the stochastic processes are given below:
+
+# %%
+
+n_samples = 100  # Num of samples
+nx, nt = 20, 10
+dx, dt = 0.05, 0.1
+
+x = np.linspace(0, (nx - 1) * dx, nx)
+t = np.linspace(0, (nt - 1) * dt, nt)
+xt_list = np.meshgrid(x, x, t, t, indexing='ij')  # R(t_1, t_2, x_1, x_2)
+
+# %% md
+
+# Defining the Autocorrelation Function.
+
+# %%
+
+R = np.exp(-(xt_list[0] - xt_list[1]) ** 2 - (xt_list[2] - xt_list[3]) ** 2)
+# R(x_1, x_2, t_1, t_2) = exp(-(x_1 - x_2) ** 2 -(t_1 - t_2) ** 2)
+
+KLE_Object = KarhunenLoeveExpansion2D(n_samples=n_samples, correlation_function=R, time_intervals=np.array([dt, dx]),
+                                      thresholds=[4, 5], random_state=128)
+
+# %% md
+
+# Simulating the samples.
+
+# %%
+
+samples = KLE_Object.samples
+
+# %% md
+
+# Plotting a sample of the stochastic field.
+
+# %%
+
+fig = plt.figure()
+plt.title('Realisation of the Karhunen Loeve Expansion for a 2D stochastic field')
+plt.imshow(samples[0, 0])
+plt.ylabel('t (Time)')
+plt.xlabel('x (Space)')
+plt.show()
+
+print('The mean of the samples is ', np.mean(samples), 'whereas the expected mean is 0.000')
+print('The variance of the samples is ', np.var(samples), 'whereas the expected variance is 1.000')

docs/source/architecture.rst

@@ -0,0 +1,169 @@
+UQpy architecture
+==================
+
+Distributions
+-----------------------
+
+The UML diagram below explains the class hierarchy of the distributions module. In case of distributions a series of
+abstract base classes are defined. :class:`.Distribution` is the basis, with classes such as :class:`.Distribution1D` \
+and :class:`.DistributionND` abstract base classes refining this hierarchy according to the distribution dimensionality.
+The 1D case if further refined with :class:`.DistributionContinuous1D` and :class:`.DistributionDiscrete1D` to take into
+account different types of random variables. Finally, the :class:`.Copula` abstract base class serves as the common
+interface for all copula implementations. Note that all the above are abstract classes and cannot be directly
+instantiated.
+
+.. image:: _static/architecture/distributions.png
+   :scale: 30 %
+   :alt:  UML diagram of distributions module.
+   :align: center
+
+Sampling
+-----------------------
+
+The sampling module contains multiple methods and the UML diagrams below are focused to one method at a time to reduce
+the complexity of the graphs. Before starting with hierarchy of each different method, we should mention that the
+classes :class:`.MonteCarloSampling`, :class:`.SimplexSampling` and :class:`.ImportanceSampling` are independent object
+and do not extend any baseclass.
+
+Starting with the Markov Chain Monte Carlo algorithms, the diagram below makes obvious the dependency between the
+different implementations. the :class:`.MCMC` abstract base class includes the common functionality between the methods,
+while the specific algorithms only override the abstract methods of that base class.
+
+
+.. image:: _static/architecture/mcmc.png
+   :scale: 30 %
+   :alt:  UML diagram of MCMC hierarchy.
+   :align: center
+
+The case of stratified sampling is more elaborate, as it includes a family of algorithms. The abstract base class
+:class:`.StratifiedSampling` serves as the interface between the concrete implementations :class:`.LatinHypercubeSampling`,
+:class:`.TrueStratifiedSampling` and :class:`.RefinedStratifiedSampling`. In the case of :class:`.LatinHypercubeSampling`
+the Strategy Design Pattern was used to extract the Latin Hypercube criteria from the sampling. Here the
+:class:`.Criterion` base class provides the interface that the specific criteria need to override. In a similar manner,
+the geometric stratification in the case of :class:`.TrueStratifiedSampling` is extracted under the
+:class:`.Strata` abstract base class. Last but not least, in the case of :class:`.RefinedStratifiedSampling`
+the different strata refinement strategies are extracted using the :class:`.Refinement` baseclass as their common interface.
+
+
+.. image:: _static/architecture/stratified_sampling.png
+   :scale: 30 %
+   :alt:  UML diagram of Stratified sampling class hierarchy.
+   :align: center
+
+In the case of :class:`.AdaptiveKriging` sampling methods, again the different learning functions are extracted into
+separate classes under the common :class:`.LearningFunction` class.
+
+.. image:: _static/architecture/adaptive_kriging_functions.png
+   :scale: 30 %
+   :alt:  UML diagram of Adaptive Kriging Hierarchy.
+   :align: center
+
+
+Transformations
+-----------------------
+
+The transformations module is one of the most simple in :py:mod:`UQpy` with three independent classes available, namely
+:class:`.Nataf`, :class:`.Correlate` and :class:`.Decorrelate`.
+
+.. image:: _static/architecture/transformations.png
+   :scale: 30 %
+   :alt:  UML diagram of Transformations module.
+   :align: center
+
+
+Stochastic Processes
+-----------------------
+
+The stochastic process module is has again simple structure with five independent classes available.
+
+.. image:: _static/architecture/stochastic_process.png
+   :scale: 30 %
+   :alt:  UML diagram of Stochastic Process module.
+   :align: center
+
+Run Model
+-----------------------
+
+In case of the RunModel module, the final algorithm to run is constructed by object composition of two different inputs.
+Initially, the type of the model to run, with :class:`.PythonModel` and :class:`.ThirdPartyModel` being the two
+available options, while the execution part is delegated to either the :class:`.SerialExecution` or :class:`.ParallelExecution`
+alternatives.
+
+.. image:: _static/architecture/run_model.png
+   :scale: 30 %
+   :alt:  UML diagram of Run Model module.
+   :align: center
+
+Inference
+-----------------------
+
+Compared to v3, the inference module has undergone a major refactoring towards v4. The initial :class:`.InferenceModel`
+class that contained all cases of computing the posterior log-likelihood is now split into three independent cases. Given
+the inference models, backward uncertainty propagation can be performed be choosing between :class:`.MLE`,
+:class:`.BayesParameterEstimation` to infer the parameter distributions of a model, or :class:`.InformationModelSelection`
+and :class:`.BayesModelSelection` to select the model that best describes the available data. In the case of
+:class:`.InformationModelSelection` the selection criteria have been extracted into separate classes under the
+:class:`.InformationCriterion` baseclass. Similarly, the evidence methods of :class:`.BayesModelSelection` are also
+parameters  that implement the abstract base class :class:`.EvidenceMethod`.
+
+.. image:: _static/architecture/inference.png
+   :scale: 30 %
+   :alt:  UML diagram of Inference module.
+   :align: center
+
+Reliability
+-----------------------
+
+The reliability module maintained the same class hierarachy as in v3, with :class:`.SubsetSimulation` being an
+independent class and :class:`.FORM` and :class:`.SORM` methods providing concrete implementations to the
+:class:`.TaylorSeries` abstract base class.
+
+.. image:: _static/architecture/reliability.png
+   :scale: 30 %
+   :alt:  UML diagram of Reliability module.
+   :align: center
+
+Surrogates
+-----------------------
+
+Another module that has extensively restructured in v4 is the surrogates. Apart from the :class:`.SROM` method which
+was retained as an independent algorithm, the previous Kriging functionality was removed. It is now replaced with
+:class:`.GaussianProcessRegression`. The functionality of the Gaussian is constructed using object composition,
+and the specific implementation of :class:`.Regression` and :class:`.Kernel` abstract base classes. An additional
+functionality of constrained surrogates is added by implementing the :class:`.ConstraintsGPR` abstract class. The
+functionality of :class:`.PolynomialChaosExpansion` was rewritten from scratch to address some performance issues of v3.
+The Strategy Design pattern was used here as well, with three abstract base classes :class:`.Polynomials`,
+:class:`.PolynomialBasis` ans :class:`.Regression` serving as the interface for the concrete classes.
+
+.. image:: _static/architecture/surrogates.png
+   :scale: 30 %
+   :alt:  UML diagram of Surrogates module.
+   :align: center
+
+Sensitivity
+-----------------------
+
+The sensitivity module has significantly benefited from the enhanced of modularity of the code introduced in v4.
+Apart from the existing independent :class:`.MorrisSensitivity` method, the :class:`.PceSensitivity` was added as an
+independent class. Finally, based on the common :class:`.Sensitivity` abstract base class, a series of new algorithms
+were introduced such as :class:`.SobolSensitivity`, :class:`.GeneralizedSobolSensitivity`, :class:`.ChatterjeeSensitivity`
+and :class:`.CramerVonMisesSensitivity`.
+
+.. image:: _static/architecture/sensitivity.png
+   :scale: 30 %
+   :alt:  UML diagram of Sensitivity module.
+   :align: center
+
+Dimension Reduction
+-----------------------
+
+The final but one of the most import modules in :py:mod:`UQpy` is dimension reduction. The :class:`.SnapshotPOD` and
+:class:`.DirectPOD` methods were retained under the :class:`.POD` abstract base class. :class:`.HigherOrderSVD` method
+was introduced as independent class, while special attention was given to Grassmann Manifolds.
+The abstract base class :class:`.GrassmannProjection` serves as an interface for different methods to project data on
+the Grassmann Manifold, with :class:`.GrassmannOperations` and :class:`.GrassmannInterpolation` support all related operations.
+
+.. image:: _static/architecture/dimension_reduction.png
+   :scale: 30 %
+   :alt:  UML diagram of Dimension Reduction module.
+   :align: center

docs/source/bibliography.bib

@@ -764,7 +764,7 @@ @article{Chatterjee
 }
 
 @misc{gamboa2020global,
-      title={Global Sensitivity Analysis: a new generation of mighty estimators based on rank statistics}, 
+      title={Global Sensitivity Analysis: a new generation of mighty estimators based on rank statistics},
       author={Fabrice Gamboa and Pierre Gremaud and Thierry Klein and Agnès Lagnoux},
       year={2020},
       eprint={2003.01772},
@@ -827,4 +827,16 @@ @article{saltelli_2002
 url = {https://www.sciencedirect.com/science/article/pii/S0010465502002801},
 author = {Andrea Saltelli},
 keywords = {Sensitivity analysis, Sensitivity measures, Sensitivity indices, Importance measures},
-}
+}
+@article{Kle2D,
+title = {Simulation of multi-dimensional random fields by Karhunen–Loève expansion},
+journal = {Computer Methods in Applied Mechanics and Engineering},
+volume = {324},
+pages = {221-247},
+year = {2017},
+issn = {0045-7825},
+doi = {https://doi.org/10.1016/j.cma.2017.05.022},
+url = {https://www.sciencedirect.com/science/article/pii/S0045782516318692},
+author = {Zhibao Zheng and Hongzhe Dai},
+keywords = {Multi-dimensional random field, Karhunen–Loève expansion, Random field simulation, Fredholm integral equation},
+}

docs/source/conf.py

@@ -97,7 +97,8 @@
         "../code/sensitivity/generalised_sobol",
         "../code/sensitivity/comparison",
         "../code/stochastic_processes/bispectral",
-        "../code/stochastic_processes/karhunen_loeve",
+        '../code/stochastic_processes/karhunen_loeve_1d',
+        '../code/stochastic_processes/karhunen_loeve_2d',
         "../code/stochastic_processes/spectral",
         "../code/stochastic_processes/translation",
         "../code/reliability/form",
@@ -136,7 +137,8 @@
         "auto_examples/sensitivity/generalised_sobol",
         "auto_examples/sensitivity/comparison",
         "auto_examples/stochastic_processes/bispectral",
-        "auto_examples/stochastic_processes/karhunen_loeve",
+        'auto_examples/stochastic_processes/karhunen_loeve_1d',
+        'auto_examples/stochastic_processes/karhunen_loeve_2d',
         "auto_examples/stochastic_processes/spectral",
         "auto_examples/stochastic_processes/translation",
         "auto_examples/reliability/form",

docs/source/index.rst

@@ -106,6 +106,7 @@ The default logging level is set to **ERROR**. The user can change the logging s
    /surrogates/index
    /transformations/index
    /utilities/index
+   architecture.rst
    paper.rst
    bibliography.rst
    news_doc

docs/source/stochastic_process/index.rst

@@ -47,5 +47,6 @@ the desired functionality. It does not require modification of any existing clas
 
     Spectral Representation Method <spectral_representation>
     Bispectral Representation Method <bispectral_representation>
-    Karhunen Loeve Expansion <karhunen_loeve>
+    Karhunen Loeve Expansion <karhunen_loeve_1d>
+    Karhunen Loeve Expansion 2D <karhunen_loeve_2d>
     Translation Processes <translation>

docs/source/stochastic_process/karhunen_loeve_2d.rst

@@ -0,0 +1,33 @@
+Karhunen Loève Expansion for Multi-Dimensional Fields
+----------------------------
+
+The Karhunen Loève Expansion expands the stochastic field as follows:
+
+.. math:: A(x, t) = \sum_{n=1}^{\infty} \sum_{k=1}^{\infty}\eta_{nk}(\theta) \sqrt{\lambda_n(x)} f_n(t, x) \sqrt{\mu_{nk}} g_{nk}(x)
+
+where :math:`\eta_{nk}(\theta)` are uncorrelated standardized normal random variables and :math:`\lambda_n(x)` and :math:`f_n(x, t)` are the eigenvalues and eigenvectors repsectively of the "quasi" one dimensional covariance function :math:`C(x, t_1, t_2)`. :math:`\mu_{nk}` and :math:`g_{nk}(x)` are the eigenvalues and eigenvectors of the derived "one" dimensional covariance function :math:`H(x_1, x_2)`. Additional details regarding the simulation formula can be found
+at :cite:`Kle2D`
+
+KarhunenLoeve2D Class
+^^^^^^^^^^^^^^^^^^^^
+
+The :class:`.KarhunenLoeve2D` class is imported using the following command:
+
+>>> from UQpy.stochastic_process.KarhunenLoeveExpansionTwoDimension2D import KarhunenLoeveExpansion
+
+Methods
+"""""""
+.. autoclass:: UQpy.stochastic_process.KarhunenLoeveExpansion2D
+    :members: run
+
+Attributes
+""""""""""
+.. autoattribute:: UQpy.stochastic_process.KarhunenLoeveExpansion2D.samples
+.. autoattribute:: UQpy.stochastic_process.KarhunenLoeveExpansion2D.xi
+
+Examples
+""""""""""
+
+.. toctree::
+
+   Karhunen Loeve Examples <../auto_examples/stochastic_processes/karhunen_loeve_2d/index>

requirements.txt

@@ -1,4 +1,4 @@
-numpy == 1.22.3
+numpy == 1.21.4
 scipy == 1.8.0
 matplotlib == 3.5.2
 scikit-learn == 1.0.2

src/UQpy/distributions/baseclass/Copula.py

@@ -2,7 +2,7 @@
     DistributionContinuous1D,
 )
 from abc import ABC
-
+from typing import Union
 
 class Copula(ABC):
 
@@ -44,7 +44,7 @@ def update_parameters(self, **kwargs: dict):
             self.parameters[key] = kwargs[key]
 
     @staticmethod
-    def check_marginals(marginals: list[DistributionContinuous1D]):
+    def check_marginals(marginals: Union[list, DistributionContinuous1D]):
         """
         Perform some checks on the marginals, raise errors if necessary.
 

src/UQpy/reliability/taylor_series/FORM.py

@@ -185,7 +185,7 @@ def run(self, seed_x: Union[list, np.ndarray] = None,
         g_record.append(0.0)
         dg_u_record = np.zeros([self.n_iterations + 1, self.dimension])
 
-        while not converged:
+        while not converged and k < self.n_iterations:
             self.logger.info("Number of iteration: %i", k)
             # FORM always starts from the standard normal space
             if k == 0:
@@ -314,8 +314,6 @@ def run(self, seed_x: Union[list, np.ndarray] = None,
 
             self.logger.info("Error: %s", error_record[-1])
 
-            if converged is True or k > self.n_iterations:
-                break
 
         if k > self.n_iterations:
             self.logger.info("UQpy: Maximum number of iterations {0} was reached before convergence."