Joint Distribution of Eigenvalues of Linear Stochastic Systems S - - PowerPoint PPT Presentation

19 April 2005

Joint Distribution of Eigenvalues

• f Linear Stochastic Systems

S Adhikari

Department of Aerospace Engineering, University of Bristol, Bristol, U.K. URL: http://www.aer.bris.ac.uk/contact/academic/adhikari/home.html

Joint Distribution of Eigenvalues – p.1/36

19 April 2005

Outline of the Presentation

Random eigenvalue problem Existing methods Exact methods Perturbation methods Asymptotic analysis of multidimensional integrals Joint moments and pdf of the natural frequencies Numerical examples & results Conclusions

Joint Distribution of Eigenvalues – p.2/36

19 April 2005

Random Eigenvalue Problem

The random eigenvalue problem of undamped or proportionally damped linear systems: K(x)φj = ω2

jM(x)φj

(1) ωj natural frequencies; φj eigenvectors; M(x) ∈ RN×N mass matrix and K(x) ∈ RN×N stiffness matrix. x ∈ Rm is random parameter vector with pdf px(x) = e−L(x) −L(x) is the log-likelihood function.

Joint Distribution of Eigenvalues – p.3/36

19 April 2005

The Objectives

The aim is to obtain the joint probability density function of the natural frequencies and the eigenvectors in this work we look at the joint statistics of the eigenvalues while several papers are available on the distribution of individual eigenvalues, only first-order perturbation results are available for the joint pdf of the eigenvalues

Joint Distribution of Eigenvalues – p.4/36

19 April 2005

Exact Joint pdf

Without any loss of generality the original eigenvalue problem can be expressed by H(x)ψj = ω2

jψj

(2) where H(x) = M−1/2(x)K(x)M−1/2(x) ∈ RN×N and ψj = M1/2φj

Joint Distribution of Eigenvalues – p.5/36

19 April 2005

Exact Joint pdf

The joint probability (following Muirhead, 1982) density function of the natural frequencies of an N-dimensional linear positive definite dynamic system is given by p Ω (ω1, ω2, · · · , ωN) = πN 2/2 Γ(N/2)

• i<j≤N
• ω2

j − ω2 i

• O(N)

pH

• ΨΩ2ΨT

(dΨ) (3) where H = M−1/2KM−1/2 & pH(H) is the pdf of H.

Joint Distribution of Eigenvalues – p.6/36

19 April 2005

Limitations of the Exact Method

the multidimensional integral over the

• rthogonal group O(N) is difficult to carry out in

practice and exact closed-form results can be derived only for few special cases the derivation of an expression of the joint pdf of the system matrix pH(H) is non-trivial even if the joint pdf of the random system parameters x is known

Joint Distribution of Eigenvalues – p.7/36

19 April 2005

Limitations of the Exact Method

even one can overcome the previous two problems, the joint pdf of the natural frequencies given by Eq. (3) is ‘too much information’ to be useful for practical problems because it is not easy to ‘visualize’ the joint pdf in the space of N natural frequencies, and the derivation of the marginal density functions of the natural frequencies from Eq. (3) is not straightforward, especially when N is large.

Joint Distribution of Eigenvalues – p.8/36

19 April 2005

Eigenvalues of GOE Matrices

Suppose the system matrix H is from a Gaussian

• rthogonal ensemble (GOE). The pdf of H:

pH(H) = exp

• −θ2Trace
• H2

+ θ1Trace (H) + θ0

• The joint pdf of the natural frequencies:

p Ω (ω1, ω2, · · · , ωN) = exp

• −

N

• j=1

θ2ω4

j − θ1ω2 j − θ0

• i<j
• ω2

j − ω2 i

• Joint Distribution of Eigenvalues – p.9/36

19 April 2005

Perturbation Method

Taylor series expansion of ωj(x) about the mean x = µ ωj(x) ≈ ωj(µ) + dT

ωj(µ) (x − µ)

+ 1 2 (x − µ)T Dωj(µ) (x − µ) Here dωj(µ) ∈ Rm and Dωj(µ) ∈ Rm×m are respec- tively the gradient vector and the Hessian matrix of ωj(x) evaluated at x = µ.

Joint Distribution of Eigenvalues – p.10/36

19 April 2005

Joint Statistics

Joint statistics of the natural frequencies can be

• btained provided it is assumed that the x is
• Gaussian. Assuming x ∼ N(µ, Σ), first few

cumulants can be obtained as κ(1,0)

jk

= E [ωj] = ωj + 1 2Trace

• DωjΣ
• ,

κ(0,1)

jk

= E [ωk] = ωk + 1 2Trace (DωkΣ) , κ(1,1)

jk

= Cov (ωj, ωk) = 1 2Trace

• DωjΣ
• (DωkΣ)
• + dT

ωjΣdωk

Joint Distribution of Eigenvalues – p.11/36

19 April 2005

Multidimensional Integrals

We want to evaluate an m-dimensional integral over the unbounded domain Rm: J =

• R

m e−f(x) dx

