ForK Library
Functions/Subroutines

kriginggradpredictGEK.f90 File Reference

Subroutine to predict mean derivative values from Gradient-Enhanced Kriging Surface (requires buildkrigingGEK to be called first) More...

Go to the source code of this file.

Functions/Subroutines

subroutine kriginggradpredictGEK (ndim, ntot, X, gtot, pts, dims, stot, H, beta, V, hyper, mtot, Xm, gtotm, ptsm, dimsm, Gm, dYm, covarflagi)
 Predicting derivative values from a gradient-enhanced Kriging model follows the same procedure as predicting derivative values from function-only Kriging models (see kriginggradpredict for details). For a gradient-enhanced model, derivative values are predicted as:

Detailed Description

Subroutine to predict mean derivative values from Gradient-Enhanced Kriging Surface (requires buildkrigingGEK to be called first)

Definition in file kriginggradpredictGEK.f90.


Function Documentation

subroutine kriginggradpredictGEK ( integer,intent(in)  ndim,
integer,intent(in)  ntot,
real(8),dimension(ndim,ntot),intent(in)  X,
integer,intent(in)  gtot,
integer,dimension(gtot),intent(in)  pts,
integer,dimension(gtot),intent(in)  dims,
integer,intent(in)  stot,
real(8),dimension(stot,ntot+gtot),intent(in)  H,
real(8),dimension(stot),intent(in)  Beta,
real(8),dimension(ntot+gtot),intent(in)  V,
real(8),dimension(ndim+3),intent(in)  hyper,
integer,intent(in)  mtot,
real(8),dimension(ndim,mtot),intent(in)  Xm,
integer,intent(in)  gtotm,
integer,dimension(gtotm),intent(in)  ptsm,
integer,dimension(gtotm),intent(in)  dimsm,
real(8),dimension(stot,gtotm),intent(in)  Gm,
real(8),dimension(gtotm),intent(out)  dYm,
integer,intent(in)  covarflagi 
)

Predicting derivative values from a gradient-enhanced Kriging model follows the same procedure as predicting derivative values from function-only Kriging models (see kriginggradpredict for details). For a gradient-enhanced model, derivative values are predicted as:

\[ \frac{\partial y(\vec{x}_{*})}{\partial x_{k}} = g_{k}(\vec{x}_{*}) \beta + \underline{l}_{*,k}^{T} \underline{K}^{-1} (\underline{Y} - \underline{H}^{T} \beta) \]

where the underline denotes the block definition required to account for derivative values. The block definitions of the training data, collocation matrix and covariance matrix are given in the documentation of buildkrigingGEK. The matrix \( \underline{l}_{*,k}^{T} \) is a block matrix whose elements represent covariance between derivative values and function values as well as derivative values with other derivative values: where the elements of \(cov( \frac{\partial y^{*}}{\partial x_{k}}, Y)\) are given by:

\[ cov(\frac{\partial y^{*}}{\partial x_{k}}, Y_{j}) = \frac{\partial}{\partial x_{k}} k(\vec{x}_{*},X_{j}) \]

and the elements of \(cov(\frac{\partial y^{*}}{\partial x_{k}}, \delta Y)\) are given as:

\[ cov(\frac{\partial y^{*}}{\partial x_{k}}, \delta Y_{j,l}) = \frac{\partial^{2}}{\partial x_{k} \partial x_{l}} k(\vec{x}_{*},X_{j}) \]

Using the processed data produced by buildkrigingGEK, derivative values can be calculated as:

\[ \frac{\partial y(\vec{x}_{*})}{\partial x_{k}} = g_{k}(\vec{x}_{*}) \beta + \underline{l}_{*,k}^{T} \underline{V} \]

Author:
Brian Lockwood Department of Mechanical Engineering University of Wyoming
Date:
May 2, 2012
Parameters:
in)ndim : The dimension of the problem
in)ntot : The number of Training points
in)X : The location of the training points (size=[ndimxntot])
in)gtot: Number of derivative values included in training data (ndim*ntot if all derivatives are included at the training points)
in)pts: List identifying what point the derivative value is enforced at (size=[gtot] with values ranging from 1 to ntot)
in)dims: List identifying the dimension the derivative is taken with respect to (size=[gtot] with values ranging from 1 to ndim)
in)stot : Number of Terms in the regression
in)H: The collocation matrix for the regression including derivative values. (size=[stotxntot+gtot])
Columns 1:ntot are the basis evaluated at the training points
Columns ntot+1:ntot+gtot are the derivative of the basis evaluated at the training points
in) beta: Regression coefficients based on the optimal estimate for the Kriging model (size=[stot])
Supplied by buildkrigingGEK subroutine
in)V: Processed Training Data (size=[ntot+gtot])
Supplied by buildkrigingGEK subroutine
in)hyper: Hyperparameters for the Kriging Model (size=[ndim+3])
Supplied by buildkrigingGEK subroutine
in)mtot : The number of test points, the places where function prediction are desired
in)Xm : The location of the test points (size=[ndimxmtot])
in)gtotm: The total number of derivatives to be predicted (ndim*mtot if the gradient at every test point is desired)
in) ptsm: A vector identifying the test point where a particular derivative is specified (size=[gtotm] with values ranging from 1 to mtot)
in) dimsm: A vector identifying the dimension the derivative is taken with respect to (size=[gtotm] with values ranging from 1 to ndim)
in)Gm : The derivative of the regression basis evaluated at the test points (size=[stotxgtotm])
out)dYm: The predicted function values (size=[gtotm])
Using the processed data V, predicting derivative values is essentially linear with respect to the number of test points so mtot can be set to one and this subroutine can be called multiple times if desired.
in)covarflagi: Flag to govern which covariance function is used

  • covarflag==1 Uses Matern function with \(\nu=3/2\)
  • covarflag==2 Uses Matern function with \(\nu=5/2\)


The parameter \(\nu\) governs the smoothness and differentiability of the covariance function.
When using gradient values, \(\nu=1/2\) is not differentiable enough so \(\nu \geq 3/2\) must be used
Must supply the same covariance flag as used in buildkrigingGEK.

Definition at line 69 of file kriginggradpredictGEK.f90.

References covars::covarflag.

Referenced by krigingwrapper().

 All Classes Namespaces Files Functions Variables