# Hyperelastic Geometrically Nonlinear Inverse 3D-FEM Truss Analyses Based on VaReS

Authors: Klaus Bernd Sautter; Kai-Uwe Bletzinger
Journal: Advances in Civil Engineering
 (2022)

---

## Abstract

Direct usage of construction plans as input for structural analyses assumes the reference configuration to match the engineering drawings. However, the built construction is typically supposed to match the construction plans after its successful erection. In that state, the structure is usually already subjected to self-weight and maybe other loadings. Consequently, an analysis approach is necessary to find the unknown reference configuration for a given, desired deformed structural shape. The standard static problem needs to be reformulated with the reference coordinates being the unknown variables. This work describes the necessary steps for geometrically and materially nonlinear truss elements based on the variation of reference strategy (VaReS) and gives a highly detailed description of all resultant system derivatives. Arbitrary hyperelastic material laws can be applied of which this work introduces the St. Venant-Kirchhoff, the Neo-Hookean, and the Ogden law. Additionally, the self-weight load case is considered, increasing the problem’s nonlinearity. Finally, two- and three-dimensional structural problems are presented to show the solution capabilities, ranging from simple 3-bar systems to larger framework bridges. While all necessary vectors and matrices are discussed and presented in great detail, a publicly available GitHub repository makes the code freely accessible as Python code.


---

## Body

## 1. Introduction

Direct usage of construction plans as input for structural analyses assumes the reference configuration to match the engineering drawings. However, the built construction is typically supposed to match the construction plans after its successful erection. In that state, the structure is usually already subjected to self-weight and maybe other loadings. Consequently, an analysis approach is necessary to find the unknown reference configuration for a given, desired deformed structural shape.In conventional structural analysis, a given configuration is typically loaded and subsequently simulated, resulting in deformations and stress states. Often, this analysis is unsuitable, and an “inverse” analysis must be performed. Especially for lightweight membrane or cable net constructions, the internal stress state is predefined, and a corresponding equilibrium shape is searched. Moreover, the desired deformed shape is given for specific structures, and an unknown reference configuration needs to be found. Standard structural analysis is not applicable in these cases, and other solution methods must be utilized. Overall, the finite element method (FEM) has proven to be one of the most effective approaches and will also be the basis for the presented method.Several approaches have been developed to find an equilibrium shape for a given stress state (including external forces). Exploiting proven dynamic simulation strategies the dynamic relaxation approach [1–3] creates a dynamic problem, based on the underlying static problem. Similar to an explicit dynamic simulation, no stiffness matrices are necessary, and a simple node-based solution loop can be used to calculate the artificially damped structure until an equilibrium shape has been found.Directly working with the static problem [4–6] has developed the force density method, whose first successful and practical application can be admired in Munich in the form of the roof of the Olympic Stadium [7]. This approach assembles nodal force equilibrium equations and keeps the ratio of internal force and element length constant, creating a linear system of equations out of a normally nonlinear system. Reference [8] proposes the SGC (statically geometrically coupled) method as an extension to the force density, considering geometrical specifications. The updated reference strategy (URS), developed by [9, 10], is a generalization of the force density method. It is derived from continuum mechanics of elastic bodies and closely follows the FEM approach, allowing its application to any structure. The URS is not restricted to one-dimensional structures compared to the force density method. Based on the URS, form-finding for hybrid structures is discussed in [11], and a description of the whole engineering design cycle is given in [12]. Further approaches are available, such as analyses based on graphic statics [13–15] and constrained form-finding using automated differentiation [16].In contrast to the discussed form-finding strategies, [17, 18] developed the variation of reference strategy (VaReS), inspired by [19]. It is further discussed in [20] and later applied to the concept of isogeometric B-Rep analysis (IBRA) in [21]. VaReS is typically employed to analyze suitable cutting patterns for prestressed membrane structures. It is based on a change of the system unknowns, called the “concept of inverse mechanics” by [17], originating from [22, 23]. A similar strategy is followed by [24] to analyze sheet metal forming. Here, the reference (undeformed configuration) represents the unknown variables instead of the structural displacements. While standard structural analyses typically result in the deformed spatial configuration with a known material configuration, the inverse analysis uses an initially given spatial (deformed) configuration to solve for a suitable unknown material (reference, undeformed) configuration. Figure 1 graphically explains this concept. The idea of “inverse engineering” was already discussed in [25–27] and successfully applied to cutting pattern generations for isotropic materials in the linear and moderately nonlinear regime.Figure 1
A standard static analysis calculates the deformation of a given reference configuration, and the displacements are the unknown variables. The inverse analysis is used to determine the unknown reference configuration for a given deformed state.While a very short introduction to a two-dimensional truss element can be found in [17], this work takes up the “concept of inverse mechanics” to find unknown reference geometries for three-dimensional truss structures. With the help of the presented approach, the engineer can determine the necessary reference construction so that the desired final aim geometry is obtained for a given load state. To allow a wide range of application cases, the base equations are geometrically and materially nonlinear. Arbitrary nonlinear hyperelastic material laws can be used within this workflow, of which the St. Venant-Kirchhoff law, the Neo-Hookean law, and the Ogden law are introduced herein. Consequently, we add the effects of self-weight inclusion on the solution procedure and present various test cases to prove the effectiveness and simulation capabilities of the proposed approach. Particular focus is put on the extensively detailed presentation of all necessary derivatives and every step towards the final solution. Additionally, all system matrices, examples, and a simulation framework are made publicly accessible via a GitHub repository [28] as Python code.The work is thematically structured as follows:(i)
Section2 introduces the general truss formulation, including material law and other assumptions. This section is necessary to later follow up with the detailed investigation of all its derivatives.Subsection2.1 formulates the nonlinear truss formulation, presenting the internal element forces, based on nonlinear stress and strain measuresSubsection2.2 discusses simplifications for linear strain and deformation regimes(ii)
Section3 subsequently formulates the global problem and the resulting system of equationsSubsection3.1 discusses the standard global force equilibrium in the reference configuration.Subsection3.2 afterwards derives the inverse problem, reformulating the force equilibrium from the previous subsection.(ii)
Section3 subsequently formulates the global problem and the resulting system of equationsSubsection3.1 discusses the standard global force equilibrium in the reference configuration.Subsection3.2 afterwards derives the inverse problem, reformulating the force equilibrium from the previous subsection.Subsection3.1 discusses the standard global force equilibrium in the reference configuration.Subsection3.2 afterwards derives the inverse problem, reformulating the force equilibrium from the previous subsection.(iii)
Section4 follows up with a variety of different test cases. Two- and three-dimensional truss structures are discussed. Investigations of convergence rates are included for selected examples.(iv)
Section5 concludes this work.


## 2. Truss Formulation

The truss element in this work is described by two nodes connected by a linear, straight line, as depicted in Figure2. It carries only normal forces and keeps a constant cross section A by neglecting the Poisson’s effect.Figure 2
Linear truss element with nodal degrees of freedom.Each node has three translational degrees of freedom (dofs), which are assembled for every element in the displacement vectorue, respectively,(1)ue=ue,1ue,2=u1v1w1u2v2w2.Additionally, we introduce the discrete reference nodal coordinates per elementXe, as well as the discrete actual nodal coordinates per element xe,(2)Xe=X1Y1Z1X2Y2Z2,(3)xe=x1y1z1x2y2z2.Using the reference coordinate differencesdX, dY, dZ, the reference length L can be calculated as(4)L=dX2+dY2+dZ2,dX=X2−X1,dY=Y2−Y1,dZ=Z2−Z1.The same operation holds for the calculation of the deformed element lengthl(5)l=dx2+dy2+dz2,dx=x2−x1,dy=y2−y1,dz=z2−z1.


### 2.1. Large Deformations—Nonlinear Truss Formulation

