A Physics-Guided Neural Network Framework for Elastic Plates: Comparison of Governing Equations-Based and Energy-Based Approaches
AA P
HYSICS -G UIDED N EURAL N ETWORK F RAMEWORK FOR E LASTIC P LATES : C
OMPARISON OF G OVERNING E QUATIONS -B ASED AND E NERGY -B ASED A PPROACHES
A P
REPRINT
Wei Li a Martin Z. Bazant b, c
Juner Zhu a, b, ∗ a Department of Mechancial Engineering, Massachusetts Institute of Technology b Department of Chemical Engineering, Massachusetts Institute of Technology c Department of Mathematics, Massachusetts Institute of TechnologyOctober 15, 2020 A BSTRACT
One of the obstacles hindering the scaling-up of the initial successes of machine learning in practicalengineering applications is the dependence of the accuracy on the size of the database that “drives” thealgorithms. Incorporating the already-known physical laws into the training process can significantlyreduce the size of the required database. In this study, we establish a neural network-based compu-tational framework to characterize the finite deformation of elastic plates, which in classic theoriesis described by the Föppl–von Kármán (FvK) equations with a set of boundary conditions (BCs).A neural network is constructed by taking the spatial coordinates as the input and the displacementfield as the output to approximate the exact solution of the FvK equations. The physical information(PDEs, BCs, and potential energies) is then incorporated into the loss function, and a pseudo datasetis sampled without knowing the exact solution to finally train the neural network. The predictionaccuracy of the modeling framework is carefully examined by applying it to four different loadingcases: in-plane tension with non-uniformly distributed stretching forces, in-plane central-hole tension,out-of-plane deflection, and buckling under compression. Two ways of formulating the loss functionare compared, one based on the PDEs and BCs, and the other based on the total potential energyof the plate. Through the comparison with the finite element simulation results, it is found thatour computational framework is capable of characterizing the elastic deformation of plates with asatisfactory accuracy. Compared with incorporating the PDEs and BCs in the loss, using the totalpotential energy is a better way in terms of training accuracy and efficiency.
Keywords
Physics-informed neural network · structural mechanics · elastic plates · Ritz method
In the past half-decade, machine learning enjoyed vast researches and achieved remarkable successes in a wide spectrumof scientific problems, including image processing [1, 2], cognitive science [3], genomics [4], drug discovery [5], andmaterial designing [6], to name a few. It has shown prominent advantages over other methods in effectively handlingcomplex natural systems with a daunting number of variables. Recently, we are witnessing a growing number ofinitial successes of machine learning, especially deep learning, in modeling complex engineered systems with a highdimensionality (the number of variables and degrees of freedom) [7, 8], for example, predicting the remaining useful lifeof a battery cell based on its partial life-cycle data [9]. In most cases, machine learning algorithms serve as a data-drivenapproach. It has been proven effective to predict the performance of a system even when the underlying physics has notbeen elucidated. However, like other statistical methods such as curve fitting and feature engineering, the accuracy of ∗ Corresponding author. Emails: [email protected] (W.L.), [email protected] (M.Z.B), [email protected] (J.Z.) a r X i v : . [ c s . C E ] O c t Physics-Guided Neural Network Framework for Elastic Plates: Comparison of Governing Equations-Based andEnergy-Based Approaches
A P
REPRINT machine learning methods highly depends on the quantity and quality of the dataset [10]. In the cases where little datais accessible or a large database is not affordable, for example at the microscales and nanoscales, machine learningmethods may lose their power, thus hindering the scaling-up of those initial successes. On the contrary, physics-basedor first-principle-based models are commonly less reliant on the size of the dataset because the governing physicallaws are elucidated by human brains beforehand and only a small amount of data is required to calibrate the unknownparameters.Bridging the gap between the data-driven approaches and the physics-based approaches creates a promising opportunityto develop novel computer methods that have the potential to unite the advantages of both approaches – characterizinghigh-dimensionality systems with a small dataset. One of these methods is often referred to as physic-guided data-drivenmethods, which aim to implement the already-known physics into the data-driven approach [8, 11, 12, 13, 14, 15, 16,17]. In this paper, we focus on machine learning algorithms, particularly artificial neural networks (ANNs). Generally,there are three key elements in a machine learning application: a dataset, a model, and a training process. In the vastopen literature, it is found that the term “physics-guided machine learning” (PGML) is being used in an unregularizedmanner. Here, we summarize the recent progress on this topic by classifying the existing studies according to the keyelements that they worked on.The first category of studies implements physics into machine learning algorithms by generating a dataset followingphysical laws [18, 19]. Instead of collecting data from the expensive and time-consuming experiments, attempts weremade to generate data through first-principle theories and physics-based simulations. In this way, the PGML algorithmcan learn from the known physics behind the man-made data. For example, Chen et al. [19] predicted the phonondensity of states of crystalline solids with unseen elements using a density functional perpetration theory-based phonondatabase to train a Euclidean neural network. Another example is that Li et al. [18] generated a large numerical datasetof lithium battery failure behavior under mechanical impact loading with a well-calibrated physics-based detailed modeland trained various machine learning algorithms to get the safety envelope of the cell.The second category of PGML studies implements physical laws by designing a physics-guided machine learningmodel [8, 20, 21]. Various types of models have been used and researched for machine learning systems for a targetregression or classification problem, for example, ANNs and supported vector machines (SVMs). Many existing studiestreated the model as a black box and empirically choose a set of parameters of the model, such as the number of nodesand hidden layers of ANNs. But it has become a consensus that understanding the physics of the problem can guidethe design of the model. E’s research team [8, 20] is a clear pioneer in this aspect. In one of their recent studies [8], aneural network was designed with several subnetworks representing the solution at different time instances to solve thehigh-dimensional partial differential equations (PDEs). Several important physical equations were successfully solvedwith their algorithm, including the nonlinear Black–Scholes equation with default risk, the Hamilton–Jacobi–Bellmanequation, and the Allen–Cahn Equation.The last category of PGML studies imposes physics into the training process of the algorithm. A typical exampleis the “physics-informed neural network” (PINN) proposed by Raissi et al. [14] PINN introduces the PDEs and theassociated initial and boundary conditions (ICs and BCs) into the loss function and solves the PDEs by minimizing it.The authors successfully applied their PINN approach to solving the Burger’s equation and the Navier-Stokes equation,which are the governing equations of a variety of flows. The same idea was adopted by Lu et al. [22] to solve a seriesof other equations such as the diffusion equation, and a library of open-source codes for solving different PDEs wascreated by the authors, named “DeepXDE”. Besides PINN, E et al. [13] proposed a deep learning Ritz method wherethe loss function is defined by the energy functional corresponding to the PDEs. It is also found that this type ofloss-function-based optimization algorithms can not only predict the performance of a system by solving the governingequations but also identify the unknown parameters in a physical law through an inverse process. Raissi et al [14]showed preliminary successes to discover the unknown parameters in the Burger’s equation using PINN. In fact, thistype of general applciations of data-driven methods does not rely on machine learning algorithms. Zhao et al. [23] usedan inverse approach to learn pattern-forming equations such as Cahn-Hilliard and Allen-Cahn from image data. Zhao’sapproach turned out to be still effective even with a very small set of images. This success was recently extended byEffendy et al [24] to analyze and design the electrochemical impedance spectroscopy of energy storage systems.Although we classify here the existing PGML studies in open literature into three categories, it is worth noting thatthere is no absolute boundary between them. It is possible to combine these approaches in one algorithm, and with therapid development of machine learning technologies, we are already witnessing an increasing number of advancedalgorithms that have all the above three merits [25, 26].There is a recent trend that the success of machine learning algorithms, which was initially achieved in the fieldslike fluid dynamics and mass and heat transfer, is now being extended into the field of applied mechanics of solids.Haghighat et al. [27] applied Raissi’s PINN framework to predicting the mechanical responses of linear elastic materials,and their predictions agreed well with the finite element simulations. Wu et al. [28] designed a recurrent neural2 Physics-Guided Neural Network Framework for Elastic Plates: Comparison of Governing Equations-Based andEnergy-Based Approaches
A P
REPRINT network-accelerated multiscale model to describe the elasto-plastic behavior of heterogeneous media subjected torandom cyclic and non-proportional loading paths. A recent study by Samaniego et al. [16] adopted an energy approach(variational method) to solve the PDEs in solid mechanics with machine learning and shows a high prediction accuracyfor the given examples. Huang et al. [29] developed a machine learning-based plasticity model that can effectivelypredict the behaviors of history-dependent materials. In our opinion, one of the challenges for the machine learningapplications in predicting the mechanical responses of solids is that most of the variables (such as stress and strain)are highly tensorial, or multi-axial. As a result, each direction has its own PDEs, leading to a large total number ofequations as well as BCs and ICs. Therefore, although the above studies all compared their predictions with othernumerical methods and showed satisfactory agreements, the most fundamental questions such as how the loss functionshould be formulated are still unsettled.The purpose of the present study is to develop a neural network framework for predicting the mechanical responsesof elastic plates using the third category of PGML strategy (implementing physics into the training process). Ourwork will be distinguished from the aforementioned existing publications in a number of ways. First, we will carefullyinvestigate the reliability of the neural network framework by solving the high-order and highly non-linear equations. Inthe classic theory, the governing equation of elastic plates are the well-known Föppl-von Kármán equations equations,which consists of two second-order PDEs and one fourth-order PDE. Second, two different approaches will be used toconstruct the loss function, one based on PDEs and BCs, and the other based on the total potential energy of the wholestructure. Third, the proposed computational framework will be applied to four different loading cases of the platesfor evaluation. Non-linearities that stem from the loading condition and the geometry of the plates will be purposelyintroduced to push the developed numerical framework to its limit. The paper will start with a brief introduction ofthe classic Kirchhoff plate theory. The theory of the neural network framework is then presented, together with acomparison with the conventional purely date-driven machine learning framework. At last, the four exemplary loadingconditions will be investigated, and some key features of the framework will be discussed.
The governing equations of an elastic plate can be obtained in two ways. One is using the condition for adjacentequilibrium, and the other is following the energy method. The former applies Newton’s laws on a control volume of thestructure by balancing all the applied and reactive forces and moments in each direction. This is quite straightforwardfor simple mechanical systems where the forces and displacements are easy to be described. For complicated systemsand structures, however, using the energy method is usually more convenient. It calculates the total potential energyas a functional and establishes the equilibrium equations (Euler-Lagrange equations of the functional) based on itsvariations. From a theoretical point of view, these two ways of obtaining the governing equations are equivalent. Here,the energy method will be used.
For a given mechanical system, infinite possible configurations are able to satisfy the geometric constraints. However,only the one that also satisfies the equilibrium condition is the true configuration. The displacement field correspondingto the true configuration is the true displacement, and the virtual displacements represent all the possible configurationsconsistent with the geometric constraints. The amounts of the virtual work is then the work done by all the forces alongwith the virtual displacements. The virtual work done by the internal stress and external body or surface force arerespectively defined as internal virtual work δ U and external virtual work δ V . Among all the admissible configurations,the one corresponding to the equilibrium configuration makes the total virtual work δ W vanish. The principle of virtualdisplacement is then stated as δ Π = δ U + δ V = 0 . (1)where δ is the variation operator. With this equation, we will derive the governing equations of the plates in thefollowing by writing down the internal and external work and making use of the principles of the variation method. We adopt the classic plate theory that is based on the following three Kirchhoff hypotheses [30]: i) straight lines normalto the plate mid-surface remain normal after deformation and thus are also called transverse normals; ii) the transversenormals remain straight after defter deformation, and iii) the thickness of the plate does not change after deformation.3 Physics-Guided Neural Network Framework for Elastic Plates: Comparison of Governing Equations-Based andEnergy-Based Approaches
A P
REPRINT 𝑥𝑦 𝑧 𝑧 𝜕𝑤𝜕𝑥(𝑥, 𝑦, 𝑧) (𝑥, 𝑦, 0)𝜃 = 𝜕𝑤𝜕𝑥 𝑥𝑧 Transverse Normals 𝑢𝑢 𝑛 𝑟 𝑠 𝑞 𝑥, 𝑦 𝜎 𝑛𝑛 𝜎 𝑛𝑠 𝜎 𝑛𝑟 Figure 1: Illustration of the reference and deformed configurations of an elastic plate and the boundary conditionsFigure 1 shows a plate of thickness h in the Cartesian coordinate ( x, y, z ) . The x - y plane coincides with the geometricmid-plane of the plate and the z -direction is taken positive downward. Without loss of generality, the three componentsof displacement field along the x , y , and z axes are noted as u , u y , and w , respectively. Based on the Kirchhoffhypothesis, the strain components can be written as (refer to Appendix A.1 for a detailed derivation), ε αβ = ε αβ + z · κ αβ , (2)where ε αβ ( α, β = 1 , are the membrane strains that represent the in-plane deformation, and κ αβ ( α, β = 1 , arethe curvatures (often known as the bending strains) that comes from the transverse bending, ε αβ = 12 ( u α,β + u β,α + w ,α w ,β ) , (3) κ αβ = − w ,αβ , (4)where the comma notation “,” indicates a derivative (for example, w ,α = ∂w/∂α ), u α ( α = x, y ) are the displacementson the mid-plane ( u α ( x α ) = u α ( x α , ). w ,α ( α = x, y ) are respectively the rotation angles of the transverse normalalong x and y axes. The internal virtual work (variation of internal energy) of a plate is defined as δ U = (cid:90) Ω (cid:90) h − h (cid:0) N αβ δ ε αβ + M αβ δ κ αβ (cid:1) d z d x d y, (5)where N αβ and M αβ are respectively the thickness-integrated forces and moments (also known as the membrane forcesand bending moments), N αβ = (cid:90) h − h σ αβ d z, (6) M αβ = (cid:90) h − h σ αβ z d z. (7)4 Physics-Guided Neural Network Framework for Elastic Plates: Comparison of Governing Equations-Based andEnergy-Based Approaches A P
REPRINT
We consider a distributed transverse ( z -direction) pressure q t on the top surface and a set of mixed traction-displacementboundary conditions. For generality, we define a local coordinate ( n, s, z ) at any point on the boundary where n and s are corresponding to the normal and tangential directions. The applied in-plane traction can thereby be described asthe normal stress ˆ σ nn , tangential stress ˆ σ ns , and the transverse shear stress ˆ σ nz . Hence, the external virtual work can becalculated as V = − (cid:90) Ω q t δ w d x d y − (cid:90) Γ σ (cid:16) (cid:98) N nn δ u n − (cid:99) M nn δ w ,s + (cid:98) N ns δ u s − (cid:99) M ns δ w ,s + (cid:98) N nz δ w (cid:17) d s, (8)where u n and u s are respectively the displacements along the boundary normal and tangential direction, u n and u s are respectively the corresponding displacements at the mid-plane, w ,n and w ,s are respectively the rotation angels ofthe transverse normal along the boundary normal and tangential directions, and the applied thickness-integrated forcesand moments are defined in the same way as in Eq. (6) and Eq. (7).To express the variation of internal energy in terms of virtual displacements, integration by parts need to be performedseveral times (details can be found in Appendix A.2). According to the principle of virtual displacements and rearrangingthe coefficients, we finally have δ U + δ V = − (cid:90) Ω (cid:110) N αβ,β δ u α + (cid:104) ( N αβ w ,β ) ,α + M αα,ββ − q t (cid:105) δ w (cid:111) d x d y + (cid:90) Γ σ (cid:104) N αβ n β δ u α + ( N αβ w ,β n α + M αβ,β n α ) δ w − M αβ n β δ w ,α + (cid:98) N nn δ u n − (cid:99) M nn ∂ δ w∂n + (cid:98) N ns δ u s − (cid:99) M ns ∂ δ w∂s + (cid:98) N nz δ w (cid:105) d s. (9)The right-hand-side of Eq. (9) consists of two integration terms. The first is over the whole domain Ω , and the second isalong the traction-based BCs Γ σ . The principle of virtual displacement δ Π = 0 implies that any small changes in thedisplacement field should not change the totally potential energy. Mathematically, both integration terms should bezero. Considering that the first term is zero, the coefficients of δ u α ( α = x, y ), and δ w must all be zero. We can get thefollowing governing equations (i.e. the Euler-Lagrange equations of the total potential energy functional), P α ≡ N αβ,β = 0 , P z ≡ ( N αβ w ,β ) ,α + M αβ,αβ − q t = 0 . (10)Eq. (10) are the general form of the governing equations of plates. Its linear elasticity special case is the well-knownFöppl–von Kármán (FvK) equations, named after August Föppl [31] and Theodore von Kármán [32]. To derive them,the constitutive equations of an isotropic elastic plate should be established, N αβ = C (cid:2) (1 − ν ) ε αβ + νε γγ δ αβ (cid:3) , (11) M αβ = D [(1 − ν ) κ αβ + νκ γγ δ αβ ] , (12)where δ αβ is the Kronecker delta ( δ αβ = 1 , if α = β and δ αβ = 0 , if α (cid:54) = β ), C is the stretching stiffness, or axialrigidity, and D is the bending stiffness, or flexural rigidity, C = Eh − ν , D = Eh − ν ) . (13)Here E is the Young’s modulus and ν is the Poisson’s ratio. Substituting Eq. (2), Eq. (3), and Eq. (11) into Eq. (10),we can get the FvK equations for the isotropic elastic plate in terms of displacement. In Eq. (10), the two in-plane equations have second-order spatial derivatives and the out-of-plane equation has thefourth-order spatial derivatives. Therefore, a total of eight boundary conditions are required. Like the derivation using5 Physics-Guided Neural Network Framework for Elastic Plates: Comparison of Governing Equations-Based andEnergy-Based Approaches
A P
REPRINT adjacent equilibrium, it is possible to identify all these BCs on each edge of the plate. But the advantage of usingthe energy method is that it can give not only the governing PDEs but also the complete description of the boundaryconditions that leads to a unique solution. By setting the term of the integration along the stress boundary (the secondintegration along Γ σ ) in Eq. (9) to zero, we can obtain the boundary conditions. The quantities with a variationare referred to as the primary variables that constitute the geometric boundary conditions and the coefficients of thevariations are referred to as the secondary variables that constitute the natural boundary conditions. We can see thereare five primary variables u x , u y , w, w ,x , and w ,y for a plate with the edges aligned with the x and y axis, whichindicates a total of ten boundary conditions (five geometric and five natural boundary conditions). This is seeminglyinconsistent with the eight boundary conditions from the order analysis of the PDEs. The reason is that there are onlyfour independent primary variables among the aforementioned five variables. Only the rotation about the normal axisis considered in the plate theory. Interested readers are referred to Appendix A.3 for the detailed derivation. After atransformation from the global Cartesian coordinate to the local coordinate of ( n, s, z ) , the second integration termbecomes (cid:90) Γ σ (cid:104)(cid:16) N nn − (cid:98) N nn (cid:17) δ u n + (cid:16) N ns − (cid:98) N ns (cid:17) δ u s + (cid:16) V n − (cid:98) V n (cid:17) δ w + (cid:16) M nn − (cid:99) M nn (cid:17) δ w ,n (cid:105) d s, (14)where V n = ( N xx w ,x + N xy w ,y ) n x + ( N yy w ,y + N xy w ,x ) n y + N xx,x n x + M yy,y n y + M xy,x n y + M xy,y n x + M ns,s . (15)It is then clear that the four primary variables are u n , u s , w , and w ,n , respectively corresponding to the in-planedisplacement in the normal direction, the in-plane displacement in the tangential direction, the out-of-plane deflection,and the rotation along the normal axis, and the four secondary variables are N nn , N ns , V n , and M nn , respectivelycorresponding to the in-plane normal force, the in-plane tangential force, the shear force, and the bending moment. The FvK equations are notoriously difficult to solve. Existing successful methods include semi-analytical solutions[33, 34, 35] and numerical simulations [34, 36, 37]. The former relies on making reasonable simplifications to makethe equations solvable, and the latter commonly makes use of finite element (FE) methods. It should be noted thatboth methods can only approximate the exact solution instead of directly solving the FvK equations. In this section,we explore the third possible method by developing a learning framework. The overall flowchart of the algorithm isillustrated in Figure 2. In general, the algorithm consists of three parts, an artificial neural network, a loss function, anda training dataset, which will be elaborated on in the following three sub-sections.
A fully connected neural network is constructed to approximate the exact solution of the displacement field. It consistsof the input layer, output layer, and hidden layers in between. Here, we take the spatial coordinates ( x, y ) as the inputsand the displacement fields ( u x , u y , w ) to be predicted as the outputs (see Figure 3a).Considering an L -layer neural network, or a ( L − -hidden layer neural network, with P k neurons in the k -th layer( P = 2 is the dimension of inputs and P L = 3 is the dimension of outputs), for the j -th neuron in the k -th layer,the output A kj is obtained by taking the weighted average of the outputs of the previous layer and then applying anactivation function, A kj = f P k − (cid:88) i =1 W ik A k − i + b kj , (16)where W ik is the weights and b kj is the bias. f ( · ) represent the nonlinear activation function. Some common choices (seeFig. 3b) are the rectified linear unit (ReLU, f ( x ) = max { x, } ), the logistic sigmoid function ( f ( x ) = 1 / (1 + e − x ) ),and the hyperbolic tangent (Tanh, f ( x ) = ( e x − e − x ) / ( e x + e − x ) ). The above equation can be written in vector andmatrix form, A k = f (cid:0) W T A k − + b k (cid:1) , (17)6 Physics-Guided Neural Network Framework for Elastic Plates: Comparison of Governing Equations-Based andEnergy-Based Approaches A P
REPRINT 𝒖 = 𝑁𝑁 𝒙; 𝜽
Neural Network Loss Function 𝜕𝑢 𝑥 /𝜕𝑥 ; 𝜕𝑢 𝑦 /𝜕𝑥 ; 𝜕𝑤/𝜕𝑥 𝜕𝑢 𝑥 /𝜕𝑦 ; 𝜕𝑢 𝑦 /𝜕𝑦 ; 𝜕𝑤/𝜕𝑦 … Deformation gradients
Training Data
ℒ = 𝑖=1𝑄
1𝑄 𝒖 𝑖 − ෝ 𝒖 𝑖
1. Purely data-driven Purely data-drivenLarge dataset [𝒙 𝑖 , ෝ𝒖 i ] observed in simulations or experiments Physics-guided
Large dataset [𝒙 𝑖 ] without outputs ℒ = 𝑖=1𝑄
1𝑄 𝑗=13 (𝒫 𝑗𝑖 ) + ℒ BCs
2. PDE-based 𝒫 𝑧 = 𝐷∆ 𝑤 − 𝛻 ∙ 𝑵𝛻𝑤 − 𝑞 𝑡 , 𝒫 𝑥,𝑦 = 𝛁 ∙ 𝑵,
3. Energy-based
ℒ = 𝑈 + 𝑉 + 𝑇
𝑈 = න Ω 𝜀 𝑖𝑗 𝜎 𝑖𝑗 ℎ𝑑𝐴, 𝑉 = − න Ω 𝑞 𝑡 𝑤 𝑑𝐴, Total Loss ℒ Minimize 𝜽 ∗ 𝑦 𝑢 𝑥 𝑢 𝑦 𝑥 𝑤 (Penalty term) Figure 2: Flow chart of the physic-guided machine learning framework compared with the purely data-driven approach 𝑦 𝑢 𝑥 𝑢 𝑦 𝑥 𝐴 𝑗𝑘 = 𝑓 𝑖=1𝑄 𝑘−1 𝑊 𝑖𝑗 𝐴 𝑖𝑘−1 − 𝑏 𝑗𝑘 Input layer 𝑏 𝑗𝑘 𝑊 𝐴 𝑗𝑘 𝑊 𝑖𝑗 𝑊 𝑄 𝑘−1 𝑗 𝐴 𝑖 𝑘−1 𝐴 𝑄 𝑘−1 𝑘−1 𝐴 Current neuronInputs from previous layerweights OutputActivation function01 𝑥 Sigmoid( 𝑥) 𝑥 ReLU( 𝑥) (a)(b) Hidden layer(s)Output layer 𝑤 bias 01 𝑥 Tanh( 𝑥) -1 Figure 3: Construction of the artificial neural network, (a) fully-connected multi-layer network, (b) activation functions7 Physics-Guided Neural Network Framework for Elastic Plates: Comparison of Governing Equations-Based andEnergy-Based Approaches
A P
REPRINT where W = [ W ij ] ∈ R M l − × M l and b k = (cid:2) b kj (cid:3) ∈ R M l are the weight matrix and bias vector, respectively. A k represents the output vector of k -th layer. The neural network can then be defined recursively as follows,input layer : A = [ x, y ] ∈ R , hidden layers : A k = f (cid:0) W T A k − + b k (cid:1) ∈ R M l , for ≤ k ≤ L − , output layer : A L = [ u, u y , w ] = W T A k − + b k ∈ R . (18)Note that the activation function is not applied for the output layer. It should also be pointed out that the lowercase w without subscripts denotes the out-of-plane displacement whereas the capital W and W ij represent the weightsof the neural network. The weights and biases are the parameters to be trained of the neural network and there is atotal of (cid:80) Li =1 P i − P i weights and (cid:80) Li =1 P i biases ( P = 2 , P L = 3 ). It has been proven that a neural network canapproximate any targeted functions with an arbitrary width or depth. This approximation ability of the neural networkmakes it possible to represent the full filed solution of the PDEs. The Loss function is the most essential part of a neural network algorithm because the already-known physical laws willbe implemented here. As illustrated in Figure 2, we consider three ways of formulating the loss function. The first ispurely data-driven and compares the prediction of the displacement field only with the data observed from experimentsor simulations. Therefore, it is not a physics-guided algorithm. We include it here for the comparison with the othertwo loss functions, which reflect the physical laws. The second one is defined on the PDEs and BCs, and the third oneis defined on the total potential energy of the plate. It is worth noting that the PDE-based is close to the concept ofPINN proposed by Raissi et al. [14] and the energy-based can be viewed as an implementation of the Ritz method [13].
Purely data-driven
Like most conventional data-driven methods, the loss function of a regression problem can be constructed by thedifference between the predictions and the true experimental or numerical observations, such as the mean square error, L Data-driven = Q (cid:88) i =1 Q (cid:104)(cid:0) u ix − ˆ u ix (cid:1) + (cid:0) u iy − ˆ u iy (cid:1) + (cid:0) w i − ˆ w i (cid:1) (cid:105) , (19)where the superscript i indicates the i -th training sample, ˆ u ix , ˆ u iy , and ˆ w i are the observed displacements in the trainingdataset, and u ix , u iy , and w i are the predicted displacements, Q is the total number of training samples, which should besufficient large to ensure an acceptable accuracy. PDE-based
For the studied plate theory, the outputs of the neural network should satisfy the governing PDEs (Eq. (10)) and theBCs (Eq. (14)). The first way to implement the plate theory into the algorithm is, therefore, to construct the loss withthe residuals of the PDEs and BCs, L PDE-based = L PDEs + λ s L BCs + λ d L BCd , (20)where L PDEs is the residual of the PDEs, L PDEs = Q P (cid:88) i =1 Q P (cid:104)(cid:0) P ix (cid:1) + (cid:0) P iy (cid:1) + (cid:0) P iz (cid:1) (cid:105) , (21) P x,y,z are the residual values of the three governing PDEs defined in Eq. (10), Q P is the number of training sampleswithin solution domain. L BCs is the residual of the stress BCs, L BCs = Q bs (cid:88) i =1 Q bs (cid:20)(cid:16) N inn − (cid:98) N inn (cid:17) + (cid:16) N ins − (cid:98) N ins (cid:17) + (cid:16) V in − ˆ V in (cid:17) + (cid:16) M inn − (cid:99) M inn (cid:17) (cid:21) , (22)8 Physics-Guided Neural Network Framework for Elastic Plates: Comparison of Governing Equations-Based andEnergy-Based Approaches A P
REPRINT L BCd is the residual of the displacement BCs, L BCd = Q bd (cid:88) i =1 Q bd (cid:104)(cid:0) u i n − ˆ u i n (cid:1) + (cid:0) u i s − ˆ u i s (cid:1) + (cid:0) w i − ˆ w i (cid:1) + (cid:0) w i,n − ˆ w i,n (cid:1) (cid:105) . (23) λ s and λ d are the weights of loss on stress boundary and displacement boundary, the hat notation indicates the prescribedsecondary and primary variables defined in Eq. (14) at the boundaries, and the corresponding variables without a hatare the predicted outputs of the neural network, Q bs and Q bd are the number of samples at the stress boundary anddisplacement boundary, respectively.By minimizing the total loss, the PDEs and boundary conditions are satisfied. Mathematically, the neural network isable to approximate the exact solution. Energy-based
The second way to implement the plate theory into the algorithm is to directly take the total potential energy as theloss and minimize the loss according to the principle of minimum potential energy. However, as E et al. [13] pointedout in their study on the Ritz method, the challenging issue is how to incorporate the displacement boundaries into thetotal potential energy because it is not automatically included. Here, we construct the loss function based on the totalpotential energy with a penalty energy term [13, 33, 34], L Energy-based = Π = U + V + T, (24)where U is the internal energy, U = (cid:90) Ω (cid:0) N αβ ε αβ + M αβ κ αβ (cid:1) d x d y. (25) V is the virtual work done by external forces on the boundary (defined as negative to be consistent with the principlevirtual displacement) V = (cid:90) Γ σ (cid:104) − (cid:98) N nn u n − (cid:98) N ns u s − (cid:98) N nz w + (cid:99) M ns w ,s + (cid:99) M nn w ,n (cid:105) d s, (26)and T is the penalty term that enforces the displacement boundary conditions, T = ε ∗ n (cid:90) Γ d | u n − ˆ u n | d s + ε ∗ s (cid:90) Γ d | u s − ˆ u s | d s + ε ∗ w (cid:90) Γ d | w − ˆ w | d s + ε ∗ w,n (cid:90) Γ d | w ,n − ˆ w ,n | d s. (27)Here, ε ∗ n , ε ∗ s , ε ∗ w , and ε ∗ w,n are four Lagrangian multipliers that are in the dimension of force ( ε ∗ n , ε ∗ s , ε ∗ w ) or moment( ε ∗ w,n ), representing the applied forces and moments on the displacement boundaries. Note that if the displacementBCs are satisfied, this penalty term vanishes. In addition, the first variation of this penalty term is zero ( δ T = 0 ), andtherefore, taking the variation of the total potential energy functional Π still returns δ U + δ V . This is consistent withEq. (1).The energy integrations can be evaluated numerically. For integration in the domain, we first uniformly sample Q P points, the integration can then be approximated by, U num = Q P (cid:88) i =1 U i δ A i ≈ Q P (cid:88) i =1 Q P U i A t , (28)where U is internal energy density U = 12 (cid:0) N αβ ε αβ + M αβ κ αβ (cid:1) . (29) U i is, therefore, the internal energy density of the i -th sample point, δ A i is the discrete area around the sampled point,which can be estimated by A t /Q P is the total area when Q P is sufficiently large.9 Physics-Guided Neural Network Framework for Elastic Plates: Comparison of Governing Equations-Based andEnergy-Based Approaches A P
REPRINT
Similarly, V num = Q bs (cid:88) i =1 V i δ l i ≈ Q bs (cid:88) i =1 Q bs V i l t , (30)where V is the external work per unit width, V = − (cid:98) N nn u n − (cid:98) N nz w − (cid:98) N ns u s + (cid:99) M ns w s + (cid:99) M nn w n . (31) δ l i is the discrete length around the sampled point, approximately l t /Q BCs , where l t is the total length of the stressboundaries. T can be numerically interpreted in the same way as L BCd in Eq. (28), T num = λ d Q bd (cid:88) i =1 Q bd (cid:104)(cid:0) u i n − ˆ u i n (cid:1) + (cid:0) u i s − ˆ u i s (cid:1) + (cid:0) w i − ˆ w i (cid:1) + (cid:0) w i,n − ˆ w i,n (cid:1) (cid:105) / . (32)It should be noted that both PDE-based and energy-based loss functions require the calculation of the partial derivativesof the outputs with respect to the inputs. Most existing machine learning platforms (e.g. Tensorflow [38] and Pytorch[39]) are already equipped with default gradient algorithms, and users can obtain the numerical gradient efficiently.Higher-order partial derivatives can also be calculated by the algorithms but take more computation time. An alternativestrategy is to introduce the derivatives into the outputs of the neural network. For example, if we set the outputs to be ( u x , u y , u x,x , u x,y , u y,x , u y,y ) , the second-order derivatives of u x and u y can be obtained by performing only the firstderivative of the outputs. In this way, we can avoid calculating high-order derivatives. The disadvantage is that it willincrease the size of the neural network and extra constraints are needed to enforce the mathematical relations among theoutputs. For instance, in the above example, the derivative of the first output u x to the first input x should be equal tothe third output u x,x . For simplicity, we adopted the default gradient algorithm in Pytorch to ensure higher accuracy ofthe calculation of the derivatives at the sacrifice of computational efficiency. As indicated by Eq. (19), the purely data-driven model loss function can be minimized only when the displacementfield (ˆ u x , ˆ u y , ˆ w ) of a sufficiently large number of sample points can be observed. It means a large experimental ornumerical database. The exact solution of (ˆ u x , ˆ u y , ˆ w ) is always preferable rather than numerical simulation results for areliable training process. Experimental data is good but usually comes with measurement uncertainties. In addition, aswe already pointed out, high-order equations such as the FvK equations are difficult to be solved analytically. Therefore,how to obtain a satisfactory dataset to train the data-driven algorithm is a big challenge even for a simple structure likean elastic plate.The training data of the physics-guided algorithms has two important features: 1) exact solutions are not required (theloss function only contains the predicted outputs), 2) the training data can be sampled both within the solution domainand at the boundaries. Theoretically, the training dataset can be sampled in any size and strategy. Due to the numericalintegration (Eq. (28) and Eq. (30)) in the energy-based loss function, we first randomly sample the data points in auniform distribution and then resample the data during each training epoch to have more accurate and consistent results.For the PDE-based loss function, a large uniform dataset was sampled. The model was trained in small batches andthere was no need to resample the dataset for each training epoch. To validate the machine learning framework that we developed, four typical loading conditions will be characterized,and the prediction will be compared with the FE simulation results.
We start with a two-dimensional case that only involves in-plane deformation under plane-stress condition, as illustratedin Figure 4a. The two lateral edges of a 20 mm ×
20 mm ( l × l ) square elastic plate are under a non-uniform stretchingforce that follows a sinuous distribution p = sin ( yπ/l ) MPa. The other two edges (upper and lower) are traction-free.10 Physics-Guided Neural Network Framework for Elastic Plates: Comparison of Governing Equations-Based andEnergy-Based Approaches
A P
REPRINT 𝑢 𝑦 = 0𝑥𝑦 Traction condition 𝑁 𝑥𝑥 = 𝑝 ∙ ℎ, 𝑁 𝑥𝑦 = 0 𝑁 𝑦𝑦 = 0, 𝑁 𝑥𝑦 = 0 Free condition 𝑥𝑦𝑢 𝑥 = 0𝑝 = sin 𝑦 𝑙 𝜋 (a) (b) Traction condition 𝑁 𝑥𝑥 = 𝑝 ∙ ℎ, 𝑁 𝑥𝑦 = 0 𝑁 𝑦𝑦 = 0, 𝑁 𝑥𝑦 = 0 𝑛 𝑥 ∙ 𝑁 𝑥𝑥 + 𝑛 𝑦 ∙ 𝑁 𝑥𝑦 = 0𝑛 𝑦 ∙ 𝑁 𝑦𝑦 + 𝑛 𝑥 ∙ 𝑁 𝑥𝑦 = 0 Free condition Free condition 𝑥𝑦 𝑥𝑦𝑢 𝑦 = 0𝑢 𝑥 = 0𝑝 (c) (d) Figure 4: (a) Loading condition of non-uniform tension of plate and (b) one-quarter equivalent model;(c) loadingcondition of uniaxial central-hole tension and (d) one-quarter equivalent model.The origin of the Cartesian coordinate is placed at the center of the square with the x axis pointing to the right. TheYoung’s modulus of the elastic plate is 70 MPa and the Poisson’s ratio is 0.3. Due to the symmetry of the geometry andloads, only one quarter of the plate is modeled (see Figure 4b). The boundary conditions at the four edges of the quartermodel are listed as following, u x | x =0 = 0 ,u y | y =0 = 0 ,N xx | x = l = p · h = sin (cid:16) yl π (cid:17) h, N xy | x = l = 0 ,N yy | y = l = 0 , N xy | y = l = 0 . (33)The two governing PDEs for this case can be written in terms of displacements, P x ≡ E − ν (cid:18) ∂ u x ∂x + 1 − ν ∂ u x ∂y + 1 + ν ∂ u y ∂x∂y (cid:19) = 0 , P y ≡ E − ν (cid:18) ∂ u y ∂y + 1 − ν ∂ u y ∂x + 1 + ν ∂ u x ∂x∂y (cid:19) = 0 . (34)Two different loss functions are defined for comparison according to Eq. (20) and Eq. (24), where we set λ s = λ d = 1 .The PDEs-based loss function is, L PDE-based = L PDEs + L BCsx + L BCsx + L BCd , (35)11 Physics-Guided Neural Network Framework for Elastic Plates: Comparison of Governing Equations-Based andEnergy-Based Approaches A P
REPRINT and L PDEs = 1 Q P Q P (cid:88) i =1 (cid:104)(cid:0) P ix (cid:1) + (cid:0) P iy (cid:1) (cid:105) , L BCsx = 1 Q sx Q sx (cid:88) i =1 (cid:26)(cid:104)(cid:0) N ixx (cid:1) − p · h (cid:105) + (cid:0) N ixy (cid:1) (cid:27) , L BCsy = 1 Q sy Q sy (cid:88) i =1 (cid:104)(cid:0) N ixx (cid:1) + (cid:0) N ixy (cid:1) (cid:105) , L BCd = 1 Q dx Q dx (cid:88) i =1 (cid:0) u ix (cid:1) + 1 Q dy Q dy (cid:88) i =1 (cid:0) u iy (cid:1) , (36)where variables with superscript i indicate value evaluated for the i -th training sample, Q P , Q sx , Q sy Q dx , and Q dy arethe total number of training samples within the domain, on the right and upper edges (stress boundary), and on the leftand lower edges (displacement boundary), respectively. The membrane forces are computed with the displacementgradients according to Eq. 11.The energy-based loss function is, L Energy-based = U + V + T = (cid:90) l/ (cid:90) l/ Eh (cid:0) − u y (cid:1) (cid:34)(cid:18) ∂u x ∂x (cid:19) + (cid:18) ∂u y ∂y (cid:19) + 2 ν ∂u x ∂x ∂u y ∂y + 1 − ν (cid:18) ∂u x ∂y + ∂u y ∂x (cid:19) (cid:35) d x d y − (cid:90) l [ N xx ( y ) · u x ( y )] x = l/ d y + (cid:90) l/ (cid:2) u x ( y ) (cid:3) x =0 h d y + (cid:90) l/ (cid:2) u y ( x ) (cid:3) y =0 h d x. (37)The fully connected neural network defined in Eq. (18) is used, where the outputs are replaced by 2D displacementfields ( u x , u y ) . We modified the outputs in the following way, (cid:2) u (cid:48) x , u (cid:48) y (cid:3) = N ( x, y ) u x = u (cid:48) x · x,u y = u (cid:48) y · y, (38)where u (cid:48) x and u (cid:48) y are the outputs of neural network N ( x, y ) , u x and u y are the modified outputs. The displacementboundary conditions can then be satisfied ( u x | x =0 = 0 , u y | y =0 = 0 ), which means that the loss term L BCd is alwayszero. This can simplify the calculation of loss function ( T ≡ is automatically satisfied).At the same time, we solved the problem using finite element method with an extremely fine mesh size of 0.1 mm(10,000 elements in total) in Abaqus/standard. Since this is a simple mechanical problem, the results were regarded asthe extract solution to evaluate the accuracy of the machine learning methods.First, we performed the purely data-driven machine learning method to be compared with our physics-guided methods.The training data were extracted from the FE simulation results with the spatial positions ( x, y ) as the inputs andthe displacement fields ( ˆ u x , ˆ u y ) as the outputs. 10,000 samples were used for training, with a learning rate of 0.001and a batch size of 64. The hyperbolic (tanh) function was used as the activation function. Figure 5a, b, and c showa comparison between the contour plots of the magnitude of the displacement obtained by FE simulation and thepredictions of the trained neural networks with two different sizes, 5-hidden layers with 5 neurons each layer and10-hidden layers with10 neurons each layer. The longitudinal membrane force fields were then obtained according to Eq.(11), as shown in Figure 5d, e, and f, respectively. The correlation coefficients between the predicted outputs and the FEresults were calculated to evaluate the accuracy of the global prediction of the displacement and membrane force fields(see Table 1). We can see that the 5-hidden layer network can well predict both the global displacement and membraneforce fields. A larger 10-hidden layer network can slightly improve the accuracy of the displacement prediction butthe accuracy of membrane force field prediction deceases significantly (see Figure 5e and f as well as the correlation12 Physics-Guided Neural Network Framework for Elastic Plates: Comparison of Governing Equations-Based andEnergy-Based Approaches A P
REPRINT
FE Simulation Purely data-driven u (mm) u (mm) u (mm) 𝑁 𝑥𝑥 (N ∙ mm) Longitudinal membrane force
Purely data-driven (a)(d) (b)(e) (c)(f)
Magnitude of the displacement 𝑁 𝑥𝑥 (N ∙ mm) 𝑁 𝑥𝑥 (N ∙ mm) Figure 5: Predicted displacement and longitudinal membrane force of FE simulation (a,d), and purely date-drivenmodels with different sizes: (b,e) 5-hidden layer (5 neurons each layer) network, (d,f) and 10-hidden layer (10 neuronseach layer) network, under non-uniform stretching load.Table 1: Correlation coefficients of the predicted displacement and longitudinal membrane force fields for non-uniformstretching Neural network Correlation coefficients u x u y N xx N yy N xy Data-driven (2 , × , (2 , × , (2 , × , (2 , × , (2 , × , A P
REPRINT
Figure 6: Predicted displacement (a) and longitudinal membrane force (b) fields of FE simulation and purely data-drivenmodels under non-uniform stretching load
PDE-based Model Energy-based Model (a)(b) (c)(d) (e)(f) u (mm) u (mm) FE Simulation
Magnitude of the displacement u (mm) 𝑁 𝑥𝑥 (N ∙ mm) 𝑁 𝑥𝑥 (N ∙ mm) 𝑁 𝑥𝑥 (N ∙ mm) Figure 7: Predicted displacement and longitudinal membrane force of FE simulation (a,b), PDE-based model (c,d), andenergy-based model (e,f) under non-uniform stretching load.14 Physics-Guided Neural Network Framework for Elastic Plates: Comparison of Governing Equations-Based andEnergy-Based Approaches
A P
REPRINT
Figure 8: Predicted displacement (a) and longitudinal membrane force (b) of FE simulation and neural network modelsalong the diagonal line of the plate under non-uniform stretching load
The second application still focuses on the in-plane loading, but we will carefully investigate the effect of geometricnonlinearity by introducing a hole at the center of the specimen. Thus, the governing equations are the same as Eq. (34).The central hole leads to the stress concentration phenomenon at the edge of the hole, as a tough loading case for themachine learning method to solve. The problem is illustrated in Figure 4c, the hole has a diameter of 5 mm in the centerof 20 mm ×
20 mm ( l × l ) square plate. A uniformly distributed stretching load (1 MPa) is applied at both ends. Theloading case is equivalent to the one-quarter model shown in Figure 4d and the boundary conditions for the 4 straightedges and 1 arch edge are listed as follows, u x | x =0 = 0 ,u y | y =0 = 0 ,N xx | x = l = p · h, N xy | x = l = 0 ,N yy | y = l = 0 , N xy | y = l = 0 ,N xx n x + N xy n y = 0 , N xy n x + N xx n y = 0 , ( at x + y = d / , (39)where n x and n y are the direction cosines of the boundary normal vector as defined in Eq. A.6.Similar to the previous case, we first performed a FE simulation as the reference and then trained a 5-hidden layerneural network (5 neurons each layer) with the displacement fields from FE simulation results as a purely data-drivencase. We then trained the two models with physics-guided loss functions, one based on PDEs and the other based onenergy.We found that all the three ANN models can well predict the field of displacement (see Figure 9). Besides, the stressconcentration phenomenon is captured by all the models, but the concentration factor is not. The predicted distributionsof the longitudinal membrane force are plotted in Figure 10. It is clear that the prediction by the energy-based modelis the closest to the FE simulation result, significantly better than the PDE-based. Regarding the purely data-drivenmodels, a larger 5-hidden layer network (10 neurons each layer) can slightly improve the accuracy (see Figure 10e).Further increasing the network to 10 hidden layers (10 neurons each) did not overfit the displacement field as it did inthe previous loading case (Figure 10e)). This implies a higher complexity due to the geometric nonlinearity and a largersize neural network is required to approximate the exact solution, therefore it is less likely to be overfitted.A quantitative comparison of the two models with physics-guided loss functions is performed by plotting the magnitudeof N xx along the edge of the central hole in the polar coordinate, as shown in Figure 11. At the same time, thecorrelation coefficients are summarized in Table 2. It becomes clearer that the energy-based model is the “winner”. Itshould also be pointed out that even the winner cannot perfectly agree with the FE result. However, for such a nonlinearproblem with stress concentration, we cannot fully trust the FE simulations as well. We expect a more persuasivecomparison with the experimental data in future studies. Here, to further validate the energy-based method, we applied15 Physics-Guided Neural Network Framework for Elastic Plates: Comparison of Governing Equations-Based andEnergy-Based Approaches A P
REPRINT
PDE-based Model Energy-based Model (a) (b) (c) u (mm) u (mm) FE Simulation
Magnitude of the displacement u (mm) (d) (e) (f) Purely data-driven
Purely data-driven
10 hidden layers - 10 neurons
Purely data-driven u (mm) u (mm) u (mm) Figure 9: Displacement field prediction of FE simulation (a), PDE-based model (b), energy-based model (c), and purelydata-driven models with three different size neural networks: (d) 5-hidden layer (5 neurons each layer) , (e) 5-hiddenlayer (10 neurons each layer), (f) 10-hidden layer (10 neurons each layer), for the central-hole tension.Table 2: Correlation coefficients of the predicted displacement and longitudinal membrane force fields for central-holetension Neural network Correlation coefficients u x u y N xx N yy N xy Data-driven (2 , × , (2 , × , (2 , × , (2 , × , (2 , × , l y /l x ) of the central hole increases. The energy-based modelcan correctly capture this trend and is also able to give reasonable predictions of the full stress field. But we can still seethe local deviations between the energy-based model and the FE simulation. The above two examples played as very strict comparisons of the purely data-driven approach and the two models withPDE-based and energy-based loss functions while using FE simulations as the reference. It has been demonstrated thatthe energy-based model can provide the predictions that are closest to the FE results. In the following two examples, wealso compared the three ANN methods, and we reached the same conclusion. Therefore, we will not cover the details ofthe comparison for conciseness. Instead, we will only focus on the energy-based method while still using FE simulationas the reference although it is not the exact solution.In this example, we consider the square plate deflection loading case under uniform out-of-plane pressure. Unlike theprevious two, this case will involve both in-plane and out-of-plane deformation. As illustrated in Figure 13a, a 10 Pauniform transverse distributed pressure is applied on a 100 mm ×
100 mm square plate whose four edges are clamped.16 Physics-Guided Neural Network Framework for Elastic Plates: Comparison of Governing Equations-Based andEnergy-Based Approaches
A P
REPRINT
PDE-based Model Energy-based Model (a)(d) (b)(e) (c)(f)
FE Simulation
Purely data-driven
Purely data-driven
10 hidden layers - 10 neurons
Purely data-driven
Longitudinal membrane force 𝑁 𝑥𝑥 (N ∙ mm) 𝑁 𝑥𝑥 (N ∙ mm) 𝑁 𝑥𝑥 (N ∙ mm) 𝑁 𝑥𝑥 (N ∙ mm) 𝑁 𝑥𝑥 (N ∙ mm) 𝑁 𝑥𝑥 (N ∙ mm) Figure 10: Longitudinal membrane force field prediction of FE simulation (a), PDE-based model (b), energy-basedmodel (c), and purely data-driven models with three different size neural networks: (d) 5-hidden layer (5 neurons eachlayer) , (e) 5-hidden layer (10 neurons each layer), (f) 10-hidden layer (10 neurons each layer), for the central-holetension Figure 11: Predicted longitudinal membrane force of neural network models around the central hole17 Physics-Guided Neural Network Framework for Elastic Plates: Comparison of Governing Equations-Based andEnergy-Based Approaches
A P
REPRINT
FE Simulation Energy-based Model (a)(b)(c) (d)(e)(f) (g)(h)(i) 𝑁 𝑥𝑥 (N ∙ mm) 𝑙 𝑥 𝑙 𝑦 𝑁 𝑥𝑥 (N ∙ mm) 𝑁 𝑥𝑥 (N ∙ mm) 𝑁 𝑥𝑥 (N ∙ mm) 𝑁 𝑥𝑥 (N ∙ mm) 𝑁 𝑥𝑥 (N ∙ mm) Figure 12: Comparison of FE simulation and energy-based neural network results for different shape and size central-hole tension: (a) l x = l y = 4 mm, (b) l x = 2 mm , l y = 4 mm, and (c) l x = 4 mm , l y = 2 mm. The membrane forcefields predicted by the FE simulation of the three cases are respectively shown in (d), (e), and (f). The predictions of theenergy-based model are respectively shown in (g), (h), and (i).The Young’s modulus and Poisson’s ratio of the plate is set as 70 MPa and 0.3, respectively. The governing PDEs arelisted in Eq. (10). The boundary conditions are u x = 0 , u y = 0 , w = 0 , ∂w∂n = 0 , ( at x = ± or y = ± . (40)FE simulations with an element size of 0.1 mm are performed in Abaqus/standard with shell element (4-node doublycurved thin shell element, with reduced integration, hourglass control, and finite membrane strains). A 5-hidden layerneural network (5 neurons each layer) is trained with the energy-based loss function and the comparison of the predicteddeflection field with the FE simulation is shown in Figure 13b-e. A quantitative comparison of the distribution along thecentral line is shown in Figure 13f. It is found that the energy-based algorithm can still provide a satisfactory predictionof the out-of-plane displacement. It should be noted that the small deviation between the energy-based method and theFE simulation cannot be fully ascribed to the computational error of the former. This is because our physics-guidedneural network framework for elastic plates is constructed to implement the classic plate theory based on Kirchhoffhypotheses. In the FE simulations, the governing equations of the shell elements are slightly different, which dependson the integration algorithms and element size. For a stricter comparison, we can either employ experimental data ordevelop finite element simulations with the same classic plate theory.18 Physics-Guided Neural Network Framework for Elastic Plates: Comparison of Governing Equations-Based andEnergy-Based Approaches A P
REPRINT
Figure 13: Loading case of plate deflection under out-of-plane pressure (a). The predicted 3D configuration and2D out-of-plane displacement field after deflection are shown in (b) and (d) for FE simulation and (c) and (d) forenergy-based model. The deflection of the central line is compared in (e)
In the last example, we investigate the buckling of the same plate as the third example under in-plane compressiveloads. Figure 14a and b show the two different boundary conditions that were studied. One has a simply-supported leftedge ( u x = 0 , w = 0 , N yy = 0 , M yy = 0 at x = − ), and the other clamped ( u x = 0 , w = 0 , N yy = 0 , ∂w/∂x = 0 at x = − ). In both cases, there is no out-of-plane load. Therefore, trivial solutions that only involve the in-planedeformation (i.e. no out-of-plane deflection, w = 0 ) exist because the trivial solutions always satisfy the out-of-planegoverning equation. The deformation of the plate follows the trivial solutions when the load is sufficiently small, but asthe load increases, there is a point where the plate will bifurcate into a more stable configuration (with lower potentialenergy) in a buckled shape. In plate theory, the first buckling mode is usually determined by seeking the lowest totalpotential energy. Therefore, we applied the energy-based model to predict it. The PDE-based loss function is notsuitable for the buckling analysis since it inevitably converges to the trivial in-plane solution.For the neural network algorithms, a 5-hidden layer neural network (5 neurons each layer) with the energy-basedloss function was constructed. The FE simulations were performed in Abaqus/standard to get the first buckling mode.Modal analysis was first conducted to obtain the different buckling modes. The first buckling mode configuration wasthen induced as the geometric imperfection with a maximum 0.01 mm transverse deviation and the model with theimperfection is used to simulate the in-plane compression with the implicit solver. The buckled configuration predictedby the neural network algorithms for the two cases are shown in Figure 14c and d, respectively. In addition, Figure 14eand f respectively compare the deflection of the central line with the FE simulations. We can see that the bulkingconfigurations under two different boundary conditions are both correctly predicted. In this study, we chose four typical loading conditions to validate the physics-guided neural network framework thatwe developed. Although they are all simple in terms of loading and boundary conditions, they together provided arather comprehensive investigation of the accuracy of different approaches. The first example involved only in-planedeformation, but a nonlinear stretching force was applied. Among the four examples, this may be the simplest taskfor modeling. However, the differences between the data-driven approach and the physics-guided approaches andbetween the PDE-based and energy-based approaches were already very clear. The second example was also in-plane,19 Physics-Guided Neural Network Framework for Elastic Plates: Comparison of Governing Equations-Based andEnergy-Based Approaches
A P
REPRINT
Figure 14: Prediction of the buckling of elastic plates with two different types of constraints: simply supported (a) andclamped (b). The buckled configurations are shown in (c) and (d). The deflections of the central line are compared withFE simulations (e) and (f).but geometrical nonlinearity was introduced by designing a circular hole at the center of the plate. Different aspectratios of the central hole were investigated to obtain a wide range of the stress concentration factor. We observed thatthere is a small difference between the FE simulation and energy-based model. However, it is promising to find that theaccuracy of the energy-based method did not decrease as the stress concentration factor increased. In other words, thismethod is stable. The third example involved a pressure in the z -direction so that the out-of-plane governing equationcould no longer be neglected. The energy-based method still provided a satisfactory prediction. The last example wasmore challenging due to instability. No out-of-plane load was applied, but out-of-plane deformation occurred throughbuckling. The energy-based method showed a great advantage over its PDE-based counterpart because the latter alwaysconverged to the trivial solution with w = 0 . Therefore, these four examples covered almost all the important aspects ofplate deformation. There are two prominent differences between the PDE-based and energy-based loss functions. The first lies in the orderof the involved partial derivatives. The energy-based loss function deals with the strain components, which are functionsof the first derivatives of the displacement field. By constructing the neural network to directly output the straincomponents, for example in the first case, we can avoid additional computational errors coming from the derivationprocess. On the contrary, the PDE-based loss function has to include the residuals in the PDEs and BCs at the sametime. Therefore, it is almost impossible to reduce the order of equations by computational treatments, and consequently,a number of derivation processes have to be performed in the algorithm, accumulating computational errors. Thesecond difference between the PDE-based and energy-based loss functions is the number of the residual terms. Thecomplete PDE-based loss function has to sum up a total number of eleven residuals coming from three equations andeight boundary conditions. However, the energy-based deals only one residual by introducing the penalty term intothe total energy. Even if considering the penalty term as an independent quantity, we still have only two residuals forsumming up. This is a big simplification. The prominent difference between the energy-based and PDE-based lossfunctions is, therefore, mainly due to the computational errors accumulated through the derivative calculations as wellas the minimization process of the loss function. The energy-based model showed clear advantages mainly because itsderivative calculation is simple and its loss function is unified and smooth.20 Physics-Guided Neural Network Framework for Elastic Plates: Comparison of Governing Equations-Based andEnergy-Based Approaches
A P
REPRINT
These two differences are computational. Here, we provide a more in-depth discussion of the physical mechanisms.As already pointed out above, the energy approach and the PDE governing equations are mathematically equivalent.The energy-based loss (Eq. (24)) and the PDE-based loss (Eq. (20) are respectively formulated following these twoapproaches and, therefore, should also be equivalent. However, the equivalence could only be achieved when the weightratios of the PDE-based loss, λ s and λ d , can be determined in advance to have the physical meanings of displacementand force, respectively. In other words, λ s and L Bcs , λ d and L Bcd should be two pairs of conjugate variables in terms ofpotential energy. However, this is impossible not only because such values of λ s and λ d are difficult to calculate but alsobecause they are usually not uniform in all boundaries. In a practical application of neural network-based computationalframeworks, λ s and λ d are usually chosen by the user empirically. Like in [13], a ratio of 500 was chosen. We have seen that for loading cases with high nonlinearity there is still a relatively large deviation of the local predictionsbetween the energy-based neural network model and the FE results. Although we have noted that FE results are notnecessarily the exact solution and that it is unfair to attribute the deviation only to the computational error of the neuralnetwork-based algorithms, it is still necessary to point out the limitations of the proposed computational framework.The first limitation is indicated by its name – the accuracy and the applicability of this physics-guided computationalframework is largely determined by the physical laws that are implemented by human brains. In our study, the physicallaws are the classic plate theory. It is based on the strong Kirchhoff hypotheses, which will lose the applicability formoderately thick plates. The framework was developed based on these hypotheses and, therefore, inherits its limitations.The second limitation stems from the fundamental and shadow neural network we have used. Its capability toapproximate a highly non-uniform displacement or strain field is limited due to its simplicity. A deep neural networkthat involves a larger number of hidden layers is likely to increase accuracy. In the present study, we focused on theimplementation of physics into machine learning algorithms and we did not try a deeper neural network due to thelimitation in computational resources. Another approach to improve the approximation ability is to modify the structureof the neural network or seek for other machine learning models. For example, Wang et al. [12] added extra connectionsbetween non-adjacent layers to improve the approximation ability of the neural network.
As mentioned in the literature survey in the introduction, many initial successes of PGML or PINN algorithms havebeen achieved in modeling the dynamics of fluids as well as mass and heat transfer. To model the mechanical responsesof solids, a special challenge is that the variables of interest (stress and strain) are highly tensorial. As a comparison, inmany cases, only the pressure of a fluid is wanted. As a result, more PDEs and BCs have to be implemented into theloss function, leading to a low accuracy of the PDE-based algorithms as mentioned in 5.3. This point becomes veryclear through our present study – the governing equations are established in three directions, and each edge of the platehas four pairs of conjugates as BCs. It is also worth noting that our present study only considered elastic plates. Toimplement the plasticity theories, even the simplest one, there will be more intermediate state variables, thus creatingmore PDEs and ordinary differential equations (ODEs) to be solved. This will be one more special challenge for theneural network-based algorithms to take.
As a preliminary study, we demonstrated that the energy-based nerual network framework can provide a satisfactoryprediction of the mechanical response of elastic plates that is as good as the FE simulation results. This seemingly easyconclusion may cause an underrating of the contribution of this work. The conclusion can be generalized – for a systemthat is governed by a large number of PDEs and BCs, if the principle of minimum potential energy is applicable, amachine learning algorithm designed to minimize the potential energy will be more effective and efficient than directlyminimizing the total sum of all the residuals stemming from the PDEs and BCs. One potential extension of our neuralnetwork framework is modeling the thermodynamics of materials, which is also based on some energy indicators.While it is true that there is no evidence to prove the current accuracy of the framework is better than FE simulations,as the complexity of the system keeps increasing, we expect that the neural network framework will show more clearadvantages over the FE simulations. One such potential application is the modeling of multiphysics and multiscalesystems, lithium-ion batteries as a typical example. It is well-known that the FE simulations often suffer from a stringentcriterion to get converged when dealing with these systems. Energy-based models are very promising to bypass theconvergence challenge and get an approximate solution for practical applications.21 Physics-Guided Neural Network Framework for Elastic Plates: Comparison of Governing Equations-Based andEnergy-Based Approaches
A P
REPRINT
In this study, we established a physics-guided neural network-based computational framework to predict the mechanicalresponses of elastic plates. The physical laws that were implemented into the algorithm were from the classic platetheory derived following the Kirchhoff hypotheses. The governing PDEs are the well-known FvK equations, whichcan be derived from the principle of virtual displacement. In our computational framework, a neural network wasconstructed to output the displacement fields (or strain fields in some cases) with the input of spatial coordinates.Three different ways of formulating the loss function were investigated. One was purely data-driven by comparingthe predicted displacement field with the observed one from tests or FE simulations. The other two were based on thephysical laws. The PDE-based loss function was the total sum of all the residuals stemming from the PDEs and BCs,and the energy-based simply used the total potential energy as its loss. The computational framework that we developedwere then applied to four different types of loading conditions, including 1) the in-plane tension with non-uniformlydistributed stretching force to study the effect of the nonlinearity from external loads, 2) the in-plane central-holetension to investigate the nonlinearity from geometric imperfections, 3) the out-of-plane deflection to examine thecapability of modeling the out-of-plane deformation, and 4) the buckling induced by uniaxial compression to validatethe algorithm on instability analysis. In all the four cases, FE simulations with an extremely fine mesh size wereperformed as references. Through these validations and comparisons, the following conclusions can be drawn.1) The energy-based neural network framework developed in this study can approximately predict the mechanicalresponse of elastic plates with a satisfactory accuracy that is close to the FE simulations in all the four examinedloading conditions.2) The model with the energy-based loss function shows a clear advantage over the one with PDE-based in boththe modeling accuracy and efficiency.3) The purely data-driven approach that was trained with observed data from simulation has a high risk to overfitthe displacement and membrane force fields.4) The main reasons for the advantage of the energy-based model are its simplicity in the unified loss functionand lower order of derivatives involved in the algorithm, thus avoiding the accumulation of computationalerrors.It is optimistically expected that our energy-based neural network framework will have a wide spectrum of applicationsin future studies. Particularly, it provides an important energy-optimization inspiration for modeling complex engineeredsystems involving multiple scales and multiple physics.
A Appendix
A.1 Strain components
Following the Kirchhoff hypotheses, the displacement field can be expressed as: u ( x α , z ) = u ( x α ) − z · w ,x , (A.1)The general three-dimensional second-order nonlinear Green strains are ε αβ = 12 ( u ,α + u ,β + u γ,α u γ,β ) , ( α, β, γ = x, y, z ) . (A.2)We consider the moderate deformation of a plate, meaning that the transverse (i.e. out-of-plane) displacement gradients u z,x = w ,x and u z,y = w ,y can be relatively large and the in-plane displacement gradients u α,β , ( α, β = x, y ) are smalldue to the large width and length. The second-order terms in the Green strains can be therefore omitted except the w ,αβ ( α, β = x, y ) .Substitute Eq. (A.1) into the above equation, the strains can subsequently be simplified to the following strains of the2D plate theory, ε xx = 12 (cid:0) u ,α + u ,β + w ,α w ,β − z · w ,αβ (cid:1) , ( α, β = x, y ) ,ε γ = 0 , ( γ = x, y, z ) . (A.3)22 Physics-Guided Neural Network Framework for Elastic Plates: Comparison of Governing Equations-Based andEnergy-Based Approaches A P
REPRINT
A.2 Integration by parts
The virtual strains are calculated from the virtual displacements according to Eq. (2) and Eq. (3). For the first term inEq. (5) , applying integration by parts we have (cid:90) Ω N αβ δ ε αβ d x d y = (cid:90) Ω N αβ (cid:0) δ u α,β + δ u β,α + δ w ,α w ,β + w ,α δ w ,β (cid:1) = (cid:90) Ω N αβ (cid:0) δ u α,β + δ w ,α w ,β (cid:1) = (cid:90) Γ (cid:0) N αβ δ u α n β + N αβ w ,β δ w γ n α (cid:1) d s − (cid:90) Ω (cid:104) N αβ,β δ u α + ( N αβ w ,β ) ,α δ w (cid:105) d x d y, (A.4)where the partial integration is applied to get the displacement variation instead of its gradient. n = n x e x + n y e y is theoutward normal on the boundary ( n x and n y are the direction cosines of the unit normal). For the second term, we have (cid:90) Ω M αβ δ κ αβ d x d y = − (cid:90) Ω M αβ δ w ,αβ d x d y = − (cid:90) Γ [ M αβ δ w ,α ] n β d s + (cid:90) Ω M αβ,β δ w ,α d x d y = − (cid:90) Γ [ M αβ δ w ,α ] n β d s + (cid:90) Γ M αβ,β δ wn α d s − (cid:90) Ω M αβ,αβ δ w d x d y, (A.5)where the integration by parts is applied twice. A.3 Boundary conditions
We perform a coordinate transformation between the global Cartesian ( x, y, z ) coordinate and the local Cartesiancoordinate ( n, s, r ) (see Fig. 1), (cid:34) xyz (cid:35) = (cid:34) cos θ − sin θ θ cos θ
00 0 1 (cid:35) (cid:34) nsr (cid:35) = (cid:34) n x n y − n y n x
00 0 1 (cid:35) (cid:34) nsr (cid:35) , (A.6)where θ is the angle between the global x axis and the local n axis along the counterclockwise direction. Thedisplacements and stresses under the two coordinates are related by (cid:34) u v w (cid:35) = (cid:34) cos θ − sin θ θ cos θ
00 0 1 (cid:35) (cid:34) u n u s w r (cid:35) , (A.7) (cid:20) σ nn σ ns (cid:21) = (cid:20) n x n y n x n y − n x n y n x n y n x − n y (cid:21) (cid:34) σ xx σ yy σ xy (cid:35) . (A.8)According to the above relations, the stress boundary integrands in Eq. (9) can be rewritten with quantities under thelocal coordinate, ( N xx n x + N xy n y ) δ u + ( N yy n y + N xy n x ) δ v = ( N xx n x + N xy n y ) ( n x δ u n − n y δ u s ) + ( N yy n y + N xy n x ) ( n y δ u n + n x δ u s )= (cid:0) N xx n x + 2 N xy n x n y + N yy n y (cid:1) δ u n + (cid:2) N yy n x n y − N xx n x n y + N xy (cid:0) n x − n y (cid:1)(cid:3) δ u s = N nn δ u n + N ns δ u s , (A.9)and − ( M yy n y + M xx n x ) ∂ δ w∂y − ( M xx n y + M xx n x ) ∂ δ w∂x = − M nn ∂ δ w∂n − M ns ∂ δ w∂s . (A.10)23 Physics-Guided Neural Network Framework for Elastic Plates: Comparison of Governing Equations-Based andEnergy-Based Approaches A P
REPRINT
To fix the inconsistency, integration by parts is applied, (cid:90) Γ σ M ns ∂ δ w∂s d s = − [ M ns δ w ] Γ σ + (cid:90) Γ σ ∂M ns ∂s δ w d s. (A.11) [ M ns δ w ] Γ σ is zero when the stress boundary is closed or M ns = 0 . Then we can get the Eq. (14). Acknowledgment
J.Z. and W.L. are grateful to the support by AVL, Hyundai, Murata, Tesla, Toyota North America, Volkswa-gen/Audi/Porsche, and other industrial partners through the MIT Industrial Battery Consortium. M.Z.B is grateful to thesupport by Toyota Research Institute through the D3BATT Center on Data-Driven-Design of Rechargeable Batteries.Thanks are also due to the MIT-Indonesia Seed Fund to support J.Z.’s postdoctoral study.
References [1] M. Egmont-Petersen, D. de Ridder, and H. Handels. “Image processing with neural networks- A review”. In:
Pattern Recognition (2002).
ISSN : 00313203.
DOI : .[2] W. Rawat and Z. Wang. “Deep convolutional neural networks for image classification: A comprehensive review”.In: Neural Computation (2017).
ISSN : 1530888X.
DOI : .[3] R. M. French. “Introduction to Neural and Cognitive Modeling”. In: Biological Psychology (2002).
ISSN :03010511.
DOI : .[4] M. W. Libbrecht and W. S. Noble. “Machine learning applications in genetics and genomics”. In: Nature ReviewsGenetics (2015).
ISSN : 14710064.
DOI : .[5] Y. C. Lo et al. “Machine learning in chemoinformatics and drug discovery”. In: Drug Discovery Today (2018).
ISSN : 18785832.
DOI : .[6] R. Ramprasad et al. “Machine learning in materials informatics: Recent applications and prospects”. In: npjComputational Materials (2017). ISSN : 20573960.
DOI : .[7] M. Alber et al. “Integrating machine learning and multiscale modeling—perspectives, challenges, and opportuni-ties in the biological, biomedical, and behavioral sciences”. In: npj Digital Medicine (2019). ISSN : 2398-6352.
DOI : .[8] J. Han, A. Jentzen, and E. Weinan. “Solving high-dimensional partial differential equations using deep learning”.In: Proceedings of the National Academy of Sciences of the United States of America
ISSN : 10916490.
DOI : .[9] K. A. Severson et al. “Data-driven prediction of battery cycle life before capacity degradation”. In: Nature Energy
ISSN : 20587546.
DOI : .[10] A. Famili et al. “Data preprocessing and intelligent data analysis”. In: Intelligent Data Analysis (1997).
ISSN :15714128.
DOI : .[11] J. Sirignano and K. Spiliopoulos. “DGM: A deep learning algorithm for solving partial differential equations”. In: Journal of Computational Physics
ISSN : 10902716.
DOI : .[12] Z. Wang and Z. Zhang. “A mesh-free method for interface problems using the deep learning approach”. In: Journal of Computational Physics
400 (2020).
ISSN : 10902716.
DOI : .[13] E. Weinan and B. Yu. “The Deep Ritz Method: A Deep Learning-Based Numerical Algorithm for SolvingVariational Problems”. In: Communications in Mathematics and Statistics
ISSN : 2194671X.
DOI : .[14] M. Raissi, P. Perdikaris, and G. E. Karniadakis. “Physics-informed neural networks: A deep learning frameworkfor solving forward and inverse problems involving nonlinear partial differential equations”. In: Journal ofComputational Physics
378 (2019), pp. 686–707.
ISSN : 10902716.
DOI : .[15] J. X. Wang, J. L. Wu, and H. Xiao. “Physics-informed machine learning approach for reconstructing Reynoldsstress modeling discrepancies based on DNS data”. In: Physical Review Fluids
ISSN :2469990X.
DOI : .[16] E. Samaniego et al. “An energy approach to the solution of partial differential equations in computationalmechanics via machine learning: Concepts, implementation and applications”. In: Computer Methods in AppliedMechanics and Engineering
362 (2020), p. 112790.
ISSN : 00457825.
DOI : .24 Physics-Guided Neural Network Framework for Elastic Plates: Comparison of Governing Equations-Based andEnergy-Based Approaches A P
REPRINT [17] A. Karpatne et al. “Theory-guided data science: A new paradigm for scientific discovery from data”. In:
IEEE Transactions on Knowledge and Data Engineering
ISSN : 10414347.
DOI : .[18] W. Li et al. “Data-Driven Safety Envelope of Lithium-Ion Batteries for Electric Vehicles”. In: Joule (2019),pp. 1–13.
ISSN : 25424351.
DOI : .[19] Z. Chen et al. Direct prediction of phonon density of states with Euclidean neural network . 2020. eprint: arXiv:2009.05163 .[20] L. Zhang et al. “DeePCG: Constructing coarse-grained models via deep neural networks”. In:
Journal of ChemicalPhysics (2018).
ISSN : 00219606.
DOI : .[21] J. Darbon, G. P. Langlois, and T. Meng. “Overcoming the curse of dimensionality for some Hamilton–Jacobipartial differential equations via neural network architectures”. In: Research in Mathematical Sciences (2020).
ISSN : 21979847.
DOI : .[22] L. Lu et al. “DeepXDE: A deep learning library for solving differential equations”. In: (2019), pp. 1–21.[23] H. Zhao et al. “Learning the Physics of Pattern Formation from Images”. In: Physical Review Letters (2020).
ISSN : 10797114.
DOI : .[24] S. Effendy, J. Song, and M. Z. Bazant. “Analysis, Design, and Generalization of Electrochemical ImpedanceSpectroscopy (EIS) Inversion Algorithms”. In: Journal of The Electrochemical Society (2020).
ISSN : 1945-7111.
DOI : .[25] W. E, J. Han, and L. Zhang. “Integrating Machine Learning with Physics-Based Modeling”. In: (2020), pp. 1–23.[26] E. Qian et al. “Lift & Learn: Physics-informed machine learning for large-scale nonlinear dynamical systems”.In: Physica D: Nonlinear Phenomena (2020).
ISSN : 01672789.
DOI : .[27] E. Haghighat et al. “A deep learning framework for solution and discovery in solid mechanics”. In: (2020).[28] L. Wu et al. “A recurrent neural network-accelerated multi-scale model for elasto-plastic heterogeneous materialssubjected to random cyclic and non-proportional loading paths”. In: Computer Methods in Applied Mechanicsand Engineering (2020).
ISSN : 00457825.
DOI : .[29] D. Huang et al. “A machine learning based plasticity model using proper orthogonal decomposition”. In: Computer Methods in Applied Mechanics and Engineering
365 (2020), p. 113008.
ISSN : 00457825.
DOI : .[30] J. N. Reddy. Theory and Analysis of Elastic Plates and Shells . CRC press, 2006.
DOI : .[31] A. Föppl. Vorlesungen über technische Mechanik . Vol. 4. BG Teubner, 1899.[32] T. V. Kármán. “Festigkeitsprobleme im Maschinenbau”. In:
Mechanik . Springer, 1907, pp. 311–385.
DOI : .[33] J. Zhu, X. Zhang, and T. Wierzbicki. “Stretch-induced wrinkling of highly orthotropic thin films”. In: InternationalJournal of Solids and Structures
ISSN : 00207683.
DOI : .[34] E. Cerda, K. Ravi-Chandar, and L. Mahadevan. “Thin films: Wrinkling of an elastic sheet under tension”. In: Nature (2002).
ISSN : 00280836.
DOI : .[35] E. Puntel, L. Deseri, and E. Fried. “Wrinkling of a stretched thin sheet”. In: Journal of Elasticity (2011).
ISSN :03743535.
DOI : .[36] V. Nayyar, K. Ravi-Chandar, and R. Huang. “Stretch-induced stress patterns and wrinkles in hyperelastic thinsheets”. In: International Journal of Solids and Structures (2011).
ISSN : 00207683.
DOI : .[37] A. A. Sipos and E. Fehér. “Disappearance of stretch-induced wrinkles of thin sheets: A study of orthotropicfilms”. In: International Journal of Solids and Structures (2016).
ISSN : 00207683.
DOI : .[38] M. Abadi et al. Tensorflow: A system for large-scale machine learning . 2016.[39] A. Paszke et al.