Deformation of planar ring mesh under static in-plane tensile loads

A ring mesh is a large-scale manufacturable structure with versatile applications in architecture and for protective systems. However, the static and dynamic numerical simulation of a large-scale ring mesh is a resource-intensive task due to the many nonlinear contact points between the individual rings. To characterize the rigid body behavior of the ring mesh, in this paper, a representative volume element is loaded under different in-plane directions. The Green-Lagrangian strain tensor components are obtained as a result. An implicit one-step algorithm is used for this dynamic relaxation issue, modified by the use of stiff springs. Static convergence positions of a representative volume element are determined in several multiaxial tensile directions. The obtained deformation parameters can be used to simulate large deformations of large-scaled ring meshes.

ring meshes, for example, six others. The number and the location of contact points between rings depend on the load on the mesh. Lenk studied the elastic and plastic behavior of prestressed ring mesh. He established a linear orthotropic material law for a representative volume element of the ring mesh, varying the diameter ratio of the tori [11].

Numerical simulation
In the case of a nonpretensioned metallic ring mesh, its deformation results mainly from the rigid body motion of the individual rings. In a rigid body simulation, all rings and contacts can be modeled and evaluated individually. This type of simulation is used for 3D printing of a ring mesh in order to position the rings in the printing volume [12]. To characterize this behavior of the ring mesh in 2D or 3D, there are currently several approaches. Kock et al. discretize the rings by their centers in the ring mesh loaded under tension and find a diamond-shaped pattern. Depending on the load, the diagonals of the virtual diamonds change. This technique has been used to simulate overthrows of ring meshes as facades over buildings [4]. For the reconstruction of historical chain mail, Wijnhoven and Moskvin have adopted it to CAD. They described the static parameters of the ring mesh under two directions of load in the plane and bending under dead weight [13]. With these data, clothing simulations can be performed with linear elastic shells or discrete mass-spring models [14]. A characterization of the rigid body behavior of interlocking tori under in-plane loads in connection with constraints was carried out by Gingras via an optimization of a position potential. Here, the distortion of the ring mesh was fixed in one direction. Then he variated the distortion in the second direction until the potential was minimized within the constraints [15].

Representative volume elements
In addition to Lenk's material law for in-plane prestressed ring meshes [11], the representative volume element (RVE) is also used in the description of other materials like composites [16]. For example, material laws for periodic masonry walls [17,18] are also investigated using the representative volume element. Here, a virtual piece is cut out of the material, which, periodically continued, results in the entire mesh again. There is no need that the surface of the RVE matches the surface of single bodies [19][20][21]. After measuring the deformation of the RVE while loading, a material law can be obtained [22]. With this material law the macroscopic behavior of a larger structure can be represented without modeling all details inside the RVE. If the local conditions like local stresses and local strains are required, numerical simulations can be performed on multiple levels. Inside of the structure, borders can be set where the macroscopic stresses and strains are used as boundary constraints for local models [23][24][25].

Research gap
The currently performed investigations of the rigid body behavior of the ring mesh without simulation of the individual rings only consider certain load directions well. In order to describe the rigid body behavior of a ring mesh under further in-plane static load directions, sections of ring mesh with periodic boundary conditions are simulated. Then, based on the relationship between load and deformation, a material law for the homogenization of the ring mesh is sought.
The paper is structured as follows: In Sect. 2, the structure of the loaded ring mesh is described. The abstraction, the boundary conditions, and the loading are discussed, as well as the fundamentals of the numerical simulation of the rigid body are given. In Sect. 3, the investigated parameters are presented, and the results for the Green-Lagrangian strain tensor components are summarized. Finally, the evaluation of the results and an outlook on further research follow in Sect. 4.

