Indirect Displacement Control
56 Direct displacement control can be applied only on structures loaded only at one point, or when the load is transmitted by a stiff platen so that all points on the loaded surface exhibit the same displacements.
57 However, this is not always the case. As an example, consider a dam loaded by hydrostatic pressure due to reservoir overflow; see Fig. 16.9. Here, the load is applied along a large portion of
the boundary, and the shape of the corresponding displacement profile is not known in advance. Another case in which direct displacement control fails is very brittle failure characterized by a loaddisplacement diagram with a snapback, Fig. 16.10.
Al  
M 
M 
M  
/ / 
/ 
/ Al/  
f 
1 
V Al  
u 
u 
u  
Load Control 
Displacement Control 
ArcLength Control 
Figure 16.10: LoadDisplacement Diagrams with Snapback
58 Advanced incrementation control techniques abandon the assumption that the values of external loads and/or displacements at supports after each incremental step are prescribed in advance. Instead, the loading program is parameterized by a scalar load multiplier.
16.4.1 Partitioning of the Displacement Corrections Adapted from (Reich 1993)
59 Restricting the applied loading to be proportional, a scalar load parameter 3 can be used to scale an arbitrary set of applied loads f .
60 The applied loads at the start of increment i are defined as the scalarvector product. (3i fe'Ti, where 3i is the load parameter at the start of increment i. 3i is zero at the start of the first increment.
61 The applied incremental loads for increment i are defined as the scalarvector product A pit, where A/% is the incremental load parameter for increment i. The updated load parameter 3i+1 at the end of increment i is
62 The incremental displacements due to the applied incremental loads are obtained using the standard modifiedNewton algorithm, as described in Zienkiewicz & Taylor (1991). The incremental displacements AuJ+1 at the end of iteration j for a. generic increment, are defined as
where AuJ are the incremental displacements at the start, of iteration j and are the incremental displacement corrections for iteration j.
63 The incremental load parameter A3j+1 at the end of iteration j is defined in an analogous manner as
where A3j is the incremental load parameter at the start of iteration j and 53j is the incremental load parameter correction for iteration j. At. the start, of the first, iteration AuJ and A/J7 are identically zero.
64 Incremental displacement corrections are determined by solving
K duj = (f3fext + Af3jfext + 5f3jTxt  Tmt) (16.48)
where K is the global stiffness matrix and
e=1 Qe are the reactions for the state of stress at the start of iteration j.
65 Defining the residual forces TZ3 at the start, of iteration j as
Equation 16.48 can be written more simply as
ee The matrixvector product. K1 fext is invariant, for the increment, and, therefore, can be treated as a. vector constant. 5ut, which Crisfield (1981) referred to as the tangent, displacements
6 The matrixvector product. K 1 TZ3 defines the displacement, corrections dui due to the residual forces
but they are obviously not invariant for the increment. The displacement corrections for iteration j are then defined as
68 Figure 16.11 shows a flowchart for an incremental nonlinear finite element program based on the modifiedNewton algorithm with indirect displacement control capabilities. The numbers in the boxes in Figure 16.11 correspond to those appearing in Figure ??.
16.4.2 ArcLength
Adapted from (Jirasek and Bazant 2001)
69 The basic idea of a flexible incrementation control technique is that the step size is specified by a constraint equation that involves the unknown displacements as well as the load multiplier. The original motivation was provided by the requirement that the size of the step measured as the geometric distance between the initial and final state in the loaddisplacement space should be equal to a prescribed constant, Fig. 16.10.
70 Despite the apparent simplicity of the condition of a constant arc length, it must be used with caution. First of all, it is important to realize that forces and displacements have completely different units, and so the purely geometrical measure of length in the loaddisplacement space does not make a good sense. It is necessary to introduce at least one scaling factor, denoted as c, that multiplies the load parameter and converts it into a quantity with the physical dimension of displacement. The length of a step during which the load parameter changes by A^ and the displacements change by Au is then defined as
71 By adjusting the scaling factor we can amplify or suppress the relative contribution of loads and displacements. One reasonable choice is derived from the condition that, the contributions should be equal as long as the response remains linear elastic, which leads to c = \/uJue where ue is the solution of Keue = f.
72 In some cases, e.g., for frame, plate, and shell models that use both translational and rotational degrees of freedom, the components of the generalized displacement, vector u do not. have the same physical dimension. It. is then necessary to apply scaling also to the vector Au.
73 Consider an incremental solution process controled by the arclength method. In a typical step number n, we start, with displacements and load parameter ^n~^ computed in the previous step, and we search for displacements u^ and load parameter ¡jSn\ The state at the end of the step must, satisfy the equations of equilibrium between the internal forces fmt(u^) and external forces fextit^) Compared to the load control or direct, displacement, control, the load parameter is an additional unknown. The corresponding additional equation is provided by the constraint that fixes the size of the step. For example, we can require that the length of the step evaluated from formula (16.55) be equal to a prescribed value, Al. We could treat the problem as a system of N/ + 1 nonlinear equations, where N/ is the number
Figure 16.12: Two points on the loaddisplacement curve satisfying the arclength constraint of unknown displacement components (degrees of freedom), and solve it by NewtonRaphson iteration. However, a more elegant and computationally more efficient procedure treats the equilibrium equations and the constraint equation to a certain extent separately. Assume for simplicity that the loading program is described by (??). The linearized equations of equilibrium in the ¿th iteration read
K(n'I1)iu(n'i) = fo + /x(ra'i_1)f  f£r1} + (1656)
where ¿u^™'^ is the unknown displacement, correction, and is the unknown correction of the load parameter. The first three terms on the righthand side are known, and the last term is an unknown scalar multiple of a. given vector f. We can therefore separately solve equations
and then express the displacement correction as
When this expression is substituted into the constraint condition,
(An(n,il) + 6{nA)^T(A(nAl) + ^(M) + ^(A//'^"1) + fi^Mf = (A I)2 (16.59)
we obtain a quadratic equation for a single unknown, This equation usually has two real roots, corresponding to the two points of the equilibrium path that have the prescribed distance from point. , ^n~see Fig. 16.12. The correct, root, is selected depending on the sense in which we march on the equilibrium path (Crisfield, M.A. 1981), and the displacement correction is determined from (16.58). After standard updates of the displacement vector and the load parameter, the iteration cycle is repeated until the convergence criteria are satisfied.
16.4.3 Relative Displacement Criterion Adapted from (Reich 1993)
74 The standard arclength control performs well if the entire structure or its large portion participates in the failure mechanism. In cases when the failure pattern is highly localized, robustness of the technique may deteriorate. The remedy is to adapt the constraint equation to the particular problem and control the incrementation process by a few carefully selected displacement components.
75 Motivation is again provided by the physical background. If the loaddisplacement diagram of a brittle structure exhibits snapback, direct displacement control applied in an experiment leads to sudden catastrophic failure. When the displacement imposed by the loading device reaches a critical value, fracture starts propagating even though the imposed displacement at the load point is kept fixed. However, opening of the crack monotonically increases during the entire failure process, and so it can be used as a control variable. If the experimental setup is arranged such that the applied force is continuously adjusted depending on the currently measured value of the crack opening, the response can be traced in a stable manner even after the point at which the loaddisplacement diagram snaps back. The same idea can be exploited by a numerical simulation. It suffices to select a suitable linear combination of displacement components that increases monotonically during the entire failure process, and to use this combination as the control variable.
76 de Borst (1985,1986) concluded that arclength methods (Riks 1979, Ramm 1981, Crisfield 1981), which were the original IDC methods, were not satisfactory for analyses involving cracking accompanied by softening.
77 The main problem with the arclength methods, when used in this context, was that the constraint involved all displacement components equally when, in fact, only a few displacement components were dominant. The dominant displacement components were typically those for nodes at or near the crack mouth. This being the case, de Borst proposed using a transformed relative displacement component between two nodes as the constraint.
78 The transformed relative displacement component can define the crack mouth opening displacement (CMOD), crack mouth sliding displacement (CMSD), or some arbitrary displacement Au between two points on a structure. The arbitrary displacement Au may correspond to a relative displacement measured during an experiment such as the relative vertical displacement between a point on the neutral axis of a 3point bend beam over a support and the bottom of the beam at midspan.
79 As it is the most general case, the relative displacement criterion will be described in terms of the arbitrary relative displacement Au. A pair of nodes, m and n, are selected to define Au, with their total displacements being (u)m and (u)„, respectively. The direction associated with Au is defined by a unit vector v. Au is thereby defined as
80 If m and n are nodes on opposite sides of a discrete crack Au = CMOD if v is normal to the crack surface and Au = CMSD if v is tangent to the crack surface. The value for Au is prescribed for an increment and the applied loads are scaled such that the total displacements at the end of each iteration reflect that value.
si Recalling that the total displacements uJ+1 at the end of iteration j are defined as ni+l = ü:/ + ößjöüT +
the load parameter correction 53j for iteration j is
All,  VT [(u7')«  (uJ')m]  VT {6uJr)n  {5uJr)m
16.4.4 IDC Methods with Approximate Line Searches
82 Employing a procedure proposed by Crisfield (1983) for use with the arclength method, the convergence of the solution alogrithm can be accelerated by performing approximate line searches; approximate line searches under fixed (i.e. nonscalable) loads are described in Section
83 This procedure requires an extra iterative loop at the beginning of the line search loop in which a combination of 83j and sk+1 satisifying the constraint conditions (i.e. Equations ?? and 16.60) is computed. As 83j is initially computed for s1 = 1.0, any change in sk requires a corresponding change in 83j for the IDC constraint to remain satisfied. Consequently, an iterative loop, in which 83j is recomputed based on the estimated value of sk+1, is required to obtain a compatible combination of 83j and sk+1.
84 After recomputing 83j, the values of g0 and gk are also recomputed using Equation 16.30 to reflect the change in the residual loads 7caused by the new value of 8ft. f^1'' is not updated to reflect the changes in sk+1 when recomputing 7, which is strictly not correct, but it does significantly reduce the number of computations without causing any difficulties (Crisfield 1983).
85 Finally, from the new values of g0 and gk, sk+1 is reestimated using Equation 16.33. The loop is terminated when sk+1 sk+1 I new I
which generally requires only a few iterations. A flow chart of this procedure is shown in Figure 16.13.
se Since the total displacements u7+1'fc are now defined as ni+U = u:/ + sk(Su{, + 8ft 8ut) (16.64)
reflecting the introduction of the relaxation parameter sk, the IDC constraint equations must be modified accordingly. 83j for the stress criterion is now defined as
I sk (8Xt)n (n)n and 83j for the relative displacement criterion is now defined as
Responses

amanda8 years ago
 Reply

sarah2 years ago
 Reply