Standard Material Point Method Research Articles

Mapped material point method for large deformation problems with sharp gradients and its application to soil‐structure interactions

AbstractThe material point method (MPM) is often applied to large deformation problems that involve sharp gradients in the solution field. Representative examples in geomechanics are interactions between soils and various “structures” such as foundations, penetrometers, and machines, where the displacement fields exhibit sharp gradients around the soil‐structure interfaces. Such sharp gradients should be captured properly in the MPM discretization to ensure that the numerical solution is sufficiently accurate. In the MPM literature, several types of locally refined discretizations have been developed and used for this purpose. However, these local refinement schemes are not only quite complicated but also restricted to certain types of basis functions or update schemes. In this work, we propose a new MPM formulation, called the mapped MPM, that can efficiently capture sharp gradients with a uniform background grid compatible with every standard MPM basis function and scheme. The mapped MPM is built on the method of auxiliary mapping that reparameterizes the given problem in a different domain whereby sharp gradients become much smoother. Because the reparameterized problem is free of undesirably sharp gradients, it can be well solved with the standard MPM ingredients including a uniform background grid. We verify and demonstrate the mapped MPM through several numerical examples, with particular attention to soil‐structure interaction problems.

Coupling transport equations to mechanics in the material point method using an approximate full capacity matrix inverse

Coupling solutions of one or more transport equations to mechanics modeling in the material point method (MPM), extends MPM to a platform for solving multi-physics problems. Unfortunately, the standard MPM approach to transport equations using a lumped transport capacity matrix is prone to large oscillations in particle values especially near interfaces, cracks, or rapid changes in transport value. Switching from lumped capacity matrix methods to an approximate full capacity matrix approach can eliminate those oscillations. Although this switch can cause unwanted numerical diffusion, that problem can be suppressed using periodic or blending options. Several examples illustrate the new methods for heat conduction, solvent transport, poroelasticity, and phase field fracture simulations. In thermal calculations, the new methods improve accuracy for tracking thermodynamics state functions.

A subdivision-stabilized B-spline mixed material point method

Subjected to external loadings, polymeric materials, e.g., biological tissues, hydrogels, and elastomers, may undergo extreme, nearly incompressible, (self-)contact deformations. For numerical modeling employing mesh-based techniques such as the finite element method (FEM), these deformations pose significant challenges due to large distortions in the deformed geometry, accuracy issues stemming from volumetric locking effects, and increased computational cost from complex contact searches. As an alternative to mesh-based methods, the material point method (MPM), a continuum-based particle technique, is gaining attention for its ability to handle extreme distortions and capture no-slip contact without added cost. For nearly incompressible material behaviors, while mixed formulations can address locking effects by treating displacements and pressure as independent fields, they can suffer from numerical instabilities close to the incompressibility limit due to the violation of the inf-sup condition, leading to inaccurate nodal pressure solutions. Here we propose an efficient and stable mixed B-spline material point method with highest achievable regularity for quasi-compressible polymeric materials. Using the two-scale relation of B-splines, we introduce a subdivision-stabilization for the two-field mixed MPM and obtain numerically stable, oscillation-free nodal solutions with equal-order interpolations with optimal regularity. Building on the Eulerian-Lagrangian nature of MPM, a previously-converged solution framework is adopted to mitigate issues related to cell-crossing and numerical fracture artifact present in standard MPM. We assess the stability and accuracy of the developed mixed MPM at large deformations for soft materials through the benchmark Cook’s membrane problem. Additionally, we test the robustness of the proposed MPM by modeling several examples, including the compression and indentation of a circular block into a quasi-compressible substrate and the twisting deformation of a rectangular block. The findings demonstrate the MPM’s capabilities for modeling practical soft material applications.

An efficient solid shell material point method for large deformation of thin structures

