Spacetime evolutive equilibrium in Nonlinear Continuum Mechanics

In Continuum Mechanics a widely cited formula, due to Nanson (Mess Math 7:182–185, 1878), contributes the relation between the product of area and unit normal vectors relevant to corresponding surfaces in two configurations of a 3D body. A geometric treatment provides equivalent expressions of Nanson formula by direct elaborations on Euler–Jacobi volume change formula. Meaning and limits of Nanson formula are underlined. In the literature this formula has been improperly assumed to be expedient for attempts of imposing equilibrium in a reference shape of the body. A critical revision shows recourse to scaling and parallel transport of surface and bulk forces is impassable in assessing alleged referential equilibrium. A variational formulation of evolutive equilibrium in spacetime is developed and finite step elastic problems guided by control algorithms in computational procedures are illustrated. Step-by-step iterative methods of solution of nonlinear structural problems and simple counterexamples should help in convincing referential equilibria are not conceivable in Nonlinear Continuum Mechanics.

However, when developing a geometric treatment, the manifold E of spacetime events, with dim(E) = 3+1 , comes readily to be regarded as the natural mathematical setting, even in the realm of classical mechanics.
Fundamental topics such as motion, velocity and the description of change of observer, can be properly investigated only in the spacetime context wherein the basic notions of time-projection t E : E → Z onto the time-line Z , time-arrows Z and dynamical trajectory T ⊂ E can be introduced, as described in the next Sect. 2.

Spacetime context
Let us here summarise basic elements of the geometric spacetime theory.
A mathematical observer defines a projection from the event manifold E onto the time-line Z , that is a surjective submersion 1 t E : E → Z which assigns a time instant t E (e) to each spacetime event e ∈ E.
Tangent vector fields (sections in geometric language) of the tangent bundle π : T E → E , are maps v : E → T E such that the composition π • v : E → E is the identity in E.
The spacetime manifold E is fibred into spatial slices S which are integral manifolds of the kernel of the differential: 2 dt Integrability is assured by Deahna-Frobenius 3 condition ddt E = 0 which is trivially fulfilled due to nilpotency of exterior differentials [27].
The time-vertical (or spatial) sub-bundle π : V E → E of the tangent bundle π : T E → E is made of those vector fields v : E → V E such that: Each observer also defines in the spacetime manifold of events a congruence of time-lines, generated by a vector field of transversal time-arrows Z : E → T E. Time-lines do intersect each space-slice just at a single event.
A direct sum decomposition T E = V E ⊕ Span(Z) holds and tangent vectors V ∈ T E may be uniquely split as sum V = v + λ · Z with v ∈ V E and λ ∈ .
The time-vertical bundle V E with dim(V E) = 3 is endowed with a field of metric tensors which are symmetric and positive covariant tensors g ∈ Cov(V E) , so that each spatial slice is a Riemann manifold. 4 The dual bundle (V E) * = (Span(Z)) • is made of those covectors in (T E) * which vanish on time-arrows Z ∈ T E. 5 The bundle (V E) * may also be identified with the factor bundle Without loss of generality the following tuning property can be assumed: The notion of time-arrow field, referred to as observer rigging in [28] and observer field in [29,30], is a discriminant which puts into evidence the superiority of a proper 4D spacetime context when compared with the standard 3D spatial context where the time plays just the role of ordering parameter. 6 Given a pair of vectors and covector frames: Movements form a commutative group under the composition rule: The spacetime velocity of motion is defined by: Taking the derivative ∂ α=0 of Eq. (10) 1 , we have: Comparing Eq. (13) with the tuning property Eq. (3), we get the decomposition in space and time components: with the spatial component qualified by time-verticality: Moreover, by definition of push: taking the derivative ∂ α=0 in Eq. (11), we infer the spacetime velocity is pushed by the motion:

