admin管理员组

文章数量:1335357

I want to calculate the eigenvalues of the (discrete) Laplace operator in 3D. For this I want to use the matlab command eigs. I have implemented the Laplace operator as a function as it is recommended in the documentation of eigs. For this it is necessary to design the function, such that it takes one array (1D) as an input and gives one array (1D) as an output. Please note, that I know, that the spectrum of the (discrete) Laplacian is known analytically. My actual problem goes beyond and I tried to minimize my actual issue to the question you read in this post. In other words: I would like to use the function together with eigs in any case.

Here I define my function:

function Lap = Laplace(f,Lx,Ly,Lz,Nx,Ny,Nz);

kx(1:Nx/2)=(0:Nx/2-1)*2*pi/Lx;            %Calculate the momentum vector in x-direction
kx(Nx/2+1:Nx)=(-Nx/2:-1)*2*pi/Lx;

ky(1:Ny/2)=(0:Ny/2-1)*2*pi/Ly;            %Calculate the momentum vector in y-direction
ky(Ny/2+1:Ny)=(-Ny/2:-1)*2*pi/Ly;

kz(1:Nz/2)=(0:Nz/2-1)*2*pi/Lz;            %Calculate the momentum vector in z-direction
kz(Nz/2+1:Nz)=(-Nz/2:-1)*2*pi/Lz;

[Kx,Ky,Kz]=meshgrid(kx,ky,kz);            %Calculate a grid in the k-space

K=sqrt(Kx.*Kx+Ky.*Ky+Kz.*Kz);             %Calculate at every grid-point the absolute value of the k-vector

f3D = reshape(f,[Nx,Ny,Nz]);              %reshape the input vector in order to multiply grid-point by grid point with K

Lap3D = ifftn(K.*K.*fftn(f3D));           %Action of the Laplacian 

Lap = Lap3D(:);                           %Reshape back to 1D-format
end

And then I run the command:

Lx=1;    %length in x-direction
Ly=1;    %length in y-direction
Lz=1;    %length in z-direction
Nx=128;  %number of points in x-direction
Ny=128;  %number of points in y-direction
Nz=128;  %number of points in z-direction

Afun= @(x) Laplace(x,Lx,Ly,Lz,Nx,Ny,Nz);

tic
eigs(Afun,Nx*Ny*Nz,12,'smallestabs','Tolerance',1e-9 ,'Display',1,'SubspaceDimension',100)
toc

As I result I should get the 12 smallest eigenvalues, but when I run the code it produces the following output:


   1.0e-05 *

   0.206138475834832
   0.208291225315224
   0.208291225315224
   0.208291225315224
   0.208291225315225
   0.210454435947029
   0.210454435947030
   0.210454435947030
   0.210489412585876
   0.210489412585876
   0.210489412585877
   0.212627347524422

This result seems wrong. It looks like the eigenvalues are shifted. I do not see my mistake. Thanks for any kind of help.

本文标签: matlabCalculating Eigenvalues of 3Dfunctions with eigsStack Overflow