Fundamentals and modeling
In this Section, the fundamentals of the investigation of the rigid body behavior of a ring mesh under uniform static in-plane loading are described. A section of the ring mesh is considered, which is subjected to periodic boundary conditions. In the following, the ring mesh will be homogenized and abstracted as a two-dimensional disk. The choice of a disk is based on the following assumptions: • The thickness of a ring mesh is significantly less than its width and length.
• The loads in thickness direction as well as the bending moments are negligible.
The individual rings, considered as rigid bodies, each have their own body-fixed coordinate system. The connection between the rigid rings and the corresponding abstracted disk can only be performed at discrete points, chosen in our case are the centers of the tori. The reference element in the present work is bordered by the center points of four numbered rings, see Fig. 1, which represent a rectangle in the stretched initial position. Note that the mass within the reference element remains constant. Under uniform loading, both the tori and the reference elements remain in a periodic pattern even when the ring mesh deforms and twists.
The displacements of the gate center points correspond to the displacements of the braid at these points. If the displacement u of the braid between these discrete points is interpolated, finite elements are obtained [26]. It should be mentioned that the choice of the initial position u = 0 is arbitrary, and the strains are thus dependent on the initial position. A finite element defined by the tori center points will be restricted to have a homogenous constant deformation. The deformation gradient is defined, see Fig. 2. Here, e x , e y , and e z are the unity vectors of the global coordinate system, b x and b y are the dimensions of the undeformed representative volume element, and x 2 and x 3 are the positions of the tori center points for body 2 and 3. It should be noted that there is no deformation in direction of e z .
With the deformation gradient F several strain tensors can be obtained. Determining the components of a strain tensor is the goal of this paper. The Green-Lagrangian strain tensor was chosen because it also describes large deformations [26,27]. Here, I is the identity matrix, and E 11 , E 12 , and E 22 are three independent components of the Green-Lagrangian strain tensor. Calculating the components of the strain tensor has the advantage of eliminating any rotations in the deformation gradient. It should be noted that only in-plane loads and in-plane deformations are obtained; therefore, all components of E based on the z-direction are equal to 0. Since the connection between the abstracted mesh and the rings occurs at the centers of the tori, the load is also applied at these centers in the rigid body model. In a two-dimensional disk, the three independent line loads s x x , s xy = s yx and s yy can be applied. For the rigid body movement, only the ratio and orientation of the main loads s 1 and s 2 are important for determining the end position. For the model, the load parameters can be reduced to s 1 =ŝ cos β and s 2 =ŝ sin β, whereŝ is the maximum load. The main loads can be related to the load components by The main load angle α as well as the stress ratio angle β = arctan2(s 2 , s 1 ) are depicted in Fig. 3. The crux is the transfer of the line loads into the discrete forces on the centers of the tori. The respective line load is assigned to the nearest center point and integrated, see Fig. 4. The forces per line load are finally summed up. Note that due to the linear projection of the line loads onto the forces at the centers, the superposition principle can be used. Please note that the forces at the ring centers depend on their current position.

Periodic constraints
A ring mesh is positioned in the coordinate system so that the interlacing direction [11] overlays the y-direction, see Fig. 5.
There are two types of periodic constraints, one based on the tori center points x i and the other one referenced on the rotation parameters p i of the rings. Together, the position of the center points x i and the rotation parameters p i build the location z i of the tori i. Rings, each with the same orientation, form a pattern of parallelograms defined by equal distance vectors, for example Fig. 6b are all parallel to each other, as well as the two blue vectors.
The equality in rotation builds the second type of periodic constraints. For the rotation parameters in Fig. 6, the constraint equations apply.
Finally, there are multiple supports required to limit the possible displacement of the ring mesh. Movement alongside the z-direction is suppressed by fixing all center points in this direction. The translation of x 1 is suppressed along the x-and y-axis, see Fig. 6a. The ring mesh should not realign itself via rotation. With a defined angle γ , the rotation can be constrained by the relation

