Stretched grid method

This article needs additional citations for verification. Please help improve this article by adding citations to reliable sources. Unsourced material may be challenged and removed. Find sources: "Stretched grid method" – news · newspapers · books · scholar · JSTOR (June 2017) (Learn how and when to remove this message)

The stretched grid method (SGM) is a numerical technique for finding approximate solutions of various mathematical and engineering problems that can be related to an elastic grid behavior. In particular, meteorologists use the stretched grid method for weather prediction^[1] and engineers use the stretched grid method to design tents and other tensile structures.

In recent decades^[when?], the finite element and boundary element methods (FEM and BEM) have become mainstays for industrial engineering design and analysis. Increasingly larger and more complex designs are being simulated using the FEM or BEM. However, some problems of FEM and BEM engineering analysis are still on the cutting edge. The first problem is that the reliability of the analysis strongly depends upon the quality of initial data generated at the pre-processing stage. It is known that automatic element mesh-generation techniques at this stage have become commonly used tools for the analysis of complex real-world models.^[2] With FEM and BEM increasing in popularity comes the incentive to improve automatic meshing algorithms. However, all of these algorithms can create distorted and even unusable grid elements. Several techniques exist that can take an existing mesh and improve its quality. For instance smoothing (also referred to as mesh refinement) is one such method, which repositions nodal locations, so as to minimize element distortion. The stretched grid method allows one to obtain pseudo-regular meshes very easily and quickly in a one-step solution.^[3]

Assume that there is an arbitrary triangle grid embedded into a plane polygonal single-coherent contour and produced by an automeshing procedure (see fig. 1). It may be assumed further that the grid considered as a physical nodal system is distorted by a number of distortions. It is supposed that the total potential energy of this system is proportional to the length of some n {\displaystyle n} -dimensional vector with all network segments as its components.

Thus, the potential energy takes the form

Π = D ∑ j = 1 n R j 2 , {\displaystyle \Pi =D\sum _{j=1}^{n}{R_{j}}^{2},}

where n {\displaystyle n} is total number of segments in the network, R j {\displaystyle R_{j}} is the length of segment number j {\displaystyle j} , and D {\displaystyle D} is an arbitrary constant.

The length of segment number j {\displaystyle j} may be expressed by two nodal coordinates as

R = ( X 12 − X 11 ) 2 + ( X 22 − X 21 ) 2 {\displaystyle R={\sqrt {(X_{12}-X_{11})^{2}+(X_{22}-X_{21})^{2}}}}

