## 4. Numerical algorithm

The dynamic nonlinear problem is:

f int ( u , u , u , t ) = f ext ( u , u , t )

u = g ( v , t )

h ( u , t ) = 0

The general form of a LMS (linear multi-step) schema for a normalised first-order ODE (ordinary differential equation) of the form:

y = f ( y )

is:

y n = h β y n + h n y

where

h β = h β 0

and

h n y = i = 0 k ( h β i y n i α i y n i )

Adapted to a second-order ODE (ordinary-differential equation) on u where:

y = [ u u ]

we obtain:

u n = u n h n u h β

u n = u n h n u h β = u n h n u h β h n u h β 2

where

h n u = i = 1 k h β i u n i α i u n i

h n u = i = 1 k h β i u n i α i u n i

At each time integration step, the problem to solve is:

f int ( u n , u n , u n , t n ) = f ext ( u n , u n , t n )

u n = g ( v n , t n )

h ( u n , t n ) = 0

The first-derivative of the potential energy of the unconstrained problem is:

u Π n ( u n , u n , u n , t n ) = f int f ext

To impose the constraint, add Lagrange multipliers to form the Lagrangian:

L n ( u n , u n , u n , t n , Λ n ) = n + Λ n T h

The first-order derivatives are:

v L n ( u n , u n , u n , t n , Λ n ) = v T g ( f int f ext + u T h Λ )

Λ L n ( u n , u n , u n , t n , Λ n ) = h

The second-order derivatives are:

vv L n = v T g ( u f int + 1 h b u f int + 1 h b 2 u f int u f ext 1 h b u f ext ) v g + vv T g ( u f int + 1 h b u f int + 1 h b 2 u f int u f ext 1 h b u f ext ) + vv T h Λ

v Λ L n = v T g u T h

Λ v L n = u h v g

Λ Λ L n = 0