FLASH4.5 API

Generated from /asc/asci2/site/flashcode/secure/release_4p5/source/physics/Eos/Eos.F90 with ROBODoc v4.99.8 on Tue Mar 05 16:16:18 2019

TABLE OF CONTENTS


[Functions] source/physics/Eos/Eos

[top][index]

NAME

  Eos

SYNOPSIS

       call Eos(integer(IN) :: mode,
                integer(IN) :: vecLen,
                real(INOUT) :: eosData(vecLen*EOS_NUM),
      optional, real(IN)    :: massFrac(vecLen*NSPECIES),
      optional, logical(IN),target :: mask(EOS_VARS+1:EOS_NUM),
      optional, integer(IN) :: vecBegin,
      optional, integer(IN) :: vecEnd    )

DESCRIPTION

  This routine applies the equation of state to thermodynamic 
  quantities at one or more grid cells.  The number of cells is 
  determined by the argument veclen.  The routine expects data packaged
  for it in the 1d array, eosData.  The data in eosData is organized as:
  1:vecLen points contain the first variable, vecLen+1:2*vecLen points 
  contain the second variable, and so on. The number and order of
  variables in the array is determined by the constants defined in Eos.h.
  
  The routine takes different quantities as givens depending on the
  value of the mode variable: if mode=MODE_DENS_TEMP, density and
  temperature are taken as given, and pressure and internal energy are generated
  as output; if mode=MODE_DENS_EI, density and internal energy are taken as
  givens, and pressure and temperature are generated as output.  If
  mode=MODE_DENS_PRES, density and pressure are taken as givens, and
  internal energy and temperature are generated as output. Note that
  internal energy is EINT_VAR, not ENER_VAR.
  
  In addition to pressure, temperature, and internal energy, which are
  always thermodynamically consistent after this call, other quantities
  such as the various thermodynamic partial derivatives can be
  calculated based on the values in the argument, mask.  mask is a
  logical array with one entry per quantity, with the order determined
  by constants defined in Eos.h (the same as those for the eosData
  argument); .true. means return the quantity, .false. means don't.

ARGUMENTS

  mode :    Selects the mode of operation of the Eos unit.
             The valid values are MODE_DENS_EI, MODE_DENS_PRES and  
             MODE_DENS_TEMP as decribed above.

  vecLen   : number of points (cells) for which the eosData array is sized.
             If vecBegin and vecEnd are not present, this is also the
             number of points (cells) for which EOS computation is to be done.

  eosData  : This array is the data structure through which variable values are 
             passed in and out of the Eos routine. The arrays is sized as 
             EOS_NUM*vecLen. EOS_NUM, and individual input and output
             Eos variables are defined in Eos.h. The array is organizes such that
             the first 1:vecLen entries represent the first Eos variable, vecLen+1:
             2*vecLen represent the second Eos variable and so on. 

  massFrac : Contains the mass fractions of the species included in
             the simulation. The array is sized as NSPECIES*vecLen.

  mask     : Mask is a logical array the size of EOS_DERIVS (number
              of partial derivatives that can be computed, defined in
              Eos.h), where each index represents a specific partial derivative
              that can be calculated by the Eos unit. A .true. value in mask 
              results in the corresponding derivative being calculated and 
              returned. It should preferably be dimensioned as
              mask(EOS_VARS+1:EOS_NUM) in the calling routine 
              to exactly match the arguments declaration in Eos Unit.
             Note that the indexing of mask does not begin at 1, but rather at one past
             the number of variables.

             An implementation that does not need derivative quantities should
             set the mask equal to .false.

  vecBegin : Index of first cell in eosData to handle.
             Can be used to limit operation to a subrange of cells, untested.
             If not present, the default is 1.
  vecEnd   : Index of last cell in eosData to handle.
             Can be used to limit operation to a subrange of cells, untested.
             If not present, the default is vecLen.

