Mechanics of moving defects in growing sheets: 3-d, small deformation theory

Growth and other dynamical processes in soft materials can create novel types of mesoscopic defects including discontinuities for the second and higher derivatives of the deformation, and terminating defects for these discontinuities. These higher-order defects move “easily", and can thus confer a great degree of flexibility to the material. We develop a general continuum mechanical framework from which we can derive the dynamics of higher order defects in a thermodynamically consistent manner. We illustrate our framework by obtaining the explicit dynamical equations for the next higher order defects in an elastic body beyond dislocations, phase boundaries, and disclinations, namely, surfaces of inflection and branch lines.


Introduction
Hyperbolic sheets abound in nature (see Fig. 1). As Margaret Wertheim writes in her delightful essay "Corals, crochet and the cosmos: how hyperbolic geometry pervades the universe" (Wertheim 2016) -We have built a world of largely straight lines -the houses we live in, the skyscrapers we work in and the streets we drive on our daily commutes. Yet outside our boxes, nature teems with frilly, crenellated forms, from the fluted surfaces of lettuces and fungi to the frilled skirts of sea slugs and the gorgeous undulations of corals.
A natural question is -Why these shapes? One suggestion is that cells in living organisms proliferate to "maximize" their number (area) subject to any applicable constraints (Wertheim 2016) and this naturally results in hyperbolic geometries. This is a "static" argument which relates the mechanisms of growth to the resulting (quasi-2D) intrinsic geometry of living organisms. In this paper, we attempt to go beyond this "static" argument and develop models, based on thermodynamic considerations, to gain a quantitative understanding of the interplay between growth, mechanics and dynamics in soft objects. These models have the potential to describe the dynamical processes that result in the observed intricate three-dimensional (i.e. extrinsic) morphologies in nature.
Particularly striking examples of dynamical behaviors in organisms with differential growth (hyperbolic geometries) occur in sea slugs (Nudibranchia) and marine flatworms (Polycladida). These marine invertebrates are found in many environments, particularly in coral reefs. While most of them crawl on the sea floor, a few are capable of free swimming (Newman 2003). They move/swim by sending waves of undulations from the front to the back along their skirts (for sea slugs) or across their entire body (for flatworms). Figure 2 shows 4 frames from a video of a free-swimming sea slug (Jones 2010). The geometry of the slug is clearly hyperbolic. It has multiple undulations and undergoes significant bending deformations in the course of one swim cycle. While it is hard to quantify the strains within the organism it is not unreasonable to consider them small in comparison to the obvious large rotations/twist of the body.
In a different context, the interplay between growth and dynamics is also relevant to the development of leaves, flowers and other plant tissues that can be modeled as thin laminae (Liang and Mahadevan 2009;Boudaoud 2010;Liang and Mahadevan 2011;Goriely 2017;Sharon and Sahaf 2018). Laboratory experiments using hydrogels (Klein et al. 2007;Kim et al. 2012;Levin et al. 2019) have led to a semi-quantitative understanding of timedependent, dissipative deformations of thin soft materials with a prescribed prestrain. In living organisms, however, the prestrain is not prescribed a priori, and how the prestrain development may be related to mechanics is not clear. Complicated physico-chemical processes are involved that need to be incorporated into mathematical models. It therefore seems reasonable to derive systematic constraints on the mathematical description based on a careful consideration of the non-standard kinematics involved and the general principles of continuum thermomechanics.
Earlier work, reviewed briefly in "Statics and equilibria of non-Euclidean elastic sheets" section, implicates higher-order defects, in contrast to disclinations and dislocations, as playing a key role in the mechanics of intrinsically hyperbolic elastic sheets (Gemmer and Venkataramani 2011;Gemmer et al. 2016; Shearman and  Venkataramani, in preparation). This points to the need for tools to describe the evolution of (terminating) discontinuities of the second-gradient of the displacement fieldwhen viewed at the macroscopic scale -for a proper description of the soft material deformations involved. It turns out that, within a continuum mechanical perspective, this fits in nicely within the question of describing the coupled mechanics of discontinuities and singularities of the elastic displacement field and its higher derivatives up to order three. This is the question that is addressed in this paper.
While, as evident from Fig. 2, it is natural, and necessary, to consider unrestricted finite deformations when dealing with soft materials, we restrict attention to 'small deformation' kinematics in this first effort due to the extra subtleties involved with higher order defect kinematics. Thus, we consider deformations of a fixed reference configuration that may or may not be stress-free. When the configurations attained by the system remain in close proximity to this fixed configuration, this is an adequate assumption. We consistently invoke Occam's razor as a guiding principle in our development -for instance, we restrict to the use of only ordinary stresses and couple stresses since forces and moments are the only agents of mechanical stimuli that we have some intuition for. Similarly, if branch point and surface defect velocities are to be the only dissipative mechanisms requiring constitutive specification without involving their spatial derivatives, then it turns out that the appropriate variables for the analysis of thermodynamics is in terms of the 'singular parts' of the first and higher order displacement gradients, instead of the more natural singular parts of the corresponding elastic distortion gradients that naturally arise in the analysis of defect kinematics. This is in sharp contrast to dislocation and g.disclination mechanics (Acharya and Fressengeas 2015) where this distinction does not arise because of the relatively lower order kinematics involved. We develop the relationship between the two types of entities in this paper.
In this article, we develop a framework from which we can derive the dynamics of higher order defects in a thermodynamically consistent manner. In particular, the framework is applicable to the motivating problems in this Introduction, namely the dynamics of Non-Euclidean (i.e. incompatible) sheets. Our framework can accommodate the specific details of the continuum mechanics of various applications. It can describe the evolution of defects in incompatible elasticity and the elastodynamics of growing bodies. It can also describe plastic deformation resulting from moving dislocations and disclinations and analogous behaviors from the motion of higher order defects. We expect that, with the appropriate choices of thermodynamic potentials and kinetic coefficients, our framework will be useful for problems governed by the interplay between growth, defects, thermodynamics and the balance laws of continuum mechanics.