Change of observer
Definition 1 (Change of frame) In mechanics a change of frame is an automorphism ζ : E → E of timebundles t E : E → Z and t ζ : E → Z over the identity id Z : Z → Z , according to the commutative diagram: From Eq. (18), commutativity of push and exterior differential gives: Newton frame changes are characterised by no change of clock rates: The spacetime motion ϕ α : T → T is pushed by the change of frame to a motion (ζ ↑ϕ) α : T ζ → T ζ , with T ζ = ζ (T ) , as depicted by the commutative diagram: and the relevant velocities are related by push: Moreover, comparing the splitting of the spacetime velocity V ζ ↑ϕ performed by the privileged observer and pushing the splitting of the spacetime velocity V ϕ in Eq. (14), we get: The computations in Eqs. (22) and (23) clarify why in the spacetime context covariance (commonly referred to as objectivity) of the spacetime velocity V ϕ holds true under any change of observer.
This conclusion is not in contrast with the lack of objectivity concerning the spatial component v ϕ of the velocity, declared in [2, p. 555].
Indeed, the transformation rule: expressing objectivity of the spatial component v ϕ of the velocity, according to Truesdell & Noll [3] and followers, implies by Eq. (23) that ζ ↑Z = Z which means there is no spatial relative velocity between the involved frames. 9 In classic treatments [3], where time is a just an ordering parameter, the time-arrow field Z is not even defined and the discussion above is simply out of reach.

Material, spatial and referential fields
Material, spatial and referential fields are defined as follows [31]: -Material vector fields are based on the trajectory and take values on the time-vertical tangent bundle to the trajectory manifold. Material tensor fields are linear maps on material vector and covector fields. -Spatial vector fields are based on the trajectory and take values on the tangent bundle to a space-slice manifold. Spatial tensor fields are linear maps involving spatial vector and covector fields. -Referential vector fields are based on a chosen local reference shape and take values on the relevant tangent bundle.
After [3,6] referential fields are often improperly referred to as material fields in the literature. No one-toone correspondence may exist between material and spatial fields, according to the clear distinction made in [31] and recalled in above. Moreover: -Spatial vectors can be parallel transported along any line c : → E drawn in spacetime endowed with a connection. The forward parallel transport is denoted by ⇑ with inverse ⇓. Invariance under parallel transport results in vanishing of the covariant derivative ∇ along the vector t = ∂ λ=0 c(λ) tangent to the line: Parallel transport to vectors based on alien manifolds outside the ambient spacetime is not feasible due to lack of a connection. -Material vectors can be transformed along the motion by push ↑ (with its inverse ↓ ) to get other material vectors, or along a placement map to get referential vectors. Invariance under push or pull results in vanishing of the Lie 10 (or convective) derivative along the motion: Parallel and convective derivatives of any tensor field are defined by Leibniz rule, due to invariance of their scalar values. 9 The relative velocity of another observer with respect to the given observer is the push ζ ↑Z of the velocity Z of the spatially standing motion. In Newton frame-changes Eq. (20) is fulfilled and the spatial component of the relative velocity is the difference ζ ↑Z − Z [32]. 10 So named after the Norwegian geometer Marius Sophus Lie (1842-1899).

