-
Notifications
You must be signed in to change notification settings - Fork 81
Commit
This commit does not belong to any branch on this repository, and may belong to a fork outside of the repository.
- Loading branch information
Showing
4 changed files
with
121 additions
and
113 deletions.
There are no files selected for viewing
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
118 changes: 118 additions & 0 deletions
118
src/UQpy/stochastic_process/KarhunenLoeveExpansion2D.py
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,118 @@ | ||
from typing import Union | ||
|
||
import numpy as np | ||
|
||
from UQpy.utilities.ValidationTypes import RandomStateType | ||
from UQpy.stochastic_process.KarhunenLoeveExpansion import KarhunenLoeveExpansion | ||
|
||
|
||
class KarhunenLoeveExpansion2D: | ||
|
||
def __init__( | ||
self, | ||
n_samples: int, | ||
correlation_function: np.ndarray, | ||
time_intervals: Union[np.ndarray, float], | ||
thresholds: Union[list, int] = None, | ||
random_state: RandomStateType = None, | ||
random_variables=None | ||
): | ||
""" | ||
A class to simulate two dimensional stochastic fields from a given auto-correlation function based on the | ||
Karhunen-Loeve Expansion | ||
:param n_samples: Number of samples of the stochastic field to be simulated. | ||
The :meth:`run` method is automatically called if `n_samples` is provided. If `n_samples` is not | ||
provided, then the :class:`.KarhunenLoeveExpansionTwoDimension` object is created but samples are not generated. | ||
:param correlation_function: The correlation function of the stochastic process of size | ||
:code:`(n_time_intervals_dim1, n_time_intervals_dim1, n_time_intervals_dim2, n_time_intervals_dim2)` | ||
:param time_intervals: The length of time discretizations. | ||
:param thresholds: The threshold number of eigenvalues to be used in the expansion for two dimensions. | ||
:param random_state: Random seed used to initialize the pseudo-random number generator. Default is :any:`None`. | ||
:param random_variables: The random variables used to generate the stochastic field. | ||
""" | ||
|
||
self.n_samples = n_samples | ||
self.correlation_function = correlation_function | ||
assert (len(self.correlation_function.shape) == 4) | ||
self.time_intervals = time_intervals | ||
self.thresholds = thresholds | ||
self.random_state = random_state | ||
|
||
if isinstance(self.random_state, int): | ||
np.random.seed(self.random_state) | ||
elif not isinstance(self.random_state, (type(None), np.random.RandomState)): | ||
raise TypeError('UQpy: random_state must be None, an int or an np.random.RandomState object.') | ||
|
||
self.samples = None | ||
self.random_variables = random_variables | ||
|
||
self._precompute_one_dimensional_correlation_function() | ||
|
||
if self.n_samples is not None: | ||
self.run(n_samples=self.n_samples, random_variables=random_variables) | ||
|
||
def _precompute_one_dimensional_correlation_function(self): | ||
self.quasi_correlation_function = np.diagonal(self.correlation_function, axis1=0, axis2=1) | ||
self.quasi_correlation_function = np.einsum('ij... -> ...ij', self.quasi_correlation_function) | ||
self.w, self.v = np.linalg.eig(self.quasi_correlation_function) | ||
if np.linalg.norm(np.imag(self.w)) > 0: | ||
print('Complex in the eigenvalues, check the positive definiteness') | ||
self.w = np.real(self.w) | ||
self.v = np.real(self.v) | ||
if self.thresholds is not None: | ||
self.w = self.w[:, :self.thresholds[1]] | ||
self.v = self.v[:, :, :self.thresholds[1]] | ||
self.one_dimensional_correlation_function = np.einsum('uvxy, uxn, vyn, un, vn -> nuv', | ||
self.correlation_function, self.v, self.v, | ||
1 / np.sqrt(self.w), 1 / np.sqrt(self.w)) | ||
|
||
def run(self, n_samples, random_variables=None): | ||
""" | ||
Execute the random sampling in the :class:`.KarhunenLoeveExpansion` class. | ||
The :meth:`run` method is the function that performs random sampling in the :class: | ||
`.KarhunenLoeveExpansionTwoDimension` class. If `n_samples` is provided when the :class: | ||
`.KarhunenLoeveExpansionTwoDimension` object is defined, the :meth:`run` method is automatically called. The | ||
user may also call the :meth:`run` method directly to generate samples. The :meth:`run`` method of the :class: | ||
`.KarhunenLoeveExpansionTwoDimension` class can be invoked many times and each time the generated samples are | ||
appended to the existing samples. | ||
:param n_samples: Number of samples of the stochastic process to be simulated. | ||
If the :meth:`run` method is invoked multiple times, the newly generated samples will be appended to the | ||
existing samples. | ||
The :meth:`run` method has no returns, although it creates and/or appends the :py:attr:`samples` attribute of | ||
the :class:`KarhunenLoeveExpansionTwoDimension` class. | ||
""" | ||
samples = np.zeros((n_samples, self.correlation_function.shape[0], self.correlation_function.shape[2])) | ||
if random_variables is None: | ||
random_variables = np.random.normal(size=[self.thresholds[1], self.thresholds[0], n_samples]) | ||
else: | ||
assert (random_variables.shape == (self.thresholds[1], self.thresholds[0], n_samples)) | ||
for i in range(self.one_dimensional_correlation_function.shape[0]): | ||
if self.thresholds is not None: | ||
print(self.w.shape) | ||
print(self.v.shape) | ||
samples += np.einsum('x, xt, nx -> nxt', np.sqrt(self.w[:, i]), self.v[:, :, i], | ||
KarhunenLoeveExpansion(n_samples=n_samples, | ||
correlation_function= | ||
self.one_dimensional_correlation_function[i], | ||
time_interval=self.time_intervals, | ||
threshold=self.thresholds[0], | ||
random_variables=random_variables[i]).samples[:, 0, :]) | ||
else: | ||
samples += np.einsum('x, xt, nx -> nxt', np.sqrt(self.w[:, i]), self.v[:, :, i], | ||
KarhunenLoeveExpansion(n_samples=n_samples, | ||
correlation_function= | ||
self.one_dimensional_correlation_function[i], | ||
time_interval=self.time_intervals, | ||
random_variables=random_variables[i]).samples[:, 0, :]) | ||
samples = np.reshape(samples, [samples.shape[0], 1, samples.shape[1], samples.shape[2]]) | ||
|
||
if self.samples is None: | ||
self.samples = samples | ||
self.random_variables = random_variables | ||
else: | ||
self.samples = np.concatenate((self.samples, samples), axis=0) | ||
self.random_variables = np.concatenate((self.random_variables, random_variables), axis=2) |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters