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