Physics-informed neural networks for modeling rateand temperature-dependent plasticity

Rajat Arora† ,†Pratik Kakkar§ ,Amit Chakraborty§ ,Biswadip Dey§

Advanced Micro Devices (AMD), Austin, TX 78735, USA.
rajat.arora9464@gmail.com§
Siemens Technology, Princeton, NJ 08536, USA.
{pratik.kakkar, amit.chakraborty, biswadip.dey}@siemens.com

Abstract
This work presents a physics-informed neural network (PINN) based framework
to model the strain-rate and temperature dependence of the deformation fields
in elastic-viscoplastic solids. To avoid unbalanced back-propagated gradientsduring training, the proposed framework uses a simple strategy with no added
computational complexity for selecting scalar weights that balance the interplay
between different terms in the physics-based loss function. In addition, we highlighta fundamental challenge involving the selection of appropriate model outputs sothat the mechanical problem can be faithfully solved using a PINN-based approach.
We demonstrate the effectiveness of this approach by studying two test problems
modeling the elastic-viscoplastic deformation in solids at different strain rates
and temperatures, respectively. Our results show that the proposed PINN-based
approach can accurately predict the spatio-temporal evolution of deformation inelastic-viscoplastic materials.

1

Introduction

Modeling the elastic-plastic response of materials using conventional numerical methods, such as
finite element method, isogeometric analysis, or mesh-free methods, has always been computationallyexpensive due to the inherent iterative nature of discretization algorithms used in such methods. Furthermore, multitude of ‘fundamentally accurate’ theories for the high-fidelity modeling of dislocation
mediated plastic deformation at different scales [1–12] or fracture modeling in materials [13–17], isbringing these numerical solvers to their limits. In this context, PINNs offer great opportunities tospeed up (nonlinear) mechanical modeling of materials.
The idea of using neural networks to learn the solution of partial differential equations (PDEs) by
minimizing a loss function, comprising the residual error of governing PDEs and its initial/boundary
conditions, has been around for some time [18, 19]. More recently, Raissi et. al [20, 21] have extended
this concept towards PINNs which can solve the forward and inverse problems involving general
nonlinear PDEs by relying on small or even zero labeled datasets. Several applications of PINNs
can be found in the literature ranging from modeling of fluid flows and Navier Stokes equations[22–27], cardiovascular systems [28, 29], and material modeling [30–35], among others. Compared
to traditional data-driven approaches for predicting path-dependent plastic behavior in metals [36–40],
PINNs can learn high-fidelity surrogate models while simultaneously reducing (or even eliminating)the need for bigger training datasets. However, developing a physics-informed neural network tomodel the spatio-temporal variation of deformation in elastic-plastic solids, along with its dependence
on strain-rate and temperature, poses several technical challenges.
Machine Learning and the Physical Sciences workshop, NeurIPS 2022.

In this work, we take a first step in highlighting these challenges and demonstrate the strengthof PINNs for modeling elastic-viscoplastic deformation in materials. In particular, we focus onpredicting the spatio-temporally varying deformation fields (displacement, stress, and plastic strain)under different strain rates (i.e., applied loading rate) and temperatures, respectively. We present
a detailed discussion on the construction of (physics-based) composite loss along with a briefsummary on ways to avoid unbalanced back-propagated (exploding) gradients during model training.
Furthermore, a strategy with no added computational complexity for choosing the scalar weights thatbalance the interplay between different terms in composite loss is also proposed. Although the current
work focuses on the scenarios with monotonic loading paths, we note that the deformation of anelastic-viscoplastic solid is a highly nonlinear function of temperature, strain rate, spatial coordinates,and strain. This real-time stress predictive capability for elastic-viscoplastic materials enjoys special
use in the design and development of energy storage devices (e.g., lithium metal solid-state-batteries).
Specifically, the study conducted here corresponds to analyzing the effect of impact (i.e. crash) andheat to the solid lithium anode in the solid state batteries.
Notation and Terminology: Vectors and tensors are represented by bold face lower- and upper-caseletters, respectively. The symbol ‘·’ denotes single contraction of adjacent indices of two tensors(i.e. a · b = ai bi or A · n = Aij nj ). The symbol ‘:’ denotes double contraction of adjacent indicesof two tensors of rank two or higher (i.e. A√: B = Aij Bij or C : A = Cijkl Akl ). The norm of asecond order tensor A is given by ||A|| = A : A. The symbols ∇ and Div denote the gradientand the divergence operators, respectively. I denotes the second order identity tensor.

2

Background on the Deformation Behavior of Elastic-viscoplastic Solids

