Spectra and Eddy Lifetime Function Learning#

Please refer to this UML diagram to see how data classes are composed to inform the CalibrationProblem classes.

Uses of this class require the OnePointSpectraDataGenerator, please refer to that documentation for specifics about data needed for spectra fitting.

CalibrationProblem#

class drdmannturb.CalibrationProblem(device: str, nn_params: NNParameters, prob_params: ProblemParameters, loss_params: LossParameters, phys_params: PhysicalParameters, logging_directory: str | None = None, output_directory: Path | str = PosixPath('/home/ai/Documents/wind/DRDMannTurb/docs/results'))[source]#

Class which manages the spectra fitting and eddy lifetime function learning based on the deep rapid distortion model developed in Keith, Khristenko, Wohlmuth (2021).

This class manages the operator regression task which characterizes the best fitting candidate in a family of nonlocal covariance kernels that are parametrized by a neural network.

This class can also be used independently of neural networks via the EddyLifetimeType used for classical spectra fitting tasks, for instance, using the EddyLifetimeType.MANN results in a fit that completely relies on the Mann eddy lifetime function. If a neural network model is used, Torch’s LBFGS optimizer is used with cosine annealing for learning rate scheduling. Parameters for these components of the training process are set in LossParameters and ProblemParameters during initialization.

After instantiating CalibrationProblem, wherein the problem and eddy lifetime function substitution type are indicated, the user will need to generate the OPS data using OnePointSpectraDataGenerator, after which the model can be fit with CalibrationProblem.calibrate.

After training, this class can be used in conjunction with the fluctuation generation utilities in this package to generate realistic turbulence fields based on the learned spectra and eddy lifetimes.

Constructor for CalibrationProblem class. As depicted in the UML diagram, this requires of 4 dataclasses.

Parameters:

device (str,) – One of the strings "cpu", "cuda", "mps" indicating the torch device to use
nn_params (NNParameters) – A NNParameters (for Neural Network) dataclass instance, which defines values of interest eg. size and depth. By default, calls constructor using default values.
prob_params (ProblemParameters) – A ProblemParameters dataclass instance, which is used to determine the conditional branching and computations required, among other things. By default, calls constructor using default values
loss_params (LossParameters) – A LossParameters dataclass instance, which defines the loss function terms and related coefficients. By default, calls constructor using the default values.
phys_params (PhysicalParameters) – A PhysicalParameters dataclass instance, which defines the physical constants governing the problem setting; note that the PhysicalParameters constructor requires three positional arguments.
output_directory (Union[Path, str], optional) – The directory to write output to; by default "./results"
logging_directory (str | None) –

calibrate(data: tuple[Iterable[float], torch.Tensor], tb_comment: str = '', optimizer_class: torch.optim.Optimizer = torch.optim.LBFGS) → Dict[str, float][source]#

Calibration method, which handles the main training loop and some data pre-processing. Currently the only supported optimizer is Torch’s LBFGS and the learning rate scheduler uses cosine annealing. Parameters for these components of the training process are set in LossParameters and ProblemParameters during initialization.

See the .print_calibrated_params() method to receive a pretty-printed output of the calibrated physical parameters.

Parameters:

data (tuple[Iterable[float], torch.Tensor]) – Tuple of data points and values, respectively, to be used in calibration.
tb_comment (str) – Filename comment used by tensorboard; useful for distinguishing between architectures and hyperparameters. Refer to tensorboard documentation for examples of use. By default, the empty string, which results in default tensorboard filenames.
optimizer_class (torch.optim.Optimizer, optional) – Choice of Torch optimizer, by default torch.optim.LBFGS

Returns:

Dict[str, float] – Physical parameters for the problem, in normal space (not in log-space, as represented internally). Presented as a dictionary in the order {L : length scale, Gamma : time scale, sigma : spectrum amplitude}.

Raises:

RuntimeError – Thrown in the case that the current loss is not finite.

Return type:

Dict[str, float]

eval(k1: torch.Tensor) → ndarray[source]#