AbstractThe standard material point method (MPM) encounters severe numerical difficulties in simulating shell structures. In order to overcome the shortcomings of locking effects, the discretization size of background grid should be small enough, usually smaller than 1/5 of the shell thickness, which however will lead to prohibitive computational cost. A novel solid shell material point method (SSMPM) is proposed to efficiently model the large deformation of thin structures. The SSMPM describes the material domain of shell structures by shell particles with hexahedral particle domains. The locking treatments of solid shell element are then introduced in SSMPM, which results in the correction of strain field throughout the shell thickness. Namely, the assumed natural strain (ANS) method is adopted to eliminate the shear locking and trapezoidal locking, while the enhanced assumed strain (EAS) method is employed to eliminate the thickness locking. With the precise description of bending modes, a single layer of particles and a coarse background grid could be used in shell structure simulations with the SSMPM, which dramatically increases the computational efficiency. A local multi‐mesh contact method is presented to naturally couple SSMPM and MPM for the contact situations of shells with other objects. Several numerical examples, including beam vibration, pinched cylinder with free edges and full hemispherical shell, are performed to verify and validate the SSMPM, which shows that the SSMPM considerably outperforms the standard MPM in these situations. A fluid–structure interaction problem and the penetration of a thin plate are investigated based on the contact method and the results are in good agreement with those in the literature.

Analysis and mitigation of spatial integration errors for the material point method

SummaryThe material point method (MPM) is a robust numerical simulation approach for continuum mechanics problems involving large material deformations coupled to changing surface topographies. These types of problems are present in many different engineering contexts, from understanding the failure processes of earthen slopes to predicting the strengths and failure mechanisms of body armor to modeling the impact forces of waves in fluid tanks. By using a set of persistent material point tracers to follow the motion and deformation of the continuum material while solving the equations of motion on a static simulation grid, the MPM avoids several shortcomings of more traditional numerical approaches including blurring of material surfaces — as in Eulerian finite element or finite volume methods (FEMs or FVMs) — and mesh tangling — as in Lagrangian FEMs. Despite its robustness, MPM is known to develop significant numerical errors: namely, (i) the particle ringing instability and (ii) solution dependent discretization and integration errors. In this work, we present an analysis of local‐in‐time, spatial integration errors in the MPM and several techniques designed to mitigate these errors. Error mitigation approaches previously described in the literature are compared to a new method we propose for problems involving very large material deformations. The proposed method is shown to offer substantial improvement over standard MPM for simulations of fluid‐like materials without requiring significant augmentation of existing MPM frameworks.

Taylor particle-in-cell transfer and kernel correction for material point method

The material point method (MPM) has been extensively used to solve problems involving large displacements and deformations. The standard MPM formulation, which adopts fluid-implicit-particle (FLIP) transfer, is less dissipative but also less stable than particle-in-cell (PIC) transfer. The affine PIC (APIC) transfer introduces the affine velocity and was developed to realize stable simulations while overcoming dissipations of the angular momentum in PIC. This paper presents Taylor-PIC (TPIC) transfer, which is a type of APIC transfer that combines the affine velocity based on the first-order Taylor series approximation and PIC transfer. TPIC is simple (only the particle velocity gradient is required) and inherits the key advantages of the original APIC, such as less dissipation and stability. Although TPIC does not conserve angular momentum, in contrast to APIC, the velocity gradient is preserved in the transfer between particles and the grid. This velocity gradient contains the angular information in its skew-symmetric component, which allows TPIC to adequately describe angular motion, similar to APIC.Furthermore, the MPM can cause stress oscillations near boundaries when boundary conditions are explicitly imposed, e.g., when the boundary grid velocity (or momentum) is set to zero. When affine-type transfers are used, this instability is exacerbated, and simulations can easily fail. Therefore, we propose a kernel correction method based on the weighted least squares for particles near boundaries. The proposed TPIC transfer and kernel correction are validated through five types of simulations, each using a different material—linear elastic, von Mises, Newtonian fluid, and Drucker–Prager. In the simulations, incorrect results due to stress oscillations are observed even for small deformations. Applying the corrected kernels successfully removes these spurious oscillations, and the results are consistent with the analytical solutions. Moreover, the numerical results confirm the accuracy and robustness of the proposed TPIC transfer.