Genesis of Nanson's formula
A configuration Ω of a body is detected by an observer as a spatial slice S ∩ T of the trajectory manifold T . The formula contributed by Nanson in [33] deals with an outer oriented material surface Σ with dim(Σ) = 2 , drawn in a configuration Ω of geometric dimension dim(Ω) = 3.
Nanson investigated about the changes of area and unit normal versor of the material surface Σ by effect of the motion ϕ : T → T along the trajectory T of dimension dim(T ) = dim(E).
Normal lines the material surface Σ embedded in Ω are generated by those vector fields n Σ : Σ → V Σ E whose associated one-forms g · n Σ : Σ → (V Σ E) * annihilate the tangent bundle T Σ : Then dim Span(n Σ ) = 3 − 2 = 1. 11 Spatial vector fields n Σ : Σ → V Σ E normal to Σ are assumed to be g-unitary, i.e. such that g(n Σ , n Σ ) = 1 and the sign of n Σ is conforming to the outer orientation of the material surface Σ. 12 Similarly, the unit spatial vector field n ϕ(Σ) : ϕ(Σ) → V ϕ(Σ) E normal to the moved material surface ϕ(Σ) , is such that g(n ϕ(Σ) , n ϕ(Σ) ) = 1 with the corresponding covector field g·n ϕ(Σ) : ϕ(Σ) → (V ϕ(Σ) E) * fulfilling the annihilation property: The area two-forms μ Σ ∈ Max(T Σ) and μ ϕ(Σ) ∈ Max(T (ϕ(Σ))) induced on Σ and on ϕ(Σ) by the volume three-form μ g ∈ Max(V E) are defined by contractions with the pertinent unit normal field: To simplify, we often omit the dot symbol (·) of linear dependence and adopt, for the time-vertical restricted tangent movement, the familiar notation: 13 The determinant J ϕ := det(F ϕ ) of this tangent movement is defined by Euler-Jacobi identity [34,35]: 14 In a frame {a} on V E as in Eq. (4), the pullback is given by: Proposition 1 (Nanson) Area-form μ Σ and unit normal n Σ to an outer oriented surface Σ ⊂ T and area-form μ ϕ(Σ) and unit normal n ϕ(Σ) to the pushed surface ϕ(Σ) fulfil on Σ the relation: 11 Span(n Σ ) is the linear subspace of V Σ E generated by n Σ ∈ V Σ E. 12 The polar of T Σ is defined by (T Σ) 13 The notation F , with no indication of underlying map, was first proposed in [3]. 14 Euler-Jacobi formula Eq. (31) holds for any motion ϕ : T → T and any form μ g ∈ Max(V E). Indeed volume forms on a manifold are proportional one another.
Proof Let {a, b} be a frame at T Σ.
Then the triplet {n Σ , a, b} is a frame on V Σ E and the triplet is a frame on V ϕ(Σ) E. The tensor product gn Σ ⊗ n Σ defined by: is the g-orthogonal projector on Span(n Σ ) and I − gn Σ ⊗ n Σ is the complementary projector on T Σ. Similarly for n ϕ(Σ) . Hence, for any v ∈ V Σ E and F ϕ v ∈ V ϕ(Σ) E : Recalling Eq. (29), rewritten in terms of tensor products, Eq. (36) yields: Each one of Eqs. (37) states no more than: The volume of a parallelepiped is height times base-area. Nanson's formula Eq. (33) is then just a rephrasing of Euler-Jacobi identity Eq. (31), with volume given by height times base-area.

Role of the cofactor map
Central in the investigation about Nanson's formula is the notion of cofactor map, as defined e.g. in [11,Eq. (1.17)] and [17, §2.11]: The tangent F ϕ : V Σ E → V ϕ(Σ) E and its g-adjoint F A ϕ : V ϕ(Σ) E → V Σ E according to the metric g ∈ Cov(V E) are related by the defining identity, for any v ∈ V Σ E and u ϕ ∈ V ϕ(Σ) E : Equation (39) is equivalent, in terms of dual maps, to: Symmetry of g means g * = g and Eq. (40) gives: It is easy to see that: stating inverse and adjoint are commuting operations.
To clarify the role of the cofactor map let us observe that the image F ϕ n Σ ∈ V ϕ(Σ) E of a normal versor n Σ ∈ V Σ E transformed by a motion ϕ : T → T , is a vector not necessarily still normal to the moved surface ϕ(Σ). 15 On the other hand, by definition of cross product a × b ∈ V Σ E : the vecor a × b ∈ V Σ E is orthogonal to the surface Σ and hence proportional to n Σ ∈ V Σ E. Moreover, for any c ϕ ∈ V ϕ(Σ) E we get the chain of equalities: Since the expression at the l.h.s. of Eq. (45) vanishes for any c ϕ ∈ V ϕ(Σ) E tangent to ϕ(Σ) , it follows the vector field: is orthogonal to the moved surface ϕ(Σ) , hence proportional to the unit normal n ϕ(Σ) , as quoted in [11,Eq. (1.18)].
A further equivalence is proven by considering, for any pair of tangent vectors a, b ∈ T x Σ , the generated parallelogram tangent to Σ at x ∈ Σ with area A Σ and the pushed parallelogram with sides F ϕ a, F ϕ b ∈ T ϕ(x) ϕ(Σ) tangent at ϕ(x) ∈ ϕ(Σ) with area A ϕ(Σ) , defined by: The pushed area is given by: Let us prove equivalence of Eqs. (47) and (52).
In particular, setting u ϕ = n ϕ(Σ) in Eq. (54), we get: This means that the factor of proportionality between the image of n Σ through the cofactor map, given by cof(F ϕ ) · n Σ , and the normal unit vector n ϕ(Σ) to ϕ(Σ) , which are parallel vectors, is equal to the area-ratio: The equivalence between the Euler-Jacobi volumetric identity Eq. (31) and Nanson formula Eq. (33) leads to the conclusion that this formula is of a purely kinematical nature, dealing with a surface Σ embedded in a 3D configuration Ω of a body in motion along the dynamical trajectory.
The unit normals to a surface pushed by the motion are spatial vectors which are not dragged by the motion.