Rigid body simulation
In a frictionless system with convex boundary conditions, the equilibrium position can be reached from any starting point. This case is assumed for a ring mesh loaded under axial or multiaxial tension. To determine the distortions of the ring mesh, the load is applied at an initial position until the ring mesh comes to rest in the equilibrium position in the damped multibody system. A linear implicit, impulse-based algorithm by Stewart and Trinkle [28,29] A stiffness matrix K = diag(k) including each contact stiffness in the rows of k is introduced to avoid that the velocity for the separation of interpenetrating tori increases proportionally when the time step h decreases. This is because some penetration is allowed to remain in an elastic contact. At least, there are the external loads f e and the vector of centrifugal forces and gyroscopic torques f g . The static equilibrium is reached when v → 0, see Eq. (6). Thus, the convergence of the system can be represented by the magnitude of the velocity or of the kinetic energy. Another way to define the convergence of the system is by dynamic relaxation [30][31][32]. In these works, a residual force is defined by the weighted norm of the static equilibrium relative to the weighted norm of the external loads. It has to be mentioned that in Eq. (9) the reference value is always fulfilled in this context. The unilateral contacts between the tori are determined via iteratively determined contact points. The maximum quantity of contact points between two rings is two, since the tori cannot be stacked. The two rings are abstracted by their center circles, and the potential contact point is centered between points on these center circles. A calculation of a contact point is done by an alternating projection algorithm [33]. In the present work the tori are reduced to their center circles. Starting from a point on the first circle, the second point results from an orthogonal projection of the first point on the second circle. This procedure is repeated alternately until the points on both circles converge. On two geometrical identic tori the potential contact point is in the middle between those two points on the two center circles, see Fig. 7.
The nonlinear influence of contacts and rotation can have a negative influence on the stability of the rigid body simulation. Therefore, in each step, the energy of the system E is calculated. Since the rigid-body simulation is dissipative in Stewart's implicit form [28], the step size h is reduced until the energy of the total system actually decreases. The individual steps for determining the convergence positions of in-plane loaded representative volume elements were described at the beginning of this Section. For a summary of the procedure, see Algorithm 1. , f g , H, B, G b , α b , N, G n , α n solve z l+1 , v l+1 , λ b , λ n update E l+1 until E l+1 ≤ E l end while Fig. 8 Geometry of a torus, center circle with radius r 1 , cross section radius r 2 , and ring surface

Results
Under multiaxial static tensile loading, a ring mesh with periodic constraints is loaded. The characteristic value for the chain mail mesh is the ratio of the two radii in the torus r 1 r 2 , see Fig. 8. The radii r 1 = 5 mm of the center circle and r 2 = 0.5 mm of the cross section are set, but the results can be transferred for all ring meshes with the radius ratio r 1 r 2 = 10.

Setup
The section of the ring mesh consists of ten rings in the x y-plane. For the simulation, the density ρ = 7.85 · 10 −6 kg/mm 3 is also specified, as well as a constant damping of D = diag(1 kg/s, 1 kg/s, 1 kg/s, 0.1 kgm/s, 0.1 kgm/s, 0.1 kgm/s) related to the individual rings to achieve the convergence position. The contact stiffness is chosen to be very high with k = 10 6 kg/s 2 , so that the rings can be considered almost rigid in the convergence position. The maximum time step is h = 1 s. The simulation step is iterated with a reduced time step, when the energy of the system in step l + 1 is higher than in step l. In each simulation, as many steps were performed until the dimensionless weighted load R (Eq. 9) fell below the value of 10 −9 .
To ensure that the main load directions do not become negative, a test series with the load angle α in the range from 0 to 45 • and β from 0 to 90 • was performed. The step size for both angles is 5 • , see Fig. 9.   In order to represent the contact points in a meaningful way, in Figs. 10, 11, and 12 the tori are only indicated by their center lines. The maximum distortion in the x-direction and at the same time the highest compression in the y-direction occur under axial tensile load in the x-direction, see Fig. 11.
Pure shear s xy = 0, s x x = s yy = 0 is not investigated because it has a negative main load component s 2 < 0. The load case with maximum shear load occurs under a uniaxial tensile load at α = 45 • , see Fig. 12. The summary of the results is shown in Fig. 13. From the final positions the components, ig. 13 Obtained Green-Lagrangian strain tensor components of the ring mesh depending on load direction, a E 11 versus β for different α, b E 12 versus β for different α and c E 22 versus β for different α of the Green-Lagrangian strain tensor can be calculated [26].
When changing the load direction, the distortions behave smoothly as long as the same rings are in contact. When new contacts are activated or rings are released, the curves of the distortions are discontinuous. In the investigated load range, the calculated Green-Lagrangian strain tensor components in a three-dimensional plot are located on a convex plane, see Fig. 14.

