Linear Algebra Provides a Basis for Elasticity without Stress or Strain

Linear algebra provides insights into the description of elasticity without stress or strain. Classical descriptions of elasticity usually begin with defining stress and strain and the constitutive equations of the material that relate these to each other. Elasticity without stress or strain begins with the positions of the points and the energy of deformation. The energy of deformation as a function of the positions of the points within the material provides the material properties for the model. A discrete or continuous model of the deformation can be constructed by minimizing the total energy of deformation. As presented, this approach is limited to hyper-elastic materials, but is appropriate for infinitesimal and finite deformations, isotropic and anisotropic materials, as well as quasi-static and dynamic responses.


Introduction
Soft materials like rubber, foam, and many biological materials can stretch far beyond the limits of infinitesimal elasticity and yet return to their original shape when forces are removed.It is useful to be able to model these deformations, but infinitesimal elasticity based on stress and strain tensors cannot be used for these large deformations.The classical equations of finite elasticity are quite difficult, requiring 30 or more tensors to explain the theory.The description of elasticity presented here presents the same equations for both infinitesimal and finite elasticity, and requires only two tensors to define the theory.
The description of elasticity presented here will be for hyper-elastic materials.A hyper-elastic material stores energy when it is deformed and returns this energy to its surroundings when it is returned to its original state.Rubber is the most common example.The energy stored in the material can be expressed as a function of the positions of the points within the material.By minimizing the stored energy, differential equations of elasticity and forces can be found.The idea of describing deformations in terms of points and forces was first used by Euler, Lagrange, and Poisson and predates Cauchy's introduction of stress and strain [1].The earlier researchers, however, did not complete the description of the general differential equations for finite elasticity.I will do that here.To do this, I will follow the notation of Spencer [2] to describe the positions of the points within the body before and after deformation.