Statics and equilibria of non-Euclidean elastic sheets
One approach to modeling the mechanics of a growing hyperelastic body, borrowed from the literature of finite elastoplasticity, is to assume a reference configuration S and a deformation y : S → R 3 along with a multiplicative decomposition of the deformation gradient F = ∇y as F = EG (or F = F e F p in the plasticity literature) where the two-point tensor G models the effect of the growth processes in the material and E is the "residual" elastic deformation (Goriely 2017). The energy of the configuration defined by y is then given by W (E) = W (FG −1 ), where W denotes a hyperelastic energy density, vanishing on SO(3) (Goriely 2017;Lewicka et al. 2014). In particular, the material is "stress-free" if FG −1 ∈ SO(3), although, for a general G, there might not be any deformation of the body y : S → R 3 whose gradient is a rotation times G. Such objects, with no stress-free configurations in R 3 , lead to incompatible elasticity.
The non-Euclidean formalism of thin sheet elasticity (Efrati et al. 2009) is a reduced dimensional description of thin elastically incompatible objects. The reference manifold , where t, the thickness, is "small" compared to the "in-plane" dimensions of the center surface ⊂ R 2 . In this setting, the effect of the growth has a reduced dimensional description as a 2-manifold ( , g, b) where g, b are symmetric (0, 2) tensors. These tensors denote, respectively, the 'target' 1st and 2nd fundamental forms of the stress-free state of the sheet, pulled back to the reference manifold (Efrati et al. 2009). This framework also applies to incompatible elasticity (Ben Amar and Goriely 2005;Lewicka et al. 2014;Bhattacharya et al. 2016) where, in general, there exists no deformation f : → R 3 realizing a surface in ambient three-dimensional space whose first and second fundamental forms match (the push-forward of ) the targets g, b (by f ), i.e., incompatible sheets have no stress-free configurations in our three dimensional space.
Assuming the Kirchhoff-Love hypothesis (Fung 1965), so that the (3D) deformation of a thin sheet is determined by the (2D) mapping y : → R 3 on the center surface.
This allows for an asymptotic expansion of the elastic energy as a sum of stretching and bending contributions (Efrati et al. 2009;Efrati et al. 2013) : where the oriented normal field N : → S 2 , also called the Gauss Normal map (Stoker 1989), is obtained from ∇y T ·N = 0. Q 3 is a non-degenerate quadratic form, on symmetric 2 × 2 matrices, that depends on the Poisson's ratio ν of the material (Efrati et al. 2009), and dA is the area element on ( , g). For various choices of g and b and boundary conditions, the energy functional (1) describes a variety of phenomena in thin sheets, including multiple-scale buckling in free sheets with 'excess length' near an edge, e.g. torn plastic or flat leaves treated with an Auxin near the edge (Sharon et al. 2002;Sharon et al. 2004). The excess length near the edge is modeled by a metric g with negative intrinsic curvature (Efrati et al. 2013).
Starting with a fully 3D elastic energy, Lewicka and Pakzad (Lewicka and Reza Pakzad 2011) have obtained a reduced dimensional model for the limit t → 0 using -convergence. In particular, they showed that for an appropriate quadratic form Q 2 . This energy has clear similarities with the energy in (1), although the details are somewhat different. Nonetheless, in either framework, the elastic energy scales like t 3 in the thin limit t → 0 if and only if there exist finite bending energy (mathematically y ∈ W 2,2 ) isometric immersions y : ( , g) → R 3 . What is the physical import of this theorem? The Gagliardo-Nirenberg-Sobolev inequality (cf. Evans (1998)[ §5.8.1]) implies that W 2,2 surfaces with finite bending content have a "nearly continuous" tangent plane and normal, in the sense they cannot oscillate "much" on small sets. More precisely, the normal map N : → S 2 is in BMO( ), and the John-Nirenberg inequality (John and Nirenberg 1961) for BMO functions rules out singularities that correspond to O(1) oscillations in the normal on arbitrarily small sets, including sharp creases (folds), cone points (disclinations) or dislocations. Indeed the energy of elastic ridges (Lobkovsky 1996;Venkataramani 2004;Conti and Maggi 2008), E t ∼ t 8/3 , and d-cones (Ben Amar and Pomeau 1997;Cerda et al. 1999;Olbermann 2016), E t ∼ t 3 log(1/t), diverge on the scale t 3 , although the limiting shapes are 'asymptotic' isometries (Vella et al. 2015;Davidovitch et al. 2019) and arguably unstretched.