Assume f(x) is smooth and at least twice differentiable The maximum contribution to this integral comes from the neighborhood where f(x) reaches its global minimum, say θ ∈ Rm

Joint Distribution of Eigenvalues – p.12/36

19 April 2005

Multidimensional Integrals

Therefore, at x = θ ∂f(x) ∂xk = 0, ∀k

• r

df (θ) = 0 Expand f(x) in a Taylor series about θ: J =

• R

m e

−

• f(θ)+ 1

2(x−θ) TDf(θ)(x−θ)+ε(x,θ)

• dx

= e−f(θ)

• R

m e− 1 2(x−θ) TDf(θ)(x−θ)−ε(x,θ) dx

Joint Distribution of Eigenvalues – p.13/36

19 April 2005

Multidimensional Integrals

The error ε(x, θ) depends on higher derivatives of f(x) at x = θ. If they are small compared to f (θ) their contribution will negligible to the value of the

• integral. So we assume that f(θ) is large so that
• 1

f (θ)D(j)(f (θ))

• → 0

for j > 2 where D(j)(f (θ)) is jth order derivative of f(x) eval- uated at x = θ. Under such assumptions ε(x, θ) → 0.

Joint Distribution of Eigenvalues – p.14/36

19 April 2005

Multidimensional Integrals

Use the coordinate transformation: ξ = (x − θ) D−1/2

f

(θ) The Jacobian: J = Df (θ)−1/2 The integral becomes: J ≈ e−f(θ)

• R

m Df (θ)−1/2 e

− 1

2

ξ

Tξ

• dξ
• r

J ≈ (2π)m/2e−f(θ) Df (θ)−1/2

Joint Distribution of Eigenvalues – p.15/36

19 April 2005

Moments of Single Eigenvalues

An arbitrary rth order moment of the natural frequencies can be obtained from µ(r)

j

= E

• ωr

j(x)

• =
• R

m ωr

j(x)px(x) dx

=

• R

m e−(L(x)−r ln ωj(x)) dx,

r = 1, 2, 3 · · · Previous result can be used by choosing f(x) = L(x) − r ln ωj(x)

Joint Distribution of Eigenvalues – p.16/36

19 April 2005

Moments of Single Eigenvalues

After some simplifications µ(r)

j

≈ (2π)m/2ωr

j(θ)e−L(θ)

• DL(θ) + 1

rdL(θ)dL(θ)T − r ωj(θ)Dωj(θ)

• −1/2

r = 1, 2, 3, · · · θ is obtained from: dωj(θ)r = ωj(θ)dL(θ)

Joint Distribution of Eigenvalues – p.17/36

19 April 2005

Maximum Entropy pdf

Constraints for u ∈ [0, ∞]: ∞ pωj(u)du = 1 ∞ urpωj(u)du = µ(r)

j ,

r = 1, 2, 3, · · · , n Maximizing Shannon’s measure of entropy S = − ∞

0 pωj(u) ln pωj(u)du, the pdf of ωj is

pωj(u) = e−{ρ0+n

i=1 ρiui} = e−ρ0e− n i=1 ρiui,

u ≥ 0

Joint Distribution of Eigenvalues – p.18/36

19 April 2005

Maximum Entropy pdf

Taking first two moments, the resulting pdf is a truncated Gaussian density function pωj(u) = 1 √ 2πσj Φ ( ωj/σj) exp

• −(u −

ωj)2 2σ2

j

• where σ2

j = µ(2) j

− ω2

j

Ensures that the probability of any natural frequencies becoming negative is zero

Joint Distribution of Eigenvalues – p.19/36

19 April 2005

Joint Moments of Two Eigenvalues

Arbitrary r − s-th order joint moment of two natural frequencies µ(rs)

jl

= E

• ωr

j(x)ωs l (x)

• =
• R

m exp {− (L(x) − r ln ωj(x) − s ln ωl(x))} dx,

r = 1, 2, 3 · · · Choose f(x) = L(x) − r ln ωj(x) − s ln ωl(x)

Joint Distribution of Eigenvalues – p.20/36

19 April 2005

Joint Moments of Two Eigenvalues

After some simplifications µ(rs)

jl

≈ (2π)m/2ωr

j(θ)ωs l (θ) exp {−L (θ)} Df (θ)−1/2

where θ is obtained from: dL(θ) = r ωj(θ)dωj(θ) + s ωl(θ)dωl(θ) and Df (θ) = DL(θ) +

r ω2

j(θ)dωj(θ)dωj(θ)T −

r ωj(θ)Dωj(θ) + s ω2

l (θ)dωl(θ)dωl(θ)T −

s ωl(θ)Dωl(θ)

Joint Distribution of Eigenvalues – p.21/36

19 April 2005

Joint Moments of Multiple Eigenvalues

We want to obtain µ(r1r2···rn)

j1j2···jn

=

• R

m

• ωr1

j1(x)ωr2 j2(x) · · · ωrn jn(x)

• px(x) dx

It can be shown that µ(r1r2···rn)

j1j2···jn

≈ (2π)m/2 ωr1

j1 (θ) ωr2 j2 (θ) · · · ωrn jn (θ)

• exp {−L (θ)} Df (θ)−1/2