It may also be supposed that coordinate vector { X } {\displaystyle \{X\}} of all nodes is associated with the non-distorted network, and the coordinate vector { X ′ } {\displaystyle \{X'\}} is associated with the distorted network. Vector { X } {\displaystyle \{X\}} may the be expressed as

{ X } = { X ′ } + { Δ X } {\displaystyle \{X\}=\{X'\}+\{\Delta X\}}

The vector { X } {\displaystyle \{X\}} determination is related to minimization of the quadratic form Π {\displaystyle \Pi } by incremental vector { Δ   X } {\displaystyle \{\Delta \ X\}} ; that is,

∂ Π ∂ Δ X k l = 0 , {\displaystyle {\frac {\partial \Pi }{\partial \Delta X_{kl}}}=0,}

where l {\displaystyle l} is the number of interior node of the area, and k {\displaystyle k} is the number of coordinate^{[clarification needed]}

After all transformations, we may write the following two independent systems of linear algebraic equations:

[ A ] { Δ X 1 } = { B 1 } {\displaystyle [A]\{\Delta X_{1}\}=\{B_{1}\}}

[ A ] { Δ X 2 } = { B 2 } , {\displaystyle [A]\{\Delta X_{2}\}=\{B_{2}\},}

where [ A ] {\displaystyle [A]} is a symmetric matrix in the banded form similar to the global stiffness matrix of a FEM assemblage, { Δ X 1 } {\displaystyle \{\Delta X_{1}\}} and { Δ X 2 } {\displaystyle \{\Delta X_{2}\}} are incremental vectors of coordinates of all nodes at axes 1 and 2, and { B 1 } {\displaystyle \{B_{1}\}} and { B 2 } {\displaystyle \{B_{2}\}} are the right part vectors that are combined by coordinates of all nodes in axes 1 and 2.

The solution of both systems, keeping all boundary nodes conservative, obtains new interior node positions corresponding to a non-distorted mesh with pseudo-regular elements. For example, Fig. 2 presents the rectangular area covered by a triangular mesh. The initial auto mesh possesses some degenerative triangles (left mesh). The final mesh (right mesh) produced by the SGM procedure is pseudo-regular without any distorted elements.

As the above systems are linear, the procedure reduces very quickly to a one-step solution. Moreover, each final interior node position meets the requirement of coordinate arithmetic mean of nodes surrounding it^{[clarification needed]} and meets the Delaunay criteria too. Therefore, the SGM has all the positive values peculiar to Laplacian smoothing and other kinds of smoothing approaches but is much easier and more reliable because of integer-valued final matrix representations. Finally, the SGM is applicable not only to 2D meshes but to 3D meshes consisting of any uniform cells as well as to mixed or transient meshes.

Mathematically, a surface bounded by a closed curve is called minimal if its area is minimal amongst all the surfaces with the same boundary. The best-known minimal surface sample is a soap film bounded by a wire frame. Usually to create a minimum surface, a fictitious constitutive law, which maintains a constant prestress, independent of any changes in strain, is used.^[4] An alternative approximate approach is based on SGM.

The idea is to approximate a surface by a triangular grid. To converge such a grid to a grid with minimum area, one should solve the same two systems described above. Increments of the third nodal coordinates may be determined by a similar system at axis 3:

[ A ] { Δ X 3 } = { B 3 } {\displaystyle [A]\{\Delta X_{3}\}=\{B_{3}\}}

By solving all three systems simultaneously, one can obtain a new grid that will approximate desired minimal surface.

As an example, a catenoid calculated by the above approach is presented in Fig 3. The radii of the bounding rings and the height of the catenoid are equal to 1. The surface area of the catenoid determined by SGM is 2.9967189; the area of the ideal catenoid is 2.992.

For structural analysis, the configuration of the structure is generally known beforehand. This is not the case for tensile structures such as tension fabric structures. Since the membrane in a tension structure possesses no flexural stiffness, its form or configuration depends upon initial prestressing and the loads that it is subjected to. Thus, the load-bearing behaviour and the shape of the membrane cannot be separated and cannot be generally described by simple geometric models only. The membrane shape, the loads on the structure, and the internal stresses interact in a nonlinear manner to satisfy the equilibrium equations.

The preliminary design of tension structures involves the determination of an initial configuration, referred to as form finding. In addition to satisfying the equilibrium conditions, the initial configuration must accommodate both architectural (aesthetic) and structural (strength and stability) requirements. Further, the requirements of space and clearance should be met, the membrane principal stresses must be tensile to avoid wrinkling, and the radii of the double-curved surface should be small enough to resist out-of-plane loads and to insure structural stability.^[5] Several variations on form-finding approaches based on FEM have been developed to assist engineers in the design of tension fabric structures. All of them are based on the same assumption as that used for analysing the behaviour of tension structures under various loads. However, it might sometimes be preferable to use minimal surfaces in the design of tension structures.

The physical meaning of SGM consists in the convergence of the energy of an arbitrary grid structure with specified constraints to a grid that minimizes the sum of the distances between arbitrary pairs of grid nodes. It allows the minimum surface energy problem to be solved. The generalized formulation of SGM presupposes the ability to apply a set of external forces and rigid or elastic constraints to grid structure nodes, in order to model various external effects. We may obtain the following expression for such a formulation:

Π = ∑ j = 1 n D j R j 2 + ∑ i = 1 3 ( ∑ k = 1 m C i k Δ X i k 2 − ∑ k = 1 m P i k Δ X i k ) {\displaystyle \Pi =\sum _{j=1}^{n}D_{j}R_{j}^{2}+\sum _{i=1}^{3}\left(\sum _{k=1}^{m}C_{ik}\Delta X_{ik}^{2}-\sum _{k=1}^{m}P_{ik}\Delta X_{ik}\right)}

where

• n {\displaystyle n} is the total number of grid segments,
• m {\displaystyle m} is the total number of nodes,
• R j {\displaystyle R_{j}} is the length of segment j {\displaystyle j} ,
• D j {\displaystyle D_{j}} is the stiffness of segment j {\displaystyle j} ,
• Δ X i k {\displaystyle \Delta X_{ik}} is the coordinate increment of node k {\displaystyle k} at axis i {\displaystyle i} ,
• C i k {\displaystyle C_{ik}} is the stiffness of an elastic constraint in node k {\displaystyle k} at axis i {\displaystyle i} , and
• P i k {\displaystyle P_{ik}} is the external force on node k {\displaystyle k} at axis i {\displaystyle i} .

Once a satisfactory shape has been found, a cutting pattern may be generated. Tension structures are highly varied in their size, curvature and material stiffness. Cutting pattern approximation is strongly related to each of these factors. It is essential for a cutting-pattern generation method to minimize the approximation error and to produce reliable plane cloth data.

The objective is to develop the shapes described by these data, as close as possible to the ideal doubly curved strips. In general, cutting-pattern generation involves two steps. First, the global surface of a tension structure is divided into individual cloths. The corresponding cutting pattern at the second step can be found by simply taking each cloth strip and unfolding it onto a planar area. In the case of the ideal doubly curved membrane surface, the subsurface cannot be simply unfolded and must be flattened. For example, in ^[6][7], SGM has been used for the flattening problem solution.

The cutting-pattern generation problem is actually subdivided into two independent formulations. These are the generation of a distortion-free plane form unfolding each cloth strip and flattening double-curved surfaces that cannot be simply unfolded. From the position of differential geometry, both formulations are the same. We may consider it as an isometric mapping of a surface onto the plane area that will be conformal and equiareal simultaneously because of invariant angles between any curves and invariance of any pieces of area. In the case of a single-curved surface that can be unfolded precisely, equiareal mapping allows one to obtain a cutting pattern for fabric structure without any distortions. The second type of surfaces can be equiareally mapped only approximately, with some distortions of linear surface elements limited by the fabric properties. Assume that two surfaces are parameterized so that their first quadratic forms may be written as follows:

I 1 = E 1 ( u , v ) d ⁡ u 2 + 2 F 1 ( u , v ) d ⁡ u d ⁡ v + G 1 ( u , v ) d ⁡ v 2 {\displaystyle I_{1}=E_{1}(u,v)\operatorname {d} u^{2}+2F_{1}(u,v)\operatorname {d} u\operatorname {d} v+G_{1}(u,v)\operatorname {d} v^{2}}

I 2 = E 2 ( u , v ) d ⁡ u 2 + 2 F 2 ( u , v ) d ⁡ u d ⁡ v + G 2 ( u , v ) d ⁡ v 2 {\displaystyle I_{2}=E_{2}(u,v)\operatorname {d} u^{2}+2F_{2}(u,v)\operatorname {d} u\operatorname {d} v+G_{2}(u,v)\operatorname {d} v^{2}}

The condition of conformal mapping for two surfaces requires that

I 2 = λ I 1 , {\displaystyle {\sqrt {I_{2}}}=\lambda {\sqrt {I_{1}}},}

where λ {\displaystyle \lambda } is the ratio of the surface distortion due to conformal mapping.

It is known that the first quadratic form reflects the distance between two surface points ( u , v ) {\displaystyle (u,v)} and ( u + d ⁡ u , v + d ⁡ v ) {\displaystyle (u+\operatorname {d} u,v+\operatorname {d} v)} . When λ {\displaystyle \lambda } is close to 1, the above equation converges to the conditions of isometric mapping and equiareal mapping respectively because of invariant angles between any curves and invariance of any pieces of area. Remembering that the first stage of form-finding is based on a triangular mesh, and using the method of weighted residuals for the description of isometric and equiareal mapping of the minimum surface onto a plane area, we may write the following function, which is defined by the sum of integrals along segments of curved triangles:

Π = D ∑ j = 1 n ∮ S j w j ( λ I 1 − I 2 ) 2 d ⁡ s , {\displaystyle \Pi =D\sum _{j=1}^{n}\oint _{S_{j}}w_{j}\left(\lambda {\sqrt {I_{1}}}-{\sqrt {I_{2}}}\right)^{2}\operatorname {d} s,}

where

• n {\displaystyle n} is the total number of grid cells,
• w j {\displaystyle w_{j}} are the weight ratios,
• Π {\displaystyle \Pi } is the total mapping residual, and
• D {\displaystyle D} is a constant that does not influence the final result and may be used as a scale ratio.

Considering further weight ratios w j = 1 {\displaystyle w_{j}=1} , we may transform the above equation into an approximate finite sum that is a combination of linear distances between nodes of the surface grid and write the basic condition of equiareal surface mapping as a minimum of following nonlinear function:

Π = D ∑ j = 1 n ∮ S j w j ( λ R j − L j ) 2 d ⁡ s , {\displaystyle \Pi =D\sum _{j=1}^{n}\oint _{S_{j}}w_{j}\left(\lambda R_{j}-L_{j}\right)^{2}\operatorname {d} s,}

where

• R j {\displaystyle R_{j}} is the initial length of linear segment j {\displaystyle j} ,
• L j {\displaystyle L_{j}} is the final length of segment number j {\displaystyle j} , and
• λ {\displaystyle \lambda } is a distortion ratio close to 1 and may be different for each segment.

The initial and final lengths of segment number j {\displaystyle j} may be expressed as usual by two nodal coordinates as

R = ( X 12 − X 11 ) 2 + ( X 22 − X 21 ) 2 + ( X 32 − X 31 ) 2 {\displaystyle R={\sqrt {(X_{12}-X_{11})^{2}+(X_{22}-X_{21})^{2}+(X_{32}-X_{31})^{2}}}}

L = ( x 12 − x 11 ) 2 + ( x 22 − x 21 ) 2 , {\displaystyle L={\sqrt {(x_{12}-x_{11})^{2}+(x_{22}-x_{21})^{2}}},}

where

• X i k {\displaystyle X_{ik}} are coordinates of nodes of the initial segment, and
• x i k {\displaystyle x_{ik}} are coordinates of nodes of the final segment.

According to the initial assumption, we can write x 32 = x 31 = 0 {\displaystyle x_{32}=x_{31}=0} for the plane surface mapping. The expression for vectors { x } {\displaystyle \{x\}} and { X } {\displaystyle \{X\}} may be written as

{ x } = { X } + { Δ X } . {\displaystyle \{x\}=\{X\}+\{\Delta X\}.}

The definition of vector { Δ X } {\displaystyle \{\Delta X\}} is as previously

∂ Π ∂ Δ X k l = 0. {\displaystyle {\frac {\partial \Pi }{\partial \Delta X_{kl}}}=0.}

After transformations, we may write the following two independent systems of nonlinear algebraic equations:

[ A ] { Δ X 1 } = { B 1 } + { Δ P 1 } {\displaystyle [A]\{\Delta X_{1}\}=\{B_{1}\}+\{\Delta P_{1}\}}

[ A ] { Δ X 2 } = { B 2 } + { Δ P 2 } , {\displaystyle [A]\{\Delta X_{2}\}=\{B_{2}\}+\{\Delta P_{2}\},}

where all the parts of the system can be expressed as previously and { Δ P 1 } {\displaystyle \{\Delta P_{1}\}} and { Δ P 2 } {\displaystyle \{\Delta P_{2}\}} are vectors of pseudo-stresses at axes 1 and 2 that have the form

{ Δ P l t } = − { ∑ j = 1 N λ R m L m ( x l m − x l t ) } , {\displaystyle \{\Delta P_{lt}\}=-\left\{\sum _{j=1}^{N}\lambda {\frac {R_{m}}{L_{m}}}(x_{lm}-x_{lt})\right\},}

where

• N {\displaystyle N} is the total number of nodes that surround node number t {\displaystyle t} , and
• l {\displaystyle l} is the number of global axes.

The above approach is another form of SGM and allows one to obtain two independent systems of nonlinear algebraic equations that can be solved by any standard iteration procedure. The less the Gaussian curvature of the surface, the higher the accuracy of the plane mapping. As a rule, the plane mapping allows one to obtain a pattern with linear dimensions 1–2% less than corresponding spatial lines of a final surface. That is why it is necessary to provide the appropriate margins while patterning.

The typical sample of cutout—also called a gore or a patch—is presented in Figs. 9, 10, and 11.

• Mesh generation
• Types of mesh