In this section, we describe the nonlinear PDEs that govern the behavior of elastic-viscoplastic solids
under loads at small deformation (see [41] for further details). In the absence of body and inertial
forces, the strong form of the mechanical equilibrium on a volumetric domain Ω can be expressed as
Div σ = 0 in Ω, withσ · n = tbc on ∂ΩN and u = ubc on ∂ΩD ,(1)

where σ and u denote the stress and displacement, respectively, n denote the unit outward normal tothe external boundary ∂Ω of domain Ω, and tbc and ubc denote the known traction and displacement
vectors on the Neumann boundary ∂ΩN and Dirichlet boundary ∂ΩD , respectively. The totalstrain tensor ϵ, which is the symmetric part of the displacement gradient, can be decomposed
into the sum of elastic and plastic strain components denoted by ϵe and ϵp , respectively, i.e.,∇u + (∇u)T = 2(ϵe + ϵp ). The stress is given by the Hooke’s law σ = C : ϵe , where C is the
fourth order elasticity tensor. Furthermore, the plastic strain evolution is governed by ′ 1/m
  m+1Qσ′||σ ||
3 2mp
Ae− Rθ,(2)ϵ̇ =2||σ ′ ||
S
where σ ′ = σ − trace(σ)I denotes the deviatoric part of the stress tensor, A is a pre-exponential
factor, Q denotes the activation energy, R is the molar gas constant, θ denotes the temperature of thedomain Ω, and m ∈ (0, 1] is a strain-rate-sensitivity parameter.
Then, by letting S denote the material strength, its dynamics can be expressed asṠ = h(S, σ),(3)with the hardening function h(S, σ) defined as [42],



a
S
S
h(S, σ) = H0 1 −ϵ̄˙p ,
sign 1 −Ss
Ss
where H0 and n are two strain-hardening parameters. Furthermore, the saturation value of S for aQ
given strain rate and temperature, i.e., Ss , is given by Ss = S∗ (e Rθ ϵ̄˙p /A)n , where S∗ and a are two
additional strain-hardening parameters.
In this work, we select the material parameters as: A = 4.25 × 104 s−1 , Q = 37 kJ/mol, m = 0.15,
S∗ = 2 MPa, S0 = 0.95 MPa, H0 = 10 MPa, a = 2, n = 0.05, E = 7810 MPa, and ν = 0.38.
These parameters have been calibrated using the experimental data from direct tension tests on
polycrystalline lithium specimens [42, 43].
2

3

Proposed Learning Framework

Since plastic strain ϵp and stress σ are related by Hooke’s law, the elastic-plastic deformation canbe uniquely characterized by the displacement vector u and any one of the two tensors, ϵp or σ,along with the internal variable S. However, as shown in Appendix B, a PINN with such choice of
outputs suffers from degraded accuracy and convergence issues. Therefore, we propose a mixedvariable formulation and use u, σ, ϵp and S as the PINN outputs. To model the effect of strain rate
and temperature on the elastic-viscoplastic behavior in two-dimensional solids, we use two separate
PINNs to predict the output variables at any given location (x1 , x2 ). To capture strain rate dependence,
the first PINN uses scalar strain Γ and strain rate Γ̂ as two more additional inputs. The second PINN
uses scalar strain Γ and temperature θ as the additional inputs to capture temperature dependence.
These PINNs are realized via a multilayer perceptron with 9 hidden layers, 120 neurons/layer, andtanh activation function. In addition, we normalize the data along each component.
Physics-informed Loss
Function: The composite
loss function is definedas the weighted sum of aphysics loss Lphy and thedata loss computed as themean-squared-error (MSE)
between the normalized
ground truth data and
PINN outputs. The physics
loss is composed of sevendifferent components i) Lpde : PDE loss,ii) LDbc : Dirichlet boundary condition loss,iii) LN bc : Neumann boundary condition loss,iv) Lic : initial conditionloss,v) Lc : constitutive loss
corresponding to the satisfaction of constitutive law,vi) Lpsr : plastic strain rate
loss corresponding to theequation governing theevolution of plastic strain,andvii) Ls :strength loss
enforcing the material
strength evolution equation.
Please see Appendix A
for further details on theconstruction of the loss
function.

Γ = 0.01

Γ = 0.02

Γ = 0.04

Γ = 0.01

Γ = 0.02

Γ = 0.04

Network Training: After initializing the networkwith Xavier initialization[44], we train it on PyTorch using Adam optimizer [45] (initial learning
rate = 10−3 ). The network(left
is trained for 8000 epochs Figure 1: Comparison of the results obtained from the PINN model−2.5 −1
block)
withthegroundtruth
referencedata(rightblock)
for
Γ̂
=
10
s
.
during which the learning
rate is varied using ReduceLROnPlateau scheduler (patience=30).
3

