Out of memory while using FFT

12 visualizaciones (últimos 30 días)
Sreenidhi Yeturi
Sreenidhi Yeturi el 18 de Abr. de 2021
Respondida: Madhav Thakker el 18 de Mayo de 2021
I'm trying to plot the field for a FDTD simulation in frequncy domain but I keep running into memory error. I've attached my code below.
Nfft = nt *8;
df = 1/(Nfft*dt);
ftEZ = abs(fft(Ez_out,Nfft));
fmax = 600e6;
numfSamps = floor(fmax/df);
fr = df*(0:numfSamps-1);
figure;
plot(fr,ftEz(1:numfSamps));
ylabel('Amplitude (V/m)');
xlabel('Frequency (Hz)');
set(gca,'FontSize',16);
hold on;
p = ftEz(1:numfSamps);
[peaks,locs] = findpeaks(p,'MinPeakDistance',width, 'MinPeakHeight',height);
f110 = 2.998e8/2/pi*sqrt(2*pi^2);
f210 = 2.998e8/2/pi*sqrt((2*pi^2)+pi^2);
nf110 = floor(f110/df);
nf210 = floor(f210/df);
fex = [f110,f210];
aex = [ftEz(nf110), ftEz(nf210)];
plot(fex,aex,'ro');
legend('FDTD sim','Exact','location','northwest');
%%%%%Error Message%%%%%
Error using fft
Requested array exceeds the maximum possible variable size.
Error in Trail_code (line 179)
ftEZ = abs(fft(Ez_out,Nfft));
More information
  2 comentarios
