v_space_procedures Submodule

Contains module procedures associated with the VSpace class


Uses


Module Functions

pure module function getVGrid(this) result(v)

Getter for vGrid

Arguments

Type IntentOptional Attributes Name
class(VSpace), intent(in) :: this

Return Value real(kind=rk), allocatable, dimension(:)

pure module function getNumH(this) result(numH)

Return total number of resolved harmonic on grid

Arguments

Type IntentOptional Attributes Name
class(VSpace), intent(in) :: this

Return Value integer(kind=ik)

pure module function getNumV(this) result(numV)

Return number of v points on grid

Arguments

Type IntentOptional Attributes Name
class(VSpace), intent(in) :: this

Return Value integer(kind=ik)

pure module function getVCellWidths(this) result(dv)

Getter for vWidths

Arguments

Type IntentOptional Attributes Name
class(VSpace), intent(in) :: this

Return Value real(kind=rk), allocatable, dimension(:)

pure module function getVLinInterp(this) result(linInterp)

Getter for linInterp

Arguments

Type IntentOptional Attributes Name
class(VSpace), intent(in) :: this

Return Value real(kind=rk), allocatable, dimension(:)

pure module function getNearestPoints(this, v) result(pair)

Return two nearest points for given velocity value v. If the first point is 0, the velocity is below the lowest velocity in the grid, and if the second point is 0 the velocity is above the greatest v in the grid.

Arguments

Type IntentOptional Attributes Name
class(VSpace), intent(in) :: this
real(kind=rk), intent(in) :: v

Return Value integer(kind=ik), dimension(2)

pure module function getContainingVCell(this, v) result(ind)

Return index of cell which containes the given velocity value v. If the returned index is 0, the point is outside of the grid.

Arguments

Type IntentOptional Attributes Name
class(VSpace), intent(in) :: this
real(kind=rk), intent(in) :: v

Return Value integer(kind=ik)

pure module function calculateMoment(this, f, h, mOrder, g, gDependsOnX) result(res)

Calculate moment of the h-th harmonic of passed local distribution function f, optionally multiplied by g(v). If gDependsOnX is .true. assumes that g is given as a strided array with stride numV, and it is used to allocate the result - this allows g to not have a halo while f does. The moment is of order mOrder.

Arguments

Type IntentOptional Attributes Name
class(VSpace), intent(in) :: this
real(kind=rk), intent(in), dimension(:) :: f
integer(kind=ik), intent(in) :: h
integer(kind=ik), intent(in) :: mOrder
real(kind=rk), intent(in), optional, dimension(:) :: g
logical, intent(in), optional :: gDependsOnX

Return Value real(kind=rk), allocatable, dimension(:)

pure module function getShkarofskyIMat(this, index) result(res)

Return sparse row data format for the lower triangular Shkarofsky I_index integral.

Arguments

Type IntentOptional Attributes Name
class(VSpace), intent(in) :: this
integer(kind=ik), intent(in) :: index

Return Value type(SparseRowData)

pure module function getShkarofskyJMat(this, index) result(res)

Return sparse row data format for the upper triangular Shkarofsky J_index integral.

Arguments

Type IntentOptional Attributes Name
class(VSpace), intent(in) :: this
integer(kind=ik), intent(in) :: index

Return Value type(SparseRowData)


Module Subroutines

pure module subroutine initVSpace(this, gridObj)

VSpace initialization routine

Arguments

Type IntentOptional Attributes Name
class(VSpace), intent(inout) :: this
type(Grid), intent(in) :: gridObj

Grid object used to initialize VSpace