Dataset Generation: To generate the ground truth data, we develop an in-house code using deal.II[46] and solve equations (1)-(3) over a 32 × 32 grid for scalar strain values up to Γ = 0.04.
Furthermore, for studying the effect of strain rate and temperature on the spatio-temporal evolutionof deformation fields in the body, the ground truth dataset is generated for strain rates Γ̂ (s−1 ):{10−4 , 10−3 , 10−2 , 10−1 } at θ = 298K and for temperatures θ (K): {298, 318, 358, 378} at Γ̂ =10−1 s−1 . The dataset is then randomly split into a 80 : 20 ratio for training and validation purposes.

4

Results & Discussion

Case I: Strain rate dependence First we compare
the predicted values of thestress, plastic strain, anddisplacement fields in thedomain with a test datasetfor Γ̂ = 10−2.5 s−1 atΓ = 0.01, 0.02, and 0.04
(Figure 1). We can notice that the predicted values have no visible artifacts and are in great agreement with the FEM reference results. Small values of the normalized rootmean-squared-error e (reported underneath the corresponding field plot) further
confirms predictions by the
PINN match the FEM reference results remarkablywell.

Γ = 0.01

Γ = 0.02

Γ = 0.04

Γ = 0.01

Γ = 0.02

Γ = 0.04

Next, we test the predictive
power of the trained PINN
for values of inputs thatlie outside the training
data range. Specifically,
we calculate the averaged
value of the normalized
root-mean-squared-error E
for 7 different strain rate
values Γ̂ (s−1 ) = {10−0.5 ,10−.75 , 10−1.5 , 10−2.5 ,10−3.5 , 10−4.25 , 10−4.5 }and multiple strain values
in the range [0, 0.08]. We
can readily notice from
Figure 3 that the error E
is very small (≈ 1%) upto
strain Γ = 0.04 whenstrain rate lies within the
training range. However, in Figure 2: Comparison of the results obtained from the PINN model (leftthe region Γ ∈ (0.04, 0.08], block) and the ground truth reference data (right block) for T = 328K.
the error E steadily increases to ≈ 10%. Also, when the strain rate is outside the training data range, the error E is large at
all strains implying that the predicted values do not match well with the FEM data.
Case II: Temperature dependence Similar to Case I, we first compare the predicted values ofthe stress, plastic strain, and displacement fields in the domain with a test dataset for T = 328K
4

at Γ = 0.01, 0.02, and 0.04 (Figure 2). We can notice that the error e has small values and thepredictions match well with the FEM reference results.
Figure 4 shows that the error E rises to ≥ 10% as
the strains go beyond the training data range. On
the other hand, when the strain is within the training range, the errors are still < 10% even when thetemperature is outside the training range.

5

Conclusion

This work demonstrates the strength of PINNs in Figure 3: Variation of error E with strain Γ forproblems dealing with the evolution of highly nonlin- different strain rates in and outside of the training
ear deformation field in elastic-viscoplastic materials range (Γ̂ ∈ [10−4 , 10−1 ]s−1 , Γ ∈ [0, 0.04]).
under monotonous loading. In particular, we trained
two specific PINN models and applied them predicting the spatio-temporally varying deformation
field in elastic-viscoplastic materials at different strain rates, and temperatures, respectively. The predicted values are in great agreement with the ground truth reference data for the test cases discussed
in this work.
This work highlighted a fundamental challenge involving selection of appropriate model outputs sothat the mechanical problem can be faithfully solved
using neural networks. We present and compare twopotential choice of outputs for the model in Appendix
B and present detailed reasoning for preferring one
choice over the other. This work also discusses theconstruction of composite loss function, comprising
the data loss component and physics-based loss components. We also use a novel physics-based strategy
for selecting the non-dimensional scalar constants Figure 4: Variation of error E with strain Γ forthat weigh each component in the physics-based loss different temperatures in and outside of the training
function without any added computational complex- data range (θ ∈ [298, 378]K, Γ ∈ [0, 0.04]).
ity. Moreover, a novel loss criterion for residual calculation corresponding to plastic strain rate equation is proposed to alleviate issues related to unbalanced
back-propagated (exploding) gradients during model training.
The real-time stress field prediction in such highly nonlinear mechanical system paves the way formany new applications, such as design and optimization of lithium ion batteries or inverse modelingproblems which were previously computationally intractable.

Broader Impact
We introduce a learning framework for predicting nonlinear spatio-temporal variation of deformation
inside elastic-viscoplastic materials under various operating conditions (e.g., temperature, strain
rate, etc.). This is an important problem in multiple engineering applications. For example, theoutcome of this work can be used to carry out fast (near real-time) simulation of the behavior of asolid Li anode under varying loading conditions. This, in turn, can improve the process of design,development, and operation of solid-state lithium-metal batteries. However, our proposed framework
is still a conceptual proposal and has a very low (around 2) Technology Readiness Level (TRL) [47].
We are yet to fully understand its limitations and failure scenarios that can significantly influence its
real-world adoption.