Following the well-known FEM approach, the total virtual workδW is expressed by the virtual internal work δWint and the virtual external work δWext(6)δW=δWint−δWext=0.This will result in a stiffness matrix and a force vector to solve the static problem. We will discuss the derivation of the global three-dimensional stiffness matrix in the following, leading to the global problem formulation. The virtual internal work is formulated by integrating over the reference volumeV=AL as the product of the reference area A and the reference length L(7)δWint=∫VSδεdV=ALSδε=ALS∂ε∂ueδue.To allow for the analysis of arbitrarily large deformations and correctly considering rigid body rotations, nonlinear strain and stress measurements are applied. WhileS represents the 1D 2 nd Piola-Kirchhoff stress, the work-conjugated 1D Green-Lagrange strain is described by ε as follows:(8)ε=12l2−L2L2.Based on the chosen material law and its associated Helmholtz free-energy functionΨ, the 2 nd Piola-Kirchhoff stress and the respective material tangent modulus C can be evaluated as follows:(9)S=∂Ψ∂ε,C=∂S∂ε=∂2Ψ∂ε2.S and C are exemplarily given for three different hyperelastic material laws in the following. They can be plugged into the derivations of this work to analyze and simulate different deformation behaviors.(1)
St. Venant-KirchhoffS=EεC=E(2)
Neo-HookeanS=Eε/2ε+1C=Eε/2ε+1(3)
Ogdenβ1,β2 are specific Ogden parameters and need to be chosen appropriately.S=E/β1−β22ε+1β1/2−1−2ε+1β2/2−1C=E/β1−β2β1−22ε+1β1/2−2−β2−22ε+1β2/2−2Specifically, the Ogden law is of major importance as it can be easily transformed into other constitutive laws, such as the St. Venant-Kirchhoff lawβ1=4,β2=2 and the Neo-Hookean law β1=2,β2=0 [29, 30]. Special attention needs to be paid when dealing with large compressive strains. As discussed in [29], the St. Venant-Kirchhoff law is not appropriate for this case anymore, and other material laws are recommended.Finally, the internal element forcesFint,e are derived from (7) [29] as follows:(10)Fint,e=ALS∂ε∂ue.


### 2.2. Small Deformations—Linear Truss Formulation

In many engineering applications, only small deformations and strains occur. For these application cases, a simplified element formulation can be applied, facilitating both the derivation as well as the implementation. The small deformations and strains allow the use of the linear elastic Hooke’s law to express the axial stressσ by the Young’s Modulus E and the linear axial strain ϵ(11)σ=Eϵ,ϵ=∂u∂X.Here,u is the axial displacement field, and X is the coordinate running along the local element axis in its reference configuration. The strains and deformations are thus restricted to be small.


#### 2.2.1. Stiffness Matrix in 1D

We start the derivation in 1D which is here expressed for a one-dimensional truss element as shown in Figure3. The two nodal dofs u1,u2 are sufficient to describe the local element displacement field.Figure 3
Linear truss element with nodal degrees of freedom, 1D.The virtual internal work is again expressed by integrating over the reference volumeV as follows:(12)δWint=∫VσδϵdV.This is used in conjunction with linear shape functions to derive the local one-dimensional element stiffness matrixKel,(13)Kel=EAL1−1−11.


#### 2.2.2. Rotation to 3D

To analyze a global three-dimensional problem, a third dof is introduced. As depicted in Figure2, three dofs at each node describe the total elemental deformation capacities and are given in equation (1). The introduction of the three-dimensional dofs allows the expression of the one-dimensional local element stiffness matrix from equation (13) to the three-dimensional space(14)Kel=EAL100−100000000000000−100100000000000000.as well as the transformation/rotation matrix T as follows:(15)T=n000000n100000n200000000n000000n100000n200,nT=dXdYdZL.The direction vectorn is used to express the orientation of the truss in the three-dimensional space, see Figure 2. Combining the previously derived element matrices, Kel can be expressed with respect to the global, three-dimensional dofs to obtain the global element stiffness matrix 23.(16)Keg=TKelTTFinally, the internal forces for the linear truss formulation are formulated as follows:(17)Fint,e=Kegue.


## 2.1. Large Deformations—Nonlinear Truss Formulation

Following the well-known FEM approach, the total virtual workδW is expressed by the virtual internal work δWint and the virtual external work δWext(6)δW=δWint−δWext=0.This will result in a stiffness matrix and a force vector to solve the static problem. We will discuss the derivation of the global three-dimensional stiffness matrix in the following, leading to the global problem formulation. The virtual internal work is formulated by integrating over the reference volumeV=AL as the product of the reference area A and the reference length L(7)δWint=∫VSδεdV=ALSδε=ALS∂ε∂ueδue.To allow for the analysis of arbitrarily large deformations and correctly considering rigid body rotations, nonlinear strain and stress measurements are applied. WhileS represents the 1D 2 nd Piola-Kirchhoff stress, the work-conjugated 1D Green-Lagrange strain is described by ε as follows:(8)ε=12l2−L2L2.Based on the chosen material law and its associated Helmholtz free-energy functionΨ, the 2 nd Piola-Kirchhoff stress and the respective material tangent modulus C can be evaluated as follows:(9)S=∂Ψ∂ε,C=∂S∂ε=∂2Ψ∂ε2.S and C are exemplarily given for three different hyperelastic material laws in the following. They can be plugged into the derivations of this work to analyze and simulate different deformation behaviors.(1)
St. Venant-KirchhoffS=EεC=E(2)
Neo-HookeanS=Eε/2ε+1C=Eε/2ε+1(3)
Ogdenβ1,β2 are specific Ogden parameters and need to be chosen appropriately.S=E/β1−β22ε+1β1/2−1−2ε+1β2/2−1C=E/β1−β2β1−22ε+1β1/2−2−β2−22ε+1β2/2−2Specifically, the Ogden law is of major importance as it can be easily transformed into other constitutive laws, such as the St. Venant-Kirchhoff lawβ1=4,β2=2 and the Neo-Hookean law β1=2,β2=0 [29, 30]. Special attention needs to be paid when dealing with large compressive strains. As discussed in [29], the St. Venant-Kirchhoff law is not appropriate for this case anymore, and other material laws are recommended.Finally, the internal element forcesFint,e are derived from (7) [29] as follows:(10)Fint,e=ALS∂ε∂ue.


## 2.2. Small Deformations—Linear Truss Formulation

In many engineering applications, only small deformations and strains occur. For these application cases, a simplified element formulation can be applied, facilitating both the derivation as well as the implementation. The small deformations and strains allow the use of the linear elastic Hooke’s law to express the axial stressσ by the Young’s Modulus E and the linear axial strain ϵ(11)σ=Eϵ,ϵ=∂u∂X.Here,u is the axial displacement field, and X is the coordinate running along the local element axis in its reference configuration. The strains and deformations are thus restricted to be small.


### 2.2.1. Stiffness Matrix in 1D

We start the derivation in 1D which is here expressed for a one-dimensional truss element as shown in Figure3. The two nodal dofs u1,u2 are sufficient to describe the local element displacement field.Figure 3
Linear truss element with nodal degrees of freedom, 1D.The virtual internal work is again expressed by integrating over the reference volumeV as follows:(12)δWint=∫VσδϵdV.This is used in conjunction with linear shape functions to derive the local one-dimensional element stiffness matrixKel,(13)Kel=EAL1−1−11.


### 2.2.2. Rotation to 3D

To analyze a global three-dimensional problem, a third dof is introduced. As depicted in Figure2, three dofs at each node describe the total elemental deformation capacities and are given in equation (1). The introduction of the three-dimensional dofs allows the expression of the one-dimensional local element stiffness matrix from equation (13) to the three-dimensional space(14)Kel=EAL100−100000000000000−100100000000000000.as well as the transformation/rotation matrix T as follows:(15)T=n000000n100000n200000000n000000n100000n200,nT=dXdYdZL.The direction vectorn is used to express the orientation of the truss in the three-dimensional space, see Figure 2. Combining the previously derived element matrices, Kel can be expressed with respect to the global, three-dimensional dofs to obtain the global element stiffness matrix 23.(16)Keg=TKelTTFinally, the internal forces for the linear truss formulation are formulated as follows:(17)Fint,e=Kegue.


