Effacer les filtres
Effacer les filtres

how plot amplitude of mode from eigenvector?

3 vues (au cours des 30 derniers jours)
Habib
Habib le 22 Nov 2015
Hi every body.
For solving an integral (this integral is Huygens Fresnel integral and it is useful for find modes of a laser resonator), I used Gauss-Legendre quadrature method and in final I got its eigenvectors and eigenvalues. but I couldn't plot the mode's amplitude and phase. If anyone has a experience about it, I hope help?
these codes are what is that I wrote.
clc
clear all
close all
a=0;
b=1;
N=200;
N=N-1;
N1=N+1; N2=N+2;
xu=linspace(a,b,N1)';
% Initial guess
y=cos((2*(0:N)'+1)*pi/(2*N+2))+(0.27/N1)*sin(pi*xu*N/N2);
% Legendre-Gauss Vandermonde Matrix
L=zeros(N1,N2);
% Derivative of LGVM
Lp=zeros(N1,N2);
% Compute the zeros of the N+1 Legendre Polynomial
% using the recursion relation and the Newton-Raphson method
y0=2;
% Iterate until new points are uniformly within epsilon of old points
while max(abs(y-y0))>eps
L(:,1)=1;
Lp(:,1)=0;
L(:,2)=y;
Lp(:,2)=1;
for k=2:N1
L(:,k+1)=( (2*k-1)*y.*L(:,k)-(k-1)*L(:,k-1) )/k;
end
Lp=(N2)*( L(:,N1)-y.*L(:,N2) )./(1-y.^2);
y0=y;
y=y0-L(:,N2)./Lp;
end
% Linear map from[-1,1] to [a,b]
x=(a*(1-y)+b*(1+y))/2;
% Compute the weights
w=(b-a)./((1-y.^2).*Lp.^2)*(N2/N1)^2;
LL=40.92;
A=1; B=LL; C=0;
n=N+1; m=n; s=n;
landa=10.6e-4;
k=2*pi/landa;
r1=linspace(a,b,n);
r2=linspace(a,b,m);
kernel1=zeros(n,m);
kernel2=zeros(n,m);
H1=zeros(n,m);
H2=zeros(n,m);
T1=zeros(1,m);
T2=zeros(1,m);
R=100*LL;
teta0=12.22e-3*pi/360;
T1(:)=exp(-2*1i*k*teta0*r1(:));
T2(:)=exp(-1i*k*r1(:).^2./R);
for s=1:n
kernel1(:,s)=-1i*(k/B)*r2(s).*besselj(0,r2(s).*r1(:)*k/B).*exp(1i*k*(A*r2(s).^2+B*r1(:).^2)/B/2);
H1(:,s)=kernel1(:,s).*T1(:);
end
for s=1:m
kernel2(:,s)=-1i*(k/B)*r1(s).*besselj(0,r1(s).*r2(:)*k/B).*exp(1i*k*(A*r1(s).^2+B*r2(:).^2)/B/2);
H2(:,s)=kernel2(:,s).*T2(:);
end
W=diag(w);
P1=H1*W;
P2=H2*W;
HH=P1*P2;
[V,D]=eig(HH);
I want to plot below pictures:

Réponses (0)

Catégories

En savoir plus sur Interpolation dans Help Center et File Exchange

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by