References
[1]

[2]

K. L. Nielsen and C. F. Niordson, “A finite strain fe-implementation of the fleck-willis gradient theory:
Rate-independent versus visco-plastic formulation,” European Journal of Mechanics / A Solids, vol. 75,
pp. 389–398, 2019.
A. L. Greer, Y. Q. Cheng, and E. Ma, “Shear bands in metallic glasses,” Materials Science and Engineering: R: Reports, vol. 74, no. 4, pp. 71–132, 2013.

5

[3]
[4]
[5]
[6]
[7]
[8]
[9]

[10]
[11]
[12]
[13]

[14]

[15]

[16]

[17]
[18]
[19]

[20]
[21]

[22]

[23]
[24]
[25]

[26]
[27]

R. Arora, “Computational approximation of mesoscale field dislocation mechanics at finite deformation,”
Ph.D. dissertation, Carnegie Mellon University, 2019.
A. Arora, R. Arora, and A. Acharya, “Mechanics of micropillar confined thin film plasticity,” arXiv
preprint arXiv:2202.06410, 2022.
R. Arora, X. Zhang, and A. Acharya, “Finite element approximation of finite deformation dislocation
mechanics,” Computer Methods in Applied Mechanics and Engineering, vol. 367, p. 113 076, 2020.
N. A. Fleck, G. M. Muller, M. F. Ashby, and J. W. Hutchinson, “Strain gradient plasticity: Theory andexperiment,” Acta Metallurgica et Materialia, vol. 42, no. 2, pp. 475–487, 1994.
R. Arora and A. Acharya, “A unification of finite deformation J2 Von-Mises plasticity and quantitative
dislocation mechanics,” Journal of the Mechanics and Physics of Solids, vol. 143, p. 104 050, 2020.
C. F. Niordson and V. Tvergaard, “A homogenized model for size-effects in porous metals,” Journal of
the Mechanics and Physics of Solids, vol. 123, pp. 222–233, 2019.
J. Lynggaard, K. L. Nielsen, and C. F. Niordson, “Finite strain analysis of size effects in wedge indentation
into a face-centered cubic (fcc) single crystal,” European Journal of Mechanics / A Solids, vol. 76, pp. 193–
207, 2019.
M. Kuroda and V. Tvergaard, “A finite deformation theory of higher-order gradient crystal plasticity,”
Journal of the Mechanics and Physics of Solids, vol. 56, no. 8, pp. 2573–2584, 2008.
L. Evers, W. Brekelmans, and M. Geers, “Non-local crystal plasticity model with intrinsic SSD and GND
effects,” Journal of the Mechanics and Physics of Solids, vol. 52, no. 10, pp. 2379–2401, 2004.
R. Arora and A. Acharya, “Dislocation pattern formation in finite deformation crystal plasticity,” International Journal of Solids and Structures, vol. 184, pp. 114–135, 2020.
C. Miehe, F. Welschinger, and M. Hofacker, “Thermodynamically consistent phase-field models offracture: Variational principles and multi-field fe implementations,” International journal for numericalmethods in engineering, vol. 83, no. 10, pp. 1273–1311, 2010.
M. J. Borden, T. J. Hughes, C. M. Landis, and C. V. Verhoosel, “A higher-order phase-field model
for brittle fracture: Formulation and analysis within the isogeometric analysis framework,” Computer
Methods in Applied Mechanics and Engineering, vol. 273, pp. 100–118, 2014.
G. Yingjun, L. Zhirong, H. Lilin, M. Hong, H. Chuanggao, and L. Kui, “Phase field crystal study of nanocrack growth and branch in materials,” Modelling and Simulation in Materials Science and Engineering,
vol. 24, no. 5, p. 055 010, 2016.
P. Areias, T. Rabczuk, and M. Msekh, “Phase-field analysis of finite-strain plates and shells including
element subdivision,” Computer Methods in Applied Mechanics and Engineering, vol. 312, pp. 322–350,
2016.
C. Kuhn and R. Müller, “A continuum phase field model for fracture,” Engineering Fracture Mechanics,
vol. 77, no. 18, pp. 3625–3634, 2010.
I. E. Lagaris, A. Likas, and D. I. Fotiadis, “Artificial neural networks for solving ordinary and partial
differential equations,” IEEE transactions on neural networks, vol. 9, no. 5, pp. 987–1000, 1998.
I. E. Lagaris, A. C. Likas, and D. G. Papageorgiou, “Neural-network methods for boundary value
problems with irregular boundaries,” IEEE Transactions on Neural Networks, vol. 11, no. 5, pp. 1041–
1049, 2000.
M. Raissi, P. Perdikaris, and G. E. Karniadakis, “Physics informed deep learning (part i): Data-driven
solutions of nonlinear partial differential equations,” arXiv preprint arXiv:1711.10561, 2017.
M. Raissi, P. Perdikaris, and G. E. Karniadakis, “Physics-informed neural networks: A deep learning
framework for solving forward and inverse problems involving nonlinear partial differential equations,”
Journal of Computational Physics, vol. 378, pp. 686–707, 2019.
L. Sun, H. Gao, S. Pan, and J.-X. Wang, “Surrogate modeling for fluid flows based on physics-constrained
deep learning without simulation data,” Computer Methods in Applied Mechanics and Engineering,
vol. 361, p. 112 732, 2020.
C. Rao, H. Sun, and Y. Liu, “Physics-informed deep learning for incompressible laminar flows,” Theoretical and Applied Mechanics Letters, vol. 10, no. 3, pp. 207–212, 2020.
T. Zhang, B. Dey, P. Kakkar, A. Dasgupta, and A. Chakraborty, “Frequency-compensated pinns for
fluid-dynamic design problems,” arXiv Preprint, vol. 2011.01456, 2020.
X. Jin, S. Cai, H. Li, and G. E. Karniadakis, “Nsfnets (navier-stokes flow nets): Physics-informed neuralnetworks for the incompressible navier-stokes equations,” Journal of Computational Physics, vol. 426,
p. 109 951, 2021.
H. Gao, L. Sun, and J.-X. Wang, “Phygeonet: Physics-informed geometry-adaptive convolutional neuralnetworks for solving parametric pdes on irregular domain,” arXiv e-prints, arXiv–2004, 2020.
T. Zhang, B. Dey, P. Kakkar, A. Dasgupta, and A. Chakraborty, “Frequency-compensated PINNs
for Fluid-dynamics Design Problems,” in NeurIPS Workshop on Machine Learning for Engineering
Modeling, Simulation, and Design, 2020.

