Arc-length method

---

Illustration of arc-length method to trace complex curve

In numerical analysis and computational mechanics, the arc-length method, also known as the Riks method, is a path-following algorithm used to trace the equilibrium path of nonlinear systems of equations. This method is particularly useful for problems exhibiting limit points, such as snap-through or snap-back, where the standard Newton's method under load control may fail. The arc-length method overcomes this difficulty by introducing an additional unknown, namely the load parameter λ, together with an arc-length constraint equation, thereby augmenting the original system of n nonlinear equations to a system of n+1 nonlinear equations. Therefore, the method is widely used in nonlinear finite element analysis, structural stability, and other problems in which the response of a nonlinear system cannot be represented by a single-valued function of the applied load.

Basic theory

Consider a nonlinear system with n displacement variables and one load parameter λ, can be cast in the form$$\mathbf{𝐑}(\mathbf{𝐮},\lambda )={\mathbf{𝐑}}^{\mathrm{i}\mathrm{n}\mathrm{t}}(\mathbf{𝐮})-\lambda {\mathbf{𝐑}}^{\mathrm{e}\mathrm{x}\mathrm{t}}={\mathbf{𝐑}}^{\mathrm{i}\mathrm{n}\mathrm{t}}(\mathbf{𝐮})-\lambda \mathbf{𝐅}=\mathbf{𝟎},$$where R(u, λ) is the residual vector, and (u, λ) represents an equilibrium state; u is the displacement vector including n displacement variables, λ is the load parameter, and R^int and R^ext are the vectors of internal and external forces; R^ext can be completely determined by the applied force and therefore is represented by F. Note that both $\mathbf{𝐮}$ and $\lambda$ are treated as unknowns in the arc-length method.

In this case, for the nonlinear system R(u, λ) with n displacement degrees of freedom, the arc-length method increases the number of unknowns by one, namely n+1 variables. The extra variable is the load parameter λ, a proportional operator to represent the applied load level. But under such situation, the above equation is not adequate to be resolved, as it is a system of n equations but with n+1 unknowns. Thus, the supplementary equation, referred to as the arc-length constraint equation, must be included$$\mathrm{\Delta }{\mathbf{𝐮}}^{\mathrm{T}}\mathrm{\Delta }\mathbf{𝐮}+\mathrm{\Phi }\mathrm{\Delta }{\lambda }^2\left({\mathbf{𝐅}}^{\mathrm{T}}\mathbf{𝐅}\right)=\mathrm{\Delta }{s}^2,$$where (Δu, Δλ) are displacement and load increments within a step, s is the arc-length parameter which parametrizes the equilibrium path , and Δs is the arc-length parameter increment, namely the step length (also the arc radius), which generally remains constant during the step; Φ is a user-defined value, where Φ = 0 represents the cylindrical constraint surfaces and Φ = 1 denotes the spherical constraint surfaces. There is another form of the arc-length constraint equation which removes the influence of applied load level$$\mathrm{\Delta }{\mathbf{𝐮}}^{\mathrm{T}}\mathrm{\Delta }\mathbf{𝐮}+\mathrm{\Phi }\mathrm{\Delta }{\lambda }^2=\mathrm{\Delta }{s}^2,$$where both forms are acceptable.

Geometrically, the arc-length constraint equation defines a hypersurface in the augmented displacement–load space, instead of the pure displacement space considered in the Newton's method. As a result, the solution of the arc-length method at each increment is obtained as the intersection of the equilibrium path and the hypersurface, yielding both the displacement vector and the corresponding load factor at the same time. As the arc-length method allows the load parameter to increase or decrease during the computation, unstable parts of an equilibrium path can be followed. Therefore, the arc-length method is able to advance along the equilibrium curve, rather than along a coordinate direction associated only with load or displacement.

Incremental-iterative procedure

This is the illustration of incremental-iterative procedure of arc-length method using Quasi-Newton's method

Arc-Length Method Code Demonstration for Tracing Equilibrium Curves with Two Limit Points

Starting from a known equilibrium point $({\mathbf{𝐮}}_p,{\lambda }_p)$, the predictor stage, which is defined as the first iteration of a step, estimates the next iteration point on the tangential direction at the known equilibrium point. Then, the corrector stage, which is defined as the following iterations except the first, continuously provides necessary corrections to the iteration point $({\mathbf{𝐮}}_p^i,{\lambda }_p^i)$, until the pre-defined convergence tolerance is achieved. A complete step consists of both the predictor stage and the corrector stage, where Newton'method is often used inside each iteration to solve the augmented nonlinear system including both the residual vector and the arc-length constraint equation. Here, one can choose either full Newton's method or Quasi-Newton's method according to different considerations and requirements.