## 2.2.1. Stiffness Matrix in 1D

We start the derivation in 1D which is here expressed for a one-dimensional truss element as shown in Figure3. The two nodal dofs u1,u2 are sufficient to describe the local element displacement field.Figure 3
Linear truss element with nodal degrees of freedom, 1D.The virtual internal work is again expressed by integrating over the reference volumeV as follows:(12)δWint=∫VσδϵdV.This is used in conjunction with linear shape functions to derive the local one-dimensional element stiffness matrixKel,(13)Kel=EAL1−1−11.


## 2.2.2. Rotation to 3D

To analyze a global three-dimensional problem, a third dof is introduced. As depicted in Figure2, three dofs at each node describe the total elemental deformation capacities and are given in equation (1). The introduction of the three-dimensional dofs allows the expression of the one-dimensional local element stiffness matrix from equation (13) to the three-dimensional space(14)Kel=EAL100−100000000000000−100100000000000000.as well as the transformation/rotation matrix T as follows:(15)T=n000000n100000n200000000n000000n100000n200,nT=dXdYdZL.The direction vectorn is used to express the orientation of the truss in the three-dimensional space, see Figure 2. Combining the previously derived element matrices, Kel can be expressed with respect to the global, three-dimensional dofs to obtain the global element stiffness matrix 23.(16)Keg=TKelTTFinally, the internal forces for the linear truss formulation are formulated as follows:(17)Fint,e=Kegue.


## 3. Global Problem Formulation

### 3.1. Standard Static Force Equilibrium

First, the standard static force equilibrium in the reference configuration is expressed to proceed with the solution process of finding the unknown reference configuration. Then, following the FEM workflow, the global element stiffness matrices of all structural elements are assembled by their dof numbering. The dof-based assembly [31] is from here on described by the ∑ symbol, using a node n or an element e loop.(18)K=∑eKeg,Fext=∑nFext,n,Fint=∑eFint,e,u=∑eue,x=∑exe,X=∑eXe.All stiffness matrices are assembled in the global system stiffness matrixK, the internal element forces are assembled in the global internal force vector Fint, and the external nodal forces are assembled in the global load vector Fext. While Fext is, in general, assumed to be constant, later in Section 3.2.3, the influence of additional self-weight is discussed. Other entities, such as the nodal discrete displacements and coordinates, are as well assembled.Finally, we can express the global force equilibrium of the whole system in equation (19). With the reference configurationXdescribing the system unknowns the expression turns out to be nonlinear, allowing us to express a force residualrequation.(19)FintX,x−Fext=rX,x.Note: For in a standard static analysisxand thus the nodal displacementsu, would be the unknown vector. The standard problem assumes a given (known) reference configuration and solves the force equilibrium in (19) for the unknown displacements.


### 3.2. Inverse Problem

For the problem at hand, however, the unknown is the reference configurationX, making it necessary to solve a system of nonlinear equations. Figure 4 demonstrates the inverse problem, depicting the known deformed configuration x and the unknown reference configuration X.Figure 4
Unknown reference configuration deforms to the desired configuration.We propose to use a Newton–Raphson solution scheme, expressing (19) in the form of a truncated Taylor’s series, with k being the current iteration number(20)rk+1=rk+∂rk∂X︷CkΔXk=0.The necessary derivative of the residual vector is introduced to be theC system matrix(21)ΔXk=Xk+1–Xk=–C–1krk.which is as well assembled by individual element e contributions(22)C=∑eCe.The following two Sections3.2.1 and 3.2.2 discuss the necessary steps to formulate Ce assembled in equation (22).


#### 3.2.1. Large Deformations—Nonlinear Truss Formulation

With respect to the global force equilibrium in (19), the element residual contribution is formulated as(23)re=Fint,e=ALS∂ε∂ue,which is derived in equation (10).The element solution matrixCe is consequently derived by differentiating the element residual re with respect to the unknown reference configuration X.(24)Ce=∂re∂Xe=∂Fint,e∂Xe,=AS∂ε∂ue∂L∂Xe︸I+AL∂ε∂ue∂S∂Xe︸II+ALS∂2ε∂ue∂Xe︸III.(24) can be split into three parts I−III and will be accordingly formulated for clarity in the following.Part I(25)∂ε∂ue=∂ε∂l∂l∂ue∂ε∂l=lL2,∂l∂ue=1ldxe,∂L∂Xe=1LdXe.(26)dxe=−dx−dy−dzdxdydz,dXe=−dX−dY−dZdXdYdzZ,Part IIPart IIIThe preceding derivations allow the expression of the element solution matrixCe for arbitrary hyperelastic material laws. Depending on the chosen material law, the 2 nd Piola-Kirchhoff stress S and the appropriate material tangent modulus C need to be calculated. For a selected choice of suitable material models, Section 2.1 provides the necessary terms.


#### 3.2.2. Small Deformations—Linear Truss Formulation

Following the internal force expression for small deformations in (17), the elemental residual is expressed as(30)re=Fint,e=Kegue.Combining equations (1-3), the nodal dofs are expressed via the difference between the reference and the actual configuration in equation (31). This step is especially important as the reference configuration is unknown in this work, see Figure 4.(31)ue=x1−X1y1−Y1z1−Z1x2−X2y2−Y2z2−Z2.The Solution MatrixCe: From (19), it can be seen that the derivative of the residual with respect to the unknown configuration results in a nonlinear problem. Following (22), the elemental solution matrix Ce is in the following derived on element level. Applying the chain rule to the internal element residual in equation (30).(32)Ce=∂re∂Xe=∂Keg∂Xeue+Keg∂ue∂Xe︸−I=∂Keg∂Xeue−Keg=K^egu^e−Keg.The new tensorsK^eg,u^e are introduced for simplification and are described in the following. This notation allows a straightforward description of the remaining operations. u^e will result in a 36 × 6 matrix, with the elemental nodal displacement vector ue moving from the very left top to the right bottom(33)u^e=ue000000ue000000ue000000ue000000ue000000ue.At the same time,K^eg can be expressed with the help of the single respective derivative of the global element stiffness matrix, with respect to each unknown reference position in space, resulting in a 6 × 36 matrix(34)K^eg=∂Keg∂X1∂Keg∂Y1∂Keg∂Z1∂Keg∂X2∂Keg∂Y2∂Keg∂Z2.The individual derivative expressions in (34) are calculated by the chain rule applied to (16) including the derivatives of the transformation matrix and the local element stiffness matrix.Derivatives ofKeg: The respective derivatives in (34) are discussed in the following, rearranging the global element stiffness matrix from (16) by merging the length terms(35)Keg=TKelTT=EAL3T¯K¯elT¯T.The•¯ matrices are obtained by multiplying the original matrices by the length L, so that the length is separated as shown in (35), rendering a constant matrix K¯el.Now, the derivative ofKeg with respect to an arbitrary reference coordinate Xi can finally be expressed as follows:(36)∂Keg∂Xi=−3EAL4∂L∂XiT¯K¯elT¯T+EAL3∂T¯∂XiK¯elT¯T+T¯K¯el∂T¯T∂Xi.First, the derivative of the length with respect to the reference coordinates can be written in the following form and will later be used again, when including self-weight(37)∂L∂Xe=∂L∂X1∂L∂Y1∂L∂Z1∂L∂X2∂L∂Y2∂L∂Z2,=−dX−dY−dZ+dX+dY+dZL.Additionally, the derivative of the transformation matrixT¯, which has the length excluded, results in the following simple matrices:(38)∂T¯∂X1=−∂T¯∂X2=−100000000000000000000−100000000000000,∂T¯∂Y1=−∂T¯∂Y2=000000−100000000000000000000−100000000,∂T¯∂Z1=−∂T¯∂Z2=000000000000−100000000000000000000−100.


