1. Jesse Capecelatro Avatar Jesse Capecelatro
  2. jCODE
  3. jCODE-public

Wiki

Clone wiki

jCODE-public / structure_jCODE

View History

Structure of the code

The CFD code contains the following directories:

  • bin : all binary (executable) files
  • examples : several examples and test cases
  • lib : all libraries generated during compilation
  • obj : the object files generated during compilation
  • mod : the Fortran modules generated during compilation
  • src : the source code itself
    • Makefiles : example Makefile.in
    • boundary : boundary conditions
    • combustion : various combustion-related modules
    • config : modules used throughout the code including solver options, simulations flags, and geometry parameters
    • controller : adjoint-based control files
    • core : the main code logic
    • functional : various cost functionals used by the adjoint solver
    • grid : grid information (coordinates, metrics, grid-related functions)
    • ibm : immersed boundary method
    • library : useful stand-alone modules
    • monitor : modules to monitor the simulation timestep (timing, stats, etc.)
    • operator : stencil operators (derivatives, dissipation, filtering, etc.)
    • particle : modules to handle particles
    • postprocess: modules to dump visualization files and statistics
    • sources : various source terms
    • state : modules to handle the state (conserved variables, dependent / transport variables, inviscid / viscous fluxes, etc.)
    • utils : useful executables

The input deck

Assigning the number of processors

use manual domain decomposition : .true.
processor decomposition : 2 2 1

If use manual domain decomposition is set to .false., the code will decompose the processors automatically and processor decomposition does not need to be specified. Otherwise, the user will need to specify the number of processors in each direction. Note that the above example processor decomposition : 2 2 1 is specifying a total of 4 processors (2 in x and 2 in y). This must be consistent with the total number of processors used to execute the code, i.e., mpirun -np 4.

Initialization parameters

simulation name : shear layer

! Grid size
nx : 768
ny : 259
Lx : 31
Ly : 20

! Stretch the grid in y
grid stretching : .true.

! Shear layer
include solenoidal perturbations : .true.
shear layer thickness            : 0.227
shear layer velocity difference  : 0.8327759
shear layer rms fraction         : 0.075
shear layer upper temperature    : 12.79
shear layer lower temperature    : 4.26
shear layer pressure             : 2.857142857142857

! Mixture
initial mixture fraction       : 1
initial fuel mass fraction     : 0.06666
initial oxidizer mass fraction : 0.233

! Initialization files
grid file to write     : grid.xyz
solution file to write : data.ic.q

The initialization parameters should only be called during the execution of init_flow. In the above shear layer example, simulation name specifies the configuration used by init_flow. The number of grid points are defined by nx ny and nz, and the domain length is defined by Lx Ly and Lz. The next several options are unique to the initialization of the shear layer. Finally, grid file to write and solution file to write set the name of the grid file and initial conditions, respectively.

Defining patches & boundary conditions

patches : 
! Name                 Type                  normDir iMin  iMax  jMin  jMax  kMin  kMax
! -------------------------------------------------------------------------------------
  bottom               SAT_FAR_FIELD          2       1    -1     1     1     1    -1
  bottomSponge         SPONGE                 2       1    -1     1     10    1    -1
  top                  SAT_FAR_FIELD         -2       1    -1    -1    -1     1    -1
  topSponge            SPONGE                -2       1    -1     250  -1     1    -1

Grid patches are used to solve specific equations on subsets of the grid. They are typically used to define boundary conditions and target/control regions for adjoint runs.

The first column is used to define a unique name to the patch. Note, if two patches have the same name the code will quit with an error message.

The second column is the patch type. In the above example, 4 patches are defined (two of type SAT_FAR_FIELD and two sponges). Available patch types are:

  • SPONGE: defined as a volume in 3D and an area in 2D, the sponge dampens the solution to a known target solution near the boundary.
  • SAT_FAR_FIELD: inflow/outflow far-field boundary condition. Must be defined as a surface in 3D or a line in 2D.
  • SAT_SLIP_WALL: no penetration slip-wall boundary condition. Must be defined as a surface in 3D or a line in 2D.
  • SAT_ISOTHERMAL_WALL: no-slip isothermal wall boundary condition. Must be defined as a surface in 3D or a line in 2D.
  • SAT_ADIABATIC_WALL: no-slip adiabatic wall boundary condition. Must be defined as a surface in 3D or a line in 2D.
  • ACTUATOR: the control region used during an adjoint run. Must be defined as a volume in 3D and an area in 2D.
  • COST_TARGET: the target region used to compute the cost functional. Must be defined as a volume in 3D and an area in 2D.
  • VISUALIZATION: used to define the region for the visualization files. If not specified, the post processing routine will output data of the entire grid size.

The third column is used to specify the normal direction. Boundary conditions are defined using a normal pointing in towards the flow, i.e., normDir=2 is used for the bottom wall, and normDir=-2 is used for the top wall. For patches without directionality (e.g., ACTUATOR, COST_TARGET, VISUALIZATION) should set normDir=0.