Deformation as a Mapping
Define the initial location of each point within a material as the vector, X , with components i X , 1, 2, 3 i = .After the deformation, each point in the material will be at some new position, x , with components i x , 1, 2, 3 i = . The components of each point after the deformation are functions of the components of the position of each point before deformation.That is ( ) ( ) , , , , .
To match the physical reality of the deformation of a material, neither inversions nor a change of dimension of the material will be allowed (i.e. a three-dimensional material cannot be turned inside out or pressed into a plane or a line).As a result of this restriction, every point in the material after the deformation will correspond to exactly one point in the material before deformation, so that the mapping from X to x is one-to-one.
Consider now a point near the point i X , i.e.
We can find the relationship between d i X and d i x by differentiating Equation (1) and using the chain rule to give d d , here and throughout the rest of this paper, the Einstein summation notation is used so that repeated indices are summed over.In general i j x X ∂ ∂ will vary from point to point within the material (see Figure 1).The 3 × 3 matrix formed by these values, is called the deformation gradient tensor by Spencer [2].The results found so far are appropriate for any coordinate system, but I will use two specific inertial coordinate systems (i.e.fixed coordinate systems where New- ton's laws apply).One I will call the observer coordinate system.This coordinate system is the one chosen to solve some problem in (e.g. a simulation or an engineering problem).The second inertial coordinate system I will call the experimental coordinate system.This coordinate system will be the one chosen by an experimentalist who wishes to measure the energy associate with a particular deformation of a particular sample.In the observer coordinate system, I will denote the components of ij  as o ij  .In the experimental coordinate system I will denote the components ij  as e ij  .I will first describe how the experimenter should measure the material properties of the material in his coordinate system.

Measuring the Energy of Deformation
It is sufficient to limit experimental deformations to homogeneous deformations of a homogeneous portion of the material in order to define the energy of deformation.If the body is anisotropic, it is necessary to define the orientation of this anisotropy.This could be determined by a visual inspection (e.g.wood grain) or by a knowledge of how the material was made (e.g.rebar in concrete).Orient the anisotropy to align with the experimental coordinate system so that in the experimental coordinate system the anisotropic coordinate system of the material is ˆi e for 1, 2, 3 i = , and ( ) ( ) ( ) ( ) ( ) ( ) ˆˆ1, 0, 0 , 0,1, 0 , 0, 0,1 e e e = = = , or A homogeneous deformation is a deformation in which every point in the body undergoes an affine transformation.Linear algebra defines an affine transformation as a mapping that transforms any point in space i X to another point in space i x by a matrix transformation ij F , followed by a translation i d .Applying this to our material body, every point in the body after the deformation, i x , corresponds to a point before the deformation, i X , as where i X = coordinates of the position of any point before the deformation, i x = coordinates of the position of this same point after the deformation, ij F = 3 × 3 matrix of values that are the same throughout the body, and i d = components of a vector that are the same throughout the body.
Taking the derivative of Equation ( 5) and comparing the result with Equation (4), we find that in our experimental coordinate system e ij ij F =  .Equation (5) indicates that to completely define any deformation in our experimental coordinate system, we need the nine components of e ij  and the three components i d .These 12 values can be found by recording the location of any four non-coplanar points in the experimental material before and after the homogeneous deformation.The coordinates of each of these points provides i x and i X in Equation (5).Since there are four points with three components each, this gives a total of 12 equations and 12 unknowns of the form of Equation (5).From these 12 equations, all 12 components of e ij  and i d can be found for any experimental deformation.
To measure the energy during experiments, record the applied forces, m F , and the corresponding displace- ments of the points where these forces are applied, d m x with m = 1 to M, where M is the number of applied forces.The energy of deformation is the total work done by these forces, i.e. the sum of the integrals from the initial to the final position of the points where the forces are applied: The energy per original volume, E, is found by dividing the total energy by the initial volume of the sample, V.
With these measurements, an experimenter can construct a table of the stored energy per unit initial volume as a function of the 12 components of e ij  and i d .To complete the description, use a linear interpolation or an equation fit to the table of data to define E as a continuous function of ( ) This completely defines the energy per unit original volume for homogeneous deformations of the material sample.

Special Experimental Cases
We have found that in general the experimentally measured energy per unit volume, E, can be a function of as many as 12 component values (nine in e ij  and three in i d ).It is a bit much to expect an experimenter to map out a 12 dimensional energy space for every material for every application, so it is useful to find some simplifications depending upon the particular material and application we are interested in.For example, for most applications, the only external body force that needs to be considered is gravity.In that case, we only need to express the energy of translation as a function of 3 û ⋅ d , where 3 û corresponds to the direction of the gravita- tional force and d the displacement of the center of mass of the material.If the experimental coordinate sys- tem is aligned so that the third component of d is parallel with the gravitational force, we only need to include energy changes due to 3 d and can omit the two variables 1 d and 2 d in our energy function, because these displacements will result in no change in the energy of the body.
A further reduction of parameters for rotations and deformations can be found if we take a Singular Value Decomposition (SVD) of e ij  .The SVD of a matrix e ij  uniquely divides the matrix into three matrices, two rotational ( 1 R and 2 R ), and one diagonal, Λ .In particular, ( ) SVD Λ This form will be useful in describing both experimental processes to measure the energy of deformation and the equations of elasticity.Since any deformation, e ij  , can be expressed uniquely as 2 1 Λ R R , then every possi- ble deformation can be considered a combination of a rotation, 1 R followed by a stretch or compression along the three fixed orthogonal coordinate axes, Λ , followed by a final rotation, 2 R .In the most general case the energy of the material can depend upon body forces from electric or magnetic forces in addition to gravity.For example, a material with an electric dipole ( p ) or a magnetic dipole ( µ ) in an electric ( ξ ) or magnetic (  ) field, will have energy ⋅ p ξ or ⋅ µ  , respectively so that the energy will be a function of 1 R , Λ , and 2 R .If the material has a charge, q, energy must include q ⋅ d ξ .In these cases the energy of deformation may be a function of all nine variables in ij F as well as the three variables in i d .If the electric and magnetic body forces are not significant in a particular application, the three components of 2 R need not be included in calculating E since in that case the energy associated with deformations are independent of body rotations after the deformation.
If the body is isotropic, then the energy will be independent of both 1 R and 2 R since rotating the body be- fore or after the deformation will produce no change in the energy of deformation.In that case only the three diagonal values of Λ , Λ ii , are needed to describe the energy associated with deformation.If in addition, the material can be considered incompressible, then the volume of the material, V, is constant and only two independent elements of Λ need to be used to describe the energy during a deformation.In addition to these simplifications, it is sufficient to measure the energy as a function of only those changes that are expected in a particular application.So for example, if infinitesimal deformations are sufficient to model the problem at hand, only one small displacement measurement in each direction is required.Alternatively, if the body is going to be used only in extension, there is no need to measure the energy associated with compressional forces.Rivlin [3] used this approach for rubber.For his application the rubber could be considered isotropic and incompressible, so he only deformed the material sample in extension along two perpendicular directions.This is sufficient to find the energy as a function of 1 Λ and 2 Λ .Finally, if a deformation contains only rigid body motions, materials are neither compressed nor extend, so that Λ 1 ii = for 1, 2, 3 i = .Thus a rigid body rotation, Equation ( 9) gives . Since rigid body rotations can be expressed in terms of a single rotation matrix, the six parameters in 1 R and 2 R reduce to only H. H. Hardy

three.
What we have found is that we can reduce the number of variables that the energy is a function of from 12 to as few as 2 in the case of an incompressible, isotropic material where we can ignore external body forces like gravity in our application.In general, however, all 12 variables may be required and it is helpful to find the most computationally efficient way to represent the energy for each application.

Some Application Issues
We have found that ( ) can be used to simplify the experimental measurements for particular cases; however, SVD is a rather computationally heavy calculation to be used during simulations.It is therefore useful to represent the deformations in a more computationally efficient manner for applications.For example, if we consider the case where the material is isotropic and there are no body forces, the energy unit volume a function of only Λ i .Thus any three independent variables spanning the same space as Λ i may be used to characterize the energy function.The values 11 Λ , 22 Λ , and 33 Λ can be rewritten [4] as In addition where a , b , and c are the column vectors of e ij  .Thus in simulations of isotropic bodies it is not necessary to compute SVD as the simulation progresses.All that is required are the components of e ij  to compute the needed three independent values.As a result it is best in this case to redefine the experimentally defined energy, ( ) , , , , , ,  , ,  E  ( ) produces R , where  is an upper triangular matrix and R is a rotation matrix.Since the energy of deformation is independent of any rotation after the deformation, the six components of  , like the three values of Λ for the isotropic case, can be used to define the energy of deformation.The Gram-Schmidt QRD algorithm is not a particularly heavy numerical calculation and can be used in applications, but an alternative is also possible.As I have noted previously [5] the components of  can be written in terms of dot and cross products of the column vectors of e ij  .Thus the energy of deformation can be written in terms of dot and cross products of the column vectors of e ij  and once the energy is expressed in these terms, it is not necessary to calculate any other tensors than e ij  in carrying out applications; however, we do have to consider the orientation of any anisotropy for engineering applications.

Coordinate Alignment
When an anisotropic material is placed in "service" it is necessary to know the initial orientation of the anisotropy.This is because the stored energy is a function of the orientation of the isotropic material relative to the observer's coordinate system.For example consider a laminate.If the laminate is oriented so that the lamina are parallel to the x-y plane in the observer's coordinate system and extended in the x direction a fixed amount there will be a change in energy of the material.However, if the same laminate is initially oriented so that the lamina are parallel to the y-z plane in the observer's coordinate system and extended in the x direction the same fixed amount there will be a different change in energy of the material.Thus the initial orientation of the lamina in the observer's coordinate system must be known in order to correctly calculate the stored energy in the material.
When the material is in "service" any deformation has occurred, the anisotropic coordinate system in which the energy measurements were made, ( )  ˆk i e , and the observer coordinate system, ( ) ˆk i  , may not align.Define a rotation matrix, 0 ij R , where 0 ij R maps the components of the observer's coordinate system, ( ) , into the components of anisotropic coordinate system before any deformation, i.e.
( ) ( ) The local deformation which is expressed in the observer's coordinate system, o ij  , can be expressed in the experimental coordinate system, e ij  , by a simple change of coordinates of o ij  [6], ( ) Note that 0 lj R is at most only a function of i X .It is not a function of i x since it is defined before any de- formation has taken place.For an inhomogeneous material, the e ij  matrix may vary from point to point, so that the energy would be a function of i X as well as  ) So that in general,

Is This Energy a Scalar?
It may seem strange that we must transfer the coordinate values of the material in the observer's coordinate system back into the experimental coordinate system in order to find the local change in energy, but this is exactly as it should be.Energy must be a scalar, which is independent of the choice of the coordinate system.That this is the case can be seen if we express the energy in terms dot products of vectors, which are independent of the coordinate choice.Note that the local anisotropy coordinate system vectors, ( )  ˆk i e , are mapped into a new set of vectors, ( )   k i e  , which are in general neither orthogonal, nor unit vectors (e.g.see Figure 1).In the observer's coordinate system this mapping is The corresponding mapping of these same vectors in the experimental coordinate system is We can choose to express the energy as a function of the invariant, ( ) ( ) ˆi j e e ⋅  , which must have the same value in both coordinate systems.To see the connection between the maps in the different coordinate systems, expand the dot product in both systems and compare the result.In the experimental coordinate system, e e e e e e δ δ In the observer's coordinate system, where ( ) Comparing Equation ( 19) to Equation (20) returns us to Equation ( 14) since ( ) ( ) ˆi j e e ⋅  is invariant and must be the same value in both coordinate systems.This is the same type of explanation that must be used in expressing the energy in terms of the components of displacement vectors, although we usually do not discuss it in these terms.Usually the energy associated with the displacement can be calculated from a formula instead of having to carry out experiments for each material.For example, a displacement in the presence of gravity changes the energy stored in the material, but it is easily expressed as m ⋅ g d where m is the mass and g is the acceleration due to gravity and no experi- ments are necessary.However, if we did carry out the experiments, we would measure the energy of deformation in terms of the components of d in the experimental coordinate system.In that case we would need to take into account any change in the components of d in the observers coordinates.For example, assume the energy stored due to gravity is expressed in the experimental coordinate system as 3 mgd , where the force of gravity in the experimental coordinate system has been chosen to align with 3 X .If we placed this material into an observer coordinate system where gravity is in the 2 X direction, we must first rotate the coordinate system to align the 2 X direction with the 3 X direction before "looking up" the corresponding energy that we stored in our energy map we compiled by measuring the energy in the experimental coordinate system.Of course if we choose to express the energy in terms of the dot product, m ⋅ g d all of this would take care of itself.A word of caution is necessary here.The vectors ( )  ˆi e and ( ) j e  are the same vectors in both the observer and the experimental coordinate systems.On the other hand, the column vectors, a , b , and c , are only de- fined in the experimental coordinate system, and therefore the transformation in Equation ( 14) must be carried out before extracting the column vectors of e ij  to calculate the energy of deformation for an anisotropic body.This is not necessary for isotropic materials because in that case the energy is independent of all rotations and applying the transformation Equation ( 14) has no effect on the final energy calculation, i.e. the Λ i 's are the same for e ij  and o ij  .

Simulations
I have now completely defined a method to measure the energy of deformation of a homogeneous body experimentally and how to place the material in a given application.Because engineering applications are often complex, it is usually necessary to put this information into a computer simulation.The simplest approach is to just "pasting" small pieces of the material together to define the complete material.This can be done by randomly positioning points in the material and use these to divide the material to be simulated into small tetrahedrons of volume Δ k V , each bounded by four points.The initial and final locations of these four points put into Equation (5) can be used to define the 12 parameters needed to calculate the energy per unit original volume for each tetrahedron.The total energy of the system is just the weighted sum of the energy per unit volume of each tetrahedron, tot , , Δ where k is summed over all the tetrahedra in the material.Apply the boundary conditions and move the internal points to produce minimum total energy, tot E , and we have a solution.I have called this this discrete region model [4].
An alternative method is to use a continuous Euler-Lagrange technique to minimize the functional, tot material , , d The result of this approach is a set of partial differential equations which can be solved by any numerical technique (e.g.finite difference, finite element, Rayleigh-Ritz, etc.).I have called this Euler-Lagrange elasticity [5].

Differential Equations
The discrete region method collapses into Euler-Lagrange elasticity if the size of each tetrahedron approaches zero as the number of tetrahedra, N, increase without bound, i.e.
To find the differential equations of elasticity we need to minimize (or find the extrema) of tot E , i.e. tot , , d 0.
. This is a classic Calculus of Variations problem with multiple variables [7].The results of this minimization are the following three Euler equations: ( ) These three equations are quite general, being appropriate for both infinitesimal and finite deformations, for isotropic and anisotropic materials, and can include surface forces, gravity, and electrical and magnetic forces.

Special Application Cases
If the material is homogeneous, def E will not be a function of i X .If gravity is the only external body force, E can be separated into the energy of deformation, def E , and energy of body forces, body E : ( ) def body , , , If only infinitesimal deformations are needed, then E can be expanded in a Taylor's expansion which yields the same differential equations Landau derived for infinitesimal deformations [8] using classical stress and strain techniques.
If time dependence is required, define the Lagrangian, , , with ρ the mass per original volume.Finding the ex- trema of results again in three Euler equations, now of the following form: ( ) ( ) These three equations are the time dependent differential equations for hyper-elasticity [9].All that is needed now is to include boundary conditions and force.
Boundary conditions consist of Neumann and Dirichlet boundary conditions.Dirichlet boundary conditions just set the positions of boundary points of the material.Neumann boundary conditions can be expressed in terms of applied forces on the surfaces of the material [5], ( ) where d j A are the components of the original surface area where the forces are applied.Equations (30) provide differential equations to set the boundary conditions if the applied surface forces are known.

Comparison to Other Elasticity Theories
The most obvious difference of this approach and classical elasticity is that in this approach there is no definition of stress or strain.Here displacements and forces are the alternatives to stress and strain.This approach also requires the definition of only one second tensor, the deformation gradient tensor.In elasticity with stress and strain more than 30 tensors have been used to describe finite elasticity [10] [11].
For classical elasticity with stress and strain, the invariants ( ) ( ) ˆî i e e ⋅ , or ( ) ( ) i i e e ⋅   , are used to describe strain in terms of T   or T  .In that case strain is second order in the displacements, d i x In the approach pre- sented in this paper, the invariants used are ( ) ( ) ˆi j e e ⋅  , and the energy of deformation is calculated from the  matrix itself, which is only first order in d i x .Infinitesimal elasticity also requires compatibility equations relat- ing stress and strain which expresses material properties as fourth order tensors.The approach given here expresses material properties for all hyper-elastic materials as a scalar, the energy per original volume, E.
Some descriptions of elasticity define dX and dx to be represented in different coordinate bases [12].In these descriptions, ij  is then called a two-point tensor.In this presentation, all vectors and tensors are expressed in the same coordinate bases, either all in the observer or all in the experiment coordinate systems.There are no two-point tensors.
In some descriptions of elasticity, only "objective tensors" are used to formulate the constitutive equations.Objective tensors are required to be independent of the motion of the material that is being deformed.That is, these tensors should be the same in both a fixed reference coordinate system and in a coordinate system that deforms with the material [11].In this paper all physical quantities are expressed only in terms of one of two coordinate systems fixed in space before any deformation takes place.In this paper, neither the observer coordinates nor the experimental coordinates deform as the material deforms.Both are inertial coordinate systems fixed in space.(The anisotropy coordinates ( )  ˆk e do deform in space into ( ) k e  after a deformation, but these are not used as bases to describe the components of vectors or tensors.)Some details may be easily confused between this approach and the Classical approach.Three examples follow:  Λ ii are not the of invariants of the right Cauchy-Green deformation tensor, which is described as T   .Instead the Λ ii values are the diagonal elements of ( )  QRD is not the same as Polar decomposition that is used in classical elasticity theory.Polar decomposition produces R , where  is a symmetric matrix, whereas QRD produces R , where  is an upper tri- angular matrix.
 The term ( ) 30) is equivalent to Cauchy stress only for infinitesimal deformations.

