chrono::fsi::SimParams Struct Reference

Description

Structure with FSI simulation parameters.

#include <ChParams.h>

Public Types

enum  fluidity_model { frictional_plasticity, Inertia_rheology, nonlocal_fluidity }
 

Public Attributes

fluid_dynamics fluid_dynamic_type
 Type of SPH mehtod (WCSPH, IISPH, or I2SPH)
 
char out_name [256]
 Name of the output directory.
 
char demo_dir [2048]
 Demo output directory.
 
int3 gridSize
 dx, dy, dz distances between particle centers.
 
Real3 worldOrigin
 Origin point.
 
Real3 cellSize
 Size of the neighbor particle searching cell.
 
uint numBodies
 Number of FSI bodies.
 
Real3 boxDims
 Dimensions of the domain. How big is the box that the domain is in.
 
Real HSML
 Interaction Radius (or h)
 
Real INVHSML
 1.0/h
 
Real MULT_INITSPACE
 Multiplier to hsml to determine the initial separation of the fluid particles and the fixed separation for the boundary particles. More...
 
Real MULT_INITSPACE_Cables
 Multiplier to hsml in cable elements.
 
Real MULT_INITSPACE_Shells
 Multiplier to hsml in shell elements.
 
int num_neighbors
 Number of neighbor particles.
 
Real epsMinMarkersDis
 epsilon mult for minimum distance between markers (d_min = eps * HSML)
 
int NUM_BOUNDARY_LAYERS
 Number of particles layers that will be used in the boundary. Default value = 3.
 
Real toleranceZone
 Helps determine the particles that are in the domain but are outside the boundaries, so they are not considered fluid particles and are dropped at the beginning of the simulation.
 
int NUM_BCE_LAYERS
 Number of fixed particle layers to rigid/flexible bodies which act as the boundaries. More...
 
Real BASEPRES
 Relative value of pressure applied to the whole domain.
 
Real LARGE_PRES
 Artificial pressure for boundary particles. More...
 
Real3 deltaPress
 Change in Pressure. More...
 
Real3 V_in
 Inlet velocity. This is needed for inlet BC.
 
Real x_in
 Inlet position. This is needed for inlet BC.
 
Real3 gravity
 Gravity. Applied to fluid, rigid and flexible.
 
Real3 bodyForce3
 Constant force applied to the fluid. More...
 
Real rho0
 Density.
 
Real invrho0
 Density's inverse.
 
Real rho_solid
 Solid Density.
 
Real volume0
 Volume.
 
Real markerMass
 marker mass
 
Real mu0
 Viscosity.
 
Real kappa
 surface tension parameter
 
Real v_Max
 Max velocity of fluid used in equation of state. Run simulation once to be able to determine it.
 
Real EPS_XSPH
 Method to modify particle velocity.
 
Real beta_shifting
 this is the beta coefficient in the shifting vector formula. See
 
Real Vis_Dam
 Viscous damping force.
 
Real dT
 Time step. More...
 
Real dT_Flex
 Setpsize for the flexible bodies dynamics.
 
Real Co_number
 Constant in CFL condition.
 
Real dT_Max
 Maximum setpsize.
 
Real out_fps
 Output frames per second.
 
Real tFinal
 Total simulation time.
 
Real timePause
 Time that we let pass before applying body forces. More...
 
Real timePauseRigidFlex
 Time before letting rigid/flex move. More...
 
Real kdT
 Implicit integration parameter. Not very important.
 
Real gammaBB
 Equation of state parameter.
 
Real3 cMin
 Lower leftmost part of the space shown in a simulation frame. More...
 
Real3 cMax
 Upper right most part of the space shown in a simulation frame. More...
 
Real3 cMinInit
 Minimum point of the fluid domain.
 
Real3 cMaxInit
 Maximum point of the fluid domain.
 
Real3 straightChannelBoundaryMin
 Origin of the coordinate system (Point (0,0,0)). More...
 
Real3 straightChannelBoundaryMax
 Upper right most part of the system, this is Point(2 * mm, 1 * mm, 3 * mm) where mm is a constant defined at the beginning of main.cpp.This is also the rightmost particle of the 3rd layer of the top boundaryparticles.Everything above this point is not considered part of thesystem.
 
Real binSize0
 Suggests the length of the bin each particle occupies. More...
 
Real3 rigidRadius
 Radius of rigid bodies.
 
int densityReinit
 Reinitialize density after densityReinit steps. More...
 
int contactBoundary
 0: straight channel, 1: serpentine
 
BceVersion bceType
 Type of boundary conditions, ADAMI or mORIGINAL.
 
bool Conservative_Form
 Whether conservative or consistent discretization should be used.
 