Nanson's formula in the literature
In the literature Nanson's formula has been variously quoted. A flawed expression was exposed in Rodney Hill early paper [36,Eq. (31)]. With the present notations, the statement therein was: with ρ and ρ ϕ volumetric mass densities in Ω and ϕ(Ω). By fixing a frame {a} in T Ω , Euler-Jacobi volumetric identity Eqs. (31)-(32) yield for the volumes: Hence, by conservation of mass (spuriously implicit in Eq. (57)): Equation (57) exposed in [36] consists then of the flawed geometric statement: In Eq. (60) area weighted unit normals are improperly related by multiplication according to volumetric ratio and push forward.
This wrong expression has to be corrected to Eq. (52) where area weighted unit normals are related by push according to the cofactor map cof(F ϕ ) := J ϕ · F −A ϕ . As a matter of fact, the flaw in Eqs. (57)-(60)  A trivial typo, resulting however in a serious flaw, is found in more recent lecture notes by Ragnar Larsson ). This mistake is there hidden by an over-simplified notation with domain and range manifolds not evidenced.

Material surfaces
Let us now consider the case of a material surface Σ , with dim(Σ) = 2 , moving along its own trajectory In this case, the surface motion ϕ Σ : T Σ → T Σ is defined only on this material trajectory, and the time-vertical restricted map: acts only on the time-vertical (2 + 1)D tangent bundle V T Σ . Therefore, the Euler-Jacobi determinant: is of order two and evaluable just on the bundle V T Σ by: By Eq. (51), this expression is equivalent to: Unlike the previous cofactor map defined in Eq. (38), the new one: doesn't act on n Σ ∈ V Σ E to Σ , which are out of its domain of definition V T Σ . The investigation can anyway still be carried out by considering just frames on Σ and their push by the motion ϕ Σ : However, Nanson's formula is no longer applicable. In this respect the remark after Eq. (1.19) in [11] is worth to be carefully read and mentioned.

