tem_variable_loadMapping_single Subroutine

public subroutine tem_variable_loadMapping_single(expectedName, conf, thandle, varDict, varSys, nComp, solverData_evalElem, ErrCode)

Loads the variable mapping from a table for single expected name. A variable mapping is used to link a user defined variable to a variable expected from, e.g., an equation system. These mappings are stored in varDict, which basically is a dictionary, whereas the key contains the name of the expected variable and the value contains the name of the user defined variable in the variable table.

Arguments

Type IntentOptional Attributes Name
character(len=*), intent(in) :: expectedName

Expected variable name from config file

type(flu_State), intent(in) :: conf

lua config file

integer, intent(in) :: thandle

table handle

type(grw_stringkeyvaluepairarray_type), intent(inout) :: varDict

The dictionary that contains the mapping between expected variables and the actual variables defined by the user.

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

Variable system to append anonymous variables to.

integer, intent(in), optional :: nComp

Number of components we expect for this variable.

type(tem_varSys_solverData_evalElem_type), intent(in), optional :: solverData_evalElem

Routine to convert point information into an element representation.

integer, intent(out), optional :: ErrCode

Error code


Source Code

  subroutine tem_variable_loadMapping_single( expectedName, conf, thandle, &
    &                                         varDict, varSys, ncomp,      &
    &                                         solverData_evalElem, ErrCode )
    ! --------------------------------------------------------------------------
    !> Expected variable name from config file
    character(len=*), intent(in) :: expectedName
    !> lua config file
    type(flu_state), intent(in) :: conf
    !> table handle
    integer, intent(in) :: thandle
    !> The dictionary that contains the mapping between expected variables
    !! and the actual variables defined by the user.
    type(grw_stringKeyValuePairArray_type), intent(inout) :: varDict
    !> Variable system to append anonymous variables to.
    type(tem_varSys_type), intent(inout) :: varSys
    !> Number of components we expect for this variable.
    integer, optional, intent(in) :: nComp

    !> Routine to convert point information into an element representation.
    type(tem_varSys_solverData_evalElem_type), &
      &  optional, intent(in) :: solverData_evalElem

    !> Error code
    integer, optional, intent(out) :: ErrCode
    ! --------------------------------------------------------------------------
    character(len=labelLen) :: varname
    character(len=labelLen) :: anony_name
    logical :: varExist
    integer :: errorCode
    type(tem_stringKeyValuePair_type) :: kvp
    type(tem_spacetime_fun_type), pointer :: stfun(:)
    real(kind=rk) :: numtest
    logical :: check_stfun
    ! --------------------------------------------------------------------------

    varExist = aot_exists( L       = conf,        &
      &                    thandle = thandle,     &
      &                    key     = expectedName )

    check_stfun = varExist

    ! if variable exist load space time function for that variable
    if (varExist) then
      write(logUnit(3),* ) 'Expected variable '  &
        & // trim(expectedName) // ' is defined'
      ! First check, whether this variable definition is a number.
      ! (They also satisfy strings).
      ! We do not accept numbers as variable names, instead this
      ! will be read as constant stfun.
      call aot_get_val( L       = conf,         &
        &               thandle = thandle,      &
        &               key     = expectedName, &
        &               val     = numtest,      &
        &               ErrCode = errorCode     )
      if (btest(errorCode, aoterr_WrongType)) then
        ! Not a number, try to interpret it as a string.
        call aot_get_val( L       = conf,         &
          &               thandle = thandle,      &
          &               key     = expectedName, &
          &               val     = varname,      &
          &               ErrCode = errorCode     )
        if (errorCode == 0) then
          ! Found a string, use it to refer to a variable.
          write(logUnit(1),*) 'Corresponding variable for ' &
            & // trim(expectedName) // ' found: ' // trim(varname)
          check_stfun = .false.
          kvp%key = expectedName
          kvp%value = varname
          call append(me = varDict, val = kvp)
        end if
      end if
      ! If the variable is not defined yet, try to get it as an anonymous
      ! variable in a space-time function.
      if (check_stfun) then
        write(logUnit(3),*) 'Did not find a variable name for ' &
          &                 // trim(expectedName)
        write(logUnit(3),*) 'Trying to load it as an anonymous space' &
          &                 // ' time function'
        allocate(stfun(1))
        call tem_load_spacetime(me      = stfun(1),           &
          &                     conf    = conf,               &
          &                     parent  = thandle,            &
          &                     key     = trim(expectedName), &
          &                     nComp   = nComp,              &
          &                     errCode = errorCode           )
        ! use global mesh for anonymous variable
        stfun(1)%subTree%useGlobalMesh = .true.
        if (errorCode == 0) then
          kvp%key = expectedName
          anony_name = tem_spacetime_hash_id(me=stfun(1), conf=conf)
          kvp%value = 'anon_' // trim(anony_name)
          call tem_varSys_append_stFun(                   &
            & varSys              = varSys,               &
            & stfun               = stfun,                &
            & varname             = kvp%value,            &
            & nComp               = nComp,                &
            & evalType            = 'firstonly_asglobal', &
            & solverData_evalElem = solverData_evalElem   )
          call append(me = varDict, val = kvp)
          write(logUnit(1),*) 'Found a spacetime function for ' &
            &                 // trim(expectedName)
        else
          write(logUnit(1),* ) 'Expected variable ' // trim(expectedName) &
            & // ' is NOT defined'
          deallocate(stfun)
        end if
      end if
    else
      write(logUnit(1),* ) 'Expected variable ' // trim(expectedName) &
        &                  // ' is NOT defined'
      errorCode = -1
    end if

    if (present(ErrCode)) ErrCode = errorCode

  end subroutine tem_variable_loadMapping_single