int gradient_type
 Type of the gradient operator.
 
int laplacian_type
 Type of the laplacian operator.
 
bool USE_NonIncrementalProjection
 Used in the I2SPH implementation.
 
bool DensityBaseProjetion
 Set true to use density based projetion scheme in ISPH solver.
 
bool USE_LinearSolver
 If a linear solver should be used to solve Ax=b, otherwise basics methods such as Jacobi-SOR are used.
 
bool Pressure_Constraint
 Whether the singularity of the pressure equation should be fixed.
 
ChFsiLinearSolver::SolverType LinearSolver
 Type of the linear solver.
 
Real Alpha
 Poisson Pressure Equation source term constant. More...
 
Real Beta
 moving exponential weighted average coefficient
 
Real LinearSolver_Abs_Tol
 Poisson Pressure Equation residual.
 
Real LinearSolver_Rel_Tol
 Poisson Pressure Equation Absolute residual.
 
int LinearSolver_Max_Iter
 Linear Solver maximum number of iteration.
 
bool Verbose_monitoring
 Poisson Pressure Equation Absolute residual.
 
Real Max_Pressure
 Max Pressure in the pressure solver.
 
PPE_SolutionType PPE_Solution_type
 MATRIX_FREE, FORM_SPARSE_MATRIX see Rakhsha et al. More...
 
Real PPE_relaxation
 PPE_relaxation.
 
bool ClampPressure
 Clamp pressure to 0 if negative, based on the IISPH paper by Ihmsen et al. (2013)
 
Real IncompressibilityFactor
 Incompressibility factor, default = 1.
 
Real Cs
 Speed of sound.
 
bool Adaptive_time_stepping
 Works only with IISPH for now, use Co_number can be used as the criteria for adaptivity. More...
 
bool Apply_BC_U
 This option lets you apply a velocity BC on the BCE markers.
 
Real L_Characteristic
 Some length characteristic for Re number computation.
 
bool non_newtonian
 Set true to model non-newtonian fluid.
 
bool granular_material
 Set true to model granular material dynamcis using I2SPH.
 
rheology rheology_model
 Model of the rheology (Inertia_rheology or nonlocal_fluidity)
 
Real ave_diam
 average particle diameter
 
Real cohesion
 c in the stress model sigma=(mu*p+c)/|D|
 
friction_law mu_of_I
 Constant I in granular material dyanmcis.
 
Real mu_max
 maximum viscosity
 
Real mu_fric_s
 friction mu_s
 
Real mu_fric_2
 mu_2 constant in mu=mu(I)
 
Real mu_I0
 Reference Inertia number.
 
Real mu_I_b
 b constant in mu=mu(I)=mu_s+b*I
 
Real Shear_Mod
 G.
 
Real HB_sr0
 Herschel–Bulkley consistency index.
 
Real HB_k
 Herschel–Bulkley consistency index.
 
Real HB_n
 Herschel–Bulkley power.
 
Real HB_tau0
 Herschel–Bulkley yeild stress.
 
bool elastic_SPH
 Handles the WCSPH solver for fluid (0) or granular (1)
 
Real E_young
 Young's modulus.
 
Real G_shear
 Shear modulus.
 
Real K_bulk
 Bulk modulus.
 
Real Nu_poisson
 Poisson’s ratio.
 
Real Ar_stress
 Artifical stress.
 
Real Ar_vis_alpha
 Artifical viscosity coefficient.
 
Real Ar_vis_beta
 Artifical viscosity coefficient.
 
Real Fri_angle
 Frictional angle of granular material.
 
Real Dil_angle
 Dilate angle of granular material.
 
Real Coh_coeff
 Cohesion coefficient.
 
Real Q_FA
 Material constants calculate from frictional angle.
 
Real Q_DA
 Material constants calculate from dilate angle.
 
Real K_FA
 Material constants calculate from frictional angle and cohesion coefficient.
 
Real C_Wi
 Threshold of the integration of the kernel function.
 
Real boxDimX
 Dimension of the space domain - X.
 
Real boxDimY
 Dimension of the space domain - Y.
 
Real boxDimZ
 Dimension of the space domain - Z.
 
Real fluidDimX
 Dimension of the fluid domain - X.
 
Real fluidDimY
 Dimension of the fluid domain - Y.
 
Real fluidDimZ
 Dimension of the fluid domain - Z.
 
Real bodyDimX
 Size of the FSI body - X.
 
Real bodyDimY
 Size of the FSI body - Y.
 
Real bodyDimZ
 Size of the FSI body - Z.
 
Real bodyRad
 Radisu of the FSI body.
 
Real bodyLength
 Length of the FSI body.
 