Joint Distribution of Eigenvalues – p.22/36

19 April 2005

Joint Moments of Multiple Eigenvalues

Here θ is obtained from dL(θ) = r1 ωj1(θ)dωj1(θ)+ r2 ωj1(θ)dωj2(θ)+· · · rn ωjn(θ)dωjn(θ) and the Hessian matrix is given by Df (θ) = DL(θ)+

jn,rn

• j = j1, j2, · · ·

r = r1, r2, · · · r ω2

j (θ)dωj(θ)dωj(θ)T −

r ωj (θ)Dωj(θ)

Joint Distribution of Eigenvalues – p.23/36

19 April 2005

Example System

Undamped three degree-of-freedom random system:

m1 m2 m3 k

4

k

5

k

1

k

3

k

2

k

6

mi = 1.0 kg for i = 1, 2, 3; ki = 1.0 N/m for i = 1, · · · , 5 and k6 = 3.0 N/m

Joint Distribution of Eigenvalues – p.24/36

19 April 2005

Example System

mi = mi (1 + ǫmxi), i = 1, 2, 3 ki = ki (1 + ǫkxi+3), i = 1, · · · , 6 Vector of random variables: x = {x1, · · · , x9}T ∈ R9 x is standard Gaussian, µ = 0 and Σ = I Strength parameters ǫm = 0.15 and ǫk = 0.20

Joint Distribution of Eigenvalues – p.25/36

19 April 2005

Computational Methods

Following four methods are compared

• 1. First-order perturbation
• 2. Second-order perturbation
• 3. Asymptotic method
• 4. Monte Carlo Simulation (15K samples) - can be

considered as benchmark. The percentage error: Error = (•) − (•)MCS (•)MCS × 100

Joint Distribution of Eigenvalues – p.26/36

19 April 2005

Scatter of the Eigenvalues

1 1.5 2 2.5 3 3.5 4 100 200 300 400 500 600 700 800 900 1000 Samples Eigenvalues, ωj (rad/s) ω1 ω2 ω3

Statistical scatter of the natural frequencies ω1 = 1, ω2 = 2, and ω3 = 3

Joint Distribution of Eigenvalues – p.27/36

19 April 2005

Error in the Mean Values

1 2 3 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

Eigenvalue number j Percentage error wrt MCS

1st−order perturbation 2nd−order perturbation Asymptotic method

Error in the mean values

Joint Distribution of Eigenvalues – p.28/36

19 April 2005

Error in Covariance Matrix

(1,1) (1,2) (1,3) (2,2) (2,3) (3,3) 2 4 6 8 10 12 14 16 Elements of the covariance matrix Percentage error wrt MCS 1st−order perturbation 2nd−order perturbation Asymptotic method

Error in the elements of the covariance matrix

Joint Distribution of Eigenvalues – p.29/36

19 April 2005

Mean and Covariance

Using the asymptotic method, the mean and covariance matrix of the natural frequencies are

• btained as

µΩ = {0.9962, 2.0102, 3.0312}T and ΣΩ =   0.5319 0.5643 0.7228 0.5643 2.5705 0.9821 0.7228 0.9821 8.7292   × 10−2 Individual pdf and joint pdf of the natural frequencies are computed using these values.

Joint Distribution of Eigenvalues – p.30/36

19 April 2005

Individual pdf

0.5 1 1.5 2 2.5 3 3.5 4 4.5 1 2 3 4 5 6 Natural frequencies (rad/s) pdf of the natural frequencies Asymptotic method Monte Carlo Simulation

Individual pdf of the natural frequencies

Joint Distribution of Eigenvalues – p.31/36

19 April 2005

Analytical Joint pdf

Joint pdf using asymptotic method

Joint Distribution of Eigenvalues – p.32/36

19 April 2005

Joint pdf from MCS

Joint pdf from Monte Carlo Simulation

Joint Distribution of Eigenvalues – p.33/36

19 April 2005

Contours of the joint pdf

0.8 1 1.2 1.4 1.6 1.8 2 2.2 2.4 2.6 1.5 2 2.5 3 3.5 4 (ω1,ω3) (ω1,ω2) (ω2,ω3) ωj (rad/s) ωk (rad/s) Asymptotic method Monte Carlo Simulation

Contours of the joint pdf

Joint Distribution of Eigenvalues – p.34/36

19 April 2005

Conclusions

Statistics of the natural frequencies of linear stochastic dynamic systems has been considered usual assumption of small randomness is not employed in this study. a general expression of the joint pdf of the natural frequencies of linear stochastic systems has been given

Joint Distribution of Eigenvalues – p.35/36

19 April 2005

Conclusions

a closed-form expression is obtained for the general order joint moments of the eigenvalues it was observed that the natural frequencies are not jointly Gaussian even they are so individually future studies will consider joint statistics of the eigenvalues and eigenvectors and dynamic response analysis using eigensolution distributions

Joint Distribution of Eigenvalues – p.36/36

References

Muirhead, R. J. (1982), Aspects of Multivariate Statistical Theory, John Wiely and Sons, New York, USA.

36-1