Two Phase MPM Solver - geomechanics/mpm GitHub Wiki
Overall Structure
The structure of the Two-phase MPM solver is similar to the other MPM schemes. The key differences of the input .JSON file are highlighted below. The two-phase solver implemented in the current code version is the single-point two-phase formulation.
Mesh
The mesh object is generally the same as the explanation given in the explicit scheme with major differences described as follow.
{
"mesh": {
"mesh": "mesh.txt",
"check_duplicates": false,
"isoparametric": false,
"entity_sets": "entity_sets.json",
"cell_type": "ED3H8",
"io_type": "Ascii3D",
"node_type": "N3D2P",
"particle_cells": "particles-cells.txt",
"boundary_conditions": {
"velocity_constraints": [
{
"nset_id": 0,
"dir": 2,
"velocity": 0.00
},
{
"nset_id": 0,
"dir": 5,
"velocity": 0.00
},
{
"nset_id": 1,
"dir": 0,
"velocity": 0.00
},
{
"nset_id": 1,
"dir": 3,
"velocity": 0.00
},
{
"nset_id": 2,
"dir": 1,
"velocity": 0.00
},
{
"nset_id": 2,
"dir": 4,
"velocity": 0.00
}
],
"pressure_constraints": [
{
"phase_id": 1,
"nset_id": 3,
"pressure": 0.00
}
]
},
"particles_pore_pressures": {
"dir_h": 0,
"dir_v": 1,
"type": "water_table",
"water_tables": [
{
"h0": 5.0,
"position": 0.0
},
{
"h0": 3.0,
"position": 1.0
}]
}
}
}
"node_type"is either"N2D2P"or"N3D2P"."velocity_constraints"follow the convention described in the explicit scheme with direction"dir"should be set as(0|1|2)for the solid phase constraints and(3|4|5)for the fluid phase constraints in 3D. Meanwhile, in 2D, they should be set as(0|1)for the solid phase and(2|3)for the fluid phase."pressure_constraints"follow the convention described in the Navier-Stokes scheme with"phase_id"to be set at 1.
Particles Pore Pressures
Optional particle pore pressure initializer can be input via a text file (in a similar format to particle stresses). The first row of the text file provides the total number of material points. Each subsequent row provides the particle id and corresponding value of pressure.
Alternatively, if the pore pressure in the entire simulation domain is isotropic, the "isotropic" option can be specified following this format:
{
"particles_pore_pressures": {
"type": "isotropic",
"values": 10000000
}
}
where the "values" are the pore pressure with positive magnitude indicating compression.
One can also define a hydrostatic pressure evolving with depth by using the option "water_table". Here, one needs to specify the following information to help the code to define the initial water table:
"dir_v"is the vertical direction (Gravity direction) of the water table."dir_h"is the horizontal direction of the water table."water_tables"is the list of reference points indicated by"position"in the horizontal axis defined in"dir_h"-direction and"h0"is the initial water height in the direction of"dir_v". Here, the datum is always measured at the zero coordinate of the"dir_v"-axis.- It is important to note that to use the
"water_table"option, one needs to define"gravity"in the"external_loading_conditions".
Particles
The particle type should be set as "particle_type": "P2D2PHASE" or "P3D2PHASE". Meanwhile, the "material_id" should be set as a list with the first entry indicating the solid-phase material and the second entry indicating the fluid-phase material. Below is the example of the "particles" object:
"particles": [
{
"generator": {
"check_duplicates": false,
"location": "particles.txt",
"io_type": "Ascii3D",
"particle_type": "P3D2PHASE",
"material_id": [
0,
1
],
"pset_id": 0,
"type": "file"
}
}
],
Materials
The "materials" object should include at least two materials that define both the solid and fluid phases. The parameters defined in the solid-phase material are assumed to follow the effective-stress analysis which solely characterizes the mechanical behavior of the solid skeleton of a porous media.
"materials": [
{
"id": 0,
"name": "solid",
"type": "LinearElastic3D",
"density": 2600.0,
"poisson_ratio": 0.3,
"youngs_modulus": 9.80E+7,
"porosity": 0.3,
"intrinsic_permeability": false,
"k_x": 0.001,
"k_y": 0.001,
"k_z": 0.001
},
{
"id": 1,
"name": "water",
"type": "Newtonian3D",
"density": 1000,
"bulk_modulus": 2.2E9,
"dynamic_viscosity": 8.9E-4
}
],
The solid-phase material should include additional information, such as:
"porosity"which indicates the initial material averaged porosity."k_x","k_y", and"k_z"which indicates the permeability in different axes."intrinsic_permeability"is a boolean that defines the unit of"k_x","k_y", and"k_z", which can be defined as intrinsic permeability (m^2) (true) or soil-mechanics permeability (m/s) (false).
Analysis
The two-phase MPM "analysis" object is described as follow:
{
"analysis": {
"type": "MPMSemiImplicitTwoPhase3D",
"scheme_settings": {
"beta": 1
},
"velocity_update": false,
"pressure_smoothing": false,
"pore_pressure_smoothing": false,
"linear_solver": {
"assembler_type": "EigenSemiImplicitTwoPhase3D",
"solver_settings": [
{
"dof": "pressure",
"solver_type": "KrylovPETSC",
"sub_solver_type": "cg",
"preconditioner_type": "jacobi",
"max_iter": 100,
"tolerance": 1E-5,
"verbosity": 0
},
{
"dof": "acceleration",
"solver_type": "KrylovPETSC",
"sub_solver_type": "bicgstab",
"preconditioner_type": "none",
"max_iter": 100,
"tolerance": 1E-5,
"verbosity": 0
}
]
},
"free_surface_detection": {
"type": "density",
"volume_tolerance": 0.25
},
"uuid": "semi_implicit_3D",
"dt": 0.0001,
"nsteps": 2001
}
}
Type
Supported types:
| Analysis | Description |
|---|---|
| MPMExplicitTwoPhase2D | explicit 2D Two-phase Single-layer MPM |
| MPMExplicitTwoPhase3D | explicit 3D Two-phase Single-layer MPM |
| MPMSemiImplicitTwoPhase2D | semi-implicit 2D Two-phase Single-layer MPM |
| MPMSemiImplicitTwoPhase3D | semi-implicit 3D Two-phase Single-layer MPM |
If the semi-implicit scheme is picked
Scheme settings
Users can also specify the "scheme_settings" option to change the default parameter considered in the simulation.
"beta"is the parameter that defines whether a full (0) or incremental (1) projection scheme is used.
Linear Solver
The "linear_solver" option is used to control the type of assembler and linear solver needed to solve the momentum equation implicitly.
"assembler_type"should be"EigenSemiImplicitTwoPhase2D"or"EigenSemiImplicitTwoPhase3D"for 2D or 3D simulation."solver_settings"contains all necessary settings to solve the linear equationAx=b. Two solvers should be defined to solve the predicted acceleration and pressure Poisson equations."dof"is the name of degrees of freedom. Here,"pressure"and"acceleration"should be defined for the semi-implicit scheme."solver_type"is the name of the linear solver library used. Options available are:
| Analysis | Description |
|---|---|
| DirectEigen | Direct solver of Eigen Library |
| IterativeEigen | Iterative solver of Eigen Library |
| KrylovPETSC | Iterative solver of PETSc Library |
"sub_solver_type"is the name of the linear solver within the selected library."preconditioner_type"is the name of used preconditioner."max_iter"is the maximum number of iterations for the iterative linear solver."tolerance"is the relative residual tolerance."abs_tolerance"is the absolute tolerance."div_tolerance"is the divergence tolerance."verbosity"is to control how much information is output in the console. The default value is 0.
Free-surface detection
The "free_surface_detection" option is used to detect nodes, particles, and cells nearby the material's free surface where zero pressure BC should be imposed.
"type"can be"density"or"geometry". Only"density"is available for parallel solvers with MPI."volume_tolerance"defines volume tolerance where all nodes should be considered as free-surface. The default value is 0.25.