EXAMPLE

 --- A single-point at a time example, does not calculate derivatives (based on Cellular Simulation)---

  #include "constants.h"   ! for MODE_DENS_TEMP
  #include "Flash.h"       ! for NSPECIES
  #include "Eos.h"         ! for EOS_VAR order

  real  :: temp_zone, rho_zone, ptot, eint, gamma
  real, dimension(EOS_NUM)  :: eosData
  real, dimension(SPECIES_BEGIN:SPECIES_END) ::  massFraction  
  integer, dimension(2,MDIM)                 :: blockRange,blockExtent


  massFraction(:) = 1.0e-12        
  massFraction(C12_SPEC) = 1.0

  .... initiale temp_zone, rho_zone

  call Grid_getBlkIndexLimits(blockId,blockRange,blockExtent)
  do k = blockRange(LOW,KAXIS), blockRange(HIGH,KAXIS)
     do j = blockRange(LOW,JAXIS),blockRange(HIGH,JAXIS)
        do i = blockRange(LOW,IAXIS),blockRange(HIGH,IAXIS)

           eosData(EOS_TEMP) = temp_zone
           eosData(EOS_DENS) = rho_zone

           call Eos(MODE_DENS_TEMP,1,eosData,massFraction)

           ptot = eosData(EOS_PRES)
           eint = eosData(EOS_EINT)
           gamma = eosData(EOS_GAMC)
           
           call Grid_putPointData(blockId,CENTER,TEMP_VAR,EXTERIOR,iPosition,temp_zone)
           call Grid_putPointData(blockId,CENTER,DENS_VAR,EXTERIOR,iPosition,rho_zone)
           call Grid_putPointData(blockId,CENTER,PRES_VAR,EXTERIOR,iPosition,ptot)
           call Grid_putPointData(blockId,CENTER,EINT_VAR,EXTERIOR,iPosition,eint)
              if you want ENER_VAR, calculate it from EINT_VAR and kinetic energy
           call Grid_putPointData(blockId,CENTER,GAMC_VAR,EXTERIOR,iPosition,gamma)
           call Grid_putPointData(blockId,CENTER,GAME_VAR,EXTERIOR,iPosition,(ptot/(etot*sim_rhoAmbient) + 1.0))

         enddo  ! end of k loop
     enddo     ! end of j loop
  enddo        ! end of i loop

 ------------------ Row at a time example, with derivates (based on Eos_unitTest) --------

  #include "constants.h"   ! for MODE_DENS_TEMP
  #include "Flash.h"       ! for NSPECIES, EOS_NUM
  #include "Eos.h"         ! for EOS_VAR order
  integer veclen, isize, jsize, ksize, i,j,k, e
  real, dimension(:), allocatable :: eosData
  real, dimension(:), allocatable :: massFrac
  logical, dimension (EOS_VARS+1:EOS_NUM) :: mask
  real, allocatable, dimension(:,:,:,:) :: derivedVariables
  integer,dimension(2,MDIM) :: blkLimits,blkLimitsGC

   ! in the Eos_unitTest, this loops over all blocks.... here is a snippet from inside
     call Grid_getBlkIndexLimits(blockID,blkLimits,blkLimitsGC)

    !  Allocate the necessary arrays for an entire block of data
    isize = (blkLimits(HIGH,IAXIS) - blkLimits(LOW,IAXIS) + 1)
    jsize = (blkLimits(HIGH,JAXIS) - blkLimits(LOW,JAXIS) + 1)
    ksize = (blkLimits(HIGH,KAXIS) - blkLimits(LOW,KAXIS) + 1)
    vecLen=isize
    allocate(derivedVariables(isize,jsize,ksize,EOS_NUM))
    allocate(eosData(vecLen*EOS_NUM))
    allocate(massFrac(vecLen*NSPECIES))
    mask = .true.

    ! indices into the first location for these variables
    pres = (EOS_PRES-1)*vecLen
    dens = (EOS_DENS-1)*vecLen
    temp = (EOS_TEMP-1)*vecLen


    call Grid_getBlkPtr(blockID,solnData)
    do k = blkLimits(LOW,KAXIS),blkLimits(HIGH,KAXIS)
        do j = blkLimits(LOW,JAXIS), blkLimits(HIGH, JAXIS)
           do i = 1,vecLen
              massFrac((i-1)*NSPECIES+1:i*NSPECIES) = &
                   solnData(SPECIES_BEGIN:SPECIES_END,ib+i-1,j,k)
           end do

           eosData(pres+1:pres+vecLen) =  solnData(PRES_VAR,ib:ie,j,k)
           eosData(dens+1:dens+vecLen) =  solnData(DENS_VAR,ib:ie,j,k)
           ! Eos Helmholtz needs a good initial estimate of temperature no matter what the mode
           eosData(temp+1:temp+vecLen) =  solnData(TEMP_VAR,ib:ie,j,k)

           call Eos(MODE_DENS_PRES,vecLen,eosData,massFrac,mask)

           do e=EOS_VARS+1,EOS_NUM
              m = (e-1)*vecLen
              derivedVariables(1:vecLen,j-NGUARD,k-NGUARD,e) =  eosData(m+1:m+vecLen)
           end do
        end do
     end do

NOTES

  NSPECIES is defined in Flash.h.

  EOS_VARS and EOS_NUM  are defined in Eos.h.
  Calling funtions should included Eos.h, in order to get the definitions of
  Eos-specific constants to be able to populate the eosData and mask arrays.
  
  MODE_DENS_TEMP, MODE_DENS_EI, and MODE_DENS_PRES are defined in constants.h.

  All routines calling this routine should include a 
  use Eos_interface 
  statement, preferable with "ONLY" attribute e.g.
  use Eos_interface, ONLY:  Eos

  When using Helmholtz, strange behaviour occurs.  See the notes in subunit
  Helmholtz/Eos.F90 directly.  For Helmholtz, when operating in MODE_DENS_EI,
  the INPUT energy is updated.
  Similarly, when operating in MODE_DENS_PRES, the INPUT pressure is updated. 
  Physicists need to be aware of this.
  This behavior can be turned off by setting the runtime parameter
  eos_forceConstantInput.

  For Gamma and Multigamma routines, the entropy and entropy derivatives 
  calculations have not been confirmed to be correct.  Use with caution.

SEE ALSO

  Eos.h    defines the variables used.
  Eos_wrapped  sets up the data structure.