Equilibrium and stress fields revisited
The formidable theoretical construction, which stands as monumental basis of Continuum Mechanics (CM), was mainly set out in the course of the XVIII century with basic contributions by Jacob, Johann and Daniel Bernoulli, Leonhard Euler, Jean-Baptiste Le Rond d'Alembert, Joseph-Louis Lagrange, and at the beginning of the XIX century with the decisive completion due to Augustin-Louis Cauchy.
The whole edifice of Continuum Mechanics is based on stating equilibrium as a variational condition concerning a force system f Ω composed of bulk actions per unit volume on the elements of a regularity partition P(Ω) of the body configuration Ω and surficial actions per unit area on its boundary ∂P(Ω).
The linear space V Ω is made of kinematic fields v ∈ V Ω each with square integrable gradient in the elements of the regularity partition P(Ω).
Firm and bilateral boundary constraints, acting between adjacent elements of the partition P(Ω) , consist of the assignment of a linear subspace L Ω ⊂ V Ω of conforming kinematic fields.
In boundary value problems the conformity space includes all kinematic fields vanishing in a boundary layer of each element of the partition P(Ω).
Frictionless constraints are characterised by vanishing of virtual power performed by their interaction field for any conforming kinematic field.
The property of equilibrium means the force system doesn't perform virtual power for any kinematic field δv : Ω → T Ω S which is rigid, i.e. such that any line drawn in the body does not tend to change its length.
The rigidity notion was introduced as far as 1586 by Simon Stevin Principle of Solidification [40]: " The state of equilibrium of a deformable body is not altered if any part of it is replaced with a rigid body of the same geometry." Peculiar to 3D Continuum Mechanics is the implicit description of the fields of virtual infinitesimal isometries δv : This notion is based on conceiving a virtual motion δϕ λ : δT → δT along a virtual trajectory δT ⊂ S Ω drawn in the spatial slice S Ω across the configuration Ω , with δϕ 0 identity in δT .
Virtual infinitesimal isometries are equivalently expressed by vanishing of the associated mixed Euler stretching tensor [41,42]: which provides the local rate of elongation of any material line drawn in S Ω through a point of evaluation in Ω , in the sense that: beingḣ := ∂ λ=0 (δϕ λ ⇓(δϕ λ ↑h)).
In Eq. (67) defining the Euler stretching tensor, ∇ is the covariant derivative associated with the spatial metric g in the Euclid space S Ω [43].
The relevant Levi-Civita connection is metric preserving, flat and symmetric. The corresponding (path independent) foreword parallel transport by translation is ⇑ , with inverse ⇓.
The linear Euler operator D : V Ω → H Ω maps David Hilbert space V Ω of regular kinematic fields into the Hilbert space H Ω of symmetric g-mixed tensor fields which are square integrable in the elements of the partition P(Ω).
The g-adjoint ∇(δv) A of ∇(δv) is defined by the condition: Symmetric and antisymmetric parts are defined by: The inner product between mixed tensors is computed, in terms of the linear invariant J 1 and in cartesian coordinates by: According to this inner product, g-symmetric and g-antisymmetric mixed tensors are orthogonal one another.
The stretching operator: D : V Ω → H Ω is linear and continuous from the Hilbert space V Ω of the vector fields in Ω to the Hilbert space H Ω .
The notion of equilibrium is expressed by vanishing of the virtual power expended by the force system f ∈ F Ω for any virtual rigid kinematic field: 16 In Eq. (72) {F Ω , V Ω } is a pair of Hilbert spaces in duality. 17 The implicit form of the rigidity constraint allows us to introduce a field of Lagrange multipliers with the mechanical meaning of Cauchy stress field T , a basic observation attributed in [2, fn,4, p. 595] to Gabrio Piola [45].
Let μ g and ∂μ g := μ g · n be the metric volume and area-forms on Ω and on its boundary ∂Ω , with n outward unit normal vector to ∂Ω.
Volterra integral transformation in terms of the exterior differential d , yields, for 3D integral domains, the expression of Gauss-Ostrogradskij divergence theorem of a smooth vector field u : Ω → T Ω : with ∂Ω , boundary of Ω , oriented consistently with the orientation of Ω.
Let Σ Ω be the linear space of mixed tensor fields T ∈ Mix(T Ω) which are g-square integrable in each element of P(Ω) with their divergence, and the linear space V Ω of vector fields δv ∈ Mix(T Ω) which are g-square integrable in each element of P(Ω) with their gradients.
Setting u = T A δv in Eq. (73) with T ∈ Σ Ω and δv ∈ V Ω we may define the divergence Div of a mixed tensor field T by a formal application of Gottfried Wilhelm Leibniz rule of calculus: This yields George Green integral formula [46,47]: The boundary integral on the l.h.s. of Eq. (75) may be rewritten as: The gradient can be split in g-symmetric and g-antisymmetric parts: with: D(δv) = sym g ∇(δv), In Eq. (79), we may assume T = sym g (T) without loss of generality, due to orthogonality of g-symmetric and g-antisymmetric mixed tensors. Then: (T, W(δv)) g = (sym g T, W(δv)) g = 0, and the inner product (T, D(δv)) g of Eq. (75) may be replaced with the product (sym g T, D(δv)) g . The representation of the equilibrium system of forces is thus replaced with an equivalent one made of the sole: -volumetric bulk force vector field −Div(sym g T) , -surficial force vector field (sym g T)n on ∂Ω.
Symmetry of T is a convenient choice 18 and not a consequence of rotatory equilibrium, as usually stated and still recently repeated [24] in the literature.
On the other hand, the assumption of a non-symmetric stress tensor is legitimate although non-convenient for description purposes and leads to representative bulk torques per unit volume.
What is not allowable, in the context of Euler implicit representation of rigidity constraint Eq. (67), is the presence of surface torques per unit boundary area.
In the regard, most attempts of modifying Euler implicit representation have clashed against redundancy of the implicit formulations, as pointed out in [49] and recently discussed in [50].
The duality Eq. (79) trivially entails the following polarity properties: For 3D bodies, the polarity Eq. (82) can be stated in the context of Functional Analysis by means of Banach closed range theorem and Korn inequality [44].
The former polarity in Eq. (82) provides the variational condition of kinematic compatibility for stretching fields ε ∈ Ker(D ) • ⊂ H Ω .
Existence and uniqueness (to within additive rigid kinematic fields in Ker(D) ) of a kinematic field v ∈ V Ω fulfilling this variational compatibility condition is proven by Eq. (82) 1 which gives: The latter polarity in Eq. (82) provides the variational condition of static equilibrium for force fields. Existence and uniqueness (to within additional self-equilibrated fields) 19 of the stress field T fulfilling the equilibrium condition f ∈ Ker(D) • follows form Eq. (82) 2 : This proves existence of a Cauchy stress field T ∈ Mix(T Ω) in equilibrium with the force system f ∈ F Ω .
By virtue of the duality in Eq. (79), from Eq. (84) we infer the virtual power principle (VPP) [50]: In terms of mass-form m = ρ μ and of natural stress Σ = ρ −1 T , with ρ scalar mass density, the VPP writes: In classical Continuum Dynamics the bulk force system: includes, á la d'Alembert [51], both the external force system b ext and the inertia term −ρ a , with a spatial acceleration field in inertial frames and ρ scalar mass density.
To eliminate constraint reactions L • Ω ⊂ F Ω from the equilibrium condition, it suffices to restrict the space of regular kinematic fields to the conformity subspace L Ω ⊂ V Ω . The restricted stretching operator D : L Ω → H Ω admits the dual equilibrium operator: defined by the variational condition: In Eq. (88) the Hilbert space L Ω , dual of L Ω , are made of active force systems, quotient of force space F Ω modulo reactive systems in L • Ω . The redundancy issue, first enlightened and investigated in [49], has been recently further addressed in [50]. Referential formulations of equilibrium will be discussed, with relevant critical observations, in Sect. 6.