A Generalized Interpolation Material Point Method for Shallow Ice Shelves. 1: Shallow Shelf Approximation and Ice Thickness Evolution.

We develop a generalized interpolation material point method (GIMPM) for the shallow shelf approximation (SSA) of ice flow. The GIMPM, which can be viewed as a particle version of the finite element method, is used here to solve the shallow shelf approximations of the momentum balance and ice thickness evolution equations. We introduce novel numerical schemes for particle splitting and integration at domain boundaries to accurately simulate the spreading of an ice shelf. The advantages of the proposed GIMPM‐SSA framework include efficient advection of history or internal state variables without diffusion errors, automated tracking of the ice front and grounding line at sub‐element scales, and a weak formulation based on well‐established conventions of the finite element method with minimal additional computational cost. We demonstrate the numerical accuracy and stability of the GIMPM using 1‐D and 2‐D benchmark examples. We also compare the accuracy of the GIMPM with the standard material point method (sMPM) and a reweighted form of the sMPM. We find that the grid‐crossing error is very severe for SSA simulations with the sMPM, whereas the GIMPM successfully mitigates this error. While the grid‐crossing error can be reasonably reduced in the sMPM by implementing a simple material point reweighting scheme, this approach it not as accurate as the GIMPM. Thus, we illustrate that the GIMPM‐SSA framework is viable for the simulation of ice sheet‐shelf evolution and enables boundary tracking and error‐free advection of history or state variables, such as ice thickness or damage.

On the dynamic assessment of the Local-Maximum Entropy Material Point Method through an Explicit Predictor–Corrector Scheme

Material Point Method (MPM) has arisen in the recent years as an alternative to Finite Element Method (FEM) under large deformations. However, the simulation of shock waves propagation in the large deformation regime is still challenging under this approach due to the incapability of the standard MPM time integration scheme to filter spurious noises. To overcome it, we propose in this paper an explicit Predictor–Corrector time integration scheme. Its powerful performance mitigates the presence of spurious oscillations with minimal dissipation in high frequency problems. Other source of numerical noise in MPM occurs when the material points cross computational grid boundaries, being caused by the lack of smoothness of the interpolation functions. This noise results in spurious local variations at the material points, where strain–stress fields are computed. This could lead to inaccurate solutions as well as aborted simulations in the worst cases. To surmount it, we propose in this work the Local Maximum-Entropy (LME) approximation schemes as a robust substitute of the traditional shape function in MPM. The performance of both improvements is validated by the high quality results of the numerical examples.

Modelling contacts with a total Lagrangian material point method

The total Lagrangian Material Point Method (TLMPM) is an efficient formulation of the well known Material Point Method (MPM). Unlike the MPM, it is not subject to cell-crossing instabilities or numerical failure, and features better quadrature. However, the TLMPM is unable to model (no-slip) contacts seamlessly as the standard updated Lagrangian MPM. This paper presents a contact algorithm for the TLMPM able to model contacts with or without friction. Its performance is tested and compared to the standard MPM using five different tests that are popular problems from the MPM literature. Even when using simple hat weighting functions, and unlike in the MPM, the stress field in the TLMPM is fairly smooth and with this algorithm, there is no gap between two surfaces in contact. Moreover, there are less energy losses, and less error in the predicted forces.

A precise critical time step formula for the explicit material point method