Branch points and lines of inflection
The preceding remark highlights the role of the regularity of isometries. Beyond the existence/non-existence of isometries, it is crucial whether a candidate isometry is in W 2,2 . This motivates the problem: We have rigorous results showing that the problem (2) is flexible and solutions are plentiful (Gemmer et al. 2016; Shearman and Venkataramani, in preparation) (with prescribed zero-traction and moment boundary conditions, i.e. for free sheets). The proof is constructive, and uses ideas from Discrete Differential Geometry DDG (Bobenko et al. 2008;Gemmer et al. 2016). This lack of uniqueness in admissible static configurations with prescribed boundary conditions underscores the necessity of a dynamical model to 'choose' between acceptable configurations and/or describe the transitions between multiple admissible states (Gemmer and Venkataramani 2013). If y : → R 3 is C 1 , the Gauss Normal map is given by N = ∂ 1 y×∂ 2 y ∂ 1 y×∂ 2 y , where ∂ i = ∂ ∂x i for (arbitrary) coordinates (x 1 , x 2 ) on . Further, if y and g are C 2 , Gauss' Theorema Egregium implies that (2) is equivalent to the Monge-Ampere Exterior differential system (EDS) (Ivey and Landsberg 2003, §6.4): where d is the area form on the sphere S 2 and κ is the Gauss curvature. Classical results in differential geometry imply that smooth solutions of (3) with κ < 0 are hyperbolic surfaces and locally saddle shaped. In contrast, the curly mustard leaf in Fig. 3b is "frilly", i.e buckled on multiple scales with a wavelength that refines ("subwrinkles") near the edge (Sharon et al. 2004). This "looks" very unlike the smooth saddle in Fig. 3a. If ⊂ R 2 is a bounded domain with a smooth boundary, and g is a smooth metric on with negative curvature, g can be extended to a smooth metricḡ on R 2 with Gauss curvature κ[ḡ] < 0 decaying (as rapidly as desired) at infinity. The existence of isometric immersions into R 3 , of smooth metrics with decaying negative curvature (Hong 1993), therefore implies that bounded smooth hyperbolic surfaces can be smoothly and isometrically embedded in R 3 . A smooth (C 2 is sufficient) hyperbolic surface cannot refine its buckling pattern and is thus "non-frilly" (Gemmer et al. 2016;Shearman and Venkataramani, in preparation). Why do we see frilly shapes in natural surfaces, as in Fig. 3b, rather than the smooth saddles of Fig. 3a?  Fig. 3 Hyperbolic surfaces in R 3 . The inscribed (geodesic) triangle in the smooth saddle has angles that sum up to less than π, illustrating the connection between the extrinsic and intrinsic geometries -Gauss' Theorema Egregium We have addressed this puzzle in recent work (Gemmer and Venkataramani 2011;Gemmer and Venkataramani 2012;Gemmer and Venkataramani 2013;Gemmer et al. 2016;Shearman and Venkataramani, in preparation) and the short answer is that, for a given metric g, the frilly surfaces, somewhat counterintuitively, can have smaller bending energy than the smooth saddle. It is true that C 2 (twice continuously differentiable) hyperbolic surfaces are saddle-like near every point. We find a topological invariant (Shearman and Venkataramani, in preparation), the index of a branch point -intimately related to the quantity α (3) n da that emerges in "Kinematics" section and the quantity of "The discontinuity of the deformation of a non simply connected domain with prescribed third 'deformation gradient" section -that distinguishes sub-wrinkled surfaces from saddles locally. With branch points, the surfaces are only C 1,1 , but gain the additional flexibility to refine their buckling pattern, while lowering their energy (Gemmer et al. 2016). This flexibility is not available to smooth saddles, and constitutes a key property of branched (sub-wrinkled) surfaces (Gemmer et al. 2016;Shearman and Venkataramani, in preparation). Figure 4 shows the construction of a non-C 2 monkey saddle (Gemmer and Venkataramani 2011). starting from the quadratic surface w = x 2 − 3y 2 . Cutting out the sector |x| ≤ √ 3|y| and then patching congruent copies of this sector by odd reflections The existence of such straight lines is a consequence of hyperboloids being doubly ruled surfaces gives a W 2,2 surface with a continuous normal vector and bounded curvature. The "defects" in this surface include the point in the middle -a branch point and the 6 rays through this pointlines of inflection, which together constitute the asymptotic skeleton of the surface (Shearman and Venkataramani, in preparation). This construction can be extended to generate C 1,1 hyperbolic surfaces with multiple distinct branch points, and an interesting question is how these defects interact with and influence each other (Gemmer et al. 2016).
Defects are of course ubiquitous in condensed matter systems. A key feature of defects in systems driven by a free energy is that the energy density typically diverges in the vicinity of a "bare" defect (and in some cases even the total energy diverges), and as a consequence, defects are always regularized, i.e. "cored" in physical systems. This is true for dislocations and disclinations in elastic objects, for creases in crumpled sheets, for defects in liquid crystals and many other types of defects. Uniquely, branch points and lines of inflection do not carry a singular energy density (Gemmer and Venkataramani 2011;Gemmer et al. 2016), and thus do not "need" a core for energetic reasons. Nonetheless, force and moment balance implies that these defects are indeed regularized into boundary layers, of width t 1/3 , mediating jumps in the normal curvature (Gemmer and Venkataramani 2012).
Branch points and lines of inflection are thus mesoscopic defects. They contain large numbers of atoms (microscopic units) and are amenable to a continuum description, but are yet much smaller than the typical size of the sheet. Arguments from energy minimization, while implying the existence of these higher order defects, do not address the question of their evolution. One has to necessarily go beyond the elastic energy (1) and incorporate dissipative effects that are crucial in determining a thermo-mechanically consistent description of the coupled evolution of the shape y : → R 3 and the internal geometry, given by the tensors g and b.

Notation
We define the notation employed in the paper in one place for convenience. When a function on a domain is discontinuous across a (non-planar) surface S, we assume that its values along any sequence of points from either side of the surface approaching any fixed point on the surface take on a unique pair of limiting values, each element in the pair corresponding to the limit from one side. The difference of these limiting values, one for each point on the surface, is defined as the jump (denoted by · ) of the function on the surface.
We think of an n th order tensor as a linear transformation between the space of vectors (in the translation space of three-dimensional Euclidean space, also 1 st -order tensors) to the space of (n − 1) th -order tensors, with its transpose defined in the natural way as being a linear transformation from the space of (n − 1) th -tensors to the space of vectors. All tensors components will be written w.r.t. the basis, (e 1 , e 2 , e 3 ) of a fixed rectangular Cartesian coordinate system and all partial derivatives, denoted often by a subscript comma, will be w.r.t coordinates of this system. The Einstein summation convention will be used unless otherwise stated. Superposed dots will represent partial derivatives w.r.t. time. If A is a p th -order tensor then the operators ∇, div, curl may be defined as Page 8 of 29 (with invariant meaning independent of the choice of coordinate system and its basis, of course). The range of all indices above is 1 to 3 and e ijk represents a component of the third-order alternating tensor.
The symbol · i represents a contraction on i indices between two tensors. For any tensor A, we define the tensor obtained by symmetrizing in the first two indices as A (s) and the one obtained by antisymmetrizing in the first two indices from the left as A (a) . We denote the deviatoric part of a second-order tensor by the superscript dev.