6

[28]

[29]
[30]

[31]

[32]
[33]

[34]

[35]
[36]

[37]
[38]

[39]

[40]

[41]
[42]
[43]
[44]

[45]
[46]
[47]
[48]
[49]
[50]
[51]

G. Kissas, Y. Yang, E. Hwuang, W. R. Witschey, J. A. Detre, and P. Perdikaris, “Machine learning in
cardiovascular flows modeling: Predicting arterial blood pressure from non-invasive 4d flow mri datausing physics-informed neural networks,” Computer Methods in Applied Mechanics and Engineering,
vol. 358, p. 112 623, 2020.
F. Sahli Costabal, Y. Yang, P. Perdikaris, D. E. Hurtado, and E. Kuhl, “Physics-informed neural networks
for cardiac activation mapping,” Frontiers in Physics, vol. 8, p. 42, 2020.
A. Frankel, K. Tachida, and R. Jones, “Prediction of the evolution of the stress field of polycrystals
undergoing elastic-plastic deformation with a hybrid neural network model,” English, Machine Learning:
Science and Technology, vol. 1, no. 3, Jul. 2020.
R. Tipireddy, P. Perdikaris, P. Stinis, and A. Tartakovsky, “A comparative study of physics-informed
neural network models for learning unknown dynamics and constitutive relations,” arXiv preprint
arXiv:1904.04058, 2019.
E. Zhang, M. Yin, and G. E. Karniadakis, “Physics-informed neural networks for nonhomogeneousmaterial identification in elasticity imaging,” arXiv preprint arXiv:2009.04525, 2020.
X. Meng and G. E. Karniadakis, “A composite neural network that learns from multi-fidelity data:
Application to function approximation and inverse pde problems,” Journal of Computational Physics,
vol. 401, p. 109 020, 2020.
Q. Zhu, Z. Liu, and J. Yan, “Machine learning for metal additive manufacturing: Predicting temperatureand melt pool fluid dynamics using physics-informed neural networks,” Computational Mechanics,
vol. 67, no. 2, pp. 619–635, 2021.
R. Arora, “Machine learning-accelerated computational solid mechanics: Application to linear elasticity,”
arXiv preprint arXiv:2112.08676, 2021.
M. Mozaffar, R. Bostanabad, W. Chen, K. Ehmann, J. Cao, and M. Bessa, “Deep learning predicts pathdependent plasticity,” Proceedings of the National Academy of Sciences, vol. 116, no. 52, pp. 26 414–
26 420, 2019.
D. W. Abueidda, S. Koric, N. A. Sobh, and H. Sehitoglu, “Deep learning for plasticity and thermoviscoplasticity,” International Journal of Plasticity, vol. 136, p. 102 852, 2021.
L. Benabou, “Development of LSTM Networks for Predicting Viscoplasticity With Effects of Deformation, Strain Rate, and Temperature History,” Journal of Applied Mechanics, vol. 88, no. 7, May2021.
D. Huang, J. N. Fuhg, C. Weißenfels, and P. Wriggers, “A machine learning based plasticity model using
proper orthogonal decomposition,” Computer Methods in Applied Mechanics and Engineering, vol. 365,
p. 113 008, 2020.
M. B. Gorji, M. Mozaffar, J. N. Heidenreich, J. Cao, and D. Mohr, “On the potential of recurrent neuralnetworks for modeling path dependent plasticity,” Journal of the Mechanics and Physics of Solids,
vol. 143, p. 103 972, 2020.
M. E. Gurtin, E. Fried, and L. Anand, The Mechanics and Thermodynamics of Continua. Cambridge
University Press, 2010.
L. Anand and S. Narayan, “An elastic-viscoplastic model for lithium,” Journal of The Electrochemical
Society, vol. 166, no. 6, A1092, 2019.
W. S. LePage et al., “Lithium mechanics: Roles of strain rate and temperature and implications forlithium metal batteries,” Journal of The Electrochemical Society, vol. 166, no. 2, A89, 2019.
X. Glorot and Y. Bengio, “Understanding the difficulty of training deep feedforward neural networks,”
in Proceedings of the thirteenth international conference on artificial intelligence and statistics, JMLR
Workshop and Conference Proceedings, 2010, pp. 249–256.
D. P. Kingma and J. Ba, “Adam: A Method for Stochastic Optimization,” ICLR, 2015.
W. Bangerth, R. Hartmann, and G. Kanschat, “deal.II – a general purpose object oriented finite element
library,” ACM Trans. Math. Softw., vol. 33, no. 4, pp. 24/1–24/27, 2007.
S. Hirshorn and S. Jefferies, “Final Report of the NASA Technology Readiness Assessment (TRA) Study
Team,” 2016.
S. Wang, Y. Teng, and P. Perdikaris, “Understanding and mitigating gradient pathologies in physicsinformed neural networks,” arXiv preprint arXiv:2001.04536, 2020.
R. Bischof and M. Kraus, “Multi-objective loss balancing for physics-informed deep learning,” arXiv
preprint arXiv:2110.09813, 2021.
C. Rao, H. Sun, and Y. Liu, “Physics-informed deep learning for computational elastodynamics withoutlabeled data,” Journal of Engineering Mechanics, vol. 147, no. 8, p. 04 021 043, 2021.
W. M. Czarnecki, S. Osindero, M. Jaderberg, G. Świrszcz, and R. Pascanu, “Sobolev training for neuralnetworks,” arXiv preprint arXiv:1706.04859, 2017.