SummaryThe material point method (MPM) combines Eulerian method and Lagrangian method and thus both Lagrangian particle position and interaction between neighboring Eulerian grid cells will affect the simulation stability. However, the original critical time step formula in the standard MPM does not reflect the effect of particle position and neighboring cell interaction on stability and overestimates the critical time step so much that the CFL number has to be very small, even smaller than 0.1, to obtain a stable solution at extreme particle positions. Therefore, in many engineering applications, the standard MPM is very expensive due to the small CFL number. In this article, the effect of particle position and neighboring cell interaction on stability of the explicit MPM is studied. An explicit critical time step formula is obtained based on the system eigenvalues in one dimension, and is then extended to two and three dimensions. For extreme deformation problems, the geometric stiffness matrix is taken into consideration which modifies the sound speed of particles in the critical time step formula. Several tests are performed to verify our formula and show a decrease in amount of time steps used for simulation with our formula comparing with the original formula.

Stabilized material point methods for coupled large deformation and fluid flow in porous materials

The material point method (MPM) has been increasingly used for the simulation of large-deformation processes in fluid-infiltrated porous materials. For undrained poromechanical problems, however, standard MPMs are numerically unstable because they use low-order interpolation functions that violate the inf–sup stability condition. In this work, we develop stabilized MPM formulations for dynamic and quasi-static poromechanics that permit the use of standard low-order interpolation functions notwithstanding the drainage condition. For the stabilization of both dynamic and quasi-static formulations, we utilize the polynomial pressure projection method whereby a stabilization term is augmented to the balance of mass. The stabilization term can be implemented with both the original and generalized interpolation material point (GIMP) methods, and it is compatible with existing time-integration methods. Here we use fully-implicit methods for both dynamic and quasi-static poromechanical problems, aided by a block-preconditioned Newton–Krylov solver. The stabilized MPMs are verified and investigated through several numerical examples under dynamic and quasi-static conditions. Results show that the proposed MPM formulations allow standard low-order interpolation functions to be used for both the solid displacement and pore pressure fields of poromechanical formulations, from undrained to drained conditions, and from dynamic to quasi-static conditions.

A Total-Lagrangian Material Point Method for solid mechanics problems involving large deformations

The material point method (MPM) has found successful applications in many engineering problems involving large displacement, large deformation and contacts. The standard MPM formulation, which adopts piece-wise linear basis functions, suffers from the so-called cell-crossing instability, low order of convergence and numerical fracture. Modifications have been made to this standard MPM to mitigate these issues: B-spline MPM (BSMPM), the generalized interpolation material point (GIMP) and convected particle domain interpolation (CPDI) all decrease cell-crossing instabilities and increase the order of convergence, but only CPDI effectively suppresses numerical fracture. However, these methods, CPDI in particular, significantly increase the method’s implementation and computational complexity. This paper presents a total Lagrangian MPM, dubbed TLMPM, that overcomes the issues of the conventional MPM while being more efficient and easier to implement than CPDI. The method is used for impact analyses of a cylinder bar made of steel and necking and fracture of cylinder alloy specimens. The numerical solutions are in satisfactory agreement with the experiment data. No numerical fracture occurred for simulations involving very large tensile deformation without special treatment such as done in the CPDI. Convergence analyses using the method of manufactured solutions show that the TLMPM is second-order accurate for problems of which boundary is axis-aligned. For the challenging generalized vortex problem, it also converges quadratically for relatively coarse meshes. Moreover, the model is able to simulate physically based fracture using continuum damage mechanics.

On the use of domain-based material point methods for problems involving large distortion

Challenging solid mechanics problems exist in areas such as geotechnical and biomedical engineering which require numerical methods that can cope with very large deformations, both stretches and torsion. One candidate for these problems is the Material Point Method (MPM), and to deal with stability issues the standard form of the MPM has been developed into new domain-based techniques which change how information is mapped between the computational mesh and the material points. The latest of these developments are the Convected Particle Domain Interpolation (CPDI) approaches. When these are demonstrated, they are typically tested on problems involving large stretch but little torsion and if these MPMs are to be useful for the challenging problems mentioned above, it is important that their capabilities and shortcomings are clear. Here we present a study of the behaviour of some of these MPMs for modelling problems involving large elasto-plastic deformation including distortion. This is carried out in a unified implicit quasi-static computational framework and finds that domain distortion with the CPDI2 approaches affects some solutions and there is a particular issue with one approach. The older CPDI1 approach and the standard MPM however produce physically realistic results. The primary aim of this paper is to raise awareness of the capabilities or otherwise of these domain-based MPMs.