Motivation for kinematics of the theory
In this section we provide some intuition on the defect kinematics we adopt for our theory of branch point singularities. This is first done by explicitly constructing a continuously differentiable deformation of a non-simply connected domain whose second derivative has a prescribed, constant jump across a planar surface in the body.
With reference to Fig. 5, we think of occupying a simply connected d = 2 or 3dimensional domain of ambient Euclidean space. Here, it may be viewed either as a right-cylinder (d = 3) or a cross-section perpendicular to its axis (d = 2). We choose a rectangular Cartesian coordinate system with the z-axis as the axis of the cylinder; e i , i = 1, 2, 3 are the unit vectors along the x, y, z directions, respectively. c is a cylindrical subset of with rectangular cross-section centered on the z-axis. The region We will refer to c as a core.
Our goal in this section is to construct a vector fieldỹ (l) piecewise-smooth, and the jump in ∇ 2ỹ(l) across S a specified constant, with the jump blowing up as l → 0 maintaining lim l→0 ∇ỹ (l) ∈ C 0 ( h ).
A necessary condition forỹ (l) ∈ C 1 ( h ) with ∇ỹ (l) piecewise-smooth is that the jump in its second derivative across S be of the form ∇ 2ỹ(l) = A ⊗ ν, where ν is the unit normal field on S (with arbitrarily chosen orientation) and A is a R n×d valued matrix field on S. Noting that l −1 ∇ 2ỹ(l) may be formally considered an approximate discrete directional derivative of ∇ 2ỹ(l) in the direction ν (if the discontinuity were ignored), we define the field with A and ν = e 2 constants, and seek to construct solutions to the equations The restriction of Z to D is used since, even though Z is (distributionally) curl-free in h (we interpret the curl of a matrix field as row-wise curls), h is not simply connected but D is and hence we are guaranteed a solution Y in D, unique up to a constant. For any such Y field, we note that curl Y = 0 in D by the symmetry in the last two entries of Z, i.e. (∇Ye l )e k − (∇Ye k )e l = (Ze l )e k − (Ze k )e l = 0. Thus W satisfying (5) can be constructed, and W is also unique in D up to a constant for a given Y.
Arbitrarily fix one of the available Y fields. Such a Y has the explicit representation for x 0 being any point on the surface S, and the line integral is along any path from x − 0 to x contained in D. Now choose any path going from x − 0 to x + 0 (see Fig. 5) contained in D with the stipulation that it go through the points x 0 ± l 2 e 2 and the segments between x 0 and x 0 ± l 2 e 2 , respectively, are parallel to e 2 . We will assume that Y (x(s)) remains constant along the path between x 0 − l 2 e 2 and x 0 + l 2 e 2 . Therefore, using the segment Hence, the jump in Y at x 0 is given by Since x 0 ∈ S and Y such that ∇Y = Z in D were chosen arbitrarily, (7) holds for all x 0 ∈ S and admissible Y in the specified class. Thus Y is unique in that class, independent of position on S, and given by the constant −A ⊗ ν. We note that Y * : h → R n×d×d may be viewed as a discontinuous function with the specification where the points x − ∈ D belong to the − side of S at x. We now evaluate the jump in the field W on S.
As already observed, for any Y satisfying (5) a W field in D can also be constructed and this will have the representation for any path linking y to x in D. We now arbitrarily fix an admissible field Y and choose the same path from Since Y remains constant at the value given by (6) along the chosen path from x 0 − l 2 e 2 to x 0 + l 2 e 2 , using (6) and (8). Now and so that and therefore (9), (10), and (6), noting Y x 0 − l 2 e 2 = Y x 0 + l 2 e 2 , imply Thus, We now define the function W * : h → R n×d as where the points x ± belong to the ± sides of S at x, respectively. W * is a continuous function on h . We now assume that the constant A is of the form A = a ⊗ ν for a ∈ R n . Then, in D, This constant is free to choose, without loss of generality (related to the choice of Y (x − 0 ), for instance), and we assume that it satisfies If not, we explicitly solve the system Solutions exist to this system (e.g. an explicit solution on star-shaped domains can be written down by using the Riemann-Graves integral operator (Edelen 1985)) that belong to C 1 ( h ). Forcing by the Dirac distribution is not necessary; functions of (x, y) with support in a cylinder contained in c satisfying Aα e 3 da = −b for any area patch A threaded by the cylinder also suffice for generating such solutions (Acharya 2001)). Then defining we note that W dx = 0 for any closed contour encircling c and that W ∈ C 0 ( h ).
for arbitrarily fixed z ∈ h and a constant p ∈ R n . Clearly,ỹ satisfies ∇ỹ = W on h andỹ ∈ C 1 ( h ).
Consider the constant vector a ∈ R n to be parametrized by the layer width l as All fields constructed with the use of A = a (l) ⊗ ν are denoted by a superscript (l). We have the following properties: lim l→0 ∇ 2ỹ(l) is bounded on S. This conclusion also holds for any value of β ≥ 0 when l > 0 is held fixed.
Remark 0.1. While the above considerations have dealt with one singular surface, the linearity of the construct on the prescribed field Z makes it clear that exactly similar arguments hold for the superposition of a set of deformations, each element of which contains a single planar surface of discontinuity of arbitrary orientation in h terminating on c . Considering y i , i = 1 to n ∈ Z + , each corresponding to a specified Z i field, the composite, superposed deformation n i=1 y i is C 1 ( h ), with generally discontinuous second derivatives on each of the S i corresponding to the specified Z i field. This corresponds to situations with a single branch point (Gemmer and Venkataramani 2011;Gemmer and Venkataramani 2012;Gemmer and Venkataramani 2013) h . This corresponds to configurations with multi ple branch-points (Gemmer et al. 2016;Shearman and Venkataramani, in preparation). Remark 0.2. For thin objects modeled by d = 2, the construction above is a representation of folds without ridges. In "Kinematics" section and "Thermodynamics" section we develop a continuum mechanical theory that encompasses the mechanics of such folds in simply connected domains within a setting that allows for deformations with less smoothness. Remark 0.3. Consider d = 2, n = 2 and b = 0, and assume that W * (x), x ∈ h , is invertible. A field y * : D → R 2 satisfying ∇y * = W * in D can be constructed that may be interpreted as a discontinuous deformation of h . Now consider the metric g := W * T W * on h . By the Nash C 1 embedding theorem, there exists a C 1 deformation z : h → R 3 with (∇z) T ∇z = g = (∇y * ) T ∇y * .
For a mechanistic interpretation, consider the configuration in R 3 defined by z( h ) as the stress-free, global reference configuration in a higher dimensional space (R 3 ) corresponding to a stressed body with a dislocation (with excluded core) in R 2 represented by h . The stress-free reference cannot be represented by a compatible mapping of h in the lowerdimensional space R 2 ; instead, one of its stress-free representations in R 2 is defined by the configuration y * ( h ). The stress-producing elastic Right-Cauchy Green tensor field is given by (W * −1 ) T W * −1 on h .