7

Appendix
A

Construction of the Physics-informed Loss Function

The development of a PINN based approach to predict the solution of a system of nonlinear PDEs canbe viewed as an optimization problem which involves solving for (W , b) that minimizes network’s
total loss. The composite loss L comprises the summation of supervised data loss Ldata and thephysics-based loss Lphy i.e. L = Ldata + Lphy . The non-dimensional supervised data loss Ldata
measures the discrepancy between the normalized ground truth data T̄ and the neural network outputs
Ȳ and is given by
NX
data N
out
X
1
(j)
(j)
Ldata =||Ȳi − T̄i ||2 ,(4)
Ndata j=1 i=1
where Ndata denotes the number of ground truth samples and Nout is the number of scalar output
variables.
To evaluate the physics-based loss Lphy , we sample a collection of randomly distributed collocation
points discretizing the normalized input space. The whole set of collocation points is denotedby P = {Pω , Pζ,D , Pζ,N , Pη } where Pω denotes the collocation points in the entire input space[−1, 1]Nin . Pζ,D and Pζ,N denote the subset of P that intersects with the ∂ΩD and ∂ΩN , respectively.
Pη denotes the subset of P that intersects with Γ = 0 (or t = 0).
To this end, we construct the physics-based loss Lphy with seven components i) PDE loss Lpde ,ii) Dirichlet boundary condition loss LDbc , iii) Neumann boundary condition loss LN bc , iv) initialcondition loss Lic , v) constitutive loss Lc corresponding to the satisfaction of constitutive law, vi)
plastic strain rate loss Lpsr corresponding to the equation governing the evolution of plastic strain,
and vii) strength loss Ls enforcing the material strength evolution equation. Each component of Lphy
is individually calculated as follows:
Lpde = ||Div σ||2Pω
LDbc = ||u − ubc ||2Pζ,D
LN bc = ||σ · n − tbc ||2Pζ,N
Lic = ||Y − I0 ||2Pη

