Belytschko T. - Introduction (779635), страница 58
Текст из файла (страница 58)
(6.3.9).6.3.3.Newton’s Method. The most widely used and most robust methodfor the solution of the nonlinear algebraic equations (6.3.9) is Newton’s method.The method is often called the Newton-Raphson method in computationalmechanics.
It is identical to the Newton method taught in introductory calculuscourses.We first illustrate the Newton method for one equation in one unknownd without a displacement boundary condition. It is then generalized to anarbitrary number of unknowns. For the case of one unknown, (6.3.9) reduces to asingle nonlinear algebraic equation()() ()sr d n +1, tn +1 = β∆tD 2 M d n +1 − d˙ n − f d n+1 , t n+1 = 06-9(6.3.11)T. Belytschko & B.
Moran, Solution Methods, December 16, 1998The solution of (6.3.11) by Newton’s method is an iterative procedure. Theiteration number is indicated by Greek subsript: dυn +1 is the vth iteration at timestep n+1; when there is no chance for confusion, the time step number will beomitted.To begin the iterative procedure, a starting value for the unknown must bechosen; usually the value of the solution d n from the last time step is used, sod0n +1 ≡ d n Taking a Taylor expansion of the residual about the current value ofthe nodal displacement, dv and setting the resulting residual equal to zero:(0 = r dυ +1,tn +1) = r(dυ ,t ) +n+1(∂r dυ ,t n+1∂d) ∆d +O ∆d 2( )(6.3.12)where∆d = dυ +1 − dυ ,((6.3.12b))()(r dυ , tn +1 = Ma( dυ ) + f int dυ , tn +1 − f ext dυ ,t n+1)(6.3.13)If the terms which are higher order in ∆d than linear are dropped, then (6.3.12)gives a linear equation for ∆d :(0 = r dυ ,tn+1)+() ∆d∂r dυ , tn +1∂d(6.3.14)Note that in the Taylor expansion, the residual is written in terms of the timet n+1 .
The time-dependence of the residual at constant nodal displacements isusually known. For example, if the tractions and body forces are given asfunctions of time, then the time dependent part of the nodal forces is known attime t n+1 at the beginning of the iterations. Therefore the residual is alwayscomputed at time t n+1 . The above is called a linear model of the nonlinearequations. The linear model is the tangent to the nonlinear residual function; theprocess of obtaining the linear model is called linearization.Equation (6.3.14) is often called a linear model of the nonlinear equations,Schnabel (?). Solving this linear model for the incremental displacements gives∆d =− ∂r( dυ ) −1r( dυ ) ∂d (6.3.15)In the Newton procedure, the solution to the nonlinear equation is obtained byiteratively solving a sequence of linear models (6.3.15). The new value for theunknown in each step of the iteration is obtained by rewriting Eq.
(6.3.12b) asdυ +1 = dυ +∆d(6.3.16)6-10T. Belytschko & B. Moran, Solution Methods, December 16, 1998The procedure is illustrated in Fig. 6.1. The process is continued until the solutionis obtained with the desired level of accuracy.rr( d )linear model(tangent)solutionddυ+1 dυ+ 2dυFig. 6.1. Linear models for a nonlinear equationr( d) = 0 .6.3.4.Newton’s Method for n Unknowns. The generalization of thisprocedure to nDOF unknowns is accomplished by replacing the above scalarequations by matrix equations.
The counterpart of Eq. (6.3.12) becomesr( dυ ) +∂r( dυ )∆d+O ∆d2 = 0∂d( )or( )ra dυ +( ) ∆d + O ∆dn DOF∂ra d υb =1∂db∑b2b =0(6.3.17)The matrix ∂r / ∂d is called the Jacobian matrix and will be denoted by A:A=∂r,∂dorAab =∂ra(6.3.18)∂dbUsing (6.3.17) and dropping terms in higher order than linear, Eqs. (6.3.16) can berewritten asr +A∆d = 0(6.3.19)6-11T. Belytschko & B.
Moran, Solution Methods, December 16, 1998which is the linear model of the nonlinear equations. The linear model is difficultto picture for problems with more than one unknown, since r( d) maps ℜn to ℜn ,Figure 6.2 shows the first component of the residual for a function of twounknowns. The linear model is a plane tangent to the nonlinear functionr1 ( d1, d2 ) . The other residual component is another nonlinear function r2 ( d1, d2 ) ,which is not drawn.r1normal Aijtangent Aijd jd2d1Figure 6.2.
Depiction of a residual component r1 as a function of d1 and d2 and the tangent plane.The increment in the nodal displacements in the Newton iterative procedure isobtained by solving (6.3.18), which gives(∆d =− A −1r d υ ,t n+1)(6.3.20)The increment in the nodal displacements is obtained from this system of linearalgebraic equations. The solution of these equations is discussed in Section X.Once the increments in nodal displacements have been obtained, the new valuesof the nodal displacements are obtained by6-12T. Belytschko & B.
Moran, Solution Methods, December 16, 1998dυ +1 = dυ + ∆d(6.3.21)The new displacement is checked for convergence, see Section 6.3.7. If theconvergence criterion is not met, a new linear model is constructed and used tofind another increment in the nodal displacements. The procedure is repeateduntil the convergence criterion is satisfied.In computational mechanics, the Jacobian is called the effective tangentstiffness matrix and the contributions of the inertial, internal and external nodalforces are linearized separately. From (6.3.9) we can writeA=∂rs∂f int ∂f ext= D2 M+−∂d β∆t∂d∂d(6.3.22)where we have used the fact that the mass matrix in a Lagrangian mesh is constantin time and (6.3.4).
The Jacobian of the internal nodal forces is called the tangentstiffness matrix and will be denoted by K int :Kintab∂f aint=∂dbKintiIjJ∂f iIint=∂u jJKint∂f int=∂d(6.3.23)The tangent stiffness matrix is shown above in three forms. The Jacobian matrixof the external nodal forces is called the load stiffness matrix and denoted byextKab=∂f aext∂dbextKiIjJ=∂f iIext∂u jJK ext =∂f ext∂d(6.3.24)The development of these matrices is the topic of linearization and istreated in Sections 6.4 and 6.5.
Using these definitions, the Jacobian matrix(6.3.22) can be written asA=sDM + Kint − K extβ∆t 2(6.3.25)This Jacobian matrix applies to both dynamic and equilibrium problems with thedynamic switch sD set by (6.3.2).The Jacobians in (6.3.23-24) can be used to relate differentials of the nodalforces to differentials of the nodal displacements bydf int = Kint dddf ext = K extdddr = Add(6.3.26)The matrices which relate finite increments of nodal displacements to incrementsof nodal forces differ from the above.
We will use a∆fint = Kint∆ ∆d∆f ext = Kext∆ ∆d6-13∆r = A∆ ∆d(6.3.27)T. Belytschko & B. Moran, Solution Methods, December 16, 1998The matrix Kint∆ is called a secant stiffness and A ∆ the secant Jacobian. Thesecant stiffnes and secant Jacobian depend on the magnitude and direction of ∆d.This can easily be seen in one dimension as illustrated in Fig.
6.3, which showssecants for various stepsizes and two directions (there are only two in a functionof a single variable). The tangent and secant Jacobians are identical only in thelimit as ∆d → 0; for finite increments, the secant stiffness in (6.3.27) differs fromthe tangent stiffness in (6.3.23).r (d )A −∆ dA 2∆dA∆∆dA (tangent)Figure 6.3. Secant Jacobians for various step sizes and directions.Conservative Problems (Stationary Points). It is useful at this point to examinethe discrete problem corresponding to the stationary principle described in Section4.9.3.
This stationary principle only applies to conservative equilibriumproblems, but it is nevertheless provides insight into the character of nonlinearproblems. An equilibrium solution is a stationary point of the potential, so byenforcing the conditions that the derivative of the potential vanish and using(4.9.29-30) and the definition of the residual (6.3.3) we have0 =r = −∂W∂d=∂W int∂d−∂W ext∂d= f int − f ext(6.3.28)A solution is a stable equilibrium solution if it corresponds to a minimum of thepotential energy. Thus stable equilibrium solutions can be found by minimizingthe potential W. The situation is depicted in Fig. 6.3, which shows the localbehavior of a potential of two generalized displacements and the contours for thispotential.
The residual is the negative of the gradient of the potential (note thesign in the above.)The linear model for (6.3.28) is (see 6.3.17-18)6-14T. Belytschko & B. Moran, Solution Methods, December 16, 1998−rν =∂r∂d−raυ =∆d = −∂2W∂d∂d∆d = A∆d where Aab =∂2 W∂2Wor A =∂da∂db∂d∂d∂ra∂2W∂2 W∆db = −∆db = Aab ∆db where Aab =∂db∂da∂db∂da∂db(6.3.29)The matrix A when it arises from the second derivatives of a potential is called aHessian matrix. It is identical to the Jacobian, soA = Kint − K ext(6.3.30)The linearized equations for a conservative system are(Kint − K ext )∆d =− rThe above are identical to Eq. (6.3.19) except that the mass matrix is omitted,since dynamic effects cannot be included in a conservative problem. However,when the problem is posed as a minimization problem, many techniques notdirectly applicable to linear models, such as the method of steepest descent, canbe applied to the problem.