#### 3.2.3. Including Self-Weight

Additional terms need to be considered when including self-weight. As the self-weight depends on the reference element length, the unknown reference configuration influences the load case and adds new nonlinear terms tore.With the help of the mass densityρ and a global gravity-direction vector g, the additional nodal forces per element are expressed as follows:(39)Fgravity,e=0.5ρALX︷massXggearth,gT=gxgygzgxgygz.Withgx,gy,gzbeing the spatial components of the normalized gravity direction vector and the gravitational accelerationgearth≈9.81m/s2. Additionally, to the previously derived element residual re, the self-weight forces need to be included(40)re=Fint,e−Fgravity,e,adding a new term to the element solution matrix Ce. The following equation describes the complete matrix, including all necessary derivatives(41)Ce=∂Fint,e∂Xe−∂Fgravity,e∂Xe.With the help of the length derivatives in equation (35), the new term is calculated as follows:(42)∂Fgravity,e∂Xe=0.5ρAggearth∂L∂Xe.The additional terms toCe in (41) hold for both the linear and the nonlinear truss element.


#### 3.2.4. Code Implementation

In order to generate lasting added value, all derivatives are made publicly available. They can be downloaded as Python code via a GitHub repository [28]. Jupyter Notebook files are included for a direct assessment of all derivatives, functions, and the proceeding examples.


## 3.1. Standard Static Force Equilibrium

First, the standard static force equilibrium in the reference configuration is expressed to proceed with the solution process of finding the unknown reference configuration. Then, following the FEM workflow, the global element stiffness matrices of all structural elements are assembled by their dof numbering. The dof-based assembly [31] is from here on described by the ∑ symbol, using a node n or an element e loop.(18)K=∑eKeg,Fext=∑nFext,n,Fint=∑eFint,e,u=∑eue,x=∑exe,X=∑eXe.All stiffness matrices are assembled in the global system stiffness matrixK, the internal element forces are assembled in the global internal force vector Fint, and the external nodal forces are assembled in the global load vector Fext. While Fext is, in general, assumed to be constant, later in Section 3.2.3, the influence of additional self-weight is discussed. Other entities, such as the nodal discrete displacements and coordinates, are as well assembled.Finally, we can express the global force equilibrium of the whole system in equation (19). With the reference configurationXdescribing the system unknowns the expression turns out to be nonlinear, allowing us to express a force residualrequation.(19)FintX,x−Fext=rX,x.Note: For in a standard static analysisxand thus the nodal displacementsu, would be the unknown vector. The standard problem assumes a given (known) reference configuration and solves the force equilibrium in (19) for the unknown displacements.


## 3.2. Inverse Problem

For the problem at hand, however, the unknown is the reference configurationX, making it necessary to solve a system of nonlinear equations. Figure 4 demonstrates the inverse problem, depicting the known deformed configuration x and the unknown reference configuration X.Figure 4
Unknown reference configuration deforms to the desired configuration.We propose to use a Newton–Raphson solution scheme, expressing (19) in the form of a truncated Taylor’s series, with k being the current iteration number(20)rk+1=rk+∂rk∂X︷CkΔXk=0.The necessary derivative of the residual vector is introduced to be theC system matrix(21)ΔXk=Xk+1–Xk=–C–1krk.which is as well assembled by individual element e contributions(22)C=∑eCe.The following two Sections3.2.1 and 3.2.2 discuss the necessary steps to formulate Ce assembled in equation (22).


### 3.2.1. Large Deformations—Nonlinear Truss Formulation

With respect to the global force equilibrium in (19), the element residual contribution is formulated as(23)re=Fint,e=ALS∂ε∂ue,which is derived in equation (10).The element solution matrixCe is consequently derived by differentiating the element residual re with respect to the unknown reference configuration X.(24)Ce=∂re∂Xe=∂Fint,e∂Xe,=AS∂ε∂ue∂L∂Xe︸I+AL∂ε∂ue∂S∂Xe︸II+ALS∂2ε∂ue∂Xe︸III.(24) can be split into three parts I−III and will be accordingly formulated for clarity in the following.Part I(25)∂ε∂ue=∂ε∂l∂l∂ue∂ε∂l=lL2,∂l∂ue=1ldxe,∂L∂Xe=1LdXe.(26)dxe=−dx−dy−dzdxdydz,dXe=−dX−dY−dZdXdYdzZ,Part IIPart IIIThe preceding derivations allow the expression of the element solution matrixCe for arbitrary hyperelastic material laws. Depending on the chosen material law, the 2 nd Piola-Kirchhoff stress S and the appropriate material tangent modulus C need to be calculated. For a selected choice of suitable material models, Section 2.1 provides the necessary terms.


### 3.2.2. Small Deformations—Linear Truss Formulation

Following the internal force expression for small deformations in (17), the elemental residual is expressed as(30)re=Fint,e=Kegue.Combining equations (1-3), the nodal dofs are expressed via the difference between the reference and the actual configuration in equation (31). This step is especially important as the reference configuration is unknown in this work, see Figure 4.(31)ue=x1−X1y1−Y1z1−Z1x2−X2y2−Y2z2−Z2.The Solution MatrixCe: From (19), it can be seen that the derivative of the residual with respect to the unknown configuration results in a nonlinear problem. Following (22), the elemental solution matrix Ce is in the following derived on element level. Applying the chain rule to the internal element residual in equation (30).(32)Ce=∂re∂Xe=∂Keg∂Xeue+Keg∂ue∂Xe︸−I=∂Keg∂Xeue−Keg=K^egu^e−Keg.The new tensorsK^eg,u^e are introduced for simplification and are described in the following. This notation allows a straightforward description of the remaining operations. u^e will result in a 36 × 6 matrix, with the elemental nodal displacement vector ue moving from the very left top to the right bottom(33)u^e=ue000000ue000000ue000000ue000000ue000000ue.At the same time,K^eg can be expressed with the help of the single respective derivative of the global element stiffness matrix, with respect to each unknown reference position in space, resulting in a 6 × 36 matrix(34)K^eg=∂Keg∂X1∂Keg∂Y1∂Keg∂Z1∂Keg∂X2∂Keg∂Y2∂Keg∂Z2.The individual derivative expressions in (34) are calculated by the chain rule applied to (16) including the derivatives of the transformation matrix and the local element stiffness matrix.Derivatives ofKeg: The respective derivatives in (34) are discussed in the following, rearranging the global element stiffness matrix from (16) by merging the length terms(35)Keg=TKelTT=EAL3T¯K¯elT¯T.The•¯ matrices are obtained by multiplying the original matrices by the length L, so that the length is separated as shown in (35), rendering a constant matrix K¯el.Now, the derivative ofKeg with respect to an arbitrary reference coordinate Xi can finally be expressed as follows:(36)∂Keg∂Xi=−3EAL4∂L∂XiT¯K¯elT¯T+EAL3∂T¯∂XiK¯elT¯T+T¯K¯el∂T¯T∂Xi.First, the derivative of the length with respect to the reference coordinates can be written in the following form and will later be used again, when including self-weight(37)∂L∂Xe=∂L∂X1∂L∂Y1∂L∂Z1∂L∂X2∂L∂Y2∂L∂Z2,=−dX−dY−dZ+dX+dY+dZL.Additionally, the derivative of the transformation matrixT¯, which has the length excluded, results in the following simple matrices:(38)∂T¯∂X1=−∂T¯∂X2=−100000000000000000000−100000000000000,∂T¯∂Y1=−∂T¯∂Y2=000000−100000000000000000000−100000000,∂T¯∂Z1=−∂T¯∂Z2=000000000000−100000000000000000000−100.