(5)
p

2

Lc = ||σ − C : (∇u − ϵ )||Ls = ||Ṡ − h(S, σ)||2qq1/m ′Q
3 Aσ ′ − Rθ
3 ||σ ||
Lpsr = MMSE ϵ̇p − 2 ||σ′ || e
2 S
In the above, I0 denotes the initial state of the system i.e. outputs at t = 0. The loss criterion MMSE
is discussed in detail below. Lphy is then given as the weighted sum of these loss components
Lphy = λ1 Lpde + λ2 Lic + λ3 Lpsr + λ4 Ls + λ5 Lc + λ6 LDbc + λ7 LN bc(6)

where λi , (i = 1...6) are the scalar weights. Please note that each of these loss components arecomputed using automatic differentiation.
Next, we briefly discuss the two main difficulties that hinder the training of DNNs for elasticviscoplastic modeling applications.
(I) The power law dependence of the equivalent plastic strain rate ϵ̇p leads to large values of L2 normof Lpsr loss (≥ O(1018 )) which causes unstable imbalance in the magnitude of the back-propagatedgradients during the training when using common loss criterions such as Mean-Squared-Error (MSE).
Therefore, in this work, we use a novel Modified Mean Squared Error (MMSE) loss criterion to
reduce the numerical stiffness associated with equation (2) and allow stable gradients to be usedduring the training
MMSE(A) = log10 (1 + ||A||).
(7)
In the above, A denotes the residual value. The loss criterion is equivalent to the Mean Squared Error(MSE) criterion when the discrepancy between the residual values are small.
(II) The relative coefficients λi (i = 1..6) for all the losses comprising Lphy play an important role in8

mitigating the gradient pathology issue during the training [48]. There are competing effects between
these different loss components which can lead to convergence issues during the minimization of thecomposite loss L (see [22, Sec. 4.1]). While the recent advances in mitigating gradient pathologies[48, 49] might improve predictive accuracy, they introduce additional computational and memory
overhead because of the calculation of an adaptive factor for each loss component. In this work, we
devise a simple strategy, with no added computational complexity, to evaluate the coefficients which
remain constant during the course of training. The strategy is outlined as follows:• The Dirichlet boundary condition and initial condition losses (Lic and LDbc ) are calculated
in a normalized manner (scaled between [−1, 1]). So, we take λ2 = λ6 = 1.
• The other loss components are nondimensionalized using appropriate scales as shown in
Table 1. µc is a constant chosen to scale quantities with units of stress. Based on theobservation that stress is often nondimensionalized by Shear Modulus µ in conventional
numerical methods, we choose µc = 0.01µ to achieve tight tolerance on the equilibrium
equation and traction boundary conditions.
• Since material strength S and µ differ by orders of magnitude, S is nondimensionalized by
S0 .
• We nondimensionalize time by using strain rate Γ̂ , since Γ̂ sets the time scale for the
problem.
• The length is nondimensionalized by the characteristic length of the domain, chosen to be
H in this work.
Loss component

Scaling

Lpde

λ1 = µHc 2

Lic

λ2 = 1Lpsr

λ3 = Γ̂122

1
2(S0 Γ̂ )

Ls

λ4 =

Lc

λ5 = µ1c 2

LDbc

λ6 = 1

LN bc

λ7 = µ1c 2

Table 1: Scaling constants for different physics-based loss components.

B

The Rationale Governing the Selection of Output Variables

This section compares the results obtained from two PINN models - (a) Model I with displacement,
stress, plastic strain, and strength (u, σ, ϵp , S) as outputs; and (b) Model II with displacement, plasticstrain, and stress (u, ϵp , S) as outputs.
For Model II, the physics-based loss Lphy is obtained from the set of equations (5) with the following
important changes: i) The stress is directly calculated from the displacements and plastic strains
which are outputs of the neural network, i.e. σ = C : (∇u − ϵp ). This implicitly leads to satisfactionof constitutive law so the loss component Lc is ignored. ii) The data loss Ldata is also modified toaccount for the current model outputs.
The study conducted here corresponds to Case I. i.e., understanding the effect of strain rate on thespatio-temporal evolution of deformation in an elastic-viscoplastic material. The learning rate forModel II is taken to be 10−4 while keeping the collocation points and all other hyperparameters thesame for both the architectures as described in Section 4.
The convergence of the training loss for both the models is presented in Fig. 5. It can be seen thatloss reaches a stagnation value of ≈ 1.9 for model II at around 3000 epochs which is approximately
hundred times larger than the converged loss value obtained for model I. We can conclude that model
I does not suffer from any such degraded accuracy or convergence issue as indicated by Figure 5. This
9