Evolutive equilibrium
Let us consider at time t ∈ Z a body in the spatial configuration Ω whose kinematics is defined by a linear space V Ω of piecewise regular spatial vector fields v Ω : The body is assumed to be subject to affine kinematic constraints described by a linear subspace L Ω ⊂ V Ω of conforming kinematic fields and by an imposed kinematic field w Ω ∈ V Ω , under the action of a force system f Ω ∈ F Ω = (V Ω ) and of the reactive system exerted by the affine constraints.
Admissible kinematic fields at Ω must accordingly belong to the affine set A Ω := w Ω + L Ω . In a time lapse α ∈ Z the movement ϕ α : T → T along the trajectory T ⊂ E drawn by the body motion in spacetime, is governed by a control system which prescribes increments of the force system Δf Ω ∈ (V Ω ) and of driven displacement field Δw Ω ∈ V Ω and also the update of the conforming kinematic subspace from the initial L Ω ⊂ V Ω to the final one L ϕ α (Ω) ⊂ V ϕ α (Ω) after the time lapse α.
Let us assume for simplicity a smooth quasi-static evolution with the material not leaving the elastic range. A initial guess on the realisation of the body configuration and of the control upgrade at the end of the incremental step can be performed by evaluating increment Δu ∈ Δw Ω + L Ω of displacement from the configuration Ω associated with the data increment {Δf Ω , Δw Ω } at the beginning of the time step.
This evaluation of the incremental displacement is based on the rate formulation of equilibrium [50,55] in which test fields are assumed to be parallel transported by the motion along the trajectory, in accord with Euler law of Dynamics.
The incremental stretch mixed tensorḊ(v ϕ , δv) is associated with the covariant tensor [55]: so that: Elasticity is expressed by the incremental constitutive relation: The elastic compliance H(Σ) is the Hesse operator of a smooth strictly convex scalar stress potential Ξ : 20 Then H(Σ) is positive definite, hence invertible, and such is the elastic stiffness K(Σ) = H(Σ) −1 so that: In a purely elastic range: with Δu ∈ Δw Ω + L Ω ⊂ V Ω , for all test fields δv ∈ L Ω . Provided the variational problem in Eq. (98) admits a displacement solution Δu , a first trial ϕ α (Ω) for the deformed configuration is available. 21 Then, with ΔΣ given by Eqs. (96)-(97), the incremental elastic response Δr ϕ α (Ω) ∈ F ϕ α (Ω) is given by: (99) 20 The fibre derivative d F is taken along spatial directions at the pertinent time instant. 21 In general a dynamical analysis is needed to ensure existence. Uniqueness breaks down when the incremental elastic response becomes singular and instability phenomena take place. Early computational investigations about stability and accuracy were provided in [56,57]. A more recent presentation can be found in [58].
Note that conservation of mass along the motion is expressed by: The response in Eq. (99) is compared with the incremental force in which δv ϕ α (Ω) is parallel transported by the movement ϕ α : T → T : where ⇑ denotes the forward distant parallel transport in the Euclid space S from the spatial configuration Ω to the displaced one ϕ α (Ω). At this stage the force increment Δf ϕ α (Ω) ∈ F ϕ α (Ω) is updated by the control system. The incremental equilibrium gap Δr ϕ α (Ω) − Δf ϕ α (Ω) is applied to the trial configuration ϕ α (Ω) corresponding to the running iteration.
The elastic equilibrium Eq. (98) with ϕ α (Ω) taking the place of Ω , updates the current guess and another iteration for the elastic incremental displacement solution is performed. The iterative loop comes to a stop provided that a suitably chosen norm of the equilibrium gap becomes smaller than a prescribed tolerance.

