2/15/2008

Introduction to
SAMCEF MECANO

1
Introduction to SAMCEF - MECANO

Outline

 Introduction
 Generalized coordinates
 Kinematic constraints
 Time integration
 Description and paramerization of finite rotations

Acknowledgements
Michel Géradin (ULg)
Doan D. Bao (GoodYear)
Alberto Cardona (Intec, Argentina)
2
Introduction to SAMCEF - MECANO

1
 2/15/2008

Our Goal
Fixeddeformable
Fixed deformable Articulated
Articulated
structures
structures rigid systems
rigid systems

?
FiniteElement
Finite Element Multibody
Multibody
Method
Method Dynamics
Dynamics

how to analyze in a general

manner flexible articulated
structures

3
Introduction to SAMCEF - MECANO

Strategic choices
Kinematic
Kinematicdescription
description
Recursive?
Absolute?
Absolute
Structural
Structural behavior
behavior
Rigid/flexible
Flexible
Rigid ? ?

Software
architecture

Motion
Motionequations
equations
formulation
formulation
Time Numerical
Symbolic
Numerical??
Timeintegration
integration
Implicit
Explicit
Implicit?

4
Introduction to SAMCEF - MECANO

2
 2/15/2008

Finite rotation Finite motion

parameterization description

Measures of deformation
and relative motion

Rigid bodies Joints Flexible elements

zero deformation relative motion elastic strains
+ external work + internal work

Mechanical
Control system elements library
FE structural
• rigid members
library substructuring
• flexible members
• rigid articulations
• flexible articulations
• active members Modelling
Step
5
Introduction to SAMCEF - MECANO

Mechanical
Control system elements library
FE structural
• rigid members
library substructuring
• flexible members
• rigid articulations
• flexible articulations
• active members

DAE system equations

• dynamic equilibrium
• kinematic constraints

Post-processing
kinematics
internal forces
element stresses Modelling
Step
6
Introduction to SAMCEF - MECANO

3
 2/15/2008

The MECANO Preprocessing

Preprocessing
••topological
topologicaldescription
description(FE)
(FE)
software ••mechanical
mechanicalproperties
properties
••response
responsespecification
specification

Initial Assembly

Static
Staticanalysis
analysis Kinetostatic
Kinetostaticanalysis
analysis Dynamic
Dynamicanalysis
analysis
equilibrium
equilibriumconfiguration
configuration trajectories
trajectories trajectories
trajectories
static
staticforces
forces quasi-static
quasi-staticforces
forces dynamic
dynamicforces
forces

Linearized
Linearizedanalysis
analysis
eigenvalues
eigenvalues
stability
stability

Preprocessing
Preprocessing
••trajectories
trajectories
••loads
loadsand
andstresses
stresses
••animation
animation 7
Introduction to SAMCEF - MECANO

The MECANO Some rigid joint elements

element library

Spherical joint
Hinge joint
Prismatic joint

Rack-and-pinion
Cylindrical joint
Gear joint

8
Introduction to SAMCEF - MECANO Kinematic constraints

4
 2/15/2008

some non-
The MECANO element library
holonomic joints

Cam contact Wheel + tyre model

Point-on-plane

Cable + pulley

9
Introduction to SAMCEF - MECANO Kinematic constraints

The basic elastic elements

The MECANO
element library

Superelement Spring-damper-stop
Bushing

Nonlinear beam

10
Introduction to SAMCEF - MECANO Kinematic constraints

5
 2/15/2008

The MECANO element library

Some nonlinear and control features

Nonlinear constraint

Hydraulic actuator

Nonlinear feedback

11
Introduction to SAMCEF - MECANO Kinematic constraints

The MECANO element library

Rail-wheel contact

12
Introduction to SAMCEF - MECANO Kinematic constraints

6
 2/15/2008

The MECANO element library

Flexible slider

13
Introduction to SAMCEF - MECANO Kinematic constraints