Conclusion
A method of describing hyper-elasticity using linear algebra has been presented that uses points within the material and forces instead of stress and strain.The theory provides a straight forward way to measure material properties and allows the inclusion of magnetic and electric fields as well as gravity.This description uses the same equations for both finite and infinitesimal deformations.Neumann boundary conditions are expressed in terms of measured forces instead of computed stresses.The result is a complete theory of hyper-elasticity which includes infinitesimal and finite deformations, isotropic and anisotropic materials, quasi-static and dynamic elastic responses.

Figure 1 .
Figure 1.In this figure,  illustrates the observer coordinate system.ij  is the mapping defined at each point in space from the X point location before the deformation to the same point x after the deformation.Any local coordinate system before the deformation, ê , is mapped into a new local coordinate system after the deformation, e  .

1 R
after the experiments are completed.Once E has been converted from a function of Λ ii to a function of the column vectors of e ij  , it is only necessary to find the components of e ij  during simulations.The ( ) SVD e  is no longer required.If the material is anisotropic, the energy per unit volume is a function of the six independent values of Λand , Three from the diagonal elements of Λ and three angles from 1 R .QR decomposition, efficient venue for calculations than does ( )SVD e ij  .(Forthis application we must use the Gram- Schmidt QRD algorithm instead of the more common Householder QRD, because the Householder algorithm permits inversions.) course any inhomogeneous material might also be made up of different materials which have different energy maps at different locations in space, i X .)Thus in the most general case,