imin, imax, jmin, jmax, kmin, kmax are used to define the extents of the grid patch. A negative value tells the code to start counting backwards from the last index, i.e., a value of -1 = nx, -2 = nx-1, etc.

! Periodicity
periodicity type in x : PLANE
periodic length in x  : 31
}}

If no patch is defined for a given direction, the flow solution will be solved at that boundary (via a one-sided finite difference) unless periodic boundary conditions are specified. in Cartesian coordinates, periodicity is defined by ##periodicity type in x : PLANE##. A periodic length must also be defined in that direction, e.g., ##periodic length in x  : 31##. In this example, there are 768 points in x, and the domain extends from 0 at grid point 1 to 31-dx (dx=Lx/nx) at grid point 768. The code assumes x=31 is located at grid point 1. 

===Running parameters for a forward (predictive) run

{{{
#!Text      
! Restart files
save interval : 1000
grid file     : grid.xyz
solution file : data.ic.q

Specifies the filenames to restart a simulation.

! Target state
use target state      : .true.
target state file     : data.ic.q

Specifies the filename for the target state (needed for the SAT_FAR_FIELD and SPONGE patches). This is typically the same as the initial condition.

! Gravity
include gravity : .true.
gravity norm    : 0 -1 0
Froude number   : 29358

By default gravity is set to .false. To use gravity you must specify gravity norm (the direction) and Froude number (the strength).

! Thermo-fluid properties
Reynolds number : 400
Prandtl number  : 0.7
include viscous terms : .true.
viscosity power law exponent : 0.666
bulk viscosity ratio         : 0.6
ratio of specific heats      : 1.4

Various properties of the flow.

! Chemistry
equation of state    : IDEAL GAS MIXTURE
combustion model     : ONE STEP
number of species    : 2
Schmidt number 1     : 0.7
Schmidt number 2     : 0.7
species 1            : H2
species 2            : O2
inert species        : N2
reference species    : AIR
stoichiometric ratio : 8
heat release         : 0.86
Zel Dovich           : 6
Damkohler number     : 500000

Chemistry parameters. In this case the equation of state is set to IDEAL GAS MIXTURE, which tells the code we have a multi-component flow. 2 species are used, each with a Schmidt number of 0.7. The species names are used to define the molecular weights for each.

! Particle parameters
include particles          : .true.
particle density           : 833.0
two way coupling           : .true.
particle collisions        : .false.
collision time             : 15.0
coefficient of restitution : 0.85
drag model                 : SCHILLER_NAUMANN

Particle properties for particle-laden flow simulations.

! Discretization scheme
default discretization scheme : SBP 3-6
curvilinear domain            : .false.

Discretization properties.

! Time stepping options
use constant CFL mode : .true.
time step size        : 1.0
cfl                   : 0.5
number of timesteps   : 100000
report interval       : 1

Time stepping options.

! Artificial dissipation
add dissipation       : .true.
dissipation amount    : 0.005
composite dissipation : .false.

Unfortunately, repeated first derivative schemes require some type of artificial dissipation to dampen high-wave number components.

! Solution limits
enable solution limits : .true.
minimum density        : 0.1
minimum temperature    : 0.05
maximum density        : 2
maximum temperature    : 4

The simulation will end if density or temperature exceed these values.

! Forcing
force momentum : y

Force the mass flow rate to be constant in direction y.

! Output
output type : ensight
output frequency : 49.99999

Output visualization files in EnSight GOLD format at a given frequency.

Running parameters for an inverse (adjoint) run

! Adjoint optimization flags
disable adjoint solver        : .false.
use continuous adjoint        : .false.
baseline prediction available : .false.
cost functional type          : TEMPERATURE
flame temperature use time ramp : .true.
flame temperature ramp peak : 1
flame temperature ramp width : .1

Specify whether or not the adjoint solver should be used, and if so whether the continuous adjoint should be solved (if set to .false. the discrete-exact adjoint will be used) and if a baseline prediction is available. In this case the cost functional type is TEMPERATURE. Refer to functional.f90 for other available cost functional types.

# Control forcing
controller type              : IGNITION_ACTUATOR
number of control parameters : 2
sensitivity parameter 1      : ENERGY
sensitivity parameter 2      : DURATION

The controller specifies what to compute the sensitivity with respect to, and updates the parameters based on the resulting gradient.

! Gradient accuracy options
check gradient accuracy           : .true.
number of control iterations      : 48
initial actuation amount          : 1e-1
restart control iteration         : 0
actuation amount geometric growth : 0.5623413251903491
minimize cost functional          : .false.
adjoint iterations                : 1000

To test the accuracy of the adjoint gradient, set check gradient accuracy to .true. In this example, a forward (predictive run) will be simulated, followed by an adjoint run (using 1000 iterations), then 48 forward runs will be simulated, each with a different value of the control parameters. A file will be generated called gradient_error.txt with the resulting gradient error.

Updated