Material point method

An MPM simulation consists of the following stages:

(Prior to the time integration phase)

1. Initialization of grid and material points.
   1. A geometry is discretized into a collection of material points, each with its own material properties and initial conditions (velocity, stress, temperature, etc.)
   2. The grid, being only used to provide a place for gradient calculations is normally made to cover an area large enough to fill the expected extent of computational domain needed for the simulation.

(During the time integration phase - explicit formulation)

1. Material point quantities are extrapolated to grid nodes.
   1. Material point mass (${\textstyle m_{mp}}$), momenta (${\textstyle {\vec {P_{mp}}}}$), stresses (${\displaystyle {\boldsymbol {\bar {\bar {\sigma }}}}_{mp}}$), and external forces (${\displaystyle {\vec {b}}}$) are extrapolated to the nodes at the corners of the cells within which the material points reside. This is most commonly done using standard linear shape functions (${\textstyle N_{nd-mp}}$), the same used in FEM.
   2. The grid use the material point values to create the masses (${\textstyle M_{node}}$), velocities (${\textstyle {\vec {V_{node}}}}$), internal and external force vectors (${\textstyle {\vec {F_{node}^{\mathsf {internal}}}}}$,${\textstyle {\vec {F_{node}^{\mathsf {external}}}}}$) for the nodes:
      $${\displaystyle M_{node}=\sum _{mp}m_{mp}~~N_{mp-nd}}$$

      
      $${\displaystyle {\vec {V_{node}}}={1 \over M_{node}}~~\sum _{mp}{\vec {P_{mp}}}~~N_{mp-nd}}$$

      
      ${\displaystyle {\vec {F_{node}^{internal}}}=\sum _{mp}~~{\bar {\bar {\sigma }}}_{mp}~~\nabla N_{mp-nd}}$
      $${\displaystyle {\vec {F_{node}^{\mathsf {external}}}}=\sum _{mp}{\vec {b}}~~N_{mp-nd}}$$

      
2. Equations of motion are solved on the grid.
   1. Newton's 2nd Law is solved to obtain the nodal acceleration (${\displaystyle {\vec {A_{node}}}}$)
      $${\displaystyle {\vec {A_{node}}}={\vec {F_{node}^{external}+{\vec {F_{node}^{internal}}} \over M_{node}}}}$$

      
   2. New nodal velocities are found (${\displaystyle {\tilde {\vec {V_{node}}}}}$).
      $${\displaystyle {\tilde {\vec {V_{node}}}}={\vec {V_{node}}}+{\vec {A_{node}}}\mathrm {d} t}$$

      
3. Derivative terms are extrapolated back to material points
   1. Material point acceleration (${\displaystyle {\vec {a_{mp}}}}$), deformation gradient (${\displaystyle {\mathcal {\bar {\bar {F_{mp}}}}}}$) (or strain rate (${\displaystyle {\bar {\bar {{\dot {\varepsilon }}_{mp}}}}}$) depending on the strain theory used) is extrapolated from the surrounding nodes using similar shape functions to before (${\displaystyle N_{nd-mp}}$).
      $${\displaystyle {\vec {a_{mp}}}=\sum _{nd}{\vec {A_{node}}}~~N_{nd-mp}}$$

      
      ${\displaystyle {\bar {\bar {{\dot {\varepsilon }}_{mp}}}}=\sum _{nd}~{1 \over 2}~~[{\vec {V_{node}}}\nabla N_{nd-mp}+(V_{node}\nabla N_{nd-mp})^{T}]}$
   2. Variables on the material points: positions, velocities, strains, stresses etc. are then updated with these rates depending on integration scheme of choice and a suitable constitutive model.
4. Resetting of grid.

   Now that the material points are fully updated at the next time step, the grid is reset to allow for the next time step to begin.

The PIC was originally conceived to solve problems in fluid dynamics, and developed by Harlow at Los Alamos National Laboratory in 1957.[1] One of the first PIC codes was the Fluid-Implicit Particle (FLIP) program, which was created by Brackbill in 1986[2] and has been constantly in development ever since. Until the 1990s, the PIC method was used principally in fluid dynamics.

Motivated by the need for better simulating penetration problems in solid dynamics, Sulsky, Chen and Schreyer started in 1993 to reformulate the PIC and develop the MPM, with funding from Sandia National Laboratories.[3] The original MPM was then further extended by Bardenhagen et al.. to include frictional contact,[4] which enabled the simulation of granular flow,[5] and by Nairn to include explicit cracks[6] and crack propagation (known as CRAMP).

Recently, an MPM implementation based on a micro-polar Cosserat continuum[7] has been used to simulate high-shear granular flow, such as silo discharge. MPM's uses were further extended into Geotechnical engineering with the recent development of a quasi-static, implicit MPM solver which provides numerically stable analyses of large-deformation problems in Soil mechanics.[8]

Annual workshops on the use of MPM are held at various locations in the United States. The Fifth MPM Workshop was held at Oregon State University, in Corvallis, OR, on April 2 and 3, 2009.

The uses of the PIC or MPM method can be divided into two broad categories: firstly, there are many applications involving fluid dynamics, plasma physics, magnetohydrodynamics, and multiphase applications. The second category of applications comprises problems in solid mechanics.

MPM is more expensive in terms of storage than other methods, as MPM makes use of mesh as well as particle data. MPM is more computationally expensive than FEM, as the grid must be reset at the end of each MPM calculation step and reinitialised at the beginning of the following step. Spurious oscillation may occur as particles cross the boundaries of the mesh in MPM, although this effect can be minimized by using generalized interpolation methods (GIMP). In MPM as in FEM, the size and orientation of the mesh can impact the results of a calculation: for example, in MPM, strain localisation is known to be particularly sensitive to mesh refinement. One stability problem in MPM that does not occur in FEM is the cell-crossing errors and null-space errors[36] because the number of integration points (material points) does not remain constant in a cell.

• Center for Simulation of Accidental Fires and Explosions – MPM code available
• NairnMPM – open source
• MPM3D - open source (MPM3D-F90) and free trial version (MPM3D)
• Taichi - Physically Based Computer Graphics Library – open source MPM code available
• Anura3D open source – software for geotechnical problems and soil-water-structure interactions by Anura3D MPM Research Community