### 3.2.3. Including Self-Weight

Additional terms need to be considered when including self-weight. As the self-weight depends on the reference element length, the unknown reference configuration influences the load case and adds new nonlinear terms tore.With the help of the mass densityρ and a global gravity-direction vector g, the additional nodal forces per element are expressed as follows:(39)Fgravity,e=0.5ρALX︷massXggearth,gT=gxgygzgxgygz.Withgx,gy,gzbeing the spatial components of the normalized gravity direction vector and the gravitational accelerationgearth≈9.81m/s2. Additionally, to the previously derived element residual re, the self-weight forces need to be included(40)re=Fint,e−Fgravity,e,adding a new term to the element solution matrix Ce. The following equation describes the complete matrix, including all necessary derivatives(41)Ce=∂Fint,e∂Xe−∂Fgravity,e∂Xe.With the help of the length derivatives in equation (35), the new term is calculated as follows:(42)∂Fgravity,e∂Xe=0.5ρAggearth∂L∂Xe.The additional terms toCe in (41) hold for both the linear and the nonlinear truss element.


### 3.2.4. Code Implementation

In order to generate lasting added value, all derivatives are made publicly available. They can be downloaded as Python code via a GitHub repository [28]. Jupyter Notebook files are included for a direct assessment of all derivatives, functions, and the proceeding examples.


## 3.2.1. Large Deformations—Nonlinear Truss Formulation

With respect to the global force equilibrium in (19), the element residual contribution is formulated as(23)re=Fint,e=ALS∂ε∂ue,which is derived in equation (10).The element solution matrixCe is consequently derived by differentiating the element residual re with respect to the unknown reference configuration X.(24)Ce=∂re∂Xe=∂Fint,e∂Xe,=AS∂ε∂ue∂L∂Xe︸I+AL∂ε∂ue∂S∂Xe︸II+ALS∂2ε∂ue∂Xe︸III.(24) can be split into three parts I−III and will be accordingly formulated for clarity in the following.Part I(25)∂ε∂ue=∂ε∂l∂l∂ue∂ε∂l=lL2,∂l∂ue=1ldxe,∂L∂Xe=1LdXe.(26)dxe=−dx−dy−dzdxdydz,dXe=−dX−dY−dZdXdYdzZ,Part IIPart IIIThe preceding derivations allow the expression of the element solution matrixCe for arbitrary hyperelastic material laws. Depending on the chosen material law, the 2 nd Piola-Kirchhoff stress S and the appropriate material tangent modulus C need to be calculated. For a selected choice of suitable material models, Section 2.1 provides the necessary terms.


## 3.2.2. Small Deformations—Linear Truss Formulation

Following the internal force expression for small deformations in (17), the elemental residual is expressed as(30)re=Fint,e=Kegue.Combining equations (1-3), the nodal dofs are expressed via the difference between the reference and the actual configuration in equation (31). This step is especially important as the reference configuration is unknown in this work, see Figure 4.(31)ue=x1−X1y1−Y1z1−Z1x2−X2y2−Y2z2−Z2.The Solution MatrixCe: From (19), it can be seen that the derivative of the residual with respect to the unknown configuration results in a nonlinear problem. Following (22), the elemental solution matrix Ce is in the following derived on element level. Applying the chain rule to the internal element residual in equation (30).(32)Ce=∂re∂Xe=∂Keg∂Xeue+Keg∂ue∂Xe︸−I=∂Keg∂Xeue−Keg=K^egu^e−Keg.The new tensorsK^eg,u^e are introduced for simplification and are described in the following. This notation allows a straightforward description of the remaining operations. u^e will result in a 36 × 6 matrix, with the elemental nodal displacement vector ue moving from the very left top to the right bottom(33)u^e=ue000000ue000000ue000000ue000000ue000000ue.At the same time,K^eg can be expressed with the help of the single respective derivative of the global element stiffness matrix, with respect to each unknown reference position in space, resulting in a 6 × 36 matrix(34)K^eg=∂Keg∂X1∂Keg∂Y1∂Keg∂Z1∂Keg∂X2∂Keg∂Y2∂Keg∂Z2.The individual derivative expressions in (34) are calculated by the chain rule applied to (16) including the derivatives of the transformation matrix and the local element stiffness matrix.Derivatives ofKeg: The respective derivatives in (34) are discussed in the following, rearranging the global element stiffness matrix from (16) by merging the length terms(35)Keg=TKelTT=EAL3T¯K¯elT¯T.The•¯ matrices are obtained by multiplying the original matrices by the length L, so that the length is separated as shown in (35), rendering a constant matrix K¯el.Now, the derivative ofKeg with respect to an arbitrary reference coordinate Xi can finally be expressed as follows:(36)∂Keg∂Xi=−3EAL4∂L∂XiT¯K¯elT¯T+EAL3∂T¯∂XiK¯elT¯T+T¯K¯el∂T¯T∂Xi.First, the derivative of the length with respect to the reference coordinates can be written in the following form and will later be used again, when including self-weight(37)∂L∂Xe=∂L∂X1∂L∂Y1∂L∂Z1∂L∂X2∂L∂Y2∂L∂Z2,=−dX−dY−dZ+dX+dY+dZL.Additionally, the derivative of the transformation matrixT¯, which has the length excluded, results in the following simple matrices:(38)∂T¯∂X1=−∂T¯∂X2=−100000000000000000000−100000000000000,∂T¯∂Y1=−∂T¯∂Y2=000000−100000000000000000000−100000000,∂T¯∂Z1=−∂T¯∂Z2=000000000000−100000000000000000000−100.


## 3.2.3. Including Self-Weight

Additional terms need to be considered when including self-weight. As the self-weight depends on the reference element length, the unknown reference configuration influences the load case and adds new nonlinear terms tore.With the help of the mass densityρ and a global gravity-direction vector g, the additional nodal forces per element are expressed as follows:(39)Fgravity,e=0.5ρALX︷massXggearth,gT=gxgygzgxgygz.Withgx,gy,gzbeing the spatial components of the normalized gravity direction vector and the gravitational accelerationgearth≈9.81m/s2. Additionally, to the previously derived element residual re, the self-weight forces need to be included(40)re=Fint,e−Fgravity,e,adding a new term to the element solution matrix Ce. The following equation describes the complete matrix, including all necessary derivatives(41)Ce=∂Fint,e∂Xe−∂Fgravity,e∂Xe.With the help of the length derivatives in equation (35), the new term is calculated as follows:(42)∂Fgravity,e∂Xe=0.5ρAggearth∂L∂Xe.The additional terms toCe in (41) hold for both the linear and the nonlinear truss element.


## 3.2.4. Code Implementation

In order to generate lasting added value, all derivatives are made publicly available. They can be downloaded as Python code via a GitHub repository [28]. Jupyter Notebook files are included for a direct assessment of all derivatives, functions, and the proceeding examples.


## 4. Examples

Several examples are presented in the following chapter. We will employ the derived system matrices in two and three dimensions. While smaller examples serve as an introduction, more complex systems are analyzed to demonstrate the simulation capabilities. For selected examples, a discussion of the solution convergence is added. While we present the nodal solutions for small examples, the solutions for large examples are only visualized. Additionally, all chosen stiffness and force values are purely academic. Forces are given in NewtonN, distances and lengths in meter m, areas in m2, and densities in kg/m3. In the respective two-dimensional examples, the third direction, here the global z-axis, is fixed on all system nodes to prevent any movement in said direction.This section includes the following examples:(i)
Section4.1 demonstrates a simple 2D 3-bar system. The nonlinear solution is compared to the linear formulation, and convergence rates are included.(ii)
Section4.2 discusses a 2D framework bridge. The nonlinear solution is compared to the linear formulation.(iii)
Section4.3 analyses a 3D 3-bar system. The nonlinear solution is compared to the linear formulation.(iv)
Section4.4 demonstrates a 3D frame. The nonlinear solution is compared to the linear formulation.(v)
Section4.5 discusses a 3D framework bridge.(vi)
Section4.6 shows a 3D 3-bar system, including self-weight. The nonlinear solution is compared to the linear formulation, and convergence rates are included.(vii)
Section4.7 discusses the influence of different material laws on a 3D system.(viii)
Section4.8 demonstrates a large structural setup with an increased number of degrees of freedom.To verify the found reference configurations, a standard FEM analysis is performed on the found reference configuration. The results coincide with the desired final configuration and are visualized in most of the example discussions in this work. These standard FEM analyses have been performed with the code available in the GitHub repository [28] and verified with the open-source multiphysics software KRATOS [32–34].All examples are accessible via a GitHub repository [28]. The interested reader is invited to download the offered files, re-run the examples herein, and analyze new structures. A tutorial is included, describing the proper usage.


