maths_library Module

Library of mathematical and numerical routines author: P J Knight, CCFE, Culham Science Centre N/A This module contains a large number of routines to enable PROCESS to perform a variety of numerical procedures, including linear algebra, zero finding, integration and minimisation. The module is an amalgamation of the contents of several different pre-existing PROCESS source files, which themselves were derived from a number of different numerical libraries including BLAS and MINPAC. http://en.wikipedia.org/wiki/Gamma_function



Contents


Functions

public function find_y_nonuniform_x(x0, x, y, n)

Routine to find y0 such that y0 = y(x0) given a set of values x(1:n), y(1:n) author: P J Knight, CCFE, Culham Science Centre x0 : input real : x value at which we want to find y x(1:n) : input real array : monotonically de/increasing x values y(1:n) : input real array : y values at each x n : input integer : size of array This routine performs a simple linear interpolation method to find the y value at x = x0. If x0 lies outside the range [x(1),x(n)], the y value at the nearest 'end' of the data is returned. None

Arguments

Type IntentOptional AttributesName
real(kind=dp), intent(in) :: x0
real(kind=dp), intent(in), dimension(n):: x
real(kind=dp), intent(in), dimension(n):: y
integer, intent(in) :: n

Return Value real(kind=dp)

public function binomial(n, k) result(coefficient)

Arguments

Type IntentOptional AttributesName
integer, intent(in) :: n
integer, intent(in) :: k

Return Value real(kind=dp)

public pure function variable_error(variable)

Arguments

Type IntentOptional AttributesName
real(kind=dp), intent(in) :: variable

Return Value logical

public pure function integer2string(value)

Arguments

Type IntentOptional AttributesName
integer, intent(in) :: value

Return Value character(len=2)

public pure function integer3string(value)

Arguments

Type IntentOptional AttributesName
integer, intent(in) :: value

Return Value character(len=3)


Subroutines

public subroutine ellipke(sqk, kk, ek)

Routine that calculates the complete elliptic integral of the first and second kinds author: P J Knight, CCFE, Culham Science Centre sqk : input real : square of the elliptic modulus kk : output real : complete elliptic integral of the first kind ek : output real : complete elliptic integral of the second kind This routine calculates the complete elliptic integral of the first and second kinds.

The method used is that described in the reference, and the code is taken from the Culham maglib library routine FN02A. Approximations for Digital Computers, C. Hastings, Princeton University Press, 1955

Arguments

Type IntentOptional AttributesName
real(kind=dp), intent(in) :: sqk
real(kind=dp), intent(out) :: kk
real(kind=dp), intent(out) :: ek

public subroutine svd(nm, m, n, a, w, matu, u, matv, v, ierr, rv1)

Singular Value Decomposition author: P J Knight, CCFE, Culham Science Centre author: B. S. Garbow, Applied Mathematics Division, Argonne National Laboratory nm : input integer : Max number of rows of arrays a, u, v; >= m,n m : input integer : Actual number of rows of arrays a, u n : input integer : Number of columns of arrays a, u, and the order of v a(nm,n) : input/output real array : On input matrix to be decomposed; on output, either unchanged or overwritten with u or v w(n) : output real array : The n (non-negative) singular values of a (the diagonal elements of s); unordered. If an error exit is made, the singular values should be correct for indices ierr+1,ierr+2,...,n. matu : input logical : Set to .true. if the u matrix in the decomposition is desired, and to .false. otherwise. u(nm,n) : output real array : The matrix u (orthogonal column vectors) of the decomposition if matu has been set to .true., otherwise u is used as a temporary array. u may coincide with a. If an error exit is made, the columns of u corresponding to indices of correct singular values should be correct. matv : input logical : Set to .true. if the v matrix in the decomposition is desired, and to .false. otherwise. v(nm,n) : output real array : The matrix v (orthogonal) of the decomposition if matv has been set to .true., otherwise v is not referenced. v may also coincide with a if u is not needed. If an error exit is made, the columns of v corresponding to indices of correct singular values should be correct. ierr : output integer : zero for normal return, or k if the k-th singular value has not been determined after 30 iterations. rv1(n) : output real array : work array This subroutine is a translation of the algol procedure SVD, Num. Math. 14, 403-420(1970) by Golub and Reinsch, Handbook for Auto. Comp., vol II - Linear Algebra, 134-151(1971).

It determines the singular value decomposition a=usvt of a real m by n rectangular matrix. Householder bidiagonalization and a variant of the QR algorithm are used. None

Arguments

Type IntentOptional AttributesName
integer, intent(in) :: nm
integer, intent(in) :: m
integer, intent(in) :: n
real(kind=dp), intent(inout), dimension(nm,n):: a
real(kind=dp), intent(out), dimension(n):: w
logical, intent(in) :: matu
real(kind=dp), intent(out), dimension(nm,n):: u
logical, intent(in) :: matv
real(kind=dp), intent(out), dimension(nm,n):: v
integer, intent(out) :: ierr
real(kind=dp), intent(out), dimension(n):: rv1

public subroutine eshellvol(rshell, rmini, rmino, zminor, drin, drout, dz, vin, vout, vtot)

