atl_varSys_module Module

This module contains types and routines to work with the ateles variable system. This variable system is build upon the treelm routines, but contains solver specific enhancements.


Uses

Used by


Contents


Derived Types

type, public :: atl_varSys_data_type

Solver-specific structure for solver and source term variables.

Components

TypeVisibilityAttributesNameInitial
type(atl_varSys_solverData_type), private, pointer:: solverData
type(tem_pointData_list_type), private :: pointData

the point_datas need to be stored levelwise

type(tem_varSys_op_data_type), private :: opData

data array for operation or derived varibales consists the index arrys for points stored in the poingtData of input variable size is number of input variables

Solver-specific structure for solver and source term variables.

Read more…

Components

TypeVisibilityAttributesNameInitial
type(atl_scheme_type), private, pointer:: scheme_listPtr(:)=> null()
type(atl_statedata_type), private, pointer:: statedata_listPtr(:)=> null()
type(atl_kerneldata_type), private, pointer:: kerneldata_listPtr(:)=> null()
type(atl_Equations_type), private, pointer:: equationPtr=> null()
type(ply_poly_project_type), private, pointer:: polyProjectPtr(:)=> null()
type(atl_cube_elem_type), private, pointer:: mesh_listPtr(:)=> null()

Levelwise mesh information used for creating chebychev coordinates.

integer, private, pointer:: poly_proj_posPtr(:)=> null()
integer, private, pointer:: levelPointer(:)=> null()
type(treelmesh_type), private, pointer:: tree=> null()

type, private :: atl_varSys_solverVar_type

Components

TypeVisibilityAttributesNameInitial
character(len=labelLen), private :: vartype
type(atl_legpolyvar_type), private :: legpolyvar

Functions

public function atl_get_new_varSys_data_ptr(solverData) result(resPtr)

Routine to get a pointer to a new instance of atl_varSys_data_type to be used as method data for a variable in the variable system.

Read more…

Arguments

TypeIntentOptionalAttributesName
type(atl_varSys_solverData_type), intent(in), target:: solverData

The prototype is used to initialize the new instance.

Return Value type(c_ptr)

Pointer to the newly created instance.


Subroutines

public subroutine atl_init_varSys_solverData(this, tree, equation, schemeList, statedataList, kerneldataList, levelPointer, meshList, polyProjPos, polyProjList)

Arguments

TypeIntentOptionalAttributesName
type(atl_varSys_solverData_type), intent(inout) :: this
type(treelmesh_type), intent(in), target:: tree
type(atl_Equations_type), intent(in), target:: equation
type(atl_scheme_type), intent(in), target:: schemeList(tree%global%minLevel:)
type(atl_statedata_type), intent(in), target:: statedataList(tree%global%minLevel:)
type(atl_kerneldata_type), intent(in), target:: kerneldataList(tree%global%minLevel:)
integer, intent(in), target:: levelPointer(:)
type(atl_cube_elem_type), intent(in), target:: meshList(tree%global%minLevel:)

Levelwise mesh information used for creating chebychev coordinates.

integer, intent(in), target:: polyProjPos(tree%global%minLevel:)
type(ply_poly_project_type), intent(in), target:: polyProjList(:)

public subroutine atl_set_stFun_getElement(solData_evalElem, fun)

Arguments

TypeIntentOptionalAttributesName
class(tem_varSys_solverData_evalElem_type), intent(in) :: solData_evalElem

Description on how to set the element retrieval function for stfuns.

type(tem_varSys_op_type), intent(inout) :: fun

Description of the method to obtain the variables, here some preset values might be stored, like the space time function to use or the required variables.

public subroutine atl_varSys_getStateForElement(fun, varSys, elempos, time, tree, nElems, nDofs, res)

To obtain values of a given variable, it is necessary to state the treeID and time at which the variable should be evaluated. The interface is nDofs values to cover the all degrees of freedoms in the element. Of course the variable system itself also needs to be passed in, to allow the computation of other derived quantities as needed. The method description itself is passed in automatically, and has not to be provided explicitly.