The discontinuity of the deformation of a non simply connected domain with prescribed third 'deformation gradient'
Consider the domain h of Fig. 5 which is rendered simply connected by a single cutsurface S which is not necessarily planar. As before, we refer to h \S =: D. We consider Z : h → R n×d×d×d as a given field for which ((Z(x)e l )e k )e j is invariant w.r.t interchanges of e j , e k , e l for any values of j, k, l ∈ {1, . . . , d}. Furthermore, we assume that Z ∈ C 0 ( h ), and curl Z = 0 in h . We are now interested in the construction of a field y : D → R n that satisfies ∇ 3 y = Z and characterizing the jump field y on S.
Define, for x ∈ h , [ ((Z(x)e l )e k )e j ] ·E I =: Z Ijkl (x), I = 1, . . . , n and j, k, l = 1, . . . , d, where E I represents an element of an orthonormal basis in R n . Z Ijkl is symmetric in the indices j, k, l. Now construct Y : which is possible since curlZ = 0 and D being simply connected. We note that due to the symmetry of Z Ijkl in j, k and the connectedness of D. Since the construction of Y allows the free specification of its value at one point of D, it can be assumed without loss of generality that Y Ijk = Y Ikj in D.
Eq. 15 and the symmetry of Z in the last two indices imply curl Y = 0 in D. Thus it is also possible to construct W : Furthermore, (16) and the symmetry of Y in its last two indices imply that a function y : D → R n can be constructed satisfying Now, because Z is curl-free in h , we have by Stokes' theorem that Z dx=: ∈ R n×d×d a constant, for the line integral over any closed loop encircling c .
By (15) and (16), this further implies that Let x 0 , x ∈ S be connected by a curve c contained in S. Consider curves c + and c − on the ± sides of S connecting x ± 0 to x ± . Then Now, due to the symmetry of in its last two indices, Ijk c k dc j ds = 1 2 d ds ( Ijk c k c j ) and the last line integral in (21) evaluates to 1 2 ( (x − x 0 )) (x − x 0 ) so that (21) implies  (22) is an extension of Weingarten's theorem (Weingarten 1901;Volterra 1907;Delphenich b) and the Weingarten theorem for g.disclinations .

Kinematics
In this section we propose the kinematics for a model of the type of discontinuities treated in "Motivation for kinematics of the theory" section, to be broadly applied to the mechanics of materials. For that purpose, it is essential to deal with simply connected, compact domains containing the said discontinuities. The excluded core regions are now included in the domain as are the excluded surfaces of discontinuity. Roughly speaking, we consider an additive split of fields into 'regular' and 'singular' parts whenever the field in question contains high magnitudes concentrated in 'thin' regions approximating smooth lower-dimensional (< d) sets; the support of the singular part of the field contains these regions of high concentration and that of the regular part contains the support of the rest of the field, including regions supporting approximate discontinuities. Importantly, both the singular and regular parts are assumed to be at least integrable functions as we want to write governing equations for these fields in the form of pde that can at least be made sense of in some weak manner. Thus, we take a somewhat microscopic point of view, assuming that discontinuities and singularities of certain fields when viewed from a macroscopic scale have a smoother definition at a microscopic scale that we describe by additional 'eigenwall' fields. We also adopt the point of view that once macroscopic theories generate discontinuities and singularities, in most circumstances additional physical insight beyond the constraints placed by the governing equations of the macroscopic theory are required to define evolution with a modicum of uniqueness. We develop such a model in the rest of the paper. We refer to a fixed reference configuration, a simply connected compact region as B. In terms of the displacement field u and the i-eigenwall fields S (i) , i ∈ {1, 2, 3}, we define the i-elastic distortions Y (i) , i ∈ {0, . . . , 4}, as (Y (0) is analogous to the field y of "Motivation for kinematics of the theory" section, Y (1) to W, Y (2) to Y, and Y (3) to Z). Thus Y (0) = u and Y (4) , the gradient of the regular part of the gradient of the 3-elastic distortion, are assumed to have no 'singular' parts. We now define the 'composite' eigenwall fields S (i) , i = 1, 2, 3, as and we note that Physical considerations related to predicting stress fields of terminating twin boundaries and the stress-free, compatible, elastic, twinning shear distortions of through-twin boundaries ) motivate the introduction of the following Stokes-Helmholtz (SH) decompositions: We will also consider exactly analogous SH decompositions for the fields Combining (25) and (27) and noting the uniqueness of the SH decomposition we have up to at most a spatially constant function of time which we will assume to be a timeindependent constant. Defining (noting that H (4) = 0), we define the i-defect density tensors for i ∈ {1, 2, 3} from (23) and (29) as using (25) and S (4) = χ (4) = 0. Since α (i) are defined locally as a curl, the local forms of the conservation laws for topological charge content, α (i) n da, of an arbitrary area patch is given bẏ where V (i) , for each i, is a vector field. V (i) is the velocity field of the i-defect density field. Combining (30) and (31), we have that for some F (i) that can be prescribed. Eqs. 24 and (32) implẏ with the last sum vanishing for i = 1. By kinematical arguments related to allowing for transverse motion of walls characterized by localized S (i) fields on surfaces, a part of where V ⊥(i) is the velocity of the i-eigenwall field. Guided by simplicity in thermodynamic arguments that precludes the appearance of (unremovable) gradients of dislocation and eigenwall velocity fields in the expression for dissipation of the body (see "Thermodynamics" section), we make the following choice In (32) and (33), incorporating (34), V (i) and V ⊥(i) are to be constitutively specified, minimally consistent with the second law of thermodynamics to be globally satisfied for all processes of any body modeled by this theory. Surfaces of displacement discontinuity (e.g. stacking faults) are not known to move transverse to themselves; moreover, such discontinuities are often not identifiable based on knowledge of only the current state (and not of the distinguished coherent reference from which displacements are measured). Hence, we will assume V ⊥(1) ≡ 0. Elastic phase boundaries, i.e. localizations of the S (1) field along surfaces are known to move transverse to themselves, and not much is known about transverse motions of surfaces of discontinuity of the second gradient of elastic distortion, i.e. surfaces of inflection. Thus, we allow V ⊥(i) , i = 2, 3 to be nonvanishing fields in general. Hence, we have the following evolution equations for the eigenwall fields:

Thermodynamics
We assume a free-energy density function of the body with the following dependencies: using (29), (24), (28), and noting that H (4) = 0 (where the argument fields of each of the functions are evaluated at (x, t) to give the value of ψ(x, t)). Roughly speaking, the dependencies of ψ * on Y (i) , α (i) , i = 1, 2, 3 are expected to be convex and those on S (i) , i = 1, 2, 3 to be multi-well, nonconvex.
The balances of linear and angular momentum are given by where ρ is the mass density, v is the material velocity vector, T is the stress, is the couple stress, and b f , K are the body force and body-couple densities per unit volume, respectively. As usual in solid mechanics, we assume balance of mass is satisfied once the deformation map at any instant is determined by evaluating the density field on the deforming body from the formula ρ = ρ 0 det(I+∇u) , where ρ 0 is the density field on the reference configuration.
The mechanical power supplied to the body is defined as (Mindlin and Tiersten 1962) using the balances of linear and angular momentum, where n is the outward unit normal to the boundary of the body, ω := 1 2 curlv = − 1 2 X · 2 is the rotation vector where := 1 2 ∇v − (∇v) T is the rotation-rate tensor, D := 1 2 ∇v + (∇v) T is the strain-rate tensor, and M := ∇ω. Denoting the mechanical dissipation, D, or the difference between the power supplied to the body and that stored in it, is given by In the following, we deduce guidelines for constitutive specification in our model that ensure that the mechanical dissipation vanishes in the absence of eigenwall and defect field evolution in any process and is positive otherwise, a minimal necessary condition for the mathematical model to be well-posed.
To facilitate the derivation of the thermodynamic driving forces for the various defect density and eigenwall fields, we will need the following auxiliary fields P (i) , i ∈ {2, 3} defined by the solutions of the following Poisson equations: which requires that the free-energy density function should satisfy the constraint (This is formally easily arranged by taking any arbitraryψ with the dependencies of (36) 3 , and defining ψ =ψ − 3 i=2 | | −1 ∂ H (i)ψ dv · i H (i) , but its physical and rigorous mathematical implications need to be understood). Defining (that exist by a unique Stokes-Helmholtz resolution of R (i) ), will aso be required in the sequel for deriving the thermodynamic driving forces.
A long computation involving (36) 3 and the kinematics of the model defined in "Kinematics" section reveals that the mechanical dissipation may be expressed in the suggestive form Thus, a set of constitutive equations, driving forces for dissipative mechanisms (denoted below by the symbol ), and some boundary conditions for the model are (it can be checked that the rhs of (52) is deviatoric). Eqs. 51-(54) along with the constitutive choices for the defect and eigenwall velocities to be in the direction of their respective driving forces, mediated by a positive, mobility/drag scalar required on dimensional grounds, ensures non-negative dissipation. Of course, other choices consistent with positive dissipation are possible as well. The boundary conditions (53)-(54) are not the most general, but a compromise between including higher order stress tensors with dubious physical meaning beyond couple stresses and simplicity in an already involved higher order theory of defects. It is clear from (51)- (52) and (37) that the governing equations lead to sixth-order pde in the displacement field u (see "Example: a model of branch-point defects in an elastic body" section below). Remark 0.9. A minimal set of field variables to be evolved in the model are (u, S (i) , i = 1, 2, 3) governed by (37) and (35), with H (1) , H (2) determined from (27)

Example: a model of branch-point defects in an elastic body
We assume the as-received body as the reference configuration with all displacements measured from it; in particular, we assume that u(x, 0) = 0. We specialize the general formalism to a specific case by making the simplest possible choice for the free energy density (36) that shows the generalization of incompatible elasticity achieved by our work: with the ansatz that V ⊥(i) = 0, i = 1, 2, 3 and V (1) = V (2) = 0, assumptions that are consistent with non-negative dissipation. Under these conditions S (i) , i = 1, 2 do not evolve and remain fixed at their values specified through initial conditions. Let S (1) (x, t) = S (1) (x, 0) = g(x) and S (2) (x, t) = S (2) (x, 0) = b(x), and we note that g sym and b are the analogs, for a 3-d body, of the freely specified, non-evolving, g and b tensors of incompatible elasticity described in "Statics and equilibria of non-Euclidean elastic sheets" section; we note that with g specified, S (2) (x, 0) can be arbitrarily specified by making a suitable choice of the field S (2) (x, 0). We then have S (3) = S (3) + ∇ S (2) =: S. Here, C is the 4 th -order tensor of elastic moduli with major and minor symmetries, c 2 , c 3 are nonnegative scalars (in place of sixth and eighth order tensors!), d 3 is a positive scalar (that could also be a positive scalar-valued function of |curlS|), and l, 3 are positive scalars. The physical dimensions of c 2 , c 3 , d 3 , l, 3 are stress.(length) 2 , stress.(length) 4 , stress, length, and stress.(length) 6 , respectively. Since the equilibria we envisage are of nominally elastic bodies that show non-trivial shapes under no applied loads, f is generally expected to be a multi-well nonconvex function with the bottom of one well at the argument 0. Thus we are looking for the mechanics of surfaces of inflection and branch line defects in bodies with an evolving stress-free reference characterized by the choices S (1) It is interesting to note that even when g = 0 and b = 0, the energy/stress-free reference for our body is never immersible in three-dimensional Euclidean space whenever S = 0, i.e. the stress-free state is necessarily incompatible or non-realizable, since it is impossible to construct a displacement field of a 3-d body with vanishing strain, i.e., (∇u) (s) = 0, whose third gradient is non-vanishing.
The balances of linear and angular momentum (37) are solved by taking a curl of (37) 2 to obtain div T (a) = 1 2 curl div dev + 1 2 curl K, that on substitution in (37) 1 leads to Constitutive Eqs. 51-(52) are used to solve for a displacement field from (58) (when the defect fields are assumed given), thus satisfying (37) 1 , and (37) 2 is then satisfied, in terms of this displacement field, by simply evaluating T a from the equation making the assumption that the constitutively undetermined tr = 0, without loss of generality.
A detail needs to be attended to in the above considerations. The thermodynamic relation (52) requires that For the constitutive choice (57), (64) is equivalent to −c 2 e ljk S (2) [jk]l + c 3 e ljk S [jk]lm,m = 0.
Since the evolution Eq. 63 for S is of the formṠ = −curl S × V for a vector field V, it is true that X · 3Ṡ (a) = −curl X · 3 S (a) × V and hence X · 3 S (a) (x, t) = 0 is consistent with the evolution for initial data X · 3 S (a) (x, 0) = 0 and we adopt this solution. Furthermore, we assume X · 3 b (a) = 0 by specification so that X · 3 S (2)(a) (x, t) = 0 and therefore (64) is satisfied. Remark 0.12. Spatial derivatives of the 3-eigenwall field serve as a source term in (62) where ν is the unit normal to a planar surface, g is a scalar-valued function of the spatial coordinate along ν given by ζ = ν · x (say a Gaussian centered at ζ = 0), and b is a constant vector, this forcing is of the form d 3 g dζ 3 b. Eq. 63 implies that there is no evolution of the eigenwall field at locations where curl S = 0, regardless of the energetic driving force there. For example, the field S(x) = g(ν · x) b ⊗ ν ⊗ ν ⊗ ν has no 'longitudinal' variation and does not evolve according to (63). However, where t is orthogonal to ν does evolve. Physically, the eigenwall field is 'dragged' by the evolution of its core. Remark 0.13. The as-received body (where u = 0) need not be in equilibrium for generally specified eigenwall initial data g, b, S(·, 0) for which the initial acceleration field can be evaluated. The class of initial data that leads to a self-equilibrated, generally stressed reference configuration may be obtained by writing g := ∇z + g * , and solving for the vector field z from (62) (and boundary conditions) with ρ := 0 for each choice of g * , b, S(·, 0) . Remark 0.14. The model with the ansatz S (1) = 0, S (2) = 0, is worthy of study on its own merits. Remark 0.15. Configurations rendering local minima of a body with the energy density (57) may be studied by an (L 2 ) gradient flow dynamics in the fields u, S (i) ; i = 1, 2, 3 starting from arbitrarily specified initial states for these variables. Remark 0.16. The governing Eq. 62 implies that, when the elastic modulus C is homogeneous and isotropic, given by C ijkl = λu k,k δ ij + μ(δ ik δ jl + δ il δ jk ), plane waves of curl u and div u are dispersive in nature, with propagation possible in any direction in space. The dilatational waves (i.e., waves of div u) with wave number |k| and direction k |k| propagate with velocity c d := ± c 3 |k| 4 + c 2 |k| 2 + (λ + 2μ) ρ while the equivoluminal waves or 'shear waves ' (i.e., vectorial waves of curl u) propagate with velocity c s := ± c 3 |k| 4 + c 2 |k| 2 + μ ρ .
Continuous dependence w.r.t initial data of the Cauchy problem for the evolution of displacement requires c 3 ≥ 0. When c 3 = 0, c 2 must be non-negative, with the requirement that μ ≥ 0 and λ + 2μ ≥ 0 if c 2 = 0. Within these parameter regimes, linear instabilities can arise for wavenumber and parameter combinations resulting in c d or c s taking complex values.

