ForK Library
Functions/Subroutines

kriginggradvarianceGEK.f90 File Reference

Subroutine to predict the variance of 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 kriginggradvarianceGEK (ndim, ntot, X, gtot, pts, dims, stot, H, hyper, mtot, Xm, gtotm, ptsm, dimsm, Gm, Var, covarflagi)
 This subroutine calculates the variance associated with derivative predictions from a gradient-enhanced Kriging model. The details of this variance prediction are given in kriginggradvariance. The form of these variance predictions are the same for a gradient-enhanced and function-only model. The variance predictions for derivative values are given as:

Detailed Description

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

Definition in file kriginggradvarianceGEK.f90.


Function Documentation

subroutine kriginggradvarianceGEK ( 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(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)  Var,
integer,intent(in)  covarflagi 
)

This subroutine calculates the variance associated with derivative predictions from a gradient-enhanced Kriging model. The details of this variance prediction are given in kriginggradvariance. The form of these variance predictions are the same for a gradient-enhanced and function-only model. The variance predictions for derivative values are given as:

\[ V\left[\frac{\partial y(\vec{x}_{*})}{\partial x_{k}} \right]= \frac{\partial^{2}}{\partial x_{k}^{2}} k(\vec{x}_{*},\vec{x}_{*})-\underline{l}_{*,k}^T \underline{K}^{-1} \underline{l}_{*,k}+ \underline{S}_{k}(\vec{x}_*) \underline{A}^{-1} \underline{S}_{k}(\vec{x}_*)^{T} \]

See the subroutines krigingfuncpredictGEK, krigingfuncvarianceGEK and kriginggradpredictGEK for the definitions of the gradient-enhanced versions of the covariance matrices, collocation matrix, regression matrix and training data. The only previously undefined term is \(\underline{S}_{k}(\vec{x}_*)\) which is given as:

\[ \underline{S}_{k}(\vec{x}_*) = g_{k}(\vec{x}_{*}) - \underline{l}_{*,k}^T \underline{K}^{-1} \underline{H} \]

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)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)Var: The predicted function values (size=[gtotm])
Because the variance prediction requires the covariance matrix of the training data, it is more expensive than the mean prediction
If possible, predict the variance associated with multiple points using a single subroutine call
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 57 of file kriginggradvarianceGEK.f90.

References covars::covarflag, covars::d2covarfunc(), choleskymod::invertsymmetric(), choleskymod::matrixmulttrans(), choleskymod::symmatrixmulttrans(), and choleskymod::symmatvec().

Referenced by krigingwrapper().

 All Classes Namespaces Files Functions Variables