Matt J
Matt J el 18 de Abr. de 2021
We don't know what the value of nt is, but isn't the cause of the problem obvious from the erorr message? Your Nfft must be absurdly big.
Sreenidhi Yeturi
Sreenidhi Yeturi el 19 de Abr. de 2021
Thanks for commenting Matt, here's the complete code.
scalefactor=11; % This scales the resolution,
c = 2.99792458e+8;
no = 376.7303134617706554679;
mu = no/c;
e0 = 1/(no*c);
epsr = 4.2;
cfln = 0.995; % Stability criterion
tw = 1e-9; % half width of pulse 10 had low resonance
to = 5*tw; % Time delay for sources
Lx = 1;
Ly = 1;
Lz = 1;
nx = 11*scalefactor;
ny = 11*scalefactor;
nz = 11*scalefactor;
dx = Lx/nx;
dy = Ly/ny;
dz = Lz/nz;
S=5.8e1;
eps = e0;
dt = cfln/(c*sqrt((1/dx)^2+(1/dy)^2+(1/dz)^2));
T = 4.5e-9*2;
simtime = T/dt;
nt = floor(simtime / dt);
nMax = 50; % Time steps
% Injection and Output point ranges
% Source Injection points
i1_src = 4;
j1_src = 4;
k1_src = 4;
i2_src = 4;
j2_src = 4;
k2_src = 6;
% Field ouput points
i1_fld = 8;
j1_fld = 8;
k1_fld = 8;
i2_fld = 8;
j2_fld = 8;
k2_fld = 9;
% Electric & Magnetic Field Coefficients
% Electric Field
cExy = dt/(dy*eps);
cExz = dt/(dz*eps);
cEyz = dt/(dz*eps);
cEyx = dt/(dx*eps);
cEzx = dt/(dx*eps);
cEzy = dt/(dy*eps);
% Magnetic Field
cHxy = dt/(mu*dy);
cHxz = dt/(mu*dz);
cHyz = dt/(mu*dz);
cHyx = dt/(mu*dx);
cHzx = dt/(mu*dx);
cHzy = dt/(mu*dy);
% PEC Boundary conditoins
% Set all tangential E fields to zero and never update them. This will keep
% them zero. Update loops only modify inner cells
Hx=zeros(nx,ny-1,nz-1);
Hy=zeros(nx-1,ny,nz-1);
Hz=zeros(nx-1,ny-1,nz);
Ex=zeros(nx-1,ny,nz);
Ey=zeros(nx,ny-1,nz);
Ez=zeros(nx,ny,nz-1);
% Ez_out = zeros(nt,1);
for n=1:nMax
%Main update loops for the Electric fields:
%Ex
for k = 2:nz-1
for j = 2:ny-1
for i = 1:nx-1 %1->2
Ex(i,j,k) = Ex(i,j,k)+cExy*((Hz(i,j,k)-Hz(i,j-1,k))-cExz*(Hy(i,j,k)-Hy(i,j,k-1)));
end
end
end
%Ey
for k = 2:nz-1
for j = 1:ny-1
for i = 2:nx-1
Ey(i,j,k) = Ey(i,j,k)+cEyz*((Hx(i,j,k)-Hx(i,j,k-1))-cEyx*(Hz(i,j,k)-Hz(i-1,j,k)));
end
end
end
%Ez
for k = 1:nz-1 %1->2
for j = 2:ny-1
for i = 2:nx-1
Ez(i,j,k) = Ez(i,j,k)+cEzx*((Hy(i,j,k)-Hy(i-1,j,k))-cEzy*(Hx(i,j,k)-Hx(i,j-1,k)));
end
end
end
%Source Injection
t = (n-0.5)*dt;
Jz =(-2/tw)*(t-to)*exp(-(t-to)^2/tw^2);
for k = k1_src:k2_src-1
for j = j1_src:j2_src
for i = i1_src:i2_src
Ez(i,j,k) = Ez(i,j,k)+dt*Jz;
end
end
end
%Main update loops for the Magnetic fields:
%Hx
for k = 1:nz-1 %1->2 all 3
for j = 1:ny-1
for i = 1:nx %-1
Hx(i,j,k) = Hx(i,j,k)-cHxy*(Ez(i,j+1,k)-Ez(i,j,k))+cHxz*(Ey(i,j,k+1)-Ey(i,j,k));
end
end
end
%Hy
for k = 1:nz-1 %1->2 all 3
for j = 1:ny %-1
for i = 1:nx-1
Hy(i,j,k) = Hy(i,j,k)-cHyz*(Ex(i,j,k+1)-Ex(i,j,k))+cHyx*(Ez(i+1,j,k)-Ez(i,j,k));
end
end
end
%Hz
for k = 1:nz %-1 %1->2 all 3
for j = 1:ny-1
for i = 1:nx-1
Hz(i,j,k) = Hz(i,j,k)-cHzx*(Ey(i+1,j,k)-Ey(i,j,k))+cHzy*(Ex(i,j+1,k)-Ex(i,j,k));
end
end
end
%%%%%%%%%% Output %%%%%%%%%%%%
Ez_out = zeros(n);
for k = k1_fld:k2_fld-1
for j = j1_fld:j2_fld
for i = i1_fld:i2_fld
Ez_out(n) = Ez_out(n) + Ez(i,j,k);
end
end
end
clc
round(n/nMax*100) %gives progress update percentage on home screen
%surf(squeeze(Ez(:,:,0.1*scalefactor))); %Z-slice with resonator in center
surf(squeeze(Ex(1*scalefactor,:,:)));
%surf(squeeze(Ey(:,0.1*scalefactor,:)));
set(gcf,'renderer','opengl')
%hold on
%whitebg('black');
grid off
set(gcf,'Position',[50 50 1500 800]);
axis([-1*scalefactor 5*scalefactor -2*scalefactor 7*scalefactor -1/scalefactor 1/scalefactor])
view([3,4,5])
Mov(n)=getframe; % comment this out for nonanimated analysis
end
movie(Mov,30,30); % comment this out for nonanimated analysis
%Post Process
Nfft = nt *8;
df = 1/(Nfft*dt);
ftEZ = abs(fft(Ez_out,Nfft));
fmax = 600e6;
numfSamps = floor(fmax/df);
fr = df*(0:numfSamps-1);
figure;
plot(fr,ftEz(1:numfSamps));
ylabel('Amplitude (V/m)');
xlabel('Frequency (Hz)');
set(gca,'FontSize',16);
hold on;
p = ftEz(1:numfSamps);
[peaks,locs] = findpeaks(p,'MinPeakDistance',width, 'MinPeakHeight',height);
f110 = 2.998e8/2/pi*sqrt(2*pi^2);
f210 = 2.998e8/2/pi*sqrt((2*pi^2)+pi^2);
nf110 = floor(f110/df);
nf210 = floor(f210/df);
fex = [f110,f210];
aex = [ftEz(nf110), ftEz(nf210)];
plot(fex,aex,'ro');
legend('FDTD sim','Exact','location','northwest');

Iniciar sesión para comentar.

Respuestas (1)

Madhav Thakker
Madhav Thakker el 18 de Mayo de 2021
Hi Sreenidhi,
The variable nt has value 3.5886e+13 which is due to which the resulting Fourier transform array exceeds the maxiumum variable size.
Hope this helps.

Categorías

Más información sobre Fourier Analysis and Filtering en Help Center y File Exchange.

Etiquetas

Community Treasure Hunt

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

Start Hunting!

Translated by