A machine-learning framework for peridynamic material models with physical constraints
aa r X i v : . [ m a t h . NA ] J a n A machine-learning framework for peridynamicmaterial models with physical constraints
Xiao Xu ∗ Marta D’Elia † John T. Foster ‡ January 5, 2021
Abstract
As a nonlocal extension of continuum mechanics, peridynamics hasbeen widely and effectively applied in different fields where discontinu-ities in the field variables arise from an initially continuous body. Animportant component of the constitutive model in peridynamics is the influence function which weights the contribution of all the interactionsover a nonlocal region surrounding a point of interest. Recent work hasshown that in solid mechanics the influence function has a strong rela-tionship with the heterogeneity of a material’s micro-structure. How-ever, determining an accurate influence function analytically from agiven micro-structure typically requires lengthy derivations and com-plex mathematical models. To avoid these complexities, the goal ofthis paper is to develop a data-driven regression algorithm to find theoptimal bond-based peridynamic model to describe the macro-scaledeformation of linear elastic medium with periodic heterogeneity. Wegenerate macro-scale deformation training data by averaging over pe-riodic micro-structure unit cells and add a physical energy constraintrepresenting the homogenized elastic modulus of the micro-structureto the regression algorithm. We demonstrate this scheme for examplesof one- and two-dimensional linear elastodynamics and show that theenergy constraint improves the accuracy of the resulting peridynamicmodel. ∗ The Oden Institute for Computational Engineering and Sciences, The University ofTexas at Austin, TX, [email protected] † Computational Science and Analysis, Sandia National Laboratories, CA,[email protected] ‡ The Oden Institute for Computational Engineering and Sciences, The University ofTexas at Austin, TX, [email protected] Introduction
First proposed by Silling [2000], peridynamic (PD) mechanics models replacethe spatial derivatives in the classical conservation of momentum equationwith an integral functional to determine the net internal force density onmaterial points during deformation. The integral equation provides conve-nience in modeling deformation problems with displacement discontinuities(e.g. cracks). The original PD formulation that is widely used in literatureis the so-called bond-based peridynamic model , which uses a pairwise forcefunctional to describe the interaction between material particles Zheng et al.[2020]. A more general theory of peridynamics for solid mechanics, called state-based peridynamics , was later proposed by Silling et al. [2007]; how-ever, the simpler bond-based peridynamic models still have a wide range ofapplications in brittle fracture Huang et al. [2015], Wang et al. [2018]. Ad-ditionally, bond-based models are useful for demonstrating advancements incomputational techniques/implementations, homogenization theory, and/ornumerical analysis where the more complicated state-based theory may beunnecessarily distracting from the central advancement of the research.The equation of motion in the bond-based PD theory is given by ρ ( x )¨ u ( x , t ) = Z H f ( η , ξ ) d ξ + b ( x , t ) , (1)where η := u ( x + ξ ) − u ( x ) ,ρ is the mass density, H is a neighborhood of x where integration of bondforces is carried out, u is the displacement vector field, b is a prescribed bodyforce density field, and f is the pairwise force function determining the forcedensity that the particle at x + ξ exerts on the particle located by the positionvector x . Regarding the constitutive model for the force function f , a scalarfunction is used to assign weights that scale the pair-wise force associatedwith each ξ . For example, in Silling and Askari [2005] the bond-based PDmodel for a linear elastic solid undergoing small deformation defined f as f ( η , ξ ) = c ( ξ ) s ξ | ξ | , (2)where c is an influence function which can be thought of as a spring constantin this setting that is, in general, unique for each ξ . The bond stretch s isdefined as s := | ξ + η || ξ | . With a special choice of c , (1) and (2) converge to the classical Cauchylinear momentum equation for a linear elastic solid with a Poisson ratio of / (in three dimensions) when the region of integration H is reduced to aninfinitesimal volume Silling and Lehoucq [2008].2n a general setting, it is clear that the influence function plays an im-portant role in determining the magnitude of the interaction between twomaterial points and the overall effect of nonlocality of the material. It iswell-established that the influence function is a key factor contributing tothe behavior of peridynamic material models, especially in the case of wavedispersion Weckner and Abeyaratne [2005] and fracture Seleson and Parks[2011]; however, there is still no general technique to systematically deter-mine the influence function for different materials or to justify its choice fora given application.A few attempts have been made; for example, Delgoshaie et al. [2015]showed that the multi-scale connectivity of natural pore networks can beused to extract nonlocal kernel functions for use in continuum nonlocaldiffusion models. D’Elia and Gunzburger [2016] used an optimal controltechnique to identify the nonlocal diffusivity parameter for both peridy-namic and fractional nonlocal models. Sridhar et al. [2018] presented a gen-eral multi-scale elastodynamic framework based on the Floquet-Bloch trans-form. Wildman [2019] used a prescribed dispersion relation to derive the cor-responding discrete micro-modulus function coefficients. More recently, Aksoylu and Gazonas[2020] developed a selection criterion for nonlocal kernels derived from Taylorexpansions that best approximate the classical linear dispersion relations forelastic solids. You et al. [2020a] introduced a nonlocal-equation-constrainedoptimization algorithm to identify the optimal influence function for wavepropagation through a one-dimensional heterogeneous bar.In an explanation on the origin of nonlocality in solid materials, Silling[2014] demonstrated that nonlocality can arise from the small-scale het-erogeneity that is excluded through an implicit or explicit homogenizationprocedure. Inspired by this idea, Xu and Foster [2020] focused on one-dimensional elastodynamics of a periodic heterogeneous bar and built atheoretical method to determine the PD influence function from the micro-structure. However, this theoretical method is effectively intractable forhigher dimensional problems. Recently, You et al. [2020b] used machinelearning to develop invertible nonlocal models from high-fidelity syntheticdata while guaranteeing the well-posedness of the learned operator. Com-bining ideas from the previous two references, this work incorporates data-driven methods together with computational homogenization techniques toderive PD influence functions based on high-fidelity data associated to aspecific material’s micro-structure.In this work, we aim to learn the discrete PD influence function fromdata, circumventing analytical complexity of conventional computational ho-mogenization theories and avoiding costly fine-scale simulations that treatthe heterogeneities in the micro-structure explicitly. The synthetic high-fidelity data set consists of microscopic displacement data of a given mate-rial obtained by solving the dynamic equations at micro-scale with a highlyresolved micro-structural finite element model (FEM). This data set is then3oarse-grained in each unit cell and it is used as training data set. Thediscrete PD influence function is the result of a data-driven regression al-gorithm; the corresponding nonlocal model serves as a surrogate model todescribe the mechanics of the material at the macro-scale. One of our maincontributions is to demonstrate that the use of a physically-justified energyconstraint imposed on the objective function reduces the amount of train-ing data needed for a highly accurate homogenized model. This constraintenforces the Hill-Mandel macro-homogeneity condition under isotropic de-formation. Our results show that this regression scheme is able to robustlylearn the discrete PD influence function and to improve the accuracy ofthe predictions compared to standard choices of the PD influence function.In fact, standard nonlocal influence functions are nonnegative; while thisguarantees well-posedness, it compromises the ability to accurately predictthe displacement Weckner and Silling [2011]. Instead, the class of kernelslearned by our algorithm is allowed to be sign-changing.The rest of the paper is organized as follows: §2 outlines the generalmulti-scale framework for a linear elasticity problem with a periodic het-erogeneous micro-structure and the learning algorithm for discrete micro-modulus functions for the macro-scale nonlocal model. In §3 and §4, thealgorithm is applied to one- and two-dimensional elastodynamics problemsand the energy constraints for both cases are introduced. §5 presents testingresults of the learned PD influence functions and illustrates the benefit ofusing the energy constraint. §6 summarizes our contributions and providesfuture research guidelines. Consider the linear elastodynamics of an open bounded domain Ω ⊂ R n ,with n being the number of space dimension. Let the elastic medium be peri-odically heterogeneous with microstructural length scale being much smallerthan the length scale of the domain Ω l ≪ L, (3)where L and l are used to denote the length scale of domain Ω and the oneof its heterogeneous micro-structures respectively. Let { r i } ni =1 be the setof lattice vectors that describe the micro-structures unit cell. Then { r i } ni =1 form a basis for R n and we can define the periodic lattice R as R = { v ∈ R n | v = n X i =1 a i r i , a i ∈ Z } , Z is the set of all integers. Naturally, every single unit cell in Ω canbe specified using a unique set of integers a = ( a , a , . . . , a n ) as T ( a ) = { v ∈ R n | v = n X i =1 b i r i , b i ∈ [ a i , a i + 1] } . We assume that, for t ∈ [0 , T ], only small deformations occur within theperiodic medium Ω under external loads and the displacement of materialpoints in the domain is governed by the linear elasticity equations. There-fore, the displacement field u satisfies the following generalized Hook’s law σ ( x , t ) = C ( x ) : 12 ( ∇ u ( x , t ) + ( ∇ u ( x , t )) ⊺ ) x ∈ Ω and 0 ≤ t ≤ T, where σ is the stress tensor field, the colon is the double contraction opera-tor , ∇ is the spatial gradient operator and C is the material’s fourth-orderelasticity tensor. Due to the periodic nature of the material, the stiffnesstensor C must satisfy C ( x + r ) = C ( x ) x , x + r ∈ Ω and r ∈ R . Similarly, the same periodicity condition holds for material’s density ρρ ( x + r ) = ρ ( x ) x , x + r ∈ Ω and r ∈ R . The equation of linear momentum conservation is ρ ( x )¨ u ( x , t ) = ∇ · σ ( x , t ) + b ( x , t ) ∀ x ∈ Ω and 0 ≤ t ≤ T, (4a) B u ( x , t ) = q ( x , t ) ∀ x ∈ ∂ Ω and 0 ≤ t ≤ T, (4b)where the superimposed double dot on u indicates double differentiation intime, ∇· is the spatial divergence operator, b is the body force density fieldover the domain Ω and B denotes a boundary operator ( B being identitymap corresponds to a Dirichlet boundary condition and B = ∂/∂n a Neu-mann boundary condition). Note that the displacement u ( x , t ) in the aboveequations represents the accurate displacement of every material point insidethe domain Ω. Thus, u ( x , t ) should be treated as the micro-scale displace-ment in the multi-scale modeling scheme. With the purpose of describingthe system by means of a courser-scale model, we introduce the macro-scaledisplacement ¯ u ( x , t ) as the average of the micro-scale displacement insideeach unit cell.To be more specific, for each x ∈ Ω, there exists a unique set of integers a x = ( a , a , . . . , a n ) such that x ∈ T ( a x ); we define the macro-displacementat x as the average micro-displacement of the unit cell T ( a x )¯ u ( x , t ) = R T ( a x ) u ( ξ , t ) d ξ R T ( a x ) d ξ . (5) i.e. A : B = A ij B ij u , to represent the macro-scale averageover the unit cell. We propose that the macro-displacement ¯ u ( x , t ) of thesolutions to the elastodynamic equation (4) with periodically oscillatingmaterial constants satisfies the PD bond-based model¯ ρ ( x )¨¯ u ( x , t ) = L ω [¯ u ]( x ) + ¯ b ( x , t ) x ∈ Ω , ≤ t ≤ T, (6a) B I ¯ u ( x , t ) = ¯ q ( x , t ) x ∈ Ω I , ≤ t ≤ T, (6b)where B I is the corresponding nonlocal interaction operator specifying avolume constraint in an appropriate nonlocal interaction domain Ω I . Forthe bond-based PD model the nonlocal operator L ω is defined as L ω [¯ u ]( x ) = Z H ω ( ξ )(¯ u ( x + ξ , t ) − ¯ u ( x , t )) d ξ , (7)where ω is the kernel function that determines the nonlocality of the operatorand H is a neighborhood of x . In the model of peridynamics, a characteristiclength scale ǫ called the medium’s horizon is often used to identify theneighborhood H as H = { x + ξ | x ∈ Ω , x + ξ ∈ Ω , | ξ | < ǫ } . As a consequence, for the well-posedness of problem (6) we require Ω I to bea layer, or collar, of thickness ǫ surrounding the domain, see, e.g., Du et al.[2012].Owing to the definition of our macro-displacement (5) , the macro-displacement ¯ u ( x ) is a constant function inside each unit cell, which leadsus to naturally discretize the nonlocal equation (6) as the summation overthe unit cells inside H ¯ ρ i ¨¯ u i ( t ) = X j ∈H n ω i,j (¯ u j ( t ) − ¯ u i ( t )) + ¯ b i ( t ) i ∈ Ω n , ≤ t ≤ T, where Ω n is the enumeration of the unit cells in Ω and H n is the numberingof the cells in H . The subscript i, j are used to identify the unit cells and, in ω i,j , they indicate the interacting cells. Since the length scale of the unit cellis much smaller than the domain, as highlighted in (3) , we simply define ω i,j as the product of the value of the influence function at the unit cell andthe volume of the unit cell itself. With this choice, in the PD literature ω i,j is often referred to as discrete micro-modulus function between unit cell i and unit cell j . Explicit definitions of the discrete micro-modulus functionare introduced later on in the paper. In order to learn the kernel function ω ( ξ ), we assume that we are given high-fidelity solutions to the elastodynamic equations (4) at micro-scale, which6an be either generated by numerically solving the equations using FEMwith a mesh that is refined enough to describe the micro-scale displacementinside the unit cell, or directly obtained from experimental data with highaccuracy.The first step is to homogenize such solutions as described in (5) . Thisresults in a collection of macro-scale data of displacement and accelerationat several time instants { ¯ u ( x , t i ) , ¨¯ u ( x , t i ) } N t i =1 . We use these data to train our macro-scale PD equation (6) and obtainan optimal surrogate model ω ∗ for the kernel function. This is achieved bysolving an optimization problem of the following form ω ∗ = argmin ω N t X i =1 (cid:13)(cid:13)(cid:13) L ω [¯ u ( x , t i )] + ¯ b ( x , t i ) − ¯ ρ ( x )¨¯ u ( x , t i ) (cid:13)(cid:13)(cid:13) Ω , (8)where k·k Ω denotes an appropriate norm over Ω. Note that the above train-ing procedure only makes use of the macro-scale dynamic data which canonly describe the deformation of the medium under certain loads and bound-ary conditions; more importantly, the minimization problem does not con-tain any physical information about the micro-structure of the medium.Thus, the accuracy of the predictions corresponding to the optimal kernelhighly depends on the type and number of training data. This fact cancompromise the generalization properties of this algorithm and, hence, thequality of the predictions. To overcome this limitation, we add physics-based information to the cost function that take into account features ofthe micro-structure and, possibly, acts as regularizers for the optimizationproblem. With this addition, we reformulate (8) as ω ∗ = argmin ω N t X i =1 (cid:13)(cid:13)(cid:13) L ω [¯ u ( x , t i )] + ¯ b ( x , t i ) − ¯ ρ ( x )¨¯ u ( x , t i ) (cid:13)(cid:13)(cid:13) Ω + N c X i =1 α i C i ( ω ) , (9)where α i is the regularization parameter, C i is the corresponding physics-based constraint and N c is the number of constraints. In the next sections,the one- and two-dimensional elastodynamics models are further specified,and we present specific constraints, such as (15) and (16) , that describedifferent features of the micro-structure of the given medium. We consider a one-dimensional elastic composite problem similar to the onediscussed in [Xu and Foster, 2020]: a composite rod made of a periodic array7 loadl / l / l / E s , ρE c , ρ Figure 1: One-dimensional composite with periodic heterogeneityof two linearly elastic, homogeneous, and isotropic constituents with perfectinterfaces. The rod is fixed on one end and external loads are applied on theother end; we also assume that no body forces are applied along the bar,i.e. b = . We define the micro-structure of the medium as the symmetricheterogeneous unit cell shown in Figure 1, where l is the length of the unitcell, the dark block represents the stiffer constituent with elastic modulus E s and density ρ and the white block represents the more compliant constituentwith elastic modulus E c and density ρ . If we use natural ordering from leftto right as the numbering of the unit cells of the composite rod, we cansimplify the discrete PD equation (7) as ρ ¨¯ u i ( t ) = ǫ/l X j = − ǫ/l ω j (¯ u i + j ( t ) − ¯ u i ( t )) . (10)Because of the symmetry of the medium’s micro-structure, we enforce thesymmetry constraint on the discrete micro-modulus function as follows: ω − i = ω i < i ≤ ǫ/l . (11)By using the discrete expression in (10) in combination with (11), the min-imization problem (8) reduces to the following linear regression problemfor the vector v ω := ( ω , ω , . . . , ω ǫ/l ), v ∗ ω = argmin v ω N t X i =1 ǫ/l X j =1 ω j (¯ u i + j ( t i ) + ¯ u i − j ( t i ) − ¯ u i ( t i )) − ¨¯ u i ( t i ) , (12)where we used the ℓ norm with respect to the values of the operator andthe acceleration in each unit cell.We next describe how to obtain a physics-based constraint that embedsthe micro-structure information in the regression problem. In [Silling and Askari,2005] the authors consider a large homogeneous body under isotropic ex-tension and derive a relationship between the PD kernel function and theclassical bulk modulus. Inspired by this work, we also assume that the8omposite rod is under constant strain, i.e. the macro-scale displacementcan be expressed as ¯ u ( x ) = sx , where s is the corresponding macro-scaleconstant strain. As a consequence, the micro-scale displacement inside eachconstituent is also linear, which implies that the macro-scale constitutiveequation is given by ¯ σ = E hom s. (13)Here, E hom is the homogenized elastic modulus of the unit cell E hom = 21 /E s + 1 /E c . In the framework of bond-based peridynamics, we can express the stress ¯ σ at point x = 0 as the sum of all the bond forces in bonds that cross x = 0[Silling et al., 2003], i.e.¯ σ = Z ∞ Z ∞ ω ( r + s )(¯ u ( r ) − ¯ u ( − s ))dsdr , = s Z ∞ ξ ω ( ξ )d ξ. (14)Combining (13) with (14) and substituting discrete micro-modulus gives E hom = ǫ/l X i =1 ( i l ) ω i . Therefore, we define the constraint as C ( v ω ) = E hom − ǫ/l X i =1 ( il ) ω i , (15)and we refer to it as an energy constraint . Note that this corresponds toenforcing the effective elastic modulus of medium’s micro-structures. We extend the previous one-dimensional learning framework to two dimen-sions. We consider a square composite thin plate with one side fixed andexternal loads applied on the other side as shown in Figure 2. Also in thiscase, we assume that no external forces are applied, i.e. b = . The com-posite thin plate has periodic micro-structure whose lattice vectors are theunit vectors along the x - and y -axis. The micro-structure unit cell of thethin plate is a square composed of two different isotropic constituents; thecompliant inclusions (represented by white block in Figure 2) is embedded in9 oadxy l l l E s , ρE c , ρ Figure 2: Two-dimensional composite with periodic heterogeneitythe stiffer continuous matrix phase (represented by dark grey block in Fig-ure 2); the inclusions also have a square shape and are placed in the middleof the unit cell. Limited by the fact that the bond-based PD model canonly be used to describe plane stress deformation of medium whose Poissonratio is / [Trageser and Seleson, 2020], we assume that the Poisson ratioof every constituent in the medium is / for consistency. Following similarprocedures as in §3, the discretization of the two-dimensional PD equationfor macro-scale displacements, the minimization problem (8) can also bereduced to a linear regression problem for the discrete micro-modulus as(12) . We choose the square neighborhood H of point ( x, y ) H = { ( x + ξ, y + η ) | − ǫ ≤ ξ ≤ ǫ, − ǫ ≤ η ≤ ǫ } . This choice is motivated by the specific micro-structure considered in thissection and, more importantly, it does not compromise the well-posednessof the problem nor the convergence to the local limit as the nonlocal neigh-borhood shrinks, see, e.g. D’Elia et al. [2020]. We use natural ordering in x direction and y direction as the two-dimensional numbering of unit cells todiscretize the PD equation for macro-scale displacements (7) as ρ ¨¯ u p,q ( t ) = ǫ/l X i,j = − ǫ/l ω i,j (¯ u p + i,q + j ( t ) − ¯ u p,q ( t )) . Since the structure of the unit cell is symmetric about i = 0, j = 0 and i = j and each constituent is isotropic, we can enforce the following symmetryconstraints on our PD kernel function ω i,j = ω − i,j = ω i, − j = ω j,i − ǫ/l ≤ i, j ≤ ǫ/l . { ω i,j } i,j = ǫ/l i,j = − ǫ/l .For the derivation of a two-dimensional energy constraint, we again con-sider the plane stress problem of the medium under isotropic extension, i.e.¯ u ( x ) = s x where s is a constant and compute the average strain energydensity of the unit cell W uc from micro-scale displacement solutions (which,in our case, are obtained via FEM simulations). Once again, in the PDframework, by using the definition of the micro-potential [Silling and Askari,2005], the strain energy density can be expressed as W uc = 12 Z H ω ( ξ ) s | ξ | ξ . Therefore, we define the energy constraint in two-dimensions as C ( ω ) = W uc − Z H ω ( ξ ) s | ξ | ξ ! , (16)which weakly prescribes the value of the average strain energy density ofthe medium under isotropic extension (i.e. the effective bulk modulus of themicro-structure). In this section we use several one- and two-dimensional numerical exam-ples to illustrate the performance of our PD kernel learning procedure anddemonstrate the effectiveness of the energy constraint.
We consider the one-dimensional problem described in §3 with geometryparameters L = 1m and l = 0 . E s = 200GPa, E c = 5GPa and ρ = 8000kg/m . We perform FEM simulations to calculatesolutions to the small-scale elastodynamic equations (4) with the followingtime-dependent boundary condition u bc ( t ) = u a t ( t − T s ) [1 − H ( t − T s )] , where u = 1 × − m, a is a scaling factor, H is the Heaviside functionand T s = 0 . , T ], for T = 10 − s. We divide this data set into atraining set, for t ∈ [0 , T t ], and a testing data set, for t ∈ ( T t , T ) and we set T t = 0 . × − s. We choose the horizon of the PD model to be ǫ = 8 l ; in11rder to avoid nonlocal boundary effects we discard the macro-scale data ofpoints whose neighborhood is not fully contained in the domain.To test the efficacy of the energy constraint, we solve the linear regressionproblems both without and with the energy constraint, i.e. problems (8) and(9) , respectively. The optimal discrete micro-modulus functions { ω ∗ i } i =8 i = − are reported in Figure 3. It is important to note that these two discretemicro-moduli generate positive definite matrices when discretizing the PDequations (6) ; this ensures the well-posedness of solutions to (6) . − − × × × numbering i ω i ( k g · m − · s − ) . . . without energy constraintwith energy constraint Figure 3: Discrete micro-modulus function ω i for the one-dimensional testcase for both the unconstrained and constrained formulations.To test the performance of the learning algorithm, we use the two discretemicro-moduli to predict the macro-scale acceleration for t ∈ ( T t , T ]. InFigure 4 we report the predicted acceleration of the unit cell at the middleof the rod as a function of time; as a reference, we also report the macro-scaleresult directly calculated from FEM simulations. The predicted accelerationusing the micro-modulus with the energy constraint shows good agreementwith the FEM solutions while the prediction using micro-modulus withoutthe energy constraint has small but noticeable differences. We further testthese two micro-moduli by predicting the deformation of the bar after time T t , i.e. solving the bond-based PD equations (6) from time T t to time T and we report the predicted macro-scale displacement of the middle unitcell in Figure 5.Due to the accumulation of error, the prediction error of the micro-modulus without the energy constraint becomes more pronounced whereasthe results corresponding to the use of the energy constraint still maintainhigh accuracy for all time.In order to better show the improvements of the energy constraint on theperformance of our algorithm, we vary the size of the training data (i.e. varythe value of T t ) and compute the relative prediction error in ℓ norm. The12 . . . . . − . × − . × − . × . × . × t (s) a cce l e r a t i o n ¨¯ u ( m / s ) . . . without constraintFEM 0 . . . . . t (s) with constraintFEM Figure 4: Testing results on the acceleration of the middle unit cell . . . . . . . . . . t (s) d i s p l a ce m e n t ¯ u ( m ) without constraintFEM 0 . . . . . t (s) with constraintFEM Figure 5: Testing results on the displacement of the middle unit cell13esults shown in Figure 6 indicate that, when the amount of training datais small, the energy constraint provides more accurate predictions, whereaswhen the amount of training data increases, the performance of the learningalgorithm without the energy constraint is similar to that with the energyconstraint. . . . . . . − . − . − . − . − . − . T t (s) t e s t e rr o r without constraintwith constraint Figure 6: Relative testing error of ¨¯ u for different size of training dataWe further validate the performance of the learning algorithm by us-ing the optimal micro-modulus function learned in the previous experi-ments to predict the deformation of the rod under different loads. Wechange the time-dependent displacement boundary condition to u bc ( t ) = u a sin(2 πt/T s )[1 − H ( t − T s )], solve the bond-based PD equations withthe optimal micro-modulus, and compare the displacement solutions withthe macro-scale displacement calculated using FEM. The comparison resultsare shown in Figure 7: we observe that the bond-based PD model success-fully describes the deformation of the periodic heterogeneous rod under adifferent loads. This shows that our algorithm generalizes well for differentloading scenarios. We consider the two-dimensional plane stress problem described in §4 withgeometry parameters L = 1m and l = 1 /
3m and material properties valuesare: E s = 200GPa, E c = 5GPa, and ρ = 8000kg/m . The time-dependentdisplacement boundary condition used for FEMsimulations is given by¯ u x ( L, y ) = u a t ( t − T s ) [1 − H ( t − T s )] 0 ≤ y ≤ L, (17)where u = 1 × − m, a is a scaling factor, H is the Heaviside func-tion and T s = 0 . T = 5 × − s and T t = 0 . × − s.14 . . . . . − . − . − . . . . . t (s) d i s p l a ce m e n t ¯ u ( m ) . without constraintFEM 0 . . . . . t (s) with constraintFEM Figure 7: Validation results on the displacement of the middle unit cellWe choose the horizon of the PD model to be ǫ = 6 l and, as done in theone-dimensional case, we do not consider macro-scale data of points whosenonlocal neighborhood is not contained in the domain.We select an appropriate weighted ℓ norm for the objective functions,perform linear regression with and without energy constraint and report theoptimal values of the discrete micro-modulus function { ω ∗ i,j } i,j =6 i,j = − at j = 0in Figure 8. Also in this case these two micro-moduli generate positive def-inite matrices when discretizing equation (6) . We test the two optimal − − − . × . × . × numbering i ω i , ( k g · m − · s − ) . . . without energy constraintwith energy constraint Figure 8: Discrete micro-modulus ω i, for two-dimensional elasticitymicro-moduli on the testing data set and report the predicted macro-scaleacceleration ¨¯ u x and ¨¯ u y of the middle unit cell as a function of time in Fig-ure 9, compared with corresponding macro-scale FEM solutions. It is evi-15 × − × − × × × × a cce l e r a t i o n ¨¯ u x ( m / s ) . . without constraintFEM with constraintFEM0 . . . . . − × − × × × t (s) a cce l e r a t i o n ¨¯ u y ( m / s ) . . without constraintFEM 0 . . . . . t (s) with constraintFEM Figure 9: Testing results on the acceleration of the middle unit cell16ent that the accuracy of the energy-constrained prediction is better whencompared to the unconstrained one, especially for the y -component of ac-celeration. This is not unexpected since the magnitude of the training datais higher in the x direction. As reported in the previous section, we alsosolve the bond-based PD equations using the two optimal micro-moduli in( T t , T ) and report the corresponding ¯ u x together with the reference FEMsolutions in Figure 10. In Figures 11 and 12 we report the relative testingerror (defined as in the previous section) in correspondence of different sizesof the training data set (i.e. different values of T t ) for ¨¯ u x and ¨¯ u y respectively.Similar considerations as in the one-dimensional case can be inferred. . . . . . − . − . . . . t (s) d i s p l a ce m e n t ¯ u x ( m ) . without constraintFEM 0 . . . . . t (s) with constraintFEM Figure 10: Testing results on the displacement of the middle unit cell . × − . × − . × − . × − − . − . − . − . . T t (s) t e s t e rr o r without constraintwith constraint Figure 11: Relative testing error of ¨¯ u x for different size of training dataFinally, we test the optimal micro-modulus functions used in the previ-17 . × − . × − . × − . × − − . − . − . − . − . − . . T t (s) t e s t e rr o r without constraintwith constraint Figure 12: Relative testing error of ¨¯ u y for different size of training data . . . . − × − − × − × − t (s) d i s p l a ce m e n t ¯ u x ( m ) . . . without constraintFEM 0 . . . . t (s) with constraintFEM Figure 13: Validation results on x displacement of the middle unit cell18us tests on another macro-scale data set generated using a different load.Specifically, we generate the training and validation data set with the fol-lowing boundary condition, which corresponds to a shear load¯ u y ( L, y ) = u a t ( t − T s ) [1 − H ( t − T s )] 0 ≤ y ≤ L. The prediction results for the x -displacement are shown in Figure 13; onceagain, the accuracy of the predictions implies that our energy-constrainedlearning algorithm generalizes well to loads that are different from the oneused for training. l l l E s , ρE c , ρ Figure 14: Two-dimensional heterogeneous unit cell − − − × × × numbering i ω i , ( k g · m − · s − ) . . . without energy constraintwith energy constraint Figure 15: Discrete micro-modulus ω i, for two-dimensional elasticityWith the purpose of proving the generality of our algorithm, we con-sider different medium of the same size as the one in Figure 2, but withdifferent micro-structure (see the illustration in Figure 14). We use thetime-dependent boundary condition reported in (17) and generate themacro-scale data which is divided into training set and testing data set with19 × − × − × × × × a cce l e r a t i o n ¨¯ u x ( m / s ) . . without constraintFEM with constraintFEM0 . . . . . . − × − × × × t (s) a cce l e r a t i o n ¨¯ u y ( m / s ) . . without constraintFEM 0 . . . . . . t (s) with constraintFEM Figure 16: Testing results on the acceleration of the middle unit cell20 = 7 . × − and T t = 2 . × − . We choose the horizon for this mediumto be ǫ = 6 l , and perform linear regression on the training set with andwithout the energy constraint. The resulting discrete micro-modulus func-tions are plotted in Figure 15. We use these two discrete micro-modulusfunctions to predict the macro-scale acceleration of the testing data set andplot the comparison results of the middle unit cell in Figure 16. The ac-curacy of our results indicates that the learning algorithm performs wellregardless of the micro-structure. The importance of the PD influence function has recently gained more at-tention from the PD research community, but what has been missing inthe literature is a systematic way for it to be determined for a given mate-rial and application setting. In this work, we use a bond-based PD modelto describe the linear elastic deformations for materials with periodic het-erogeneity. We used a highly-resolved micro-structural FEM to solve theclassical elastodynamic equations for a short time periods and upscaled thesolutions to an averaged macro-scale deformation. We used these solutionsas training data in a machine-learning framework to identify the optimaldiscrete micro-modulus for the PD model. In the regression algorithm, anenergy constraint that represents the average elastic modulus of the micro-structure was added to the objective function. The testing results indicatedthat the homogenized macro-scale deformation can be predicted by the re-sulting micro-modulus and the energy constraint helps in constructing thePD model with better accuracy using less data.In the interest of simplicity, this work only focuses on the application ofbond-based PD models and an energy constraint corresponding to isotropicdeformation, which reduces the regression algorithm to simple linear regres-sion. Future work should include extensions of the algorithm to more com-plex state-based PD models and constructing additional mathematically-and physically-justified constraints, e.g. constraints that can guarantee thewell-posedness of resulting PD model and/or constraints that contain moredetailed information of the medium’s micro-structure.
This work was partially supported by the Sandia National Laboratories(SNL) Laboratory-directed Research and Development program and by theU.S. Department of Energy, Office of Advanced Scientific Computing Re-search under the Collaboratory on Mathematics and Physics-Informed Learn-ing Machines for Multiscale and Multiphysics Problems (PhILMs) project.SNL is a multimission laboratory managed and operated by National Tech-21ology and Engineering Solutions of Sandia, LLC., a wholly owned sub-sidiary of Honeywell International, Inc., for the U.S. Department of Energy’sNational Nuclear Security Administration under contract DE-NA-0003525.This paper (SAND2021-0028) describes objective technical results and anal-ysis. Any subjective views or opinions that might be expressed in this paperdo not necessarily represent the views of the U.S. Department of Energy orthe United States Government.
References
B. Aksoylu and G. A. Gazonas. On the choice of kernel function in non-local wave propagation.
Journal of Peridynamics and Nonlocal Mod-eling , 2(4):379–400, 2020. doi: 10.1007/s42102-020-00034-x. URL https://doi.org/10.1007/s42102-020-00034-x .A. H. Delgoshaie, D. W. Meyer, P. Jenny, and H. A. Tchelepi. Non-localformulation for multiscale flow in porous media.
Journal of Hydrology ,531:649–654, 2015.M. D’Elia and M. Gunzburger. Identification of the diffusion parameter innonlocal steady diffusion problems.
Applied Mathematics and Optimiza-tion , 73:227–249, 2016.M. D’Elia, M. Gunzburger, and C. Vollman. A cookbook for finite elementmethods for nonlocal problems, including quadrature rule choices and theuse of approximate neighborhoods. arXiv:2005.10775, 2020.Q. Du, M. Gunzburger, R. Lehoucq, and K. Zhou. Analysis and approx-imation of nonlocal diffusion problems with volume constraints.
SIAMReview , 54(4):667–696, 2012.D. Huang, G. Lu, and P. Qiao. An improved peridynamic approach forquasi-static elastic deformation and brittle fracture analysis.
InternationalJournal of Mechanical Sciences , 94:111–122, 2015.P. Seleson and M. Parks. On the role of the influence function in the peri-dynamic theory.
International Journal of Multiscale Computational En-gineering , 9(6):689–706, 2011.S. A. Silling. Reformulation of elasticity theory for discontinuities and long-range forces.
Journal of the Mechanics and Physics of Solids , 48(1):175–209, 2000.S. A. Silling. Origin and effect of nonlocality in a composite.
Journal ofMechanics of Materials and Structures , 9(2):245–258, 2014.22. A. Silling and E. Askari. A meshfree method based on the peridynamicmodel of solid mechanics.
Computers & structures , 83(17):1526–1535,2005.S. A. Silling and R. B. Lehoucq. Convergence of peridynamics to classicalelasticity theory.
Journal of Elasticity , 93(1):13, 2008.S. A. Silling, M. Zimmermann, and R. Abeyaratne. Deformation of a peri-dynamic bar.
Journal of Elasticity , 73(1-3):173–190, 2003.S. A. Silling, M. Epton, O. Weckner, J. Xu, and E. Askari. Peridynamicstates and constitutive modeling.
Journal of Elasticity , 88(2):151–184,2007.A. Sridhar, V. G. Kouznetsova, and M. G. Geers. A general multiscaleframework for the emergent effective elastodynamics of metamaterials.
Journal of the Mechanics and Physics of Solids , 111:414–433, 2018.J. Trageser and P. Seleson. Bond-based peridynamics: A tale of two poisson’sratios.
Journal of Peridynamics and Nonlocal Modeling , 2(3):278–288,2020.Y. Wang, X. Zhou, Y. Wang, and Y. Shou. A 3-d conjugated bond-pair-based peridynamic formulation for initiation and propagation of cracks inbrittle solids.
International Journal of Solids and Structures , 134:89–115,2018.O. Weckner and R. Abeyaratne. The effect of long-range forces on thedynamics of a bar.
Journal of the Mechanics and Physics of Solids , 53(3):705–728, 2005.O. Weckner and S. A. Silling. Determination of nonlocal constitutive equa-tions from phonon dispersion relations.
International Journal for Multi-scale Computational Engineering , 9(6), 2011.R. A. Wildman. Discrete micromodulus functions for reducing wave disper-sion in linearized peridynamics.
Journal of Peridynamics and NonlocalModeling , (1):56–73, 2019.X. Xu and J. T. Foster. Deriving peridynamic influence functions for one-dimensional elastic materials with periodic microstructure. arXiv preprintarXiv:2003.05520 , 2020.H. You, Y. Yu, S. Silling, and M. D’Elia. Data-driven learning ofnonlocal models: from high-fidelity simulations to constitutive laws.arXiv:2012.04157, 2020a. 23. You, Y. Yu, N. Trask, M. Gulian, and M. D’Elia. Data-driven learning ofrobust nonlocal physics from high-fidelity synthetic data. arXiv preprintarXiv:2005.10076 , 2020b.G. Zheng, G. Shen, Y. Xia, and P. Hu. A bond-based peridynamic modelconsidering effects of particle rotation and shear influence coefficient.