### 4.1. 2D 3-Bar System

The first two-dimensional system describes a 3-bar structure as visualized in Figure5. The Young’s modulus E is chosen to be1e6N/m2, and the cross section isA=1e–2m2, while the force vector on the free middle note is 4e2,−8e2N. The St. Venant-Kirchhoff material law is used for the nonlinear truss formulation, and the standard, small deformation Hooke’s law is for the linear truss formulation.Figure 5
Unknown and known configuration of a 3-bar structure. The red arrows demonstrate the force vector.In Figure5, the gray line presents the known deformed configuration, while the solid black and blue lines are the results of the inverse analyses herein. The different resulting reference configurations for the linear and the nonlinear truss formulation can be observed. While the linear truss formulation calculates the force equilibrium in the reference configuration and is restricted to small deformations, the nonlinear truss formulation formulates its force equilibrium in the deformed configuration and can represent large deformations. Additionally, the dotted yellow line describes the result of a standard FEM analysis on the found, previously unknown reference configurations. It can be observed that the solution of the standard FEM analysis coincides with the given desired aim configuration of the inverse analysis.The found reference configuration of the middle free node is given in the following Table1.Table 1
Solution of the initially unknown and known configuration of the 3-bar structure in Figure5.
XZxzLinear1.00939197743751.18585025374761.0−1.0Nonlinear1.02734126841671.33049320124031.0−1.0Additionally, Figure6 presents the convergence of the problem residual with respect to the Newton-Raphson iterations k. The L − 2 norm •2 is chosen to depict a convergence rate that is at least quadratic.Figure 6
Convergence of the 3-bar system.


### 4.2. 2D Framework Bridge

A truss framework bridge is modeled using multiple structural elements in two dimensions. The structure in Figure7 has a Young’s modulus E=3e5N/m2, a cross sectionA=1e–2m2, and two vertical forces in y-direction of magnitude 1e2N. The forces are intentionally chosen to be arbitrary and nonsymmetric to test the simulation capabilities. The St. Venant-Kirchhoff material law is used for the nonlinear truss formulation, and the standard, small deformation Hooke’s law is for the linear truss formulation. Figure 7 nicely shows the solution for the unknown reference configuration and the result of a proceeding standard FEM analysis. The found reference shapes deform nicely into the desired aim configuration. The difference between the linear and nonlinear truss formulation again can be observed. Figure 8 presents the course of the residual over iterations steps k. The gradient for rk2<1e−1 demonstrates the correct implementation.Figure 7
Unknown and known configuration of a two-dimensional bridge structure. The red arrows demonstrate the force vector.Figure 8
Convergence of the 2D framework bridge system.


### 4.3. 3D 3-Bar System

To begin with 3D systems, a 3-bar structure is investigated. The structure has a Young’s modulusE=1e6N/m2, a cross sectionA=1e–2m2, and a force vector of 1e3,0,2e3N is applied to the free node. The St. Venant-Kirchhoff material law is used for the nonlinear truss formulation, and the standard, small deformation Hooke’s law is for the linear truss formulation. Solutions for the unknown reference configuration are visualized in Figure 9 and given in Tables 2 and 3.Figure 9
Unknown and known configuration of a 3-bar structure. The red arrow demonstrates the force vector1e3,0,2e3N.Table 2
Solution of the initially unknown configuration and aim configuration of the 3-bar structure in Figure9 with linear truss elements.
XYZ0.4775519165819412−0.480257269602899331.4935108461479483xyz1.00.02.0Table 3
Solution of the initially unknown configuration and aim configuration of the 3-bar structure in Figure9 with nonlinear truss elements.
XYZ0.6553972454362691−0.4246380164345481.735432507307086xyz1.00.02.0The found reference configuration of the middle free node is given in the following.Again, Figure10 presents the convergence for both the linear and the nonlinear truss elements.Figure 10
Convergence of the 3D 3-bar system.


### 4.4. 3D Frame

Figure11 depicts the next three-dimensional setup. The frame is supposed to deflect to a regular rectangular shape, as demonstrated by the gray dotted lines. The system uses a Young’s modulus E=1e6N/m2 and a cross sectionA=1e–2m2, and the nodal force vector is 300,300,300N. The St. Venant-Kirchhoff material law is used for the nonlinear truss formulation, and the standard, small deformation Hooke’s law is for the linear truss formulation.Figure 11
Unknown and known configuration of a three-dimensional frame structure. The red arrow demonstrates the force vector300,300,300N.The respective convergence rates are given in Figure12.Figure 12
Convergence of the 3D frame system.


### 4.5. 3D Framework Bridge

Subsequently, we present the example of a 3D framework bridge structure, as shown in Figure13. While a single force loads the structure nonsymmetric with 0,−800,0N, the system’s Young’s modulus is E=1e7N/m2, and the elements are assigned a cross sectionA=1e–2m2. The St. Venant-Kirchhoff material law is used for the nonlinear truss formulation, and the standard, small deformation Hooke’s law is for the linear truss formulation. As the three dimensionality makes a proper visualization difficult, the found reference configuration is compared to the desired aim configuration in multiple spatial projections. Figure 13 depicts a three-dimensional side view.Figure 13
Unknown and known configuration of a 3D bridge structure. The red arrow demonstrates the force vector0,−800,0N.Figures14 and 15 show a projection in the x−z plane for different levels of z. The vertical point load of −800 N in y-direction is visualized in Figure 15.Figure 14
Bridge configurations inx−y plane for z=0.Figure 15
Bridge configurations inx−y plane for z=1. The red arrow describes the vertical point load.In the following Figures16 and 17, the projection on the x−z plane is shown for different y levels.Figure 16
Bridge configurations inx−z plane for y=0.Figure 17
Bridge configurations inx−z plane for y=1.The convergence rate for the 3D bridge system is visualized in Figure18. While the linear truss formulation needs several iterations to reach rk2<1e−1, it thereafter shows proper convergence.Figure 18
Convergence of the 3D bridge system.


### 4.6. Self-Weight

