Non Linear Eigenvalue problem

I have been trying to solve a Non linear Eigenvalue problem using fsolve and Newton's iteration method and have not been successful. The matrix which I am looking to solve:

A=[2*w -300 0 0;sin(w/2) cos(w/2) -sin(w/2) cos(w/2); 2*cos(w/2) -2*sin(w/2) -cos(w/2) sin(w/2); 0 0 cos(w) -sin(w)]; (Found in a paper, using it as a practice case)

My Newton method code:

%Define intial values and tolerances for the variable
w0=0.1;
tol=2;
maxiter=1000;
w=w0;
wold=w0;
lambda=0.1;

%Start Iteration
for i=1:maxiter
    %Define A and B
A=[2*w -300 0 0;sin(w/2) cos(w/2) -sin(w/2) cos(w/2); 2*cos(w/2) -2*sin(w/2) -cos(w/2) sin(w/2);
0 0 cos(w) -sin(w)];

B=[-2 0 0 0;-0.5*cos(w/2) 0.5*sin(w/2) 0.5*cos(w/2) -0.5*sin(w/2); sin(w/2) cos(w/2) -0.5*sin(w/2) -0.5*cos(w/2);
    0 0 sin(w) cos(w)];
C=inv(B);

%Find Eigen value for the intermediate step
beta=eig(C*A);
epsilon=min(abs(beta));

%Update the variable
w=w0+epsilon;
err=abs(epsilon);
wold=w;
if(err<tol)
    break;
end
end

Fsolve code

function fval=fun4evp(w)
A=[2*w -300 0 0;sin(w/2) cos(w/2) -sin(w/2) cos(w/2); 2*cos(w/2) -2*sin(w/2) -cos(w/2) sin(w/2);
0 0 cos(w) -sin(w)];
fval=det(A);
end

wsol=fsolve(@(w)fun4evp,0.1);

Thanks

[Gyan Swarup Nag]

clear all;

% One root at a time

% Number of iteration

k=0; % Number of iteration

ks=300;

epsi=1*0.01;

w0=-1

while 1e-14<abs(epsi)

k=k+1;

A=[2*w0 -ks 0 0;

sin(w0/2) cos(w0/2) -sin(w0/2) -cos(w0/2);

2*cos(w0/2) -2*sin(w0/2) -cos(w0/2) sin(w0/2);

0 0 cos(w0) -sin(w0)];

B=-[2 0 0 0;

cos(w0/2)/2 -sin(w0/2)/2 -cos(w0/2)/2 sin(w0/2)/2;

-sin(w0/2) -cos(w0/2) sin(w0/2)/2 cos(w0/2)/2;

0 0 -sin(w0) -cos(w0)];

e = eig(A, B);

[dum,jj]=min(abs(e));

epsi=e(jj)

w0=w0+epsi

end

[w0 k]

Answer 1

[Matt J]

 wsol=fsolve(@fun4evp,0.1),

[Matt J]

As for "Newton's iterations", I don't see how

w=w0+min(abs(beta));

is a meaningful update. For one thing, it doesn't utilize wold. If you really meant to have

w=wold+min(abs(beta))

it still doesn't look right, because w will always be increasing with iterations. What if the algorithm needs to take a step backwards?

[Abhishek Gaikwad]

!!!That's a bad mistake.Thanks. What about the Newton's method I am not sure about setting the tolerance. Currently, it is equal to the eigenvalue of C*A matrix which never goes to q very low value so i am able to solve the equation but the solution is incorrect.

[Abhishek Gaikwad]

The idea was to use the taylor series expansion and expand the matrix in the vicinity of w0 giving rise to epsilon which I found equal to the eigenvalue of the equation (A-lamba*B).

Answer 2

[Christine Tobler]

Sorry for answering so long after the initial post. One option you could consider for solving a nonlinear eigenvalue problem is the Rational Krylov Toolbox for MATLAB: http://guettel.com/rktoolbox/examples/html/example_nlep.html.

(edit to fix the link)