Although the original arc-length method was independently developed by Wempner and Riks in the 1970s, the modern incremental-iterative procedures were subsequently developed by Crisfield into the forms widely used today in the commercial finite element software. By re-formatting the augmented nonlinear system, Crisfield contributed to the standardization of the original arc-length method for coding, making it possible to be implemented into the standard finite element method.

By using the Taylor series expansion, the iteration point of the $(i+1)$-th iteration of the p-th step $({\mathbf{𝐮}}_p^{i+1},{\lambda }_p^{i+1})$ can be expressed in terms of the last known iteration point $({\mathbf{𝐮}}_p^i,{\lambda }_p^i)$$$\mathbf{𝐑}({\mathbf{𝐮}}_p^{i+1},{\lambda }_p^{i+1})\approx {\mathbf{𝐑}}^{\mathrm{i}\mathrm{n}\mathrm{t}}\left({\mathbf{𝐮}}_p^i\right)+{\left[{\mathbf{𝐊}}_{\mathrm{T}}\right]}_{{\mathbf{𝐮}}_p^i}\delta {\mathbf{𝐮}}_p^{i+1}-({\lambda }_p^i+\delta {\lambda }_p^{i+1})\mathbf{𝐅}=\mathbf{𝟎},$$where ${\left[{\mathbf{𝐊}}_{\mathrm{T}}\right]}_{{\mathbf{𝐮}}_p^i}$ denotes the tangent stiffness matrix evaluated at ${\mathbf{𝐮}}_p^i$, and $\delta {\mathbf{𝐮}}_p^{i+1}$ and $\delta {\lambda }_p^{i+1}$ are the displacement and load corrections to be determined at this iteration.

According to Crisfield's suggestion, the above equation is re-written as$$\delta {\mathbf{𝐮}}_p^{i+1}=-{\left[{\mathbf{𝐊}}_{\mathrm{T}}\right]}_{{\mathbf{𝐮}}_p^i}^{-1}[{\mathbf{𝐑}}^{\mathrm{i}\mathrm{n}\mathrm{t}}\left({\mathbf{𝐮}}_p^i\right)-{\lambda }_p^i\mathbf{𝐅}]+\delta {\lambda }_p^{i+1}{\left[{\mathbf{𝐊}}_{\mathrm{T}}\right]}_{{\mathbf{𝐮}}_p^i}^{-1}\mathbf{𝐅}$$$$=\delta {\mathbf{𝐮}}_R^i+\delta {\lambda }_p^{i+1}\delta {\mathbf{𝐮}}_F^i$$with$$\delta {\mathbf{𝐮}}_R^i=-{\left[{\mathbf{𝐊}}_{\mathrm{T}}\right]}_{{\mathbf{𝐮}}_p^i}^{-1}[{\mathbf{𝐑}}^{\mathrm{i}\mathrm{n}\mathrm{t}}\left({\mathbf{𝐮}}_p^i\right)-{\lambda }_p^i\mathbf{𝐅}]=-{\left[{\mathbf{𝐊}}_{\mathrm{T}}\right]}_{{\mathbf{𝐮}}_p^i}^{-1}{\mathbf{𝐑}}_p^i$$$$\delta {\mathbf{𝐮}}_F^i={\left[{\mathbf{𝐊}}_{\mathrm{T}}\right]}_{{\mathbf{𝐮}}_p^i}^{-1}\mathbf{𝐅},$$where $\delta {\mathbf{𝐮}}_R^i$ is the residual displacement increment used to correct the displacement discrepancy induced by the unbalanced forces, $\delta {\mathbf{𝐮}}_F^i$ is the reference load displacement increment which defines the tangential direction of the equilibrium path in the displacement space , and ${\mathbf{𝐑}}_p^i$ is the residual vector at the $i$-th iteration of the $p$-th load step. All of these quantities can be computed from known information.

