Two Phase MPM Solver

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.

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
				}
			]
		}
	}

}

• "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 and (3|4|5) for the fluid phase 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

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) or soil-mechanics permeability (m/s).

Analysis

The two-phase MPM "analysis" object are 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

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 equation Ax=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.