Arguments

TypeIntentOptionalAttributesName
class(tem_varSys_op_type), intent(in) :: fun

Description of the method to obtain the variables, here some preset values might be stored, like the space time function to use or the required variables.

type(tem_varSys_type), intent(in) :: varSys

The variable system to obtain the variable from.

integer, intent(in) :: elempos(:)

Position of the TreeID of the element to get the variable for in the global treeID list.

type(tem_time_type), intent(in) :: time

Point in time at which to evaluate the variable.

type(treelmesh_type), intent(in) :: tree

global treelm mesh info

integer, intent(in) :: nElems

Number of values to obtain for this variable (vectorized access).

integer, intent(in) :: nDofs

Number of degrees of freedom within an element.

real(kind=rk), intent(out) :: res(:)

Resulting values for the requested variable.

Linearized array dimension: (n requested entries) x (nComponents of this variable) x (nDegrees of freedom) Access: (iElem-1)fun%nComponentsnDofs + (iDof-1)*fun%nComponents + iComp

public subroutine atl_varSys_getStateForPoint(fun, varSys, point, time, tree, nPnts, res)

Interface description for a variable access method (single point).

Read more…

Arguments

TypeIntentOptionalAttributesName
class(tem_varSys_op_type), intent(in) :: fun

Description of the method to obtain the variables, here some preset values might be stored, like the space time function to use or the required variables.

type(tem_varSys_type), intent(in) :: varSys

The variable system to obtain the variable from.

real(kind=rk), intent(in) :: point(:,:)

Three-dimensional coordinates at which the variable should be evaluated. Only useful for variables provided as space-time functions.

type(tem_time_type), intent(in) :: time

Point in time at which to evaluate the variable.

type(treelmesh_type), intent(in) :: tree

global treelm mesh info

integer, intent(in) :: nPnts

Number of values to obtain for this variable (vectorized access).

real(kind=rk), intent(out) :: res(:)

Resulting values for the requested variable.

Dimension: n requested entries x nComponents of this variable Access: (iElem-1)*fun%nComponents + iComp

public subroutine atl_varSys_setupStateIndices(fun, varSys, point, offset_bit, iLevel, tree, nPnts, idx)

This routine takes points coordinates, stores them in the method_data and return indices where points are located in the growing array of points or values ( sometimes we do not need to store the points ) It is need to setup points for every variable. Points will be provided by boundaries or sources depends on what uses the variable. This points do not change with time . This indices will be stored in corresponding boundary or source to evaluate a value on that point later using tem_varSys_proc_getValOfIndex.

Arguments

TypeIntentOptionalAttributesName
class(tem_varSys_op_type), intent(in) :: fun

Description of the method to obtain the variables, for this routine we need the location where to store the points.

type(tem_varSys_type), intent(in) :: varSys

The variable system to obtain the variable from.

real(kind=rk), intent(in) :: point(:,:)

arrays of points for which the indices are returned

character, intent(in), optional :: offset_bit(:)

Offset bit encoded as character for every point. If not present default is to center i.e offset_bit = achar(1+4+16)

integer, intent(in) :: iLevel

the point data need to be loaded levelwise, we need the current level

type(treelmesh_type), intent(in) :: tree

global treelm mesh info

integer, intent(in) :: nPnts

number of points

integer, intent(out) :: idx(:)

public subroutine atl_varSys_getStateValOfIndex(fun, varSys, time, iLevel, idx, idxLen, nVals, res)

Routine for gettint the actual value for a given array of indices. The indices belong to the grwarray of points storing levelwise in Pointdata%pntLvl(iLevel). Hence this routines takes the indeices as input, can refer to the pointData and evaluate the variable and returns the values

Read more…

Arguments

TypeIntentOptionalAttributesName
class(tem_varSys_op_type), intent(in) :: fun

Description of the method to obtain the variables,

type(tem_varSys_type), intent(in) :: varSys

