Lagrangian Crumpling Equations

0.16 MB·English

by  Mark A. Peterson

Lagrangian Crumpling Equations
9
0
0 Mark A. Peterson
2
Mount Holyoke College

n
a January 27, 2009
J

7
2
Abstract

]
t A concise method for following the evolving geometry of a moving
f
o surface using Lagrangian coordinates is described. All computations can
s be done in the fixed geometry of the initial surface despite the evolving
. complexity of the moving surface. The method is applied to three prob-
t
a lems in nonlinear elasticity: the bulging of a thin plate under pressure
m
(theoriginal motivationforF¨oppl-vonKarmantheory),thebucklingofa
- spherical shell underpressure, and the phenomenon of capillary wrinkles
d
induced by surface tension in a thin film. In this last problem the inclu-
n
sion ofagravitational potentialenergyterm inthetotalenergyimproves
o
the agreement with experiment.
c
[

1 1 Introduction
v
6
The elasticity theory of thin shells is largely differential geometry by another
6
name. Inthis paper I describe a method for followingthe differentialgeometric
1
4 dataofasurfaceasitmoves,andillustrateitsapplicationtonon-linearelasticity
. theory. Theequationsofthemethodarecompletelygeneralforsmoothsurfaces,
1
0 and so could in principle describe the complex motions of crumpling up to the
9 formation of singularities.
0 Problemsinvolvingelasticmembraneshavebeenapproachedinseveralways,
v: including numericalsimulationby triangulatedsurfaces,using a polyhedralap-
i proximation to differential geometry [1]. Another approach has been to use
X
differential geometry and scaling laws to understand the line and point singu-
r
larities of crumpled surfaces analytically [2, 3, 4, 5, 6, 7], and numerically[8].
a
The method of this paper generalizes familiar methods of mechanical engineer-
ing for the non-linear elasticity theory of thin shells [9, 10, 11] in going beyond
second order, and in treating initially curved surfaces in a unified way.
Section 2 establishes notation for the differential geometry of a moving sur-
face and shows how to use Lagrangian coordinates to simplify its description,
themainideaofthispaper. Section3summarizestheobservationsoftheprevi-
ous section in a system of differential equations for the evolving surface and its
strains. Section 4 compares this approachto F¨oppl-von Karman(FvK) theory,
1
andsolvesthe motivating problemfor thattheory,the bulging ofathin rectan-
gular plate subject to pressure, by integrating the evolution equations forward
in time. Section 5 uses second order expansions of the crumpling equations to
describe the buckling of a sphere under pressure. Section 6 uses the insights of
CerdaandMahadevan[12]togiveamoredetaileddescriptionofaphenomenon
recentlydiscussedin[13],capillarywrinklesinducedbysurfacetensioninathin
film. A previously unnoticed discrepancy with experiment is partially resolved
with the inclusion of the gravitationalpotential energy of the system.
2 Geometrical Methods
In terms of smooth coordinates (x1,x2,x3) in space one can describe the de-
formation of a material object by the trajectories of its constituent particles,
solutions of equations of motion
dxi
=Vi(x,t) (1)
dt
where
V =Vi∂ (2)
i
isthe vectorfieldgeneratingthe flow,andtisaparameteralongthe flow. Inte-
gratingthe systemforwardtot=1,onecanalsothink ofVi asadisplacement,
a slight abuse of notation that should be clear from context. Metric relations
among the particles are given by
ds2 =g dxidxj (3)
ij
where g is a Riemannian metric tensor, perhaps, but not necessarily, the Eu-
ij
clidean metric.
I coordinatize the material object by Lagrangiancoordinates, convected by
the flow, i.e., every material point keeps the same coordinates that it had orig-
inally. In this case the changing metric relationship of material points, namely
the changein the expressionEq.(3), is due entirely to the changein the metric
components g , because dxi, which for this purpose simply assigns to a line
ij
segment the coordinate difference of its endpoints, is invariant. The rate of
changeas a consequenceofthis deformationinthe components ofthe metric g,
or of any second rank tensor G, expressed in convectedcoordinates, is givenby
the Lie derivative[14][15]
£ G(∂ ,∂ )=VG +G([∂ ,V],∂ )+G(∂ ,[∂ ,V]) (4)
V j k jk j k j k
Here [ , ] is the Lie bracket of vector fields. It is more common to express
objectslikethis,derivativesoftensorswhicharethemselvestensors,intermsof
the covariant derivative with respect to the metric connection, and to employ
the conventions of raising and lowering indices with g and its matrix inverse
ij
gij, such that, for example the covector with components
V =g Vj (5)
i ij
2
has covariantderivative with respect to xk (denoted V ), in terms of the ordi-
i;k
nary partial derivative (denoted V ) given by
i,k
V =V +Γj V (6)
i;k i,k ik j
where the coefficients of connection Γ are
1
Γj = gjm(g g g ) (7)
ik 2 ik,m− mi,k− km,i
It is straightforwardto verify for any second rank tensor G that
µν
£ G(∂ ,∂ )=gij(V G +V G +V G ) (8)
V k ℓ j kl;i j;k iℓ j;ℓ ki
In particular, if G is the metric tensor g, which is a covariant constant, we
recover the well known result
£ g(∂ ,∂ )=V +V =2U (9)
V k ℓ ℓ;k k;ℓ kℓ
where U is the rate of strain tensor of the flow V (or the first order strain of
the displacement V). Nothing said above was specific to three dimensions, and
therefore every statement can be interpreted as referring to a surface with a
Riemannian structure if the indices take only two values and not three. From
now on I shall use Latin indices for tensors in three-space, and Greek indices
for tensors on a surface.
Now consider a smooth material surface M, so thin that one may regard it
as 2-dimensional,andlet (x2,x3) be coordinatesin this surface, while x1 =z is
displacementalongthenormaltothesurface,withthepositivedirectionchosen
conventionally, such that the surface M is z = 0. Such a coordinate system
exists for a neighborhood of M such that z <1/C, where C is the supremum
| |
over M of both principal curvatures in absolute value. The metric tensor in
these coordinates takes the form
1 0
g = . (10)
0 g +2zh +z2k
µν µν µν
(cid:18) (cid:19)
The tensor g , with Greek indices taking values (2,3), is the first fundamental
µν
form of M, h is the second fundamental form, and k = hλh is the third
µν µν µ λν
fundamental form. All these tensors are associated with the surface M, and
not with the ambient space. They do not depend on z, i.e. all z dependence in
Eq. (10) is explicit. The plus sign on the middle term is a conventional choice.
On a sphere, for example, one could take the positive direction for z to be the
outernormaldirection,andtheprincipalcurvaturesofthespheretobepositive.
Now let a vector field (a,Vµ) be prescribed on M with normal component
a(x2,x3) and tangential components Vµ(x2,x3), and extend it to a neighbor-
hood of M as
W =a∂ +Vµ∂ zG′µνa ∂ , (11)
z µ ,µ ν
−
whereinitiallythetensorG′µν =gµν. Inashorttime ∆t,the flowgeneratedby
the velocity field W changes the metric tensor components by approximately
∆g =∆t£ g (12)
W
3
The tensor g+∆g regarded as a tensor on 3-space expresses the ambient Eu-
clideangeometryinLagrangiancoordinates. Ifg+∆gisrestrictedtothesurface
z = 0 and indices (2,3), one has the slightly altered first fundamental form of
M
G =(g +∆g ) , (13)
µν µν µν |z=0
expressingthe non-Euclideangeometry ofthe slightly alteredM induced by its
embedding in the ambient Euclidean space. The term linear in z in Eq. (11)
was chosen to maintain the block diagonal form of Eq. (10) to first order in z.
Therefore, taking the z-derivative, one has the slightly altered second funda-
mental form of M,
∂
H = (g +∆g ) , (14)
µν µν µν
∂z
(cid:20) (cid:21)|z=0
The third fundamentalform couldnotbe computed in this way,but it is deter-
mined by H ,
µν
K =H Hλ . (15)
µν µλ ν
I now imagine taking a sequence of such small steps, and I will continue to
denote by G andH the evolvingfirstandsecondfundamentalforms giving
µν µν
the Riemannian structure on M induced by the embedding in Euclidean space.
I will not make use of this Riemannian structure for computations, however.
There is anothernaturalRiemannianstructureonM, namely thatgivenby
theoriginal,undeformedfirstfundamentalformg ,togetherwithitsassociated
µν
connection, etc., which I shall continue to use, being careful not to give it erro-
neous interpretations. This Riemannian structure, unlike G , has no obvious
µν
geometrical meaning on the deformed surface, but it is still useful in a formal
way. Another possible interpretation, deliberately suppressing the geometrical
meaningofG ,is toimagineasurfacethatis notdeformedbythe flowW but
µν
carries tensor fields G and H , initially coinciding with g and h , that
µν µν µν µν
aredeformed by W. That these happen to be the first andsecondfundamental
formsofanevolvingsurfaceisforgotten. Inthispicturetheundeformedg has
µν
an obvious geometrical meaning as the metric on the underlying, unchanging
surface which is the arena for the evolving G and H .
µν µν
In Eq. (11) I introduced the tensor G′µν, initially gµν. More generally G′µν
istheinverseofG asamatrix. ItisatensorfieldonM,butitisnotobtained
µν
from G by raising indices. Raising indices is an operation accomplished by
µν
gµν, my chosen Riemannian structure, not by G′µν. The prime on G′ is a
reminder that it is not some version of the tensor G.
I have shown how G and H change, to first order, under a deformation
µν µν
(a,Vµ) of M, assumed now always to be extended off M as in Eq (11). In
turn, (a,Vµ) might evolve so as to reduce at each step a free energy functional
depending on G and H . In this way I will arrive at crumpling equations,
µν µν
a system of differential equations for (a,Vµ), G , and H , describing the
µν µν
evolution of M. Before considering the equation for (a,Vµ), though, there is
another issue to consider.
4
This formulation leaves implicit what the evolving surface actually looks
like, since mere knowledge of G and H is not a convenient description of
µν µν
M. To keeptrackofthe positions ofpoints onthe surface,one shouldintegrate
Eq.(1)usingcomponentsofW(0,xµ)=(a,Vµ)withrespecttofixedCartesian
coordinate axes. Let XA(x1,x2,x3,t) be a Cartesian coordinate function in
space. It is time independent in the physical sense, but its functional form
depends on time because the xi evolve in time. The 1-form dXA = XA dxj
,j
assigns the XA component WA to the vector W. This 1-form evolves in time
at the rate given by the Lie derivative
£ dXA(∂ ) = W dXA(∂ )+dXA([∂ ,W]) (16)
W i i i
= WjXA +Wj XA (17)
,ij ,i ,j
= (XA Wj) (18)
,j ,i
Thus UA = UiXA, the XA-component of any vector field U = Ui∂ at time t,
i i
can be found using XA(x1,x2,x3,t) solving
i
∂XA ∂(XAWj)
i = j (19)
∂t ∂xi
with appropriate initial conditions. By the definition of the coordinate x1 =z,
the Cartesian coordinate XA is an affine linear function of z. It is essential
thereforetoexpandXAWj onlytofirstorderinz inEq.(19). Tobecompletely
j
explicit, XA is independent of z and we can represent
1
XA =YA(x2,x3)+zZA(x2,x3). (20)
µ µ µ
Then Eq. (19) says
∂X1A = ZAVµ YAa G′µν (21)
∂t µ − µ ,ν
∂YA
µ = (XAa+XAVν) (22)
∂t 1 ν ,µ
∂ZA
µ = (ZAVν YAa G′νλ) (23)
∂t ν − ν ,λ ,µ
The linear approximation I have made in the neighborhood of M obscures the
fact that if W were made to carry affine normal lines to affine normal lines
exactly, as one could always require by a suitable nonlinear extension W of
(a,Vµ) off M, then XAWj would be exactly an affine linear function of z
j
without approximation. The evolution of M is the same for any extension,
however, so what looks like a linear approximation in the method is actually
exact.
As a special case, I describe motion at constant velocity, i.e., ∂WA/∂t = 0
for each component A. Then
∂ XAWj ∂Wj
0= j =(XAWk) Wj +XA (24)
∂t k ,j j ∂t
(cid:0) (cid:1)
5
Thus the components of W must evolve according to
∂Wk
= Xk(XAWj) Wi (25)
∂t − A ,j ,i
Here Xk is the inverseofXA, consideredas a matrix. Eq.(25)for straightline
A j
motion is recognizable as
∂Wk
+Wj Wk =0 (26)
j
∂t ∇
where isthecovariantderivativewithrespecttothemetricconnectionofthe
k
∇
Euclidean metric in 3-space expressed in the evolving Lagrangian coordinates.
I emphasize that I have chosen, however, not to use the evolving geometry but
rather the fixed initial geometry of M for all computations, a great simplifica-
tion.
3 Evolution Equations
By the arguments of the previous section the surface M evolves according to
∂G
κλ = VµG +Vµ G +Vµ G +2aH (27)
∂t κλ;µ ;κ µλ ;λ µκ κλ
∂Hκλ = aK a + 1a G′µν( G +G +G )
κλ ,λ;κ ,µ κλ;ν νλ;κ νκ;λ
∂t − 2 −
+VµH +Vµ H +Vµ H (28)
κλ;µ ;κ µλ ;λ µκ
(29)
Using these relations one can find how other geometric quantities change, for
example the area element √Gdx2dx3, involving the determinant of the first
fundamental form
G=G G G G (30)
22 33 23 32
−
The result is
∂√G =(Vµ√G) +aG′µνH √G (31)
,µ µν
∂t
Integrating one finds √G and hence dilation strain. The strain tensor
1
(G g ) (32)
µν µν
2 −
can be found by integrating Eq. (27). A natural definition for nonlinear shear
strain S is
µν
∂S 1 ∂G 1 ∂√G
µν µν
= G . (33)
µν
∂t 2 ∂t − √G ∂t !
The subtracted term removes the contribution of dilation strain. S is not
µν
traceless, in general, beyond first order.
6
4 Comparison with F¨oppl-von Karman Approach
A simple example illustrates the use of this formalism and points out its rela-
tionshiptoF¨oppl-vonKarman(FvK)theory[9]. FvKconsiderstheequilibrium
state of a thin membrane subject to external forces and boundary conditions.
Since the metric strain within a membrane is typically small, even for large
normal displacements, it makes sense to continue to use linear stress-strain re-
lationships. The strain may, however, be a nonlinear function of displacement,
and therefore displacement may be nonlinearly related to stress. FvK thus
produces nonlinear equations for the equilibrium shape of an elastic membrane
subject to external stress.
Historically this idea was implemented by expanding the strain tensor to
first orderin tangentialdisplacement but secondorder in normaldisplacement.
I derive the FvK strain by solving the evolution equations to first order in Vµ
and secondorder in a, continuing to use the notation of previous sections, with
the initial velocity vector
W(0) =a∂ +V(0)µ∂ zG′µνa ∂ (34)
z µ ,µ ν
−
of Eq. (11). I am using the superscript(0) to indicate the initial value, which is
also the zeroth approximationfor an iterative solution. Other initial values are
g = G(0) = δ and h = H(0) = 0. I use Picard’s method to generate the
µν µν µν µν µν
solution to the differential system Eqs. (19), (25, (27), and (28) iteratively as a
power series in t, taking M to be the Euclidean plane with the usual Cartesian
coordinates. In this case there is no distinction between indices up and indices
down,andcovariantderivativesareordinarypartialderivatives. Iteratingonce,
and ignoring quadratic terms except in a gives
G(1) = δ +t(V(0)+V(0)) (35)
κλ κλ κ,λ λ,κ
H(1) = ta (36)
κλ − ,λκ
V(1) = taa (37)
µ ,µ
Iterating a second time, still ignoring quadratic terms except in a, gives
G(2) =δ +t(V(0)+V(0)) 2ta +t2a a (38)
κλ κλ κ,λ λ,κ − ,κλ ,κ ,λ
Finally, evaluating at t=1, gives the FvK metric strain
1 1
(G(2) δ )= (V(0)+V(0)+a a ) a (39)
2 κλ − κλ 2 κ,λ λ,κ ,κ ,λ − ,κλ
ThisisthecomputationalstartingpointforFvKtheory. Therestofthattheory
follows from minimizing the elastic energy, expressed as a quadratic functional
of this strain and the first order bending strain H(1), to find the equilibrium
µν
shape.
The approachofthis paper isto developthe nonlinearstrainasthe solution
toadifferentialsystem. FromthatpointofviewthederivationofEq.(39)isnot
7
very natural, since to obtain it one must artificially impose the condition that
the trajectories of the particles are straight lines, a condition that introduces,
viaEq.(37),asecondordercorrectionintothestrainthatisnecessarytoobtain
Eq. (39). Although one can certainly parameterize the possible final shapes of
M by displacementofparticles alongstraightlines, it is a differentthing to say
that particles actually move along straight lines. FvK theory does not claim
this, and in that sense it is not a dynamical theory. A dynamical theory would
determine the evolution of the velocity vector (a,Vµ) by some local physical
law, replacing Eq. (25) in the differential system. It would be a simpler theory,
both conceptually and computationally, in that solving it would only require
integrating a differential system forward in time. I will do the obvious thing
and choose W to reduce the elastic energy at each step, seeking the minimum.
A typical phenomenological elastic energy functional is
E =E +E +E (40)
d s c
where
2
Λ √G
E = 1 √gdx2dx3 (41)
d
2 ZM √g − !
E = µ SκλS √gdx2dx3 (42)
s κλ
ZM
E = κ G′µνH gµνh 2√gdx2dx3 (43)
c µν µν
2 −
ZM(cid:16) (cid:17)
andwhereΛ,µ,κarethe2Dcompressionmodulus,shearmodulus,andbending
modulus respectively. The area element involves √g, not √G, because the
energy due to metric strain is better understood to be per unit mass, not per
unitarea,andthe massisconvectedwiththematerialcoordinates. Thesystem
will move,if possible, to lower its energy, so one must compute the variationof
E with respect to a small normal displacement δa and tangential displacement
δVµ
δE δE
δE = δVµ+ δa √gdx2dx3 (44)
δVµ δa
ZM(cid:20) (cid:21)
The work done on M in deforming it represents energy given up by some other
partof the system, so this workshould be added with a minus signto the total
change in energy. Work done by pressure P in a small normal deformation δa,
for instance, is
W =P δa√Gdx2dx3 (45)
ZM
where now one must use the physical area element √Gdx2dx3 on M. A small
displacement in the direction opposite to this “gradient,” i.e.
δE √G
a = L +P (46)
a
−δa √g!
δE
Vµ = L (47)
V −δVµ
(cid:18) (cid:19)
8
willlowerthe energyandmovethe systemtowarda localminimum. The linear
operators L and L include a projection onto the space of admissible veloc-
a V
ity vector fields. They must define positive semi-definite quadratic forms with
respect to the inner product givenby integrationoverM. Apart from these re-
quirements,they will varywith the application. This is just the familiarnotion
of conjugate gradient. One could also think of L and L together as defin-
a V
ing a generalized mobility tensor, because it transforms generalized force into
velocity. If one only wants to know the final state, one could try to choose L
a
and L so as to reach equilibrium in the most efficient way. In any case, the
V
dynamics of the system is not completely determined by the elastic energies,
andadditionalphysicalconsiderationsmustbeaddedto completethetheoryin
a specific application.
Eqs. (46) and (47), together with the evolution equations of Section 3, are
what I mean by Lagrangian crumpling equations. The original problem ad-
dressed by FvK theory, the bulging of a square plate fixed on the boundary
and subject to pressure,can be solvedstraightforwardlyin this way. Represent
all geometric data by discretization on a square grid of points of the original
square. Spectral methods (fast Fourier transform with anti-aliasing) make the
computation efficient, and the gradient flow converges quickly to a solution.
5 Buckling of a sphere under pressure
I consider an elastic spherical shell subject to pressure P, described by the
phenomenological energies of Eqs. (41), (42), (43), and (45). For small enough
pressurethesphereisuniformlycompressed,butaspressureincreasesitbuckles.
I will describe the buckling by using expansions of strain to second order in
displacement where necessary, not the FvK expansion, but the “dynamic” one
ofthis paper,foundby solvingthe crumplingequationsiteratively. Itturns out
that the expansion must include more terms than FvK.
For a sphere of radius R, in terms of spherical polar coordinates (θ,φ),
g = diag(R2,R2sin2θ) (48)
µν
h = g /R (49)
µν µν
k = g /R2 (50)
µν µν
TakingR=1,andregardingallquantities nowasdimensionless,the perturbed
geometric quantities in a general displacement (a,Vµ) are
G = g +V +V +2ag (51)
µν µν µ;ν ν;µ µν
H = g +ag a +V +V (52)
µν µν µν ,µ;ν µ;ν ν;µ
−
1
√G = √g[1+(Vµ +2a)+ (Vµ Vν +VµVν
;µ 2 ;µ ;ν ;ν;µ
+ 4aVµ +2Vµa aaµ +2a2)] (53)
;µ ,µ− ;µ
Theareaelement√Ghadtobe foundtosecondorderindisplacement. Tofirst
9
order in displacement the shear strain in the sphere is
1
S = (V +V Vλ g ) (54)
µν 2 µ;ν ν;µ− ;λ µν
Parameterize the displacement by coefficients (a ,b ,c ), such that
ℓm ℓm ℓm
a = a Y (55)
ℓm ℓm
ℓm
X
Vµ = gµν b Y +ǫµν c Y (56)
ℓm ℓm,ν ℓm ℓm,ν
ℓm ℓm
X X
where the Y are spherical harmonics and ǫ = ǫ =sinθ, ǫ =ǫ =0 is
ℓm 32 23 22 33
−
the antisymmetric tensor. Then for example the change in the mean curvature
of the perturbed sphere is
δH =GµνH gµνh = [ℓ(ℓ+1) 2]a Y (57)
µν µν ℓm ℓm
− −
ℓm
X
so that the curvature energy is
κ
E = [ℓ(ℓ+1) 2]2 a 2 (58)
c ℓm
2 − | |
ℓm
X
Itvanishes for ℓ=1, asit must by Galileaninvariance,andit is independent of
the tangential displacement Vµ. The other energy expressions are
W = 4πPa Y P ℓ(ℓ+1)a b +2P a 2
00 00 ℓm ℓm ℓm
− | |
ℓm ℓm
X X
1
+ Pa Y [ 2ℓ(ℓ+1)a b +ℓ(ℓ+1)a 2+2a 2] (59)
00 00 ℓm ℓm ℓm ℓm
2 − | | | |
ℓm
X
Λ
E = [ ℓ(ℓ+1)b +2a ]2
d ℓm ℓm
2 −
ℓm
X
+ Λa Y [ 2ℓ(ℓ+1)a b +ℓ(ℓ+1)a 2+2a 2] (60)
00 00 ℓm ℓm ℓm ℓm
− | | | |
ℓm
X
µ
E = ℓ(ℓ+1)[ℓ(ℓ+1) 2](b 2+ c 2) (61)
s ℓm ℓm
2 − | | | |
ℓm
X
These expansions have been carried out to second order in all coefficients, but
they anticipate that a is the same order as a 2 for ℓ > 1, so that some
00 ℓm
| |
terms quadratic in a appear to be third order. I also anticipate that the first
00
response to pressure is a uniform compression
P
a (62)
00
| |∼ Λ
sothatconsistencyrequiresP <<Λ. Nowseektheminimumofthetotalenergy
E =W +E +E +E (63)
tot d s c
10