In the following, the inclusion of the self-weight is discussed with respect to the problem solving procedure. The same 3-bar system as in Section 4.3 is used withE=1e6N/m2, A=1e–2m2, g=0,1,0, and ρ=1e3kg/m3. The St. Venant-Kirchhoff material law is used for the nonlinear truss formulation, and the standard, small deformation Hooke’s law is for the linear truss formulation. The theoretical background and the respective symbols are discussed in Section 3.2.3. Figure 19 visualizes the problem setup and the found, initially unknown reference configuration in the following.Figure 19
Unknown and known configuration of a three-dimensional truss structure, subjected to self-weight. The red arrow indicates the normalized gravity direction vectorg, acting in the y-direction.After solving for the unknown reference configuration, the following discrete nodal position is found. Table4 presents the result for the linear truss formulation, while Table 5 relates to the nonlinear formulation.Table 4
Solution of the initially unknown configuration and aim configuration of the 3-bar structure in Figure19 with linear truss elements.
XYZ0.874387790220−0.26275037784601.9263479975443xyz1.00.02.0Table 5
Solution of the initially unknown configuration and aim configuration of the 3-bar structure in Figure19 with nonlinear truss elements.
XYZ0.8431479069299−0.34373540177141.9298317925520xyz1.00.02.0While the inclusion of the self-weight leads to additional forces, it influences the setup of the solution matrixCe, too. The necessary changes to Ce are discussed in Section 3.2.3, especially in equation (39). While, in some cases, it is also possible to find a solution for the self-weight problem employing the original matrix from equations (27) or (32), a lower convergence rate is expected. This behavior is visualized in Figure 20, where the same problem is solved twice. Ones, using the original matrix from equation (30) for linear trusses and equation (29) for nonlinear trusses (here named “incomplete,” due to its lack of necessary terms) and again with the complete matrix from equation (41). The latter equation includes all derivatives of the new residual vector, containing the self-weight loads dependent on the unknown reference configuration. The large differences in convergence rates can be seen in Figure 20, proving the correct derivation and implementation of the necessary additional derivative terms.Figure 20
Convergence of the nonlinear self-weight problem, employing either the complete solution matrixCe from equation 57 or the incomplete matrix from equation 32, respectively, 43. Both converge to the approximately same solution but show different convergence rates.


### 4.7. Different Hyperelastic Material Laws

To discuss the influence of the three different hyperelastic material laws introduced in Section 2.1, a 2D truss system is investigated, as depicted in Figure21. The Young’s modulus is always E=1e6N/m2, and the cross section isA=1e–2m2. While the St. Venant-Kirchhoff law and the Neo-Hookean law do not need further material parameters, β1=10,β2=0 is assigned to the Ogden law. Only the nonlinear truss formulation can be applied in conjunction with the nonlinear hyperelastic material laws. Figure 21 shows the different reference configurations for each of the material law. All three cases yield the desired aim configuration when simulated with a standard FEM analysis.Figure 21
Unknown and known configuration of a 2D nonlinear truss structure. The three different hyperelastic material laws are compared: St. Venant-Kirchhoff, Neo-Hookean, Ogden.Figure22 visualizes the convergence rate for each of the simulations. At least quadratic convergence is observed for all three cases.Figure 22
Convergence of the 2D nonlinear truss structure for different hyperelastic material laws.The solution for the unknown reference configuration yields the following y-coordinate for the middle node, given in Table6.Table 6
Unknown reference configuration and desired aim configuration for various hyperelastic material laws.
St. Venant-KirchhoffNeo-HookeanOgdenβ1=10,β2=0Y2.02.02.0Y1.31134341281.16952007371.4055031762


### 4.8. TUM Logo—Large Number of Degrees of Freedom

Finally, to conclude, we present an example with a larger number of dofs. The logo of the Technical University of Munich (TUM) is taken as an inspiration for the following example. As depicted in Figure23, the whole structure is discretized with nonlinear truss elements, using the St. Venant-Kirchhoff material law. The Young’s modulus is chosen to be E=2.069e11N/m2, and the cross section of each truss isA=1e–2m2. Supports are introduced on the left lower side as shown in Figure 23, and the right lower side is loaded with a y-load of 5e5N on 11 nodes. As all displacements in the global z-direction are fixed, the total number of 3241 nodes yields ≈7e3 dofs.Figure 23
Known aim configuration. The whole model consists of single nonlinear truss elements, as depicted in the blue box. As all displacements out-of-plane are constrained, a total number of 3241 nodes yield≈7e3 dofs.Figure24 finally visualizes the found reference configuration and the desired aim configuration, which coincides with the deformed reference structure, obtained by a standard FEM analysis on the found reference configuration.Figure 24
Unknown reference configuration and desired aim configuration for the nonlinear truss TUM logo. Boundary conditions are given in Figure23.


## 4.1. 2D 3-Bar System

The first two-dimensional system describes a 3-bar structure as visualized in Figure5. The Young’s modulus E is chosen to be1e6N/m2, and the cross section isA=1e–2m2, while the force vector on the free middle note is 4e2,−8e2N. The St. Venant-Kirchhoff material law is used for the nonlinear truss formulation, and the standard, small deformation Hooke’s law is for the linear truss formulation.Figure 5
Unknown and known configuration of a 3-bar structure. The red arrows demonstrate the force vector.In Figure5, the gray line presents the known deformed configuration, while the solid black and blue lines are the results of the inverse analyses herein. The different resulting reference configurations for the linear and the nonlinear truss formulation can be observed. While the linear truss formulation calculates the force equilibrium in the reference configuration and is restricted to small deformations, the nonlinear truss formulation formulates its force equilibrium in the deformed configuration and can represent large deformations. Additionally, the dotted yellow line describes the result of a standard FEM analysis on the found, previously unknown reference configurations. It can be observed that the solution of the standard FEM analysis coincides with the given desired aim configuration of the inverse analysis.The found reference configuration of the middle free node is given in the following Table1.Table 1
Solution of the initially unknown and known configuration of the 3-bar structure in Figure5.
XZxzLinear1.00939197743751.18585025374761.0−1.0Nonlinear1.02734126841671.33049320124031.0−1.0Additionally, Figure6 presents the convergence of the problem residual with respect to the Newton-Raphson iterations k. The L − 2 norm •2 is chosen to depict a convergence rate that is at least quadratic.Figure 6
Convergence of the 3-bar system.


## 4.2. 2D Framework Bridge

A truss framework bridge is modeled using multiple structural elements in two dimensions. The structure in Figure7 has a Young’s modulus E=3e5N/m2, a cross sectionA=1e–2m2, and two vertical forces in y-direction of magnitude 1e2N. The forces are intentionally chosen to be arbitrary and nonsymmetric to test the simulation capabilities. The St. Venant-Kirchhoff material law is used for the nonlinear truss formulation, and the standard, small deformation Hooke’s law is for the linear truss formulation. Figure 7 nicely shows the solution for the unknown reference configuration and the result of a proceeding standard FEM analysis. The found reference shapes deform nicely into the desired aim configuration. The difference between the linear and nonlinear truss formulation again can be observed. Figure 8 presents the course of the residual over iterations steps k. The gradient for rk2<1e−1 demonstrates the correct implementation.Figure 7
Unknown and known configuration of a two-dimensional bridge structure. The red arrows demonstrate the force vector.Figure 8
Convergence of the 2D framework bridge system.


## 4.3. 3D 3-Bar System

To begin with 3D systems, a 3-bar structure is investigated. The structure has a Young’s modulusE=1e6N/m2, a cross sectionA=1e–2m2, and a force vector of 1e3,0,2e3N is applied to the free node. The St. Venant-Kirchhoff material law is used for the nonlinear truss formulation, and the standard, small deformation Hooke’s law is for the linear truss formulation. Solutions for the unknown reference configuration are visualized in Figure 9 and given in Tables 2 and 3.Figure 9
Unknown and known configuration of a 3-bar structure. The red arrow demonstrates the force vector1e3,0,2e3N.Table 2
Solution of the initially unknown configuration and aim configuration of the 3-bar structure in Figure9 with linear truss elements.
XYZ0.4775519165819412−0.480257269602899331.4935108461479483xyz1.00.02.0Table 3
Solution of the initially unknown configuration and aim configuration of the 3-bar structure in Figure9 with nonlinear truss elements.
XYZ0.6553972454362691−0.4246380164345481.735432507307086xyz1.00.02.0The found reference configuration of the middle free node is given in the following.Again, Figure10 presents the convergence for both the linear and the nonlinear truss elements.Figure 10
Convergence of the 3D 3-bar system.


## 4.4. 3D Frame