Uniqueness of the displacement field and boundary conditions
Our model encompasses a model of third-order elasticity in the absence of dissipative defect evolution, and involves the thermodynamically motivated higher-order boundary conditions (53)-(54). Here, we use a uniqueness argument (in a putative smooth class of solutions) to deduce a full set of boundary conditions for the problem (62) when the S field is assumed specified. We abstract the results of the exercise in this special case related to the 'quadratic' energy (57) to identify a likely set of sufficiently general boundary conditions for the determination of the displacement field for processes consistent with the general constitutive statement (36). Consider two solutions u (1) and u (2) of (62) corresponding to identical S, K, b f , g, b fields. Denote the difference displacement as u := u (1) − u (2) and its velocity v =u. Then u satisfies and taking the inner-product of the difference velocity with the equation and integrating in space, we have Let us now assume that both u (1) and u (2) satisfy (53)-(54) consistent with (57). Then the last line of (65) vanishes due to the boundary condition (54) and the line before that due to (53).
Let the stress field arising from (u (i) , S), i = 1 2, be T (i) = T (s)(i) +T (a)(i) , in accord with (59), (60), and (61). Then the third line from the bottom of (65) may be interpreted as il v i n l da and if we now additionally require that solutions satisfy specified tractions and velocities (or displacements) on mutually complementary parts of the boundary of the body for all times, then this term vanishes. Consequently, we are left with d dt and if u (1) and u (2) both satisfy specified initial conditions on the displacement and velocity fields, then the bracketed quantity, an integral of sums of squares (in fact, the potential and kinetic energies of the body subjected to the difference displacement) vanishes at all times. This proves that the difference velocity vanishes point-wise, and the initial condition on the difference displacement implies that u (1) = u (2) for all (x, t). Obviously, the dynamic problem allows the prediction of unique rigid motions. In statics, i.e., when the inertia term is absent, one takes the inner product of the governing equation for the difference displacement with the difference displacement, and obtains, for the same boundary conditions (except only the displacement can now be specified on the part of the boundary complementary to where tractions are specified), All integrands are non-negative implying that the strain, or the symmetrized displacement gradient, vanishes (recall the minor symmetries of C) which, by compatibility, further implies that the displacement field is unique if a displacement boundary condition is specified and otherwise it is unique up to an infinitesimally rigid deformation. Thus, the higher order boundary conditions (53)-(54), along with classical displacement and traction boundary conditions may be expected to define a well-set problem (for the displacement field) in the case of the general constitutive Eq. 36. Of course, the traction now involves a stress tensor that has an antisymmetric part, and is constitutively dependent on higher order displacement gradients.

A 'plate' idealization
For simplicity we consider g = 0 and b = 0. Let the reference B be a plate of thickness 2t, i.e., B = {(x 1 , x 2 , x 3 )|(x 1 , x 2 , 0) ∈ B 2 , x 3 ∈[ −t, +t] }, where B 2 is a flat 2-dimensional simply connected domain. Defining the through-the-thickness average of a function as and the notation we now seek the governing equations for u and S, under the ansatz that S = S and ρ = ρ, i.e., S and ρ do not vary through the thickness of the plate, and K = b = 0. It is also assumed that a component of S vanishes if any of its last three indices takes the value 3. We use the notation that all lowercase Greek indices vary from 1 to 2 while lowercase Latin indices span from 1 to 3. While not essential, the assumptions l = 2t, c 2 = Et 2 and c 3 = Et 4 , where E is the Young's modulus of the material can be made to draw an analogy with classical plate theory (the curvature-related elastic energy term in the thickness-integrated expression of (57) would then be proportional to t 3 ). For 0 < t 1, whenever S = 0, there is energy and stress in the body, possibly small, with the corresponding thickness-integrated 'elastic' energy of the plate (arising from the first three terms in (57)), alternatively the 'plate elastic energy' , scales as ∼ t 5 , assuming energy is minimized, there are no external forcing or constraints, and 3 > 0 to rule out any possibility of a singular energy. Our governing Eqs. 62 or (66) do not require that energy be minimized, so that scaling of the thickness-integrated elastic energy w.r.t t as t → 0 in the model can well contain lower order bending O(t 3 ) , and even stretching (O(t)), contributions.
Applying the averaging operator to (62) and noting that we obtain Similarly, In Eq. 66, the terms beyond the first line represent forcings in the transverse direction to the plate and need to be specified (it would be physically legitimate to assume many of these terms to vanish); the third line of (67) has similar meaning and needs specification. The functions u, S represent the fundamental fields of the plate theory, governed by (66)-(67). Evaluating T (a) from (59) in terms of (u, S) solving (66)-(67) and K would imply the satisfaction of balance of angular momentum (i.e., moment balance) in the throughthe-thickness averaged sense. Remark 0.17. We note that non-evolving and non-vanishing S (1) = g, S (2) = b 'target' composite eigenwall fields can be included in the considerations of this plate idealization with only slight increase of tedium in bookkeeping.
Within the context of energy minimization and for t > 0, if curl curl S (1)(s) T = 0, i.e. S (1)(s) satisfies the St.-Venant compatibility condition, then an infinitesimal isometry exists (the reference configuration is assumed to be simply-connected) and the plate elastic energy scales as ∼ t 3 or of smaller magnitude; if S (1)(s) is not compatible, then the energy has to scale as ∼ t. We note that when S (1)(s) is compatible, unless S (2) = ∇ 2 v, where v is s.t. (∇v) (s) = S (1)(s) so that ∇ 2 v is unique, the plate elastic energy is going to scale as ∼ t 3 . The requirement S (2) = ∇ 2 v is non-generic for a freely-specifiable S (2) field that, however, is satisfied by the choice S (1) = 0, S (2) = 0. Thus, in most circumstances the plate energy is expected to scale as ∼ t 3 , if the plate energy is minimized.