After that, by substituting the above equation into the arc-length constraint equation, the quadratic equation with $\delta {\lambda }_p^{i+1}$ as the only variable can be derived as$$a(\delta {\lambda }_p^{i+1}{)}^2+2b(\delta {\lambda }_p^{i+1})+c=0,$$where$$a={\left[\delta {\mathbf{𝐮}}_F^i\right]}^{\mathrm{T}}\delta {\mathbf{𝐮}}_F^i+\mathrm{\Phi }{\mathbf{𝐅}}^{\mathrm{T}}\mathbf{𝐅}$$$$b={\left[{\mathbf{𝐯}}^i\right]}^{\mathrm{T}}\delta {\mathbf{𝐮}}_F^i+\mathrm{\Phi }{\lambda }_p^i{\mathbf{𝐅}}^{\mathrm{T}}\mathbf{𝐅}$$$$c={\left[{\mathbf{𝐯}}^i\right]}^{\mathrm{T}}{\mathbf{𝐯}}^i+\mathrm{\Phi }({\lambda }_p^i{)}^2{\mathbf{𝐅}}^{\mathrm{T}}\mathbf{𝐅}-\mathrm{\Delta }{s}^2$$with$${\mathbf{𝐯}}^i=\mathrm{\Delta }{\mathbf{𝐮}}_p^i+\delta {\mathbf{𝐮}}_R^i.$$The roots of the quadratic equation can be solved using the root-finding formula$$\delta {\lambda }_p^{i+1}=\frac{-b\pm \sqrt{b^2-ac}}{a},\text{when }{b}^2-ac>0.$$This equation will produce two roots, the correct one with forward search direction can be selected out based on multiple approach. After that, by substituting the $\delta {\lambda }_p^{i+1}$ back, the displacement and load corrections for this iterations are determined.;

The process is repeated for successive steps, producing a sequence of equilibrium points$$({\mathbf{𝐮}}_1,{\lambda }_1),({\mathbf{𝐮}}_2,{\lambda }_2),({\mathbf{𝐮}}_3,{\lambda }_3),\dots$$on the nonlinear equilibrium path.

Limitations

Although the arc-length method significantly improves the ability to trace nonlinear equilibrium paths beyond limit points, it does not eliminate all numerical difficulties. Its performance still depends on the nonlinear characteristics of the problem, the predictor-corrector strategy, and the numerical implementation. Several issues remain important in practical applications:

• Initialization: At the initial step, no previous information is available. Therefore, the initial step length and the first predictor direction must be specified by the user.
• Step length increment: The step length affects both convergence and computational efficiency. Excessively large increments may cause divergence, whereas too small increments increase the computational cost.

• Root selection: The arc-length constraint generally produces two roots. Selecting the correct root corresponding to the forward equilibrium path is essential to avoid tracing back or convergence issue.

• Bifurcation: The arc-length method does not automatically switch to secondary equilibrium branches at bifurcation points. Additional branch-switching techniques are generally required.

• Convergence: As with other Newton-type methods, convergence depends on the tangent stiffness matrix and the predictor. Severe geometric or material nonlinearities, contact, and softening behavior may reduce convergence robustness.

Variants and extensions

Numerous variants and extensions of the arc-length method have been developed to improve its robustness, convergence, and computational efficiency. Most of these developments are intended to address the numerical difficulties encountered in practical nonlinear analysis, including root selection, adaptive step-length control, and branch switching at bifurcation points.

Various root-selection criteria have been proposed for arc-length method with arc-length constraint equation to ensure that the solution proceeds along the desired equilibrium path. Representative approaches include criteria based on the sign of the tangent stiffness determinant, incremental work, stiffness parameters,secant path methods.

Adaptive adjustment of the step length has also been extensively investigated. A common strategy is to multiply the previously converged step length by a scaling factor, which is usually determined according to the number of iterations required in the previous increment. If the iteration count exceeds a prescribed limit, the step length is typically reduced, often by one-half, to improve convergence.

To trace multiple equilibrium branches, the arc-length method is frequently combined with branch-switching techniques based on eigenvalue analysis or perturbation methods. These approaches enable the solution to transition from the primary equilibrium path to secondary branches at bifurcation points, thereby extending the applicability of the method to more complicated post-buckling problems.

Recent developments have further incorporated some advanced techniques to deal with problems involving severe geometric and material nonlinearities, contact, damage, and fracture. Despite the diversity of modern implementations, most arc-length methods share the same underlying continuation principle and differ primarily in the formulation of the constraint equation, predictor strategy, and step length control.

Applications

The arc-length method is used in a wide range of nonlinear mechanical problems, including:

• post-buckling analysis of plates and shells;
• snap-through analysis of shallow arches and shell structures;
• nonlinear analysis with material non-linearity
• fracture and damage simulations in which the load-displacement response contains snapping behavior.

In many commercial finite element software, arc-length-type procedures are often provided under names such as Static Riks, arc-length method or arc-length control.

See also

• Newton's method
• Nonlinear finite element analysis
• Structural stability
• Numerical continuation
• Post-buckling analysis

External links

• Arc-length method implementation on MATLAB File Exchange
• Non-Linear Finite Element Analysis of Solids and Structures

This article "Arc-length method" is from Wikipedia. The list of its authors can be seen in its historical and/or the page Edithistory:Arc-length method. Articles copied from Draft Namespace on Wikipedia could be seen on the Draft Namespace of Wikipedia and not main one.

---