Real bodyIniPosX
 Initial position of the FSI body - X.
 
Real bodyIniPosY
 Initial position of the FSI body - Y.
 
Real bodyIniPosZ
 Initial position of the FSI body - Z.
 
Real bodyIniVelX
 Initial velocity of the FSI body - X.
 
Real bodyIniVelY
 Initial velocity of the FSI body - Y.
 
Real bodyIniVelZ
 Initial velocity of the FSI body - Z.
 
Real bodyIniAngVel
 Initial angular velocity of the FSI body.
 
Real bodyMass
 Mass of the FSI body.
 
Real bodyDensity
 Density of the FSI body.
 

Member Data Documentation

◆ Adaptive_time_stepping

bool chrono::fsi::SimParams::Adaptive_time_stepping

Works only with IISPH for now, use Co_number can be used as the criteria for adaptivity.

dT_Max is set either from the Co_number, or the time step the is required for outputting data

◆ Alpha

Real chrono::fsi::SimParams::Alpha

Poisson Pressure Equation source term constant.

This is used for controlling the noise in the FS forces, F_new=beta*F_old+(1-beta)*F_t beta=1 highly damped force, beta=0 : no modificaiton to force

◆ binSize0

Real chrono::fsi::SimParams::binSize0

Suggests the length of the bin each particle occupies.

Normally this would be 2*hsml since hsml is the radius of the particle, but when we have periodic boundary condition varies a little from 2 hsml.This may change slightly due to the location of the periodic BC.

◆ bodyForce3

Real3 chrono::fsi::SimParams::bodyForce3

Constant force applied to the fluid.

Flexible and rigid bodies are not affected by this force directly, but instead they are affected indirectly through the fluid.

◆ cMax

Real3 chrono::fsi::SimParams::cMax

Upper right most part of the space shown in a simulation frame.

This point is usually outside thetolerance zone.

◆ cMin

Real3 chrono::fsi::SimParams::cMin

Lower leftmost part of the space shown in a simulation frame.

This point is usually outside the tolerance zone.

◆ deltaPress

Real3 chrono::fsi::SimParams::deltaPress

Change in Pressure.

This is needed for periodic BC. The change in pressure of a particle when it moves from end boundary to beginning.

◆ densityReinit

int chrono::fsi::SimParams::densityReinit

Reinitialize density after densityReinit steps.

Note that doing this more frequently helps in getting more accurate incompressible fluid, but more stable solution is obtained for larger densityReinit

◆ dT

Real chrono::fsi::SimParams::dT

Time step.

Depending on the model this will vary and the only way to determine what time step to use is to run simulations multiple time and find which one is the largest dT that produces a stable simulation.

◆ LARGE_PRES

Real chrono::fsi::SimParams::LARGE_PRES

Artificial pressure for boundary particles.

Make sure fluid particles do not go through the boundaries.Note that if time step is not small enough particles near the boundaries might build up huge pressures and will make the simulation unstable.

◆ MULT_INITSPACE

Real chrono::fsi::SimParams::MULT_INITSPACE

Multiplier to hsml to determine the initial separation of the fluid particles and the fixed separation for the boundary particles.

This means that the separation will always be a multiple of hsml. Default value = 1.0.

◆ NUM_BCE_LAYERS

int chrono::fsi::SimParams::NUM_BCE_LAYERS

Number of fixed particle layers to rigid/flexible bodies which act as the boundaries.

Default value = 2.

◆ PPE_Solution_type

PPE_SolutionType chrono::fsi::SimParams::PPE_Solution_type

MATRIX_FREE, FORM_SPARSE_MATRIX see Rakhsha et al.

2018 paper for details omega relaxation in the pressure equation, something less than 0.5 is necessary for Jacobi solver

◆ straightChannelBoundaryMin

Real3 chrono::fsi::SimParams::straightChannelBoundaryMin

Origin of the coordinate system (Point (0,0,0)).

In this case (straigh channel) this point is where the leftmost particle of the 3rd layer of boundary particles is located. Everything below this point is outside the tolerance zone, this means that everything below is not considered part of the system

◆ timePause

Real chrono::fsi::SimParams::timePause

Time that we let pass before applying body forces.

This is done to allow the particles to stabilize first.Run the fluid only during this time, with dTm = 0.1 * dT

◆ timePauseRigidFlex

Real chrono::fsi::SimParams::timePauseRigidFlex

Time before letting rigid/flex move.

Keep the rigid and flex stationary during this time (timePause + timePauseRigidFlex) until the fluid is fully developed


The documentation for this struct was generated from the following file:
  • /builds/uwsbel/chrono/src/chrono_fsi/physics/ChParams.h