B-spline based boundary conditions in the material point method

The material point method is an increasingly popular method for tackling solid mechanics problems involving large deformations. However, there are issues associated with applying boundary conditions in the method and, to date, no general approach for imposing both Neumann and Dirichlet boundary conditions has been proposed. In this paper, a new B-spline based boundary method is developed as a complete methodology for boundary representation, boundary tracking and boundary condition imposition in the standard material point method. The B-spline interpolation technique is employed to form continuous boundaries which are independent of the background mesh. Dirichlet boundary conditions are enforced by combining the B-spline boundaries with an implicit boundary method. Neumann boundary conditions are included by direct integration of surface tractions along the B-spline boundary. This general boundary method not only widens the problems that can be analysed by all variants of the material point method, when implemented using an implicit solver, but is also applicable to other embedded and non-matching mesh approaches. Although the Dirichlet boundary conditions are restricted to implicit methods, boundary representation, tracking and Neumann boundary condition enforcement can be applied to explicit and implicit methods.

Conservative Taylor least squares reconstruction with application to material point methods

SummaryWithin the standard material point method (MPM), the spatial errors are partially caused by the direct mapping of material‐point data to the background grid. In order to reduce these errors, we introduced a novel technique that combines the least squares method with the Taylor basis functions, called the Taylor least squares (TLS), to reconstruct functions from scattered data while preserving their integrals. The TLS technique locally approximates quantities of interest such as stress and density, and when used with a suitable quadrature rule, it conserves the total mass and linear momentum after transferring the material‐point information to the grid. The integration of the technique into MPM, dual domain MPM, and B‐spline MPM significantly improves the results of these methods. For the considered examples, the TLS function reconstruction technique resembles the approximation properties of highly accurate spline reconstruction while preserving the physical properties of the standard algorithm.

Hysteretic behavior using the explicit material point method

The material point method (MPM) is an advancement of particle in cell method, in which Lagrangian bodies are discretized by a number of material points that hold all the properties and the state of the material. All internal variables, stress, strain, velocity, etc., which specify the current state, and are required to advance the solution, are stored in the material points. A background grid is employed to solve the governing equations by interpolating the material point data to the grid. The derived momentum conservation equations are solved at the grid nodes and information is transferred back to the material points and the background grid is reset, ready to handle the next iteration. In this work, the standard explicit MPM is extended to account for smooth elastoplastic material behavior with mixed isotropic and kinematic hardening and stiffness and strength degradation. The strains are decomposed into an elastic and an inelastic part according to the strain decomposition rule. To account for the different phases during elastic loading or unloading and smoothening the transition from the elastic to inelastic regime, two Heaviside-type functions are introduced. These act as switches and incorporate the yield function and the hardening laws to control the whole cyclic behavior. A single expression is thus established for the plastic multiplier for the whole range of stresses. This overpasses the need for a piecewise approach and a demanding bookkeeping mechanism especially when multilinear models are concerned that account for stiffness and strength degradation. The final form of the constitutive stress rate–strain rate relation incorporates the tangent modulus of elasticity, which now includes the Heaviside functions and gathers all the governing behavior, facilitating considerably the simulation of nonlinear response in the MPM framework. Numerical results are presented that validate the proposed formulation in the context of the MPM in comparison with finite element method and experimental results.

Smoothing algorithm for stabilization of the material point method for fluid–solid interaction problems