Discussion
Starting from the work of the brothers Cosserat (1909) (as presented in Truesdell and Toupin (1960)), through those of Toupin (1964), Green andRivlin (1964), Mindlin (1962;1964), on to that of Fleck et al. (1994), Fleck and Hutchinson (2001), Hutchinson (2012) and of Gurtin (2002), Gurtin and Anand (2009), higher order theories of continuum mechanics have made an appearance off and on and have been noted for their intricacy and elegance, but always, arguably, with the nagging question of the physical justification (in their details 1 ) in view of their added complexity. Our work aims to provides a concrete, tangible, and compelling justification -that the precise treatment of defects in the deformation and its higher order gradients is the raison d'être for higher order theory in continuum mechanics.
Our work is in the context of non-Euclidean elastic sheets with negative in-plane Gauss curvature. These objects are ubiquitous in nature and they display varied and intricate multi-scale behaviors (Sharon et al. 2002;Audoly and Boudaoud 2003;Klein et al. 2007;Kim et al. 2012;Gemmer and Venkataramani 2013). Their elastic behavior is significantly different from that of elastic plates or spherical shells (Gemmer et al. 2016;Shearman and Venkataramani, in preparation). In particular, they have "large" continuous families of low-energy states obtained from piecewise isometries, with each piece possessing additional "bending" degrees of freedom. Thin hyperbolic free sheets are thus easily deformed by weak stresses and their morphology is strongly dependent on the dynamics of the growth/swelling processes, material imperfections, or other weak external forces. This naturally motivates the need for tools to describe singularities/defects in these sheets, their interactions and the resulting dynamics.
Mesoscopic defects in hyperbolic sheets, associated with their "soft" modes of deformation, include lines of inflection that terminate at branch points (Gemmer and Venkataramani 2012;Shearman and Venkataramani, in preparation). These are higher-order defects (termination of jumps in curvature) unlike the more common types of defects, disclinations and dislocations. Irreversible effects in the dynamics of disclinations and dislocations are associated with (macroscopic) plastic behaviors -stress-free large deformations, internal stresses, and microstructure -in solids. A natural question therefore is -what are the macroscopic manifestations of moving lines/surfaces of inflection and branch points/lines?
In this work we have begun to address this question in the context of 'small deformations' from a (potentially stressed, when occupied) reference configuration. A detailed analysis and characterization of the kinematics of branch point defects and the discontinuities in the deformation that they induce is achieved. This analysis, in its essence, is a non-trivial adaptation and extension of the ideas of Weingarten (1901) and Volterra (1907), from the dawn of elastic defect theory, to a context not restricted within the kinematics of only strain (the symmetrized gradient of the displacement, as well as its nonlinear analog) and its incompatibilities, and shows the natural way forward for deducing the constraints on possible jumps in deformation, i.e. global constraints, for locally compatible higher order deformation gradients, albeit on domains with the simplest non-trivial topology 2 . We then develop a thermodynamically consistent theory for the dissipative dynamics of such defects in a nominally elastic solid, allowing for their interaction with dislocation, g.disclination, grain, and phase boundary defects. The constitutive guidance provided by this thermodynamic argument ensures that the model is equipped with an energy (in)equality, a crucial necessary condition for its physical and mathematical well-posedness. The analysis uncovers the non-Newtonian, energetic driving forces on these defects that couple their dynamics and mutual interactions to applied loads and the deformation of the body 3 . Evolution of the defect fields subject to such driving forces necessarily reduces the system free-energy by design, within an overall dynamics that accounts for material inertia and is not restricted to its free-energy decreasing with time (depending on the external driving). As an example, we explicitly demonstrate the full set of governing equations for the case of branch point defects in an elastic material and develop a 'plate' theory idealization for it. The development of the finite deformation version of the model poses no conceptual or technical barriers 4 based on our prior work in g.disclination mechanics (Acharya and Fressengeas 2015), but this same work makes it clear that the bookkeeping tasks in pushing through the analysis are going to be formidable.
We observe in passing that while we have been interested in developing a theory for branch point/line defects and lines/surfaces of inflection, i.e. a theory for the discontinuities and singularities of the deformation and its gradients up to order three, the analysis makes it clear that the mathematical/continuum mechanical formalism extends to describing the discontinuities and singularities of any finite integer order gradient of the deformation, while including only stresses and couple stresses. As already observed in Acharya and Fressengeas (2015), using the Second Law in global form is crucial for this, albeit at the expense of the application of limited (but adequate, as we show in "Uniqueness of the displacement field and boundary conditions" section) higher-order boundary conditions about which not much is physically known anyway.
As a final comment, we note that a geometric model of growth mechanics, based on Riemannian geometry and including evolution, has been proposed in Yavari (2010). The viewpoint is different from ours and, in particular, the mechanics of incompatibility based on a Riemannian metric cannot describe (without non-trivial extension) the 'softer' branch point defects we focus on. We expect that one can recast our continuum mechanical kinematic constructs within a differential geometric structure involving the specification of a moving frame, and higher-order constructs based on such a field, thereby making connections with the "geometric" viewpoint of growth mechanics.