Figure11 depicts the next three-dimensional setup. The frame is supposed to deflect to a regular rectangular shape, as demonstrated by the gray dotted lines. The system uses a Young’s modulus E=1e6N/m2 and a cross sectionA=1e–2m2, and the nodal force vector is 300,300,300N. The St. Venant-Kirchhoff material law is used for the nonlinear truss formulation, and the standard, small deformation Hooke’s law is for the linear truss formulation.Figure 11
Unknown and known configuration of a three-dimensional frame structure. The red arrow demonstrates the force vector300,300,300N.The respective convergence rates are given in Figure12.Figure 12
Convergence of the 3D frame system.


## 4.5. 3D Framework Bridge

Subsequently, we present the example of a 3D framework bridge structure, as shown in Figure13. While a single force loads the structure nonsymmetric with 0,−800,0N, the system’s Young’s modulus is E=1e7N/m2, and the elements are assigned a cross sectionA=1e–2m2. The St. Venant-Kirchhoff material law is used for the nonlinear truss formulation, and the standard, small deformation Hooke’s law is for the linear truss formulation. As the three dimensionality makes a proper visualization difficult, the found reference configuration is compared to the desired aim configuration in multiple spatial projections. Figure 13 depicts a three-dimensional side view.Figure 13
Unknown and known configuration of a 3D bridge structure. The red arrow demonstrates the force vector0,−800,0N.Figures14 and 15 show a projection in the x−z plane for different levels of z. The vertical point load of −800 N in y-direction is visualized in Figure 15.Figure 14
Bridge configurations inx−y plane for z=0.Figure 15
Bridge configurations inx−y plane for z=1. The red arrow describes the vertical point load.In the following Figures16 and 17, the projection on the x−z plane is shown for different y levels.Figure 16
Bridge configurations inx−z plane for y=0.Figure 17
Bridge configurations inx−z plane for y=1.The convergence rate for the 3D bridge system is visualized in Figure18. While the linear truss formulation needs several iterations to reach rk2<1e−1, it thereafter shows proper convergence.Figure 18
Convergence of the 3D bridge system.


## 4.6. Self-Weight

In the following, the inclusion of the self-weight is discussed with respect to the problem solving procedure. The same 3-bar system as in Section 4.3 is used withE=1e6N/m2, A=1e–2m2, g=0,1,0, and ρ=1e3kg/m3. The St. Venant-Kirchhoff material law is used for the nonlinear truss formulation, and the standard, small deformation Hooke’s law is for the linear truss formulation. The theoretical background and the respective symbols are discussed in Section 3.2.3. Figure 19 visualizes the problem setup and the found, initially unknown reference configuration in the following.Figure 19
Unknown and known configuration of a three-dimensional truss structure, subjected to self-weight. The red arrow indicates the normalized gravity direction vectorg, acting in the y-direction.After solving for the unknown reference configuration, the following discrete nodal position is found. Table4 presents the result for the linear truss formulation, while Table 5 relates to the nonlinear formulation.Table 4
Solution of the initially unknown configuration and aim configuration of the 3-bar structure in Figure19 with linear truss elements.
XYZ0.874387790220−0.26275037784601.9263479975443xyz1.00.02.0Table 5
Solution of the initially unknown configuration and aim configuration of the 3-bar structure in Figure19 with nonlinear truss elements.
XYZ0.8431479069299−0.34373540177141.9298317925520xyz1.00.02.0While the inclusion of the self-weight leads to additional forces, it influences the setup of the solution matrixCe, too. The necessary changes to Ce are discussed in Section 3.2.3, especially in equation (39). While, in some cases, it is also possible to find a solution for the self-weight problem employing the original matrix from equations (27) or (32), a lower convergence rate is expected. This behavior is visualized in Figure 20, where the same problem is solved twice. Ones, using the original matrix from equation (30) for linear trusses and equation (29) for nonlinear trusses (here named “incomplete,” due to its lack of necessary terms) and again with the complete matrix from equation (41). The latter equation includes all derivatives of the new residual vector, containing the self-weight loads dependent on the unknown reference configuration. The large differences in convergence rates can be seen in Figure 20, proving the correct derivation and implementation of the necessary additional derivative terms.Figure 20
Convergence of the nonlinear self-weight problem, employing either the complete solution matrixCe from equation 57 or the incomplete matrix from equation 32, respectively, 43. Both converge to the approximately same solution but show different convergence rates.


## 4.7. Different Hyperelastic Material Laws

To discuss the influence of the three different hyperelastic material laws introduced in Section 2.1, a 2D truss system is investigated, as depicted in Figure21. The Young’s modulus is always E=1e6N/m2, and the cross section isA=1e–2m2. While the St. Venant-Kirchhoff law and the Neo-Hookean law do not need further material parameters, β1=10,β2=0 is assigned to the Ogden law. Only the nonlinear truss formulation can be applied in conjunction with the nonlinear hyperelastic material laws. Figure 21 shows the different reference configurations for each of the material law. All three cases yield the desired aim configuration when simulated with a standard FEM analysis.Figure 21
Unknown and known configuration of a 2D nonlinear truss structure. The three different hyperelastic material laws are compared: St. Venant-Kirchhoff, Neo-Hookean, Ogden.Figure22 visualizes the convergence rate for each of the simulations. At least quadratic convergence is observed for all three cases.Figure 22
Convergence of the 2D nonlinear truss structure for different hyperelastic material laws.The solution for the unknown reference configuration yields the following y-coordinate for the middle node, given in Table6.Table 6
Unknown reference configuration and desired aim configuration for various hyperelastic material laws.
St. Venant-KirchhoffNeo-HookeanOgdenβ1=10,β2=0Y2.02.02.0Y1.31134341281.16952007371.4055031762


## 4.8. TUM Logo—Large Number of Degrees of Freedom

Finally, to conclude, we present an example with a larger number of dofs. The logo of the Technical University of Munich (TUM) is taken as an inspiration for the following example. As depicted in Figure23, the whole structure is discretized with nonlinear truss elements, using the St. Venant-Kirchhoff material law. The Young’s modulus is chosen to be E=2.069e11N/m2, and the cross section of each truss isA=1e–2m2. Supports are introduced on the left lower side as shown in Figure 23, and the right lower side is loaded with a y-load of 5e5N on 11 nodes. As all displacements in the global z-direction are fixed, the total number of 3241 nodes yields ≈7e3 dofs.Figure 23
Known aim configuration. The whole model consists of single nonlinear truss elements, as depicted in the blue box. As all displacements out-of-plane are constrained, a total number of 3241 nodes yield≈7e3 dofs.Figure24 finally visualizes the found reference configuration and the desired aim configuration, which coincides with the deformed reference structure, obtained by a standard FEM analysis on the found reference configuration.Figure 24
Unknown reference configuration and desired aim configuration for the nonlinear truss TUM logo. Boundary conditions are given in Figure23.


## 5. Conclusion

In contrast to standard static analyses of structures, searching for a suitable reference geometry is not as straightforward. A reformulation of the structural force equilibrium equation is necessary. We are changing the unknowns from the structural displacements to the unknown reference configuration. This change of unknowns yields a different solution matrix, compared to the stiffness matrix of standard static analyses. Including load cases, which depend on the reference configuration (such as self-weight), adds additional nonlinear terms to the problem statement. This work presents the described problem reformulation in terms of unknowns for geometrically nonlinear truss elements combined with nonlinear hyperelastic material laws. Additionally, simplifications to linear truss formulations are given, offering a wide range of detailed derivations. To solve the newly obtained nonlinear system of equations, we propose the gradient-based Newton-Raphson solution scheme. Consequently, the preceding sections discuss and present arising system derivatives in great detail. Every part of each necessary solution matrix is derived and consequently given. Additionally, we published everything in a freely accessible GitHub repository [28]. While all system matrices are available, each presented example is added to the repository and thus freely available.The derivation and implementation are checked by various two- and three-dimensional test cases, ranging from small systems to larger, more complex structures. The solution scheme herein, including all system matrices, can now be used to analyze workshop geometries for truss structures.



---

*Source: 3573608-2022-12-01.xml*