Figure 5: Comparison of training history for Model I and II differing only in the model outputs.
result is an extension of the similar observation for the purely linear elastic calculations presented in[50] to the general elastic-plastic modeling case discussed here.
While the exact reasons for such a behavior are still unclear, we highlight the main differences between
the two models. First, the stress calculated in model II is sensitive to the noise in the gradients of u.
Second, we note that highest order of the spatial derivatives occurring in the composite loss function
is one and two for models I and II, respectively. Moreover, in elastic/elastic-plastic deformations theorder of displacement field magnitudes in the x1 and x2 direction can be vastly different becauseof the loading setup and Poisson’s effect. We believe that these factors combine together to giverise to convergence issue and degraded accuracy when using model II. The use of improved training
technique [51], which also approximates target derivatives along with target values, may alleviatethese issues for model II but that may involve added computational complexity and remains the
subject of future investigation.

C

Ablation Studies on the Effect of Network Depth and WidthFirst, we perform a (non-exhaustive) parametric study to identify a suitable number of hidden layers
Nl and number of neurons per layer Nn needed to model the deformation field with an acceptableaccuracy. We train 6 neural networks with the following architectures: i) 5 × 100, ii) 6 × 100, iii)

Figure 6: Comparison of the composite loss L associated with different PINN architectures (Nl × Nn ) for
Case I. Nl and Nn denotes the number of layers and number of nodels per layer, repectively.

6 × 120, iv) 7 × 120, v) 8 × 120, and vi) 9 × 120. Figure 6 presents the training history for eachof these architectures. As expected, we see a merit in increasing both Nn and Nl initially but the
final value of the composite loss stops improving when the number of layers are increased from 7
to 9 keeping Nn fixed at 120. These three network architectures (7 × 120, 8 × 120, and 9 × 120)
reduce the nondimensional composite loss by almost five orders of magnitude (from 102 to ∼ 10−3 ).
The values of the corresponding validation losses are monitored to notice any overfitting issues. We
use the neural network with architecture 9 × 120 for generating results associated with Case I in
Section 4.
10

Similar to Case I, we first conduct a study to gain insight into the effect of Nl and Nn on the composite
loss L and train the six aforementioned neural network architectures (i.e., 5 × 100, 6 × 100, 6 × 120,
7 × 120, 8 × 120, 9 × 120). Figure 7 presents the training history for each of these architectures
which shows similar trend as in Figure 6. Therefore, we use the neural network with architecture
9 × 120 for generating results associated with Case II in Section 4.

Figure 7: Comparison of the composite loss L associated with different PINN architectures (Nl × Nn ) for
Case II.

11

Checklist
1. For all authors...
(a) Do the main claims made in the abstract and introduction accurately reflect the paper’s
contributions and scope? [Yes](b) Did you describe the limitations of your work? [Yes](c) Did you discuss any potential negative societal impacts of your work? [N/A](d) Have you read the ethics review guidelines and ensured that your paper conforms tothem? [Yes]2. If you are including theoretical results...
(a) Did you state the full set of assumptions of all theoretical results? [N/A](b) Did you include complete proofs of all theoretical results? [N/A]
3. If you ran experiments...
(a) Did you include the code, data, and instructions needed to reproduce the main experimental results (either in the supplemental material or as a URL)? [No](b) Did you specify all the training details (e.g., data splits, hyperparameters, how they
were chosen)? [Yes] See section 2 and 3.
(c) Did you report error bars (e.g., with respect to the random seed after running experiments multiple times)? [No](d) Did you include the total amount of compute and the type of resources used (e.g., typeof GPUs, internal cluster, or cloud provider)? [No]4. If you are using existing assets (e.g., code, data, models) or curating/releasing new assets...
(a) If your work uses existing assets, did you cite the creators? [Yes](b) Did you mention the license of the assets? [No](c) Did you include any new assets either in the supplemental material or as a URL? [No](d) Did you discuss whether and how consent was obtained from people whose data you’reusing/curating? [N/A](e) Did you discuss whether the data you are using/curating contains personally identifiable
information or offensive content? [N/A]
5. If you used crowdsourcing or conducted research with human subjects...
(a) Did you include the full text of instructions given to participants and screenshots, if
applicable? [N/A](b) Did you describe any potential participant risks, with links to Institutional Review
Board (IRB) approvals, if applicable? [N/A](c) Did you include the estimated hourly wage paid to participants and the total amountspent on participant compensation? [N/A]12