Eulerian Objects

  • Mohamed Houssem Kasbaoui

Eulerian Objects

The derived types eulerian_obj_r and eulerian_obj_i define real and integer Eulerian data on LEAP blocks. Eulerian data represents fields . These include velocity fields , pressure field , and scalar fields .

Staggering

Eulerian objects may have a staggering value 0, 1, 2, or 3, depending on where the data is defined in a computational cell. Cell-centered quantities have staggering =0. Face-centered quantities have a staggering =1, 2, or 3, if the data is shifted from the cell center to the -, -, or -face. See the figure below for illustration:

use leapKinds
use leapBlock
use leapParallel
use leapEulerian
type(block_obj)      :: block
type(parallel_obj)   :: parallel
type(eulerian_obj_r) :: phi

! Perform initializations
call parallel%Initialize()
call block%Initialize(ngc=1,parallel)
! ... (setup the block)

! Define phi as a cell-centered (stag=0) eulerian object
! on our given block - This will allocate the internal
! arrays
call phi%Initialize(name="Phi", block, parallel, stag=0)

associate (xm => block%xm, ym => block%ym, zm => block%zm)
  do k=block%lo(3),block%hi(3)
    do j=block%lo(2),block%hi(2)
      do i=block%lo(1),block%hi(1)
        phi%cell(i,j,k) = cos(xm(i))*cos(ym(j))*cos(zm(k))
      end do
    end do
  end do
end associate

call phi%UpdateGhostCells()

Overloaded operators

To make working with eulerian objects a little easier, LEAP has a limited number of overloaded operators. For example, addition ofEulerian objects, which goes as

do k=block%lo(3),block%hi(3)
  do j=block%lo(2),block%hi(2)
    do i=block%lo(1),block%hi(1)
      C%cell(i,j,k) = A%cell(i,j,k) + B%cell(i,j,k)
    end do
  end do
end do

call C%UpdateGhostCells()

can be replaced with

C = A+B

call C%UpdateGhostCells()

The same applies for subtraction. For multiplication by scalar, the scalar must appear on the right side:

! This overloaded scalar multiplication
B = A*0.5_wp
call B%UpdateGhostCells()

! is equivalent to this
do k=block%lo(3),block%hi(3)
  do j=block%lo(2),block%hi(2)
    do i=block%lo(1),block%hi(1)
      B%cell(i,j,k) = A%cell(i,j,k) * 0.5_wp
    end do
  end do
end do

call B%UpdateGhostCells()

Ghost cells treatment

There are two main operations involving ghostcells: UpdateGhostCells and AddUpGhostCells. The first operation will ensure that ghostcells contain updated information from neighboring blocks. The following schematic illustrates this procedure in a simplified example with 1 layer of ghost cells (ngc=1).

Note that the ghost cells left and right to the domain boundaries are not updated, unless if the domain is periodic, then the procedure is similar to what is shown with the exchange occuring between the two blocks at both ends of the domain.

The AddUpGhostCells operation take values in the ghostcells and adds them to neighboring blocks, then finishes up with an UpdateGhostCells. With the domain decomposition above, the sequence of operations goes like this:

Eulerian sets

The eulerian_set derived type is essentially a container to store and facilitate handling multiple eulerian objects. Here is an example where two eulerian objects of type real and integer are added to the same eulerian set:

use leapKinds
use leapBlock
use leapParallel
use leapEulerian
type(block_obj)      :: block
type(parallel_obj)   :: parallel
type(eulerian_set)   :: fields
type(eulerian_obj_r) :: A
type(eulerian_obj_i) :: B

! Perform initializations
call parallel%Initialize()
call block%Initialize(ngc=1,parallel)
! ... (setup the block)

! Initialize the eulerian set
call fields%Initialize(block, parallel)

! Add object A with staggering 1 (x1-face centered)
! this will also initialize and allocate A
call fields%Add(name='A', stag=1, A)

! Add object B with staggering 0 (cell centered)
! this will also initialize and allocate B
call fields%Add(name='B', stag=0, B)

Once done using the eulerian objects contained in an eulerian set, finalizing the latter will also finalize (and deallocate) the eulerian objects.

! ... (do some work)

! Once done, finalize the set to also deallocate all associated data
! (including the objects A and B)
call fields%Finalize()

Writing/Reading Eulerian Objects

All I/O operations are handeled by the eulerian_set derived type. LEAP supports reading/writing with the H5HUT format. Write operations with the SILO format are also supported.

To read from file follow this example:

use leapKinds
use leapBlock
use leapParallel
use leapEulerian
type(block_obj)      :: block
type(parallel_obj)   :: parallel
type(eulerian_set)   :: fields
type(eulerian_obj_r) :: A
type(eulerian_obj_i) :: B
integer              :: iter
real(wp)             :: time

! Perform initializations
call parallel%Initialize()
call block%Initialize(ngc=1,parallel)
! ... (setup the block)

! Initialize the eulerian set
call fields%Initialize(block, parallel)

! Add objects to this set
call fields%Add(name='A', stag=1, A)
call fields%Add(name='B', stag=0, B)

! Specify the name of the file to read
call fields%SetReadFileName('fields.h5')

! Read only 'A', the iteration and time values
call fields%Read(iter, time,[character(len=str8):: 'A'])

! Read 'A' and 'B', the iteration and time values
call fields%Read(iter, time,[character(len=str8):: 'A', 'B'])

! Read all objects previously added to this set, the iteration and time values
call fields%Read(iter, time)

To write data to disk, follow this example:

use leapKinds
use leapBlock
use leapParallel
use leapEulerian
type(block_obj)      :: block
type(parallel_obj)   :: parallel
type(eulerian_set)   :: fields
type(eulerian_obj_r) :: A
type(eulerian_obj_i) :: B
integer              :: iter
real(wp)             :: time

iter = 0
time = 1.0_wp
! Perform initializations
call parallel%Initialize()
call block%Initialize(ngc=1,parallel)
! ... (setup the block)

! Initialize the eulerian set
call fields%Initialize(block, parallel)

! Add objects to this set
call fields%Add(name='A', stag=1, A)
call fields%Add(name='B', stag=0, B)

! Specify the name of the file to read
call fields%SetWriteFileName('fields.h5')

! Change the overwrite attribute from true (default) to false
call fields%SetOverwrite(.false.)

! Write nly 'A', the iteration and time values
call fields%Write(iter, time,[character(len=str8):: 'A'])

! wWrite 'A' and 'B', the iteration and time values
call fields%Write(iter, time,[character(len=str8):: 'A', 'B'])

! Writeall objects previously added to this set, the iteration and time values
call fields%Write(iter, time)