Thread Subject: Help!!! Eigenvalue program for frequency dependent matrix

I am trying find the natural frequency of a vibration problem

The standard eigenvalue problem looks like

Kx=lMx

for which there is a simple eig(K,M) command that does the
trick.

The trouble is that the K matrix is frequency dependent in
this case. So, each element of the K matrix is a function of
frequency. How can I approach this type of a problem in matlab?

Help will be much appreciated!!

Thanks,
YM

"Yogesh " <yahoo@trashmail.net> wrote in message
<g51ddf$pm4$1@fred.mathworks.com>...
> I am trying find the natural frequency of a vibration problem
>
> The standard eigenvalue problem looks like
>
> Kx=lMx
>
> for which there is a simple eig(K,M) command that does the
> trick.
>
> The trouble is that the K matrix is frequency dependent in
> this case. So, each element of the K matrix is a function of
> frequency. How can I approach this type of a problem in
matlab?
>
> Help will be much appreciated!!
>
> Thanks,
> YM

Hi

Your problem may be rewritten into the form

    K(omega).X = M.X.S

where
K(omega) is a stiffness matrix of a system dependent on an
operational frequency (speed),
M is a mass matrix of the system,
X is a modal matrix of the system composed of modal vectors,
S is a diagonal spectral matrix with eigenvalues on diagonal.

To each operational frequency omega belong matrices X and S,
which are different from those belonging to other omega. The
system is "evolutionary".

You wanted to solve the eigenvalue problem for omega as a
function of S. It is possible to solve it in iterations for
a single eigenvalue s_k, if you reguire that, say, omega=
real(s_k). Then, you may use the MATLAB function fzero
applied to the difference omega-real(s_k).

I hope that it helps.

Mira

"Miroslav Balda" <balda.nospam@cdm.it.cas.cz> wrote in
message <g51igi$8g1$1@fred.mathworks.com>...
:
SNIP
:

> You wanted to solve the eigenvalue problem for omega as a
> function of S. It is possible to solve it in iterations for
> a single eigenvalue s_k, if you reguire that, say, omega=
> real(s_k). Then, you may use the MATLAB function fzero
> applied to the difference omega-real(s_k).

Sorry, there should be omega = imag(s_k).

Mira

"Miroslav Balda" <balda.nospam@cdm.it.cas.cz> wrote in
message <g51igi$8g1$1@fred.mathworks.com>...
:
SNIP
:

> You wanted to solve the eigenvalue problem for omega as a
> function of S. It is possible to solve it in iterations for
> a single eigenvalue s_k, if you reguire that, say, omega=
> real(s_k). Then, you may use the MATLAB function fzero
> applied to the difference omega-real(s_k).

Sorry, there should be omega = imag(s_k).

Mira

"Miroslav Balda" <miroslav.nospam@balda.cz> wrote in message
<g53alh$2tc$1@fred.mathworks.com>...
> "Miroslav Balda" <balda.nospam@cdm.it.cas.cz> wrote in
> message <g51igi$8g1$1@fred.mathworks.com>...
> :
> SNIP
> :
>
> > You wanted to solve the eigenvalue problem for omega as a
> > function of S. It is possible to solve it in iterations for
> > a single eigenvalue s_k, if you reguire that, say, omega=
> > real(s_k). Then, you may use the MATLAB function fzero
> > applied to the difference omega-real(s_k).
>
> Sorry, there should be omega = imag(s_k).
>
> Mira
>

Hi Mira,

Thanks for your response. I had implemented something
similar before posting it on the newsgroup.

Here is what I was doing:

1. Assume a certain omega.
2. Find K for this omega
3. Find eigenvalues
4. Check is the eigenvalue is the same as the omega assumed
in step-1.
5. If yes, this omega is an eigenvalue. Store the freuqency
and eigenvector.

There are "numerical" problems with this approach though. I
know that the frequency will lie between 0-10,000 Hz. Now,
your ability to predict the eigenvalue depends upon how good
your approximation in (1) is. This, in my opinion, is not
the most elegant approach.

I tried vector iteration approaches but they can only give
the first and last eigenvalues.

I appreciate your inputs and thank you again.

On 11 Jul, 23:08, "Yogesh " <ya...@trashmail.net> wrote:
> "Miroslav Balda" <miroslav.nos...@balda.cz> wrote in message
> Here is what I was doing:
>
> 1. Assume a certain omega.
> 2. Find K for this omega

Depending on exactly what you do, this will be the wrong answer.
If this is a wave propagation problem (the terminology suggests
it is) then the eigenvalues you are looking for are pairs of
frequencies and wavenumbers (w,k) that make the boundary
conditions in the physical environment match.

For instance, if you work with acoustic waves in a perfect
underwater waveguide, the reflection coefficient at the
water surface is -1 and at the bottom +1.

So you want to find the pair(s) (w,k) that meets the
boundary conditions *both* at the surface and at the
bottom. Depending on exactly what you look for, there may
be any number of solutions.

Check out the book 'Computational Ocean Acoustics' by
Jensen & al to find a number of ways to find eigenvalues
in wave propagation problems.

Rune

"Yogesh " <yahoo@trashmail.net> wrote in message
<g58i3i$8f6$1@fred.mathworks.com>...

Sorry that I am replaying so late. I was off the Internet.
My reply concernes discrete mechanical systems, say
multibody systems.
 :
 SNIP
 :
> Thanks for your response. I had implemented something
> similar before posting it on the newsgroup.
>
> Here is what I was doing:
>
> 1. Assume a certain omega.
> 2. Find K for this omega
> 3. Find eigenvalues
> 4. Check is the eigenvalue is the same as the omega assumed
> in step-1.
> 5. If yes, this omega is an eigenvalue. Store the freuqency
> and eigenvector.
>

In principal, the procedure is in order. However, I have
also found problems with convergence. This has been a reasom
why I aplied the function LMFnlsq (FEX Id 17524). The best
way how to show it is a solved example:
%Main script:
% Yogesh Balda, 2008-07-14
%%%%%%%%%%%%
omega = 300; % Initial guess for omega near to first
eigenvalue
omega = LMFnlsq('KomgM',omega,'Display',1)

% Example of a cost function:
function ds = KomgM(omg)
% ds difference between eigenvalue and omega
% Mass matrix:
M = diag([5, 32.12, 32.12, 5]); % [kg]
% Stiffness matrix:
K = [0.4397, -0.6143, 0.4095, -0.0683 % [N/m]
     -0.6143, 1.6381, -1.4334, 0.4095
      0.4095, -1.4334, 1.6381, -0.6143
     -0.0683, 0.4095, -0.6143, 0.4397]*1e7;
K = K-10*omg*diag(ones(4,1)); Frequeny dependent K
s = sqrt(eig(K,M));
omega = s(k);

And results for the second natural frequency:
>> Yogesh
***************************************************************************
  itr nfJ SUM(r^2) x dx l
      lc
***************************************************************************
   0 1 1.4677e+004 3.0000e+002 0.0000e+000
0.0000e+000 1.0000e+000
   1 2 1.2516e-006 4.2099e+002 -1.2099e+002
0.0000e+000 1.0000e+000
   2 3 6.6392e-017 4.2099e+002 1.1173e-003
0.0000e+000 1.0000e+000
 
   3 4 6.6392e-017 4.2099e+002 2.3314e-008
0.0000e+000 1.0000e+000
omega =
  420.9890

It is obvious that the solution has been found in 4
iterations and that square of delta s(2)-omega dropped down
by 21 orders! It is necessary to stress out that the code is
not optimized.

Hope it helps

Mira