Convergence study
In the convergence study, two aspects of the evaluation of the load cases are investigated. First, the simulation time, and second, the path to convergence for two exemplary load cases are discussed.
For each simulation of the 190 load cases, two characteristic values were documented, which provide information about the duration. First there is the number of steps n and second there are the cumulated time steps n l=1 h l , where h l is the time step size of the step number l. The number of steps n indicates how many time steps were calculated until convergence was reached.
In Fig. 15, the number of time steps n is depicted versus the load case number. It can be seen that most of the load cases require between 150 and 300 time steps. There also exist several downward outliers with only There are two possible reasons for the large difference between the minimum and maximum values of n and n l=1 h l . Collisions between noninterlinking rings slow down the system abruptly and shorten the path to convergence. An increase in simulation time can occur due to a larger displacement, especially when the forces The deviations from the linear trend in Fig. 16 for high values of n show that the step size h was reduced to lower values for several load cases. The reason for this may be a high nonlinearity in the occurring contacts.
Finally, two exemplary load cases are investigated in more detail. These are the load case number 1 with α = β = 0 • with a small number of steps and the load case number 190 with α = 45 • and β = 90 • with many steps. Figure 17 shows the time evolution of the components of the strain tensor and of the convergence parameters for the two load cases. In these two uniaxial load cases the system converges when noninterlinking rings get into contact.
In both cases, the components of the Green-Lagrangian strain tensor converge to their final values. The two convergence parameters R and E kin behave similarly and go toward zero. The reason for this is a high damping in the system. It can be seen that the convergence parameters slowly decrease, followed by a larger decrease at approx. 30 s (Fig. 17a, c) and at approx. 1200 s (Fig. 17b, d). This kink in the curves results from noninterlinking rings colliding and abruptly slowing the system to convergence.
It has to be mentioned that for times up to approx. 400 s oscillations in R and E kin occur (Fig. 17d) having no significant impact on the strain tensor components. The system in Fig. 17c has a higher E kin than in Fig. 17d and converges within a much shorter time.
It can be seen that, in the method described in the present paper, the simulation time can vary drastically for different load cases.

Discussion
With the described method, convergence positions in the rigid body model can be achieved, and subsequently two-dimensional distortions of the four-in-one linked ring mesh can be calculated. The advantage of the periodic constraints is the reduction in the number of rings required in this method compared to a homogeneous loaded large-scaled ring mesh simulation.
The calculated distortions can be applied to all ring meshes with the same ratio r 1 /r 2 . For other ratios the procedure can be performed analogically. Only tensile loads were investigated because otherwise more rings get into contact. Then the reference volume element requires more rings, and other effects can negatively influence the results. These are, for example, buckling or an ambiguous convergence position. In addition, no bending of the braid is considered with these simulations.
The novelty of the present approach compared to other works available in the literature is that Green-Lagrangian strain tensor components can be simply obtained for the different load directions.
The present method gives some extensions and advantages compared to other approaches. However, for determining the strain parameters there are some challenges left: first, there are sometimes high simulation times for specific load cases. Second, the system parameters, i.e., the damping and the convergence criterion, have to be chosen carefully.
The distortions can be used in future work for macroscopic deformation simulation of a ring mesh for static and dynamic loading. The continuum described by the macroscopic distortions can be simulated by using numerical methods. For example, discrete elements or finite elements can be defined. The obtained distortions under load can serve as the starting point for both, elastic and plastic calculations.
In future research, the obtained rigid body behavior will be incorporated into a material law. The goal is to homogenize the planar ring mesh in static and dynamic simulations with no need to simulate all tori individually for large-scale ring meshes.
Funding Open Access funding enabled and organized by Projekt DEAL.

Declarations
Conflict of interest The authors declare that there is no conflict of interest.
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/.