Calls the calibrated model on \(k_1\). This can be done after training or after loading trained model parameters from file.

Parameters:: k1 (torch.Tensor) – Tensor of \(k_1\) data
Returns:: np.ndarray – Evaluation of the model represented in a Numpy array (CPU bound)
Return type:: ndarray

eval_grad(k1: torch.Tensor)[source]#

Evaluates gradient of \(k_1\) via Autograd

Parameters:: k1 (torch.Tensor) – Tensor of :math:\(k_1\) data
Returns:: np.ndarray – Numpy array of resultant gradient (CPU bound)

eval_trainable_norm(ord: float | str | None = 'fro')[source]#

Evaluates the magnitude (or other norm) of the: trainable parameters in the model.

Note

The EddyLifetimeType must be set to one of TAUNET or CUSTOMMLP, which involve a network surrogate for the eddy lifetime.

Parameters:: ord (Optional[Union[float, str]]) – The order of the norm approximation, follows torch.norm conventions.
Raises:: ValueError – If the OPS was not initialized to one of TAUNET, CUSTOMMLP.

format_input(k1: torch.Tensor) → torch.Tensor[source]#

Wrapper around clone and cast \(k_1\) to torch.float64

Parameters:: k1 (torch.Tensor) – Tensor of \(k_1\)
Returns:: torch.Tensor – Copy of :math:`k_1 casted to doubles
Return type:: torch.Tensor

format_output(out: torch.Tensor) → ndarray[source]#

Wrapper around torch’s out.cpu().numpy(). Returns a CPU tensor.

Parameters:: out (torch.Tensor) – Tensor to be brought to CPU and converted to an np.ndarray
Returns:: np.ndarray – Numpy array of the input tensor
Return type:: ndarray

init_device(device: str)[source]#

Initializes the device (CPU or GPU) on which computation is performed.

Parameters:: device (str) – string following PyTorch conventions – “cuda” or “cpu”

initialize_parameters_with_noise()[source]#: Simple routine to introduce additive white noise to the OPS parameters.

log_dimensional_scales() → None[source]#

Sets the quantities for non-dimensionalization in log-space.

Note

The first 3 parameters of self.parameters() are exactly

LengthScale
TimeScale
Spectrum Amplitude

Return type:: None

num_trainable_params() → int[source]#

Computes the number of trainable network parameters

in the underlying model.

The EddyLifetimeType must be set to one of the following, which involve a network surrogate for the eddy lifetime:

TAUNET

CUSTOMMLP

Returns:: int – The number of trainable network parameters in the underlying model.
Raises:: ValueError – If the OPS was not initialized to one of TAUNET, CUSTOMMLP
Return type:: int

property parameters: ndarray#

Returns all parameters of the One Point Spectra surrogate model as a single vector.

Note

The first 3 parameters of self.parameters() are exactly #. LengthScale

TimeScale
Spectrum Amplitude

Returns:: np.ndarray – Single vector containing all model parameters on the CPU, which can be loaded into an object with the same architecture with the parameters setter method. This automatically offloads any model parameters that were on the GPU, if any.

plot(Data: tuple[Iterable[float], torch.Tensor] | None = None, model_vals: torch.Tensor | None = None, plot_tau: bool = True, save: bool = False, save_dir: Path | str | None = None, save_filename: str = '')[source]#

Plotting method which visualizes the spectra fit as well as the learned eddy lifetime function, if plot_tau=True. By default, this operates on the data used in the fitting, but an alternative \(k_1\) domain can be provided and the trained model can be re-evaluated.

Parameters:

Data (tuple[Iterable[float], torch.Tensor], optional) – Tuple of data points and corresponding values, by default None
model_vals (torch.Tensor, optional) – Evaluation of the OPS on the data, by default None in which case Data must provided (since the function will call OPS on the provided Data)
plot_tau (bool, optional) – Indicates whether to plot the learned eddy lifetime function or not, by default True
save (bool, optional) – Whether to save the resulting figure, by default False
save_dir (Optional[Union[Path, str]], optional) – Directory to save to, which is created safely if not already present. By default, this is the current working directory.
save_filename (str, optional) – Filename to save the final figure to, by default drdmannturb_final_spectra_fit.png

Raises:

ValueError – Must either provide k1space or re-use what was used for model calibration; thrown in the case neither is specified.
ValueError – Must either provide data points or re-use what was used for model calibration; thrown in the case neither is specified.
ValueError – Thrown in the case that save is true but neither the save_dir or output_directory are provided.

plot_losses(run_number: int)[source]#

A wrapper method around the plot_loss_logs helper, which plots out the loss function terms, multiplied by their associated hyperparameters

Parameters:: run_number (int) – The number of the run in the logging directory to plot out. This is 0-indexed.

print_calibrated_params()[source]#

Prints out the optimal calibrated parameters L, Gamma, sigma, which are stored in a fitted CalibrationProblem object under the calibrated_params dictionary.

Raises:: ValueError – Must call .calibrate() method to compute a fit to physical parameters first.

save_model(save_dir: str | None = None)[source]#

Saves model with current weights, model configuration, and training histories to file. The written filename is of the form save_dir/<EddyLifetimeType>_<DataType>.pkl

This routine stores

NNParameters

ProblemParameters

PhysicalParameters

LossParameters

Optimized Parameters (self.parameters field)

Parameters:: save_dir (Optional[str], optional) – Directory to save to, by default None; defaults to provided output_dir field for object.
Raises:: ValueError – No output_directory provided during object initialization and no save_dir provided for this method call.

OnePointSpectraDataGenerator#

class drdmannturb.OnePointSpectraDataGenerator(zref: float, ustar: float = 1.0, data_points: Iterable[Tuple[torch.tensor, float]] | None = None, data_type: DataType = DataType.KAIMAL, k1_data_points: torch.Tensor | None = None, spectra_values: torch.Tensor | None = None, spectra_file: Path | str | None = None, seed: int = 3)[source]#

The one point spectra data generator, which evaluates one of a few spectral tensor models across a grid of \(k_1\) wavevector points which is used as data in the fitting done in the OnePointSpectra class.

The type of spectral tensor is determined by the DataType argument, which determines one of the following models:

DataType.KAIMAL, which is the Kaimal spectra.
DataType.VK, which is the von Karman spectra model.
DataType.CUSTOM, usually used for data that is processed from real-world data. The spectra values are to be provided as the spectra_values field, or else to be loaded from a provided spectra_file. The result is that the provided data are matched on the wavevector domain.
DataType.AUTO, which generates a filtered version of provided spectra data. The filtering is based on differential evolution to perform a non-linear fit onto functions of the following form:

\begin{align} & \frac{k_1 F_{11}\left(k_1 z\right)}{u_*^2}=J_1(f):=\frac{a_1 f}{(1+b_1 f)^{c_1}} \\ & \frac{k_1 F_{22}\left(k_1 z\right)}{u_*^2}=J_2(f):=\frac{a_2 f}{(1+b_2 f)^{c_2}} \\ & \frac{k_1 F_{33}\left(k_1 z\right)}{u_*^2}=J_3(f):=\frac{a_3 f}{1+ b_3 f^{ c_3}} \\ & -\frac{k_1 F_{13}\left(k_1 z\right)}{u_*^2}=J_4(f):=\frac{a_4 f}{(1+ b_4 f)^{c_4}}, \end{align}

with \(F_{12}=F_{23}=0\). Here, \(f = (2\pi)^{-1} k_1 z\). In the above, the \(a_i, b_i, c_i\) are free parameters which are optimized by differential evolution. The result is a spectra model that is similar in form to the Kaimal spectra and which filters/smooths the spectra data from the real world and eases fitting by DRD models. This option is highly suggested in cases where spectra data have large deviations.

Note

The one-point spectra for \(F_{13}\) are NOT to be pre-multiplied with a negative, this data generator automatically performs this step both when using DataType.CUSTOM and DataType.AUTO.

Parameters:

zref (float) – Reference altitude value
ustar (float) – Friction velocity, by default 1.0.
data_points (Iterable[Tuple[torch.tensor, float]], optional) – Observed spectra data points at each of the \(k_1\) coordinates, paired with the associated reference height (typically kept at 1, but may depend on applications).
data_type (DataType, optional) – Indicates the data format to generate and operate with, by default DataType.KAIMAL
k1_data_points (Optional[Any], optional) – Wavevector domain of \(k_1\), by default None. This is only to be used when the AUTO tag is chosen to define the domain over which a non-linear regression is to be computed. See the interpolation module for examples of unifying different \(k_1\) domains.
spectra_file (Optional[Path], optional) – If using DataType.CUSTOM or DataType.AUTO, this is used to indicate the data file (a .dat) to read from. Since it is not used by others, it is by default None
spectra_values (torch.Tensor | None) –
seed (int) –

Raises:

ValueError – In the case that DataType.CUSTOM is indicated, but no spectra_file is provided
ValueError – In the case that DataType.AUTO data is indicated, but no k1_data_pts is provided
ValueError – Did not provide DataPoints during initialization for DataType method requiring spectra data.

eval_Kaimal(k1: float) → torch.Tensor[source]#

Evaluates the one-point spectra as proposed by Kaimal et al in 1972. Clasically motivated by measurements taken over a flat homogeneous terrain in Kansas, the one-point spectra were proposed as

\begin{align} & \frac{k_1 F_{11}\left(k_1 z\right)}{u_*^2}=J_1(f):=\frac{52.5 f}{(1+33 f)^{5 / 3}} \\ & \frac{k_1 F_{22}\left(k_1 z\right)}{u_*^2}=J_2(f):=\frac{8.5 f}{(1+9.5 f)^{5 / 3}} \\ & \frac{k_1 F_{33}\left(k_1 z\right)}{u_*^2}=J_3(f):=\frac{1.05 f}{1+5.3 f^{5 / 3}} \\ & -\frac{k_1 F_{13}\left(k_1 z\right)}{u_*^2}=J_4(f):=\frac{7 f}{(1+9.6 f)^{12 / 5}}, \end{align}

with \(F_{12}=F_{23}=0\). Here, \(f = (2\pi)^{-1} k_1 z\). This method returns a \(3\times 3\) matrix whose entries are determined by the above equations.

Parameters:: k1 (torch.Tensor) – First dimension of the wavevector \(k_1\), this is the domain over which the associated spectra are defined.
Returns:: torch.Tensor – The \(3 \times 3\) matrix with entries determined by the Kaimal one-point spectra.
Return type:: torch.Tensor

eval_VK(k1: float) → torch.Tensor[source]#

Evaluation of von Karman spectral tensor given in the form

\[\Phi_{i j}^{\mathrm{VK}}(\boldsymbol{k})=\frac{E(k)}{4 \pi k^2}\left(\delta_{i j}-\frac{k_i k_j}{k^2}\right)\]

which utilizes the energy spectrum function

\[E(k)=c_0^2 \varepsilon^{2 / 3} k^{-5 / 3}\left(\frac{k L}{\left(1+(k L)^2\right)^{1 / 2}}\right)^{17 / 3},\]

where \(\varepsilon\) is the viscous dissipation of the turbulent kinetic energy, \(L\) is the length scale parameter and \(c_0^2 \approx 1.7\) is an empirical constant. Here, the physical constants are taken to be :math:L = 0.59.

Parameters:: k1 (torch.Tensor) – First dimension of the wavevector \(k_1\), this is a single coordinate of the domain over which the associated spectra are defined.
Returns:: torch.Tensor – The \(3 \times 3\) matrix with entries determined by the von Karman spectra.
Return type:: torch.Tensor

generate_Data(DataPoints: Iterable[Tuple[torch.tensor, float]]) → tuple[torch.Tensor, torch.Tensor][source]#

Generates a single spectral tensor from provided data or from a surrogate model. The resulting tensor is of shape (number of \(k_1\) points):math:times 3 times 3, ie the result consists of the spectral tensor evaluated across the provided range of \(k_1\) points. The spectra model is set during object instantiation.

Note

The DataType.CUSTOM type results in replication of the provided spectra data.

Parameters:: DataPoints (Iterable[Tuple[torch.tensor, float]]) – Observed spectra data points at each of the \(k_1\) coordinates, paired with the associated reference height (typically kept at 1, but may depend on applications).
Returns:: tuple[torch.Tensor, torch.Tensor] – Evaluated spectral tensor on each of the provided grid points depending on the DataType selected.
Raises:: ValueError – DataType set such that an iterable set of spectra values is required. This is for any DataType other than CUSTOM and AUTO.
Return type:: tuple[torch.Tensor, torch.Tensor]

plot(x_interp: ndarray | None = None, spectra_full: ndarray | None = None, x_coords_full: ndarray | None = None)[source]#

Utility for plotting current spectra data over the wavevector domain. Note that if the datatype is chosen to be DataType.AUTO, the interpolation coordinates in the frequency space must be provided again to determine the domain over which to plot filtered spectra.

Parameters:

x_interp (Optional[np.ndarray], optional) – Interpolation domain in normal space (not log-space), if used DataType.AUTO, by default None
spectra_full (Optional[np.ndarray], optional) – Full tensor of spectra data, by default None
x_coords_full (Optional[np.ndarray], optional) – Full \(k_1\) log-space coordinates, by default None

Raises:

ValueError – Provide interpolation domain (normal space) for filtered plots.
ValueError – Provide original spectra and associated domains as a single stack of np.arrays, even if all x coordinates are matching.

OnePointSpectra#

class drdmannturb.OnePointSpectra(*args: Any, **kwargs: Any)[source]#

One point spectra implementation, including set-up of eddy lifetime function approximation with DRD models, or several classical eddy lifetime functions.

Initialization of the one point spectra class. This requires the type of eddy lifetime function to use, the power spectra type (currently only the von Karman spectra is implemented), the neural network parameters to use if a DRD model is selected, and whether or not to learn \(\nu\) in

\[\tau(\boldsymbol{k})=\frac{T|\boldsymbol{a}|^{\nu-\frac{2}{3}}}{\left(1+|\boldsymbol{a}|^2\right)^{\nu / 2}}, \quad \boldsymbol{a}=\boldsymbol{a}(\boldsymbol{k}).\]

Here,

\[\boldsymbol{a}(\boldsymbol{k}):=\operatorname{abs}(\boldsymbol{k})+\mathrm{NN}(\operatorname{abs}(\boldsymbol{k}))\]

if a neural network is used to learn the eddy lifetime function. For a discussion of the details and training, refer to the original DRD paper by Keith et al.

Non-neural network eddy lifetime functions are provided as well, specifically the Mann model. The default power spectra used is due to von Karman.

Parameters:

type_eddy_lifetime (EddyLifetimeType, optional) – Type of eddy lifetime function \(\tau\) to use.
physical_params (PhysicalParameters,) – Object specifying physical parameters of the problem.
type_power_spectra (PowerSpectraType, optional) – Type of power spectra function to use, by default PowerSpectraType.RDT, the only one currently implemented.
nn_parameters (NNParameters, optional) – Dataclass containing neural network initialization if a neural network is used to approximate the eddy lifetime function, by default None.
learn_nu (bool, optional) – Whether or not to learn \(\nu\), by default False.

Raises:

ValueError – “Selected neural network-based eddy lifetime function approximation but did not specify neural network parameters. Pass a constructed NNParameters object.”

EddyLifetime(k: torch.Tensor | None = None) → torch.Tensor#

Evaluation of eddy lifetime function \(\tau\) constructed during object initialization. This may be the Mann model or a DRD neural network that learns \(\tau\).

Parameters:: k (Optional[torch.Tensor], optional) – Wavevector domain on which to evaluate the eddy lifetime function, by default None, which defaults to grids in logspace(-3, 3).
Returns:: torch.Tensor – Evaluation of current eddylifetime function on provided domain.
Raises:: Exception – Did not specify an admissible EddyLifetime model. Refer to the EddyLifetimeType documentation.
Return type:: torch.Tensor

InitialGuess_EddyLifetime()#

Initial guess at the eddy lifetime function which the DRD model uses in learning the \(\tau\) eddy lifetime function. By default, this is just the \(0\) function, but later functionality may allow this to be dynamically set.

Returns:: float – Initial guess evaluation, presently, the \(0\) function.

PowerSpectra()#

Calls the RDT Power Spectra model with current approximation of the eddy lifetime function and the energy spectrum.

Returns:: torch.Tensor – RDT power spectra evaluation.
Raises:: Exception – In the case that the Power Spectra is not RDT and therefore incorrect.

exp_scales() → tuple[float, float, float][source]#

Exponentiates the length, time, and spectrum amplitude scales,

Note

The first 3 parameters of self.parameters() are exactly

LengthScale
TimeScale
SpectrumAmplitude

Returns:: tuple[float, float, float] – Scalar values for each of the length, time, and magnitude scales, in that order.
Return type:: tuple[float, float, float]

forward(k1_input: torch.Tensor) → torch.Tensor[source]#

Evaluation of one point spectra

\[\widetilde{J}_i(f ; \boldsymbol{\theta})=C k_1 \widetilde{F}_{i i}\left(k_1 z ; \boldsymbol{\theta}\right),\]

defined by equations 6 (a-d) in the original DRD paper. Here,

\[\widetilde{F}_{i j}\left(k_1 ; \boldsymbol{\theta}\right)=\int_{-\infty}^{\infty} \int_{-\infty}^{\infty} \Phi_{i j}^{\mathrm{DRD}}(\boldsymbol{k}, \boldsymbol{\theta}) \mathrm{d} k_2 \mathrm{~d} k_3.\]

Parameters:: k1_input (torch.Tensor) – Discrete \(k_1\) wavevector domain.
Returns:: torch.Tensor – Network output
Return type:: torch.Tensor

get_div(Phi: torch.Tensor) → torch.Tensor#

Evaluates the divergence of an evaluated spectral tensor. This is evaluated simply as \(\textbf{k} \cdot \Phi_{\textbf{k}}\) and normalized by the trace.

Parameters:: Phi (torch.Tensor) – Discrete evaluated spectral tensor \(\Phi(\textbf{k}, \tau)\), which may or may not depend on the eddy lifetime function. For instance, if the von Karman model is used, no \(\tau\) dependence is present.
Returns:: torch.Tensor – Divergence of the spectral tensor in the frequency domain.
Return type:: torch.Tensor

init_mann_approximation()[source]#

Initializes Mann eddy lifetime function approximation by performing a linear regression in log-log space on a given wave space and the true output of

\[\frac{x^{-\frac{2}{3}}}{\sqrt{{ }_2 F_1\left(1 / 3,17 / 6 ; 4 / 3 ;-x^{-2}\right)}}\]

This operation is performed once on the CPU.

quad23(f: torch.Tensor) → torch.Tensor#

Computes an approximation of the integral of the discretized function \(f\) in the dimensions defined by \(k_2\) and \(k_3\) using the trapezoidal rule:

\[\int \int f(k_1 ; \mathbf{\theta}) d k_2 dk_3.\]

Parameters:: f (torch.Tensor) – Function evaluation (tensor) to integrate over the frequency domain constructed during initialization.
Returns:: torch.Tensor – Evaluated double integral.
Return type:: torch.Tensor

set_scales(LengthScale: float64, TimeScale: float64, Magnitude: float64)[source]#

Sets scalar values for values used in non-dimensionalization.

Parameters:

LengthScale (np.float64) – Length scale.
TimeScale (np.float64) – Time scale.
Magnitude (np.float64) – Spectrum amplitude magnitude.