Phenomena involving general solid–water interactions such as flows with debris are challenging to model numerically because they are not easily represented using solid- or fluid-oriented methods. The material point method (MPM) provides a unified multi-material interaction platform potentially capable of modeling complex solid–water flow phenomena. However, it is necessary to address volumetric locking for (nearly) incompressible materials when modeling fluids, while also stabilizing integration errors that arise in standard MPM. This paper examines these challenges in depth, and presents a flux-based smoothing algorithm designed to address integration-error-induced destabilization via controlled strain energy dissipation. The effectiveness of the algorithm is demonstrated with two simple but fundamental fluid/solid problems, and with an application to a complex solid–water dynamic interaction problem. Results show the flux-based smoothing algorithm is capable of stabilizing the side-effects of numerical integration errors, while at the same time remaining inactive if there is no integration-error-induced oscillation. Based on this study, the flux-based smoothing algorithm is suggested as a stabilization scheme for MPM when using constant-interpolated hybrid elements.

IGIMP: An implicit generalised interpolation material point method for large deformations

The Material Point Method (MPM) uses a combined Eulerian-Lagrangian approach to solve problems involving large deformations. A problem domain is discretised as material points which are advected on a background grid. Problems are encountered with the original MPM when material points cross between grid cells, and this has been tackled by the development of the Generalised Interpolation MPM, where material points’ domains of influence extend beyond the currently occupied grid cell. In this paper, the Generalised Interpolation Material Point (GIMP) Method has been implemented implicitly in a manner that allows a global stiffness matrix to be constructed similar to that in the Finite Element Method (FEM) by combining contributions from individual elements on the background grid. An updated Lagrangian finite deformation framework has been used to ensure non-linear behaviour within each of the loadsteps. The weighting functions used for this which make the GIMP method different to standard MPM are presented and the implementation is explained. Specific details on computing the deformation gradient to be consistent with the updated Lagrangian framework and the updating of the material point influence domains are outlined, both of which are currently unclear in the published literature. It is then shown through numerical examples that for both small and large deformation elastic and elasto-plastic problems, the implicit GIMP method agrees well with analytical solutions and exhibits convergence properties between that of linear and quadratic FEM.

Application of a Second-order Implicit Material Point Method

In this work, we present the application of a newly developed implicit second-order Material Point Method (MPM) on offshore geotechnical applications. The presented second-order MPM uses a special set of piecewise quadratic shape functions to circumvent the well-known issue of producing zero nodal mass contributions. To mitigate the effect of the standard MPM to produce highly oscillating stresses across cell interfaces, we have carried over our ideas obtained from the derivation of the second-order MPM to the Dual Domain Material Point (DDMP) Method, too. The resulting second-order DDMP Method produces a smoother stress distribution across the entire computational domain while being able to profit from the improved convergence rates of second-order finite elements. In a numerical example from geotechnical engineering applications, we illustrate the practical application of our enhanced Material Point and DDMP Methods by simulating a cone penetration.

Contact algorithms for the material point method in impact and penetration simulation

AbstractThe inherent no‐slip contact constraint in the standard material point method (MPM) creates a greater penetration resistance. Therefore, the standard MPM was not able to treat the problems involving impact and penetration very well. To overcome these deficiencies, two contact methods for MPM are presented and implemented in our 3D explicit MPM code, MPM3D. In MPM, the impenetrability condition may not satisfied on the redefined regular grid at the beginning of each time step, even if it has been imposed on the deformed grid at the end of last time step. The impenetrability condition between bodies is only imposed on the deformed grid in the first contact method, while it is imposed both on the deformed grid and redefined regular grid in the second contact method. Furthermore, three methods are proposed for impact and penetration simulation to determine the surface normal vectors that satisfy the collinearity conditions at the contact surface. The contact algorithms are verified by modeling the collision of two elastic rings and sphere rolling problems, and then applied to the simulation of penetration of steel ball and perforation of thick plate with a particle failure model. In the simulation of elastic ring collision, the first contact algorithm introduces significant disturbance into the total energy, but the second contact algorithm can obtain the stable solution by using much larger time step. It seems that both contact algorithms give good results for other problems, such as the sphere rolling and the projectile penetration. Copyright © 2010 John Wiley & Sons, Ltd.