The variable system to obtain the variable from.

type(tem_time_type), intent(in) :: time

Point in time at which to evaluate the variable.

integer, intent(in) :: iLevel

Level on which values are requested

integer, intent(in) :: idx(:)

Index of points in the growing array and variable val array to return. Size: n

integer, intent(in), optional :: idxLen(:)

With idx as start index in contiguous memory, idxLength defines length of each contiguous memory Size: nVals

integer, intent(in) :: nVals

Number of values to obtain for this variable (vectorized access).

real(kind=rk), intent(out) :: res(:)

Resulting values for the requested variable.

Dimension: n requested entries x nComponents of this variable Access: (iElem-1)*fun%nComponents + iComp

public subroutine atl_create_fortranVar(me)

This routine creates hard coded variables which are required by ateles.

Read more…

Arguments

TypeIntentOptionalAttributesName
type(tem_variable_type), intent(out), allocatable:: me(:)

Contains list of hard coded variables

public subroutine atl_varSys_load_user(L, parent, specifics, appender, iError)

Method to load user defined variables from the configuration.

Read more…

Arguments

TypeIntentOptionalAttributesName
type(flu_State) :: L

Lua script to load the variable data from.

integer, intent(in) :: parent

Parent table in the Lua script to read the variable from.

type(c_ptr), intent(out) :: specifics

Data to read for the solver specific variable

procedure(tem_append_solverVar_method), pointer:: appender

Function pointer to use for appending the solver variable.

integer, intent(out) :: iError

Indication whether the attempted reading was successful.

private subroutine atl_generic_fromNodal_getElement(fun, varSys, elempos, time, tree, nElems, nDofs, res)

Routine to obtain a modal representation for a variable, which is only available in nodal form, like space-time functions.

Read more…

Arguments

TypeIntentOptionalAttributesName
class(tem_varSys_op_type), intent(in) :: fun

Description of the method to obtain the variables, here some preset values might be stored, like the space time function to use or the required variables.

type(tem_varSys_type), intent(in) :: varSys

The variable system to obtain the variable from.

integer, intent(in) :: elempos(:)

Position of element in tree%treeID to get the variable for.

type(tem_time_type), intent(in) :: time

Point in time at which to evaluate the variable.

type(treelmesh_type), intent(in) :: tree

global treelm mesh info

integer, intent(in) :: nElems

Number of elements to obtain for this variable (vectorized access).

integer, intent(in) :: nDofs

Number of degrees of freedom within an element.

real(kind=rk), intent(out) :: res(:)

Resulting values for the requested variable.

Linearized array dimension: (nComponents of resulting variable) x (nDegrees of freedom) x (nElems) Access: (iElem-1)fun%nComponentsnDofs + (iDof-1)*fun%nComponents + iComp

private subroutine atl_generic_fromConst_getElement(fun, varSys, elempos, time, tree, nElems, nDofs, res)

Routine to obtain a modal representation for a variable, which is only available in nodal form, like space-time functions.

Arguments

TypeIntentOptionalAttributesName
class(tem_varSys_op_type), intent(in) :: fun

Description of the method to obtain the variables, here some preset values might be stored, like the space time function to use or the required variables.

type(tem_varSys_type), intent(in) :: varSys

The variable system to obtain the variable from.

integer, intent(in) :: elempos(:)

Position of element in tree%treeID to get the variable for.

type(tem_time_type), intent(in) :: time

Point in time at which to evaluate the variable.

type(treelmesh_type), intent(in) :: tree

global treelm mesh info

integer, intent(in) :: nElems

Number of elements to obtain for this variable (vectorized access).

integer, intent(in) :: nDofs

Number of degrees of freedom within an element.

real(kind=rk), intent(out) :: res(:)

Resulting values for the requested variable.

Linearized array dimension: (nComponents of resulting variable) x (nDegrees of freedom) x (nElems) Access: (iElem-1)fun%nComponentsnDofs + (iDof-1)*fun%nComponents + iComp