Routine to calculate the inboard, outboard and total volumes of a toroidal shell comprising two elliptical sections author: P J Knight, CCFE, Culham Science Centre rshell : input real : major radius of centre of both ellipses (m) rmini : input real : horizontal distance from rshell to outer edge of inboard elliptical shell (m) rmino : input real : horizontal distance from rshell to inner edge of outboard elliptical shell (m) zminor : input real : vertical internal half-height of shell (m) drin : input real : horiz. thickness of inboard shell at midplane (m) drout : input real : horiz. thickness of outboard shell at midplane (m) dz : input real : vertical thickness of shell at top/bottom (m) vin : output real : volume of inboard section (m3) vout : output real : volume of outboard section (m3) vtot : output real : total volume of shell (m3) This routine calculates the volume of the inboard and outboard sections of a toroidal shell defined by two co-centred semi-ellipses. Each section's internal and external surfaces are in turn defined by two semi-ellipses. The volumes of each section are calculated as the difference in those of the volumes of revolution enclosed by their inner and outer surfaces.

See also eshellarea Internal CCFE note T&M/PKNIGHT/PROCESS/009, P J Knight: Surface Area and Volume Calculations for Toroidal Shells

Arguments

Type IntentOptional AttributesName
real(kind=dp), intent(in) :: rshell
real(kind=dp), intent(in) :: rmini
real(kind=dp), intent(in) :: rmino
real(kind=dp), intent(in) :: zminor
real(kind=dp), intent(in) :: drin
real(kind=dp), intent(in) :: drout
real(kind=dp), intent(in) :: dz
real(kind=dp), intent(out) :: vin
real(kind=dp), intent(out) :: vout
real(kind=dp), intent(out) :: vtot

public subroutine dshellvol(rmajor, rminor, zminor, drin, drout, dz, vin, vout, vtot)

Routine to calculate the inboard, outboard and total volumes of a D-shaped toroidal shell author: P J Knight, CCFE, Culham Science Centre rmajor : input real : major radius to outer point of inboard straight section of shell (m) rminor : input real : horizontal internal width of shell (m) zminor : input real : vertical internal half-height of shell (m) drin : input real : horiz. thickness of inboard shell at midplane (m) drout : input real : horiz. thickness of outboard shell at midplane (m) dz : input real : vertical thickness of shell at top/bottom (m) vin : output real : volume of inboard straight section (m3) vout : output real : volume of outboard curved section (m3) vtot : output real : total volume of shell (m3) This routine calculates the volume of the inboard and outboard sections of a D-shaped toroidal shell defined by the above input parameters. The inboard section is assumed to be a cylinder of uniform thickness. The outboard section's internal and external surfaces are defined by two semi-ellipses, centred on the outer edge of the inboard section; its volume is calculated as the difference in those of the volumes of revolution enclosed by the two surfaces.

See also dshellarea Internal CCFE note T&M/PKNIGHT/PROCESS/009, P J Knight: Surface Area and Volume Calculations for Toroidal Shells

Arguments

Type IntentOptional AttributesName
real(kind=dp), intent(in) :: rmajor
real(kind=dp), intent(in) :: rminor
real(kind=dp), intent(in) :: zminor
real(kind=dp), intent(in) :: drin
real(kind=dp), intent(in) :: drout
real(kind=dp), intent(in) :: dz
real(kind=dp), intent(out) :: vin
real(kind=dp), intent(out) :: vout
real(kind=dp), intent(out) :: vtot

public subroutine dshellarea(rmajor, rminor, zminor, ain, aout, atot)

Routine to calculate the inboard, outboard and total surface areas of a D-shaped toroidal shell author: P J Knight, CCFE, Culham Science Centre rmajor : input real : major radius of inboard straight section (m) rminor : input real : horizontal width of shell (m) zminor : input real : vertical half-height of shell (m) ain : output real : surface area of inboard straight section (m3) aout : output real : surface area of outboard curved section (m3) atot : output real : total surface area of shell (m3) This routine calculates the surface area of the inboard and outboard sections of a D-shaped toroidal shell defined by the above input parameters. The inboard section is assumed to be a cylinder. The outboard section is defined by a semi-ellipse, centred on the major radius of the inboard section.

See also dshellvol Internal CCFE note T&M/PKNIGHT/PROCESS/009, P J Knight: Surface Area and Volume Calculations for Toroidal Shells

Arguments

Type IntentOptional AttributesName
real(kind=dp), intent(in) :: rmajor
real(kind=dp), intent(in) :: rminor
real(kind=dp), intent(in) :: zminor
real(kind=dp), intent(out) :: ain
real(kind=dp), intent(out) :: aout
real(kind=dp), intent(out) :: atot

public subroutine eshellarea(rshell, rmini, rmino, zminor, ain, aout, atot)

Routine to calculate the inboard, outboard and total surface areas of a toroidal shell comprising two elliptical sections author: P J Knight, CCFE, Culham Science Centre rshell : input real : major radius of centre of both ellipses (m) rmini : input real : horizontal distance from rshell to inboard elliptical shell (m) rmino : input real : horizontal distance from rshell to outboard elliptical shell (m) zminor : input real : vertical internal half-height of shell (m) ain : output real : surface area of inboard section (m3) aout : output real : surface area of outboard section (m3) atot : output real : total surface area of shell (m3) This routine calculates the surface area of the inboard and outboard sections of a toroidal shell defined by two co-centred semi-ellipses.

See also eshellvol Internal CCFE note T&M/PKNIGHT/PROCESS/009, P J Knight: Surface Area and Volume Calculations for Toroidal Shells

Arguments

Type IntentOptional AttributesName
real(kind=dp), intent(in) :: rshell
real(kind=dp), intent(in) :: rmini
real(kind=dp), intent(in) :: rmino
real(kind=dp), intent(in) :: zminor
real(kind=dp), intent(out) :: ain
real(kind=dp), intent(out) :: aout
real(kind=dp), intent(out) :: atot