Referential equilibria: critical notes
In application-oriented treatments, a referential configuration Ω ref of a 3D body is a domain of simple shape in 3 (simplexes or parallelepides) which are put in diffeomorphic correspondence with similar elements of an approximate finite partition of the actual body configuration Ω under investigation.
The one-to-one differentiable placement map p : Ω ref → Ω does not admits an explicit inverse map, as a rule. This is the case, for instance, in the highly successful Finite Element Method (FEM) developed since the sixties of the past century.
In theoretical treatments exposed in the literature, the situation is different and the referential configuration is considered to be occupied by the whole body, or by an element of it, at some instant of time in the relevant ambient Euclid space S.
There is clear evidence that equilibrium cannot be investigated in referential configurations other than isometric or uniformly scaled avatars considered in structural calculations, an opinion shared by any welltrained structural engineer.
In spite of this, the confusion-illusion induced by scholars who considered referential equilibrium as a valid convenient alternative, quickly became viral and is nowadays widespread in the relevant literature, as witnessed by the list in Sect. Therein external forces acting on the actual configuration Ω of the body are assumed to consist of boundary forces per unit area s ∈ T ∂Ω S and of bulk forces per unit volume b ∈ T Ω S.
According to Eq. (10.5) 1 of [11], these force densities are (tacitly) parallel transported from the actual configuration Ω to a reference configuration along the map p −1 : Ω → Ω ref and (explicitly) scaled in such a way that the force resultant is left invariant.
The portrayed surficial and volumetric forces in Eq. (10.5) of [11]: are, in our notation, given by: With a common abuse in Euclid spacetime, the foreword parallel transport ⇑ and its backward inverse ⇓ are usually assumed to be the identity and therefore ignored. In [6, Ch.IX §27] Gurtin says: "It is not convenient to work with T (Cauchy stress) since the deformed configuration is not known in advance." 22 In [11], after Eqs. (10.5), Podio-Guidugli says: "With these definitions we can portray the surface-force and volume-force fields s, b in the reference shape." The portrait adopted in [11] and here reported in Eq. (104) is motivated by the change of integration domain from Ω ref to Ω and from ∂Ω ref to ∂Ω and by the corresponding changes of volumetric and surficial Jacobi determinants.
Resultant force and resultant moment are evaluated in the actual configuration Ω by the integrals: with r : S Ω → Ω position vector field with respect to a given pole O ∈ S Ω . In the reference manifold Ω ref a formulation in a convected frame gives: Indeed, as remarked in [3,], according to Antonio Signorini [59,60] the vanishing of Eq. (106) 2 is to be considered as a compatibility condition to be verified a posteriori, that is once the placement map will eventually be available for the evaluation of the deformed configuration. No strategy to get this goal was, however, envisaged.
In the treatment originally exposed in [3,] the integration domains Ω ref and Ω = p(Ω ref ) were improperly denoted by a same symbol B and in the treatments [6,11] the position vector field p −1 • r : S Ω → Ω ref is deceptively denoted by the same symbol as the position vector r : S Ω → Ω in the actual configuration Ω.
As described in Sect. 5 dealing with evolutive equilibrium, in a nonlinear structural analysis at each instant of time the control system assigns loads and kinematic constraints acting on the actual configuration of the body at that time instant. Familiar examples are a bicycle being driven along a country path or a car moving along city streets.
A shape-independent system of loads is therefore conceivable only in the trivial case of spatially standing motions. 23 Remark 2 In [3,Eq. (44.20) p. 128] an attempt is suggested to recover a referential equilibrium condition in the standard form: 22 That the deformed configuration be not known in advance is an essential feature of structural analysis where the deformed configuration in an incremental time step, or rather a suitable approximation of an element of it, is the objective of the computation. 23 In the footnote on p. 40 of [11, Ch.II §12]: "While we here confine attention to shape-independent systems of loads, to describe realistic body-environment interactions it is essential to study loads that depend on shape." This statement shows the author of [11] realised that something was not going the right way. The passage from the compatibility condition Eq. (106) 2 to the vanishing of the resultant moment in Eq. (107) 2 , was attempted by appealing to a theorem due to Daniel Augusto da Silva [61], exposed in [3,Eq. (44.18) The whole procedure is of no interest in Structural Mechanics where kinematic constraints are usually considered acting on the body and further conditions of admissibility, in addition to equilibrium, are imposed on the stress field, as usually made for instance in constitutive relations of perfect plasticity.
For completeness sake, let us reproduce hereafter the treatment usually described in the literature. In the sequel the tangent placement is denoted by: and Eq. (38) is rewritten as: Putting Σ = ∂Ω ref in Nanson formula Eq. (52) 1 we get: Equation (79) provides for the traction the formula s := T · n ∂Ω in terms of the Cauchy stress T and of the normal versor n ∂Ω to ∂Ω.
Then Nanson formula Eq. (110) leads to the evaluation: Moreover, the portrait proposed in Eq. (10.5) of [11] and reproduced here in Eq. (104) 1 with an explicit parallel transport: when substituted in Eq. (111), gives: Under the special assumption that the reference placement belongs to the dynamical trajectory, 25 the spatial metric tensor g : V E → (V E) * is well defined therein and the following analysis may be applied.
This reshaping was aimed to preserve translational equilibria in rewriting the equilibrium conditions pertaining to the deformed configuration Ω in terms of fields defined in the initial configuration Ω ref . As observed above, the argument is not applicable when dealing with rotatory equilibria where the placement map p : Ω ref → Ω plays an essential role. As

Concluding remarks
In computational contexts, local reference manifolds are alien geometric elements conceived to perform constitutive evaluations by means of forth and back local transformations of material fields between actual and reference shapes, with the choice of a local reference element dictated by convenience.
Physical evidence shows that a material body is a time-parametrised family of configurations, each made of simultaneous events dragged by the motion along the trajectory and characterised by mass invariance.
Adoption of this physically based definition eliminates the need for the introduction of an alien (that is, not included in the trajectory) body configuration B , whose points are candidates for being called particles, and the need for endowing it with specific topological properties. Rather, a particle is what experiments reveal: -a one-parameter family (a time-parametrised curve) of motion-related material events along the trajectory.
Only invertible smooth correspondences between referential elements and elements of a partition of the current configuration are usually available.
Push-pull transformation of material tensor fields to referential ones, as specified in Sect. 2.3, can be performed at corresponding points and may be convenient for numerical integration in computational procedures, such as the finite element method (FEM) and similar ones.
Parallel transport of spatial fields is not applicable when the reference configuration is an alien manifold, a convenient finite subdomain with simple shape in , 2 or 3 . This is the case occurring in engineering methods of analysis and their applications. Simple direct questions arise: "What structural engineer would ever carry out the calculation of an arch stone bridge by writing the equilibrium equations on a straight beam of reference?". 28 "How to compare equilibrium features of a compressed column under gravitational loads with the ones of a referential cantilever (with a 90 degrees rotated axis) under parallel transported loads?".
Equilibrium conditions ought to be imposed in variational form on the current configuration. The geometry of this configuration plays a decisive role in detecting uniqueness, bifurcation and dynamical instability (e.g. snap-through) phenomena.
Equilibrium in Continuum Mechanics, with a discussion about symmetry of Cauchy stress, has been revisited in Sect. 4.
Incremental equilibrium is based on formulation of rate equilibrium illustrated in geometric terms in [50,55].
The investigation is motivated by and aimed at improving theoretical formulations in Nonlinear Structural Mechanics.
The target is a critical analysis of previous treatments in order to promote an elimination of still persisting gaps and contradictions between formal presentations and applicative methodologies devised to render directly available to structural engineers effective computational tools of nonlinear analysis.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.
Funding Open access funding provided by Università degli Studi di Napoli Federico II within the CRUI-CARE Agreement.

Conflict of interest
The authors have no competing interests to declare that are relevant to the content of this article.