Algorithm Template for Nonlinear Dynamic Solver using the Newmark & Newton Loops

Nonlinear dynamic analysis plays an important role in seismic performance evaluation of structures. The solution of nonlinear dynamic can be solved by using Newmarks's method. This article provides a 1-DOF cookbook for review purposes and to provide basic knowledge for multi-DOF problems. The algorithm consists of Newmark loop and Newton loop as

Time Integration (Newmark)

└── for each time step

└── Newton Iteration

└── Stress Return Mapping

Fig. 1 Program structure

The algorithm can be expanded below

Newmark time step t+1

1) Linear predictor: Keff * Δu = ΔFeff

2) Update u, v, a

3) Newton correction if nonlinear:

r = fs + c*v + m*a + m*ag

Keff_t = Kt + a0*m + a1*c

Keff_t * δu = -r

return mapping → fs, Kt

4) Commit state and move to t+2

By separating the Newmark and Newton loop, the solution of complex nonlinear dynamic analysis shall become obvious and easy to manage.

For a special and most common case (average acceleration)

beta = 1/4, gamma = 1/2

a0 = 4/dt^2

a1 = 2/dt

a2 = 4/dt

a3 = a4 = 1

a5 =0

For each time step n → n+1


    Known:

        un, vn, an, state_n



    Step A: Linear Newmark predictor

        Use tangent/stiffness from previous step

        Keff = Kt_n + a0*m + a1*c



        ΔFeff = F(n+1)-F(n)

                + m*(a2*vn + a3*an)

                + c*(a4*vn + a5*an)



        Solve:

            Keff * Δu = ΔFeff



        Update:

            u = un + Δu

            a = a0*Δu - a2*vn - a3*an

            v = a1*Δu - a4*vn - a5*an



    Step B: Newton correction

        repeat until convergence:



            Return mapping:

                [fs, Kt, state_trial] = material_update(u,state_n)



            Residual:

                r = fs + c*v + m*a + m*ag(n+1)



            Effective tangent:

                Keff_t = Kt + a0*m + a1*c



            Solve:

                δu = -r / Keff_t



            Update:

                u = u + δu

                a = a + a0*δu

                v = v + a1*δu



    Step C: Commit

        U(n+1)=u

        V(n+1)=v

        A(n+1)=a

        state = state_trial