The MECANO element library

Symbolic element representation

14
Introduction to SAMCEF - MECANO Kinematic constraints

7
 2/15/2008

Domains of application

• mechanical engineering
• automotive industry
• sport industry
• aeronautics
• space structures
• biomechanics
• ...

15
Introduction to SAMCEF - MECANO

Adequacy of Mecano
Dynamic car behavior • association of kinematic
joints with finite elements
• open software through user
element capability

M. Ebalard, J. Mercier (PSA Citroën)

Utilisation du logiciel Mecano pour
le comportement routier à PSA
16
Introduction to SAMCEF - MECANO Examples of application

8
 2/15/2008

Landing gear analysis

Research conducted in framework of ELGAR consortium

17
Introduction to SAMCEF - MECANO Examples of application

Landing gear analysis

Functional phases
• deployment/retraction
• ground impact
• rolling
• breaking
• taxiing

Model components
• mechanism
• damper (internal
hydraulic system)
• tyre: contact,
deformation, friction

18
Introduction to SAMCEF - MECANO Examples of application

9
 2/15/2008

MEA antenna: result animation

19
Introduction to SAMCEF - MECANO Examples of application

Simulation of Cluster satellite boom deployment

20
Introduction to SAMCEF - MECANO Examples of application

10
 2/15/2008

Generalized coordinates for mechanism analysis

• Generalized co-ordinates
description of various types
„

the 4-bar mechanism example

„

• The finite element method for articulated systems

kinematics
„

elasto-dynamics
„

 differential-algebraic equations of motion

 numerical integration by Newmark algorithm
 the double pendulum example

21
Introduction to SAMCEF - MECANO Generalized coordinates

Introduction
• Classical approach: minimal number of DOF.
• Lagrangian coordinates: best suited to robotics
• open-tree structure + associated graph
• loop closure conditions
• lack of generality, difficult for flexible systems
• Cartesian coordinates: large scale simulation packages
• inherent topology description
• large sets of differential-algebraic equations (DAE)
• still difficult for flexible systems
• Nonlinear Finite element coordinates:
• inherent topology description
• simplification of kinematic constraints

• geometrically nonlinear effects naturally included

22
Introduction to SAMCEF - MECANO Generalized coordinates

11
 2/15/2008

The 4-
4-bar mechanism example

• Simplest planar, closed-loop mechanism

• one single variable: crank angle β

Gruebler’s formula m : number of joints

N: number of bodies
p : number of DOF removed by joint i
m
X
F = 3(N − 1) − pi F=1
i=1

23
Introduction to SAMCEF - MECANO Generalized coordinates

Minimal coordinates Expression in terms of minimal

coordinate β

θj = f (β) j = 2, . . . 4

From geometry

d2 + s2 + 2r` cos β
cos θ2 =
2ds − r2 − `2

(` − r cos β) cos θ3 + r sin β sin θ3 = s − d cos θ2

c
sin(θ3 + φ) =
ρ
θ1 = 2π − β − θ2 − θ3

24
Introduction to SAMCEF - MECANO Generalized coordinates

12
 2/15/2008

Minimal coordinates: kinematic transformer

Input-output relationship

Complex mechanisms as assembly of

transformers

25
Introduction to SAMCEF - MECANO Generalized coordinates

Lagrangian coordinates

3 generalized coordinates

q T = [ θ1 θ2 θ3 ]

2 closure conditions

`1 cos θ1 + `2 cos(θ1 + θ2 ) + `3 cos(θ1 + θ2 + θ3 ) − `4 = 0

`1 sin θ1 + `2 sin(θ1 + θ2 ) + `3 sin(θ1 + θ2 + θ3 ) = 0

Equations of motion: 3 differential

equations linked by 2 constraints
(DAE)

26
Introduction to SAMCEF - MECANO Generalized coordinates

13
 2/15/2008

Cartesian coordinates

9 generalized coordinates q T = [ x1 y1 θ1 x2 y2 θ2 x3 y3 θ3 ]

Equations of motion: 9
differential equations linked
by 8 constraints (DAE)

8 kinematic constraints
`1
x1 − 2 cos θ1 = 0 y1 − `21 sin θ1 = 0
`1
x1 + 2 cos θ1 − x2 + `22 cos θ2 = 0 y1 + `21 sin θ1 − y2 + `22 sin θ2 = 0
`2
x2 + 2 cos θ1 − x2 + `23 cos θ3 = 0 y2 + `22 sin θ2 − y3 + `23 sin θ3 = 0
`3
x3 + 2 cos θ3 = `4 y3 + `23 sin θ3 = 0

27
Introduction to SAMCEF - MECANO Generalized coordinates

Finite Element coordinates 13 generalized coordinates:

• 12 finite element DOF q
q = [ x1 y1 . . . x6 y6 ]T • 1 driving DOF β

Strong form of kinematic constraints • 3 zero strain constraints

• 4 boundary nodal constraints 1 `2i − `2i,0
²i =
x1 = 0, y1 = 0 2 `2i,0
x6 = 0, y6 = `4
²1 = 0 ²2 = 0 ²3 = 0
• 4 assembly nodal constraints
• 1 driving constraints (implicit)
x2 = x3 y2 = y3
x4 = x5 y4 = y5 y2 cos β − x2 sin β = 0
ΦD (q) = 0
28
Introduction to SAMCEF - MECANO Generalized coordinates

14
 2/15/2008

Finite element description of a mechanical system

q = DOF at structural level

q e = DOF at element level

Le = DOF localization operator(boolean)

λ = lagrangian multipliers

• Boolean constraints: localization of DOF

• Implicit constraints: algebraic treatment

q e = Le q
System topology results implicitly
from Boolean assembly and Φ(qe , q̇ e , t) = 0
constraint description

29
Introduction to SAMCEF - MECANO Kinematic constraints

Classification of kinematic constraints

• Holonomic constraints ex: prismatic joint

• Restrict the number of DOF

of the system
Φ(q, t) = 0
• Non-holonomic constraints ex: rolling motion

• Restrict the system behavior

Ψ(q, q̇, t) = 0

• Unilateral constraints Ψ(q, q̇, t) ≥ 0

30
Introduction to SAMCEF - MECANO Kinematic constraints

15
 2/15/2008

Only 3 lower pairs

Classical mechanical pairs • surface contact
• motion reversibility

• All other
mechanical pairs
are higher pairs planar
revolute prismatic

screw spherical
cylindrical

Other modes of classification:

• number of DOF (1 for all lower pairs)
• mode of closure
31
Introduction to SAMCEF - MECANO

The algebraic constrained problem

⎧ Possible methods of solution
⎪
⎪ min F(q)
⎪
⎪ q • constraint elimination
⎪
⎪
⎨ • Lagrange multipliers
subject to • penalty function
⎪
⎪ • augmented Lagrangian
⎪
⎪
⎪
⎪ • perturbed Lagrangian
⎩
Φ(q) = 0
Jacobian matrix

Formulated in terms of ∂Φm

B such that Bmj =
∂qj
Residual vector
Hessian matrix

∂F ∂2 F
r(q) = H such that Hij =
∂q ∂qi ∂qj

32
Introduction to SAMCEF - MECANO

16
 2/15/2008

Constraint elimination method

• First variation of constraints
δΦ = Bδq = 0
• Partitioning into independent and
dependent variables

B I δq I + B D δqD = 0
• Elimination of dependent variables
⎡ ⎤
I
reduction
δq = RδqI R=⎣ ⎦
matrix
−B −1
D BI
• Projected residual equation

∂F • exact verification of constraints

RT = RT r(q) = 0 • complex and computationally
∂q
expensive

33
Introduction to SAMCEF - MECANO

Augmented Lagrangian method

Augmented functional
p
Fp∗ (q, λ) = F(q) + (kλT Φ + ΦT Φ)
2 p penalty coefficient
k scaling factor
Second-order solution method: by Newton-Raphson

⎡ ⎤⎡ ⎤ ⎡ ⎤
q = q? + ∆q H + pBT B kB T ∆q −r? − BT (kλ? + pΦ)
⎣ ⎦⎣ ⎦=⎣ ⎦
λ = λ? + ∆λ
kB 0 ∆λ kΦ

• provides exact solution regardless of penalty and scaling factors

• quadratic term pBT B reinforces positive character of H

• scaling factor required for pivoting

34
Introduction to SAMCEF - MECANO

17
 2/15/2008

Constrained dynamic problems

Z t2
Hamilton’s
principle
[δ(L − kλT Φ − pΦT Φ) + δW]dt = 0
t1

Matrix form of
constrained motion M q̈ + BT (pΦ + kλ) = g(q, q̇, t)
equations Φ(q, t) = 0

• Differential - Algebraic nature (DAE)

• vanishing of penalty term at equilibrium → exact
• scaling of constraints for numerical conditioning
• multipliers λ = forces needed to close constraints
• all remaining forces (internal, external, complementary inertia) in RHS
35
Introduction to SAMCEF - MECANO Kinematic constraints

Constrained equations of motion

M q̈ + BT (pΦ + kλ) = g(q, q̇, t)

Φ(q, t) = 0

Methods of solution

• Constraint reduction
• Constraint regularization → m+n ODE + constraint stabilization
• second-order solution → linearization

But: specific problem of numerical

stability of time integration

36
Introduction to SAMCEF - MECANO Kinematic constraints

18
 2/15/2008

Linearization of motion equations

¸ ¸ ¸ ¸ ¸ ¸ ¸
M 0 ∆q̈ Ct 0 ∆q̇ Kt kBT ∆q r?
+ + = + O(∆2 )
0 0 ∆λ̈ 0 0 ∆λ̇ 0 kB ∆λ −Φ?

Residual vector of equilibrium r = g(q, q̇, t) − M q̈ − B T (pΦ + kλ)

Tangent damping matrix Tangent stiffness matrix

∂g ∂g ∂(M q̈) ∂(B T (pΦ + kλ))

Ct = − Kt = −
∂q
+
∂q
+
∂q
∂ q̇
approx.
• plays central role for
convergence in iteration ∂g ∂(B T λ)
• quality of approximation needed Kt ' − + pB T B + k
∂q ∂q
depends on penalty factor

37
Introduction to SAMCEF - MECANO Kinematic constraints

Integration of motion equations: Implicit Method

Motivations
filtering of high frequencies (elasticity)
• independence of algorithm stability on time step size
• one-step method simplicity of software
• experience in structural dynamics.

Newmark algorithm

Application to multibody systems

• appropriate integration of rotation motion
• stability in presence of constraints (DAE)
38
Introduction to SAMCEF - MECANO

19
 2/15/2008

Implicit integration
Differential - algebraic equations

r(q, q̇, q̈, λ) = M q̈ + BT λ − g(q, q̇, t) = 0

Φ(q, t) = 0

Correction of solution at time t Linearized matrices

Tangent damping
q= q∗ + ∆q ∂r ∂g
(q∗ , q̇ ∗ , q̈∗ , λ∗ ) q̇ = q̇∗ + ∆q̇ CT = =−
∂ q̇ ∂ q̇
q̈ = q̈∗ + ∆q̈
Tangent stiffness
λ= λ∗ + ∆λ
∂r ∂g ∂
Linearized equations of motion KT = =− + (M q̈ ∗ + BT λ∗ )
∂q ∂q ∂q

¸ ¸ ¸ ¸ ¸ ¸ ¸
M 0 ∆q̈ CT 0 ∆q̇ KT BT ∆q r∗
+ + =
0 0 ∆λ̈ 0 0 ∆λ̇ B 0 ∆λ −Φ∗

39
Introduction to SAMCEF - MECANO

Newmark formula

• Interpolation of displacements and velocities

q̇ n+1 = q̇ n + (1 − γ)hq̈ n + γhq̈ n+1

q n+1 = q n + hq̇ n + ( 12 − β)h2 q̈ n + +βh2 q̈ n+1

• Free parameters: average constant acceleration

1 1
γ= and β=
2 4

• Introduction of numerical damping

1 1 1
γ= +α and β= (γ + )2 α>0
2 4 2

40
Introduction to SAMCEF - MECANO

20
 2/15/2008

• Prediction at time t Integration procedure

q̈0n+1 =0
q̇0n+1 = q̇ n + (1 − γ)hq̈n
q0n+1 = q n + hq̇n + ( 12 − β)h2 q̈n
• Correction
q̈k+1 k k
n+1 = q̈ n+1 + ∆q̈ = q̈ n+1 +
1
β h2 ∆q
γ
q̇k+1 k k
n+1 = q̇ n+1 + ∆q̇ = q̇ n+1 + β h ∆q
k k
qn+1 = q n+1 + ∆q

• Correction equation
¸ ¸ ¸
ST BT ∆q r∗
=
B 0 ∆λ −Φ∗
Iteration γ 1
ST = KT + CT + M
matrix βh βh2
• Convergence check
krk < ² and kΦk < η

41
Introduction to SAMCEF - MECANO

The double pendulum • Description by Cartesian coordinates

(redundant)
q T = [ x1 y1 θ1 x2 y2 θ2 ]

• Kinetic and potential energies

m1 = m2 = 1.
l1 = l2 = 1., K = 12 m1 (ẋ21 + ẏ12 ) + 12 m2 (ẋ22 + ẏ22 )
(θ1 = π/2, θ2 = 0 P = m1 gy1 + m2 gy2
• Kinematic constraints
⎡ ⎤
x1 − `1 sin θ1
⎢ y1 + `1 cos θ1 ⎥
• Constraint matrix Φ=⎢ ⎥
⎣ x2 − x1 − `2 sin θ2 ⎦
⎡ ⎤
1 0 −`1 cos θ1 0 0 0 y2 − y1 + `2 cos θ2
⎢ 0 1 −`1 sin θ1 0 0 0 ⎥
B=⎢
⎣ −1
⎥
0 0 1 0 −`2 cos θ2 ⎦
0 −1 0 0 1 −`2 sin θ2
• External force vector and mass matrix are constant
• Tangent stiffness contains only geometric terms
K T = diag [ 0 0 `1 (λ1 sin θ1 − λ2 cos θ1 ) 0 0 `2 (λ3 sin θ2 − λ4 cos θ2 ) ]

42
Introduction to SAMCEF - MECANO

21
 2/15/2008

Newmark integration : undamped (γ = 12 , β = 14 )

• normal evolution of angular

displacements and velocities
• but: numerical instability
detected after 2 sec.

43
Introduction to SAMCEF - MECANO

Newmark integration : undamped (γ = 12 , β = 14 )

Weak instability: generated by

kinematic constraints

Numerical damping

44
Introduction to SAMCEF - MECANO

22
 2/15/2008

Newmark integration( with damping

1
(γ = 2 + α, β = 14 ( 12 + γ)2 , α = 0.015)

Normal aspect of response over 5 sec. period, but ...

45
Introduction to SAMCEF - MECANO

Newmark integration( with damping

1
(γ = 2 + α, β = 14 ( 12 + γ)2 , α = 0.015)

Unacceptable amount of energy loss

Search for better compromise

between accuracy and stability
46
Introduction to SAMCEF - MECANO

23
 2/15/2008

Kinematics of finite motion

• Spherical motion
 Eigenvalue analysis of rotation operator
 Euler theorem
 Explicit expressions of rotation operator
 Expression in terms of linear invariants
 Exponential map
• General motion of rigid body
• Velocity analysis of spherical motion
• Explicit expression of angular velocities
 Rotation parameterization
 Cartesian rotation vector

47
Introduction to SAMCEF - MECANO Finite Motion

Spherical motion Spatial coordinates

Current position:
material configuration x x = [x1 x2 x3 ]T

material coordinates
Initial position: material
configuration X X = [X1 X2 X3 ]T

Base vectors of
current configuration [~e1 ~e2 ~e3 ]

Spherical motion such that:

• length of position vector unaffected by
pure rotation
Base vectors of • relative angle between 2 directions remain
initial configuration constant
~1 E
[E ~2 E
~ 3] Linear transformation

x = RX with RT = R−1 det(R) = 1

Proper orthogonal
48
Introduction to SAMCEF - MECANO Finite Motion

24
 2/15/2008

Euler theorem

If a rigid body undergoes a motion

leaving fixed one of its points, then
a set of points of the body lying on
a line that passes through that
point, remains fixed as well.

49
Introduction to SAMCEF - MECANO Finite Motion

Explicit expressions of rotation operator

Orthonomality implies 6 constraints

R = [ r1 r2 r3 ] rTi rj = δij

3 independent rotation parameters R = R(α1 , α2 , α3 )

Outer product expression Inner product expression

⎡ ⎤
X E T1 e1 E T1 e2 E T1 e3
R= ei E Ti ⎣
R = E T2 e1 E T2 e2 E T2 e3 ⎦
i
E T3 e1 E T3 e2 E T3 e3

50
Introduction to SAMCEF - MECANO Finite Motion

25
 2/15/2008

Expression in terms of linear invariants

Rotation operator

R = I cos φ +(1 −cos φ)nnT + ñsinφ

Invariants

tr(R) = λ1 +λ2 +λ3 = 1 +2cosφ

vect(R) = nsinφ

51
Introduction to SAMCEF - MECANO Finite Motion

Exponential map

• Differentiate with respect to φ

x = RX
• verify that
dR T
R = ñ
dφ
• differential motion

dx
− ñx = 0 with x(0) = X
dφ
• integration

x = exp(ñφ) X R = exp(ñφ)

52
Introduction to SAMCEF - MECANO Finite Motion

26
 2/15/2008

General motion of rigid body

Combination of
translation and rotation

xP = x0 + RX P

spatial motion material

of origin coordinates of
point P

53
Introduction to SAMCEF - MECANO Finite Motion

Velocity analysis of spherical motion

Spatial description Material description
Differentiate

xP = RX P with RT = R−1

spatial velocities material velocities

vP = ẋP = ṘX P V P = RT vP
skew-symmetry

= ṘRT + (ṘRT )T = 0

vP = ω̃xP with ω̃ = ṘRT V P = Ω̃X P with Ω̃ = RT Ṙ

spatial angular velocities material angular velocities

vect(ṘRTT)
ωω== vect(ṘR ) Ω T
vect(RTṘ)
Ω==vect(R Ṙ)
54
Introduction to SAMCEF - MECANO Finite Motion

27
 2/15/2008

Cartesian rotation vector

defined as Ψ = nφ
Rotation operator

Trigonometric form Exponential map

sin kΨ k e 1 − cos kΨ k e e
R=I+ Ψ + ΨΨ e
R = exp(Ψ)
kΨ k kΨ k2

Angular velocities

Ω = T (Ψ)Ψ̇ ω = T T (Ψ)Ψ̇ Limit properties

Tangent operator lim T (Ψ) = lim R(Ψ) = I

k Ψk→0 k Ψk→0
µ ¶ µ ¶ e e
cos kΨk − 1 e sin kΨk Ψ Ψ
T (Ψ) = I + Ψ+ 1−
kΨk2 kΨk kΨk2

55
Introduction to SAMCEF - MECANO Finite Motion

28