gen_treelm_cube_global Subroutine

private subroutine gen_treelm_cube_global(me, conf, thandle, myPart, nParts, comm)

Generate the header for the simple full cube mesh.

Arguments

Type IntentOptional Attributes Name
type(tem_global_type), intent(out) :: me

Structure to load the mesh to

type(flu_State) :: conf

Directory containing the mesh informations

integer, intent(in) :: thandle

Handle for the table to read the description of the mesh from.

integer, intent(in) :: myPart

Partition to use on the calling process (= MPI Rank in comm)

integer, intent(in) :: nParts

Number of partitions, the mesh is partitioned into (= Number of MPI processes in comm).

integer, intent(in) :: comm

MPI Communicator to use


Calls

proc~~gen_treelm_cube_global~~CallsGraph proc~gen_treelm_cube_global gen_treelm_cube_global proc~aot_get_val~2 aot_get_val proc~gen_treelm_cube_global->proc~aot_get_val~2

Called by

proc~~gen_treelm_cube_global~~CalledByGraph proc~gen_treelm_cube_global gen_treelm_cube_global proc~tem_global_mesh_internal tem_global_mesh_internal proc~tem_global_mesh_internal->proc~gen_treelm_cube_global proc~tem_global_mesh_read tem_global_mesh_read proc~tem_global_mesh_read->proc~tem_global_mesh_internal

Contents


Source Code

  subroutine gen_treelm_cube_global( me, conf, thandle, myPart, nParts, comm )
    ! -------------------------------------------------------------------- !
    !> Structure to load the mesh to
    type(tem_global_type), intent(out) :: me
    !> Directory containing the mesh informations
    type(flu_State) :: conf
    !> Handle for the table to read the description
    !! of the mesh from.
    integer, intent(in) :: thandle
    !> Partition to use on the calling process (= MPI Rank in comm)
    integer, intent(in) :: myPart
    !> Number of partitions, the mesh is partitioned into (= Number of MPI
    !! processes in comm).
    integer, intent(in) :: nParts
    !> MPI Communicator to use
    integer, intent(in) :: comm
    ! -------------------------------------------------------------------- !
    integer :: iError
    integer :: orig_err(3)
    ! -------------------------------------------------------------------- !

    write(logUnit(1),*) 'Creating HEADER for a full cubical mesh ' &
      &                 // 'without boundaries'

    me%nParts = nParts
    me%myPart = myPart
    me%comm = comm

    ! Get the origin of the cube:
    call aot_get_val( L       = conf,                    &
      &               thandle = thandle,                 &
      &               key     = 'origin',                &
      &               val     = me%origin,               &
      &               ErrCode = orig_err,                &
      &               default = [0.0_rk, 0.0_rk, 0.0_rk] )

    ! Get the length of the cube:
    call aot_get_val( L       = conf,                  &
      &               thandle = thandle,               &
      &               val     = me%BoundingCubeLength, &
      &               ErrCode = iError,                &
      &               key     = 'length',              &
      &               default = 1.0_rk                 )

    ! Get the refinement level:
    call aot_get_val( L       = conf,             &
      &               thandle = thandle,          &
      &               val     = me%minlevel,      &
      &               ErrCode = iError,           &
      &               key     = 'refinementLevel' )

    me%label = 'Generic_Cube'
    me%predefined = 'cube'

    if (me%minlevel == 0) then
      me%BoundingCubeLength = me%BoundingCubeLength*2
      me%minlevel = 1
      me%label = 'Generic_Single'
      me%predefined = 'single'
    end if

    me%maxLevel = me%minLevel

    write(me%comment,'(a15,i7,a16,i2,a1)') &
      &  'Generated with ', nParts, ' parts on Level ', me%minlevel, '.'
    me%dirname = './'

    ! No properties in this mesh
    me%nProperties = 0
    if (associated(me%Property)) deallocate(me%property)
    allocate(me%Property(me%nProperties))

  end subroutine gen_treelm_cube_global