Why are values in my loop filled with NaN?

8 Ansichten (letzte 30 Tage)
Westin Messer
Westin Messer am 27 Nov. 2018
Bearbeitet: Stephen23 am 27 Nov. 2018
Hello,
Can anyone tell why 'u' and 'v' inside my loop are all being filled with NaN? I cant figure out what I'm doing wrong. Thanks!
close all;
clear all;
a = 0;
b = 20;
nel = 200; % number of elements
h = (b-a)/nel; % step size
nv = nel+1;% number of vertices
x = a:h:b; % mesh
dt = 0.001; % time steps
tend = 20;
e=ones(nv, 1);
%Discretization of the Laplace operator
Delta=1/(h^2)*spdiags([e -2*e e], -1:1, nv, nv);
Delta(1,2) = 2/h^2;
Delta(end,end-1) = -2/h^2;
% parameters of the model
eps = 0.02;
alpha = 0.1;
% initial guess
u = 0.5 * (x < 3);
v = 0 * x;
counter = 1;
% explicit Euler scheme
for t=dt:dt:tend
u = u + dt/eps*(u .*(1-u).*(u-alpha) - v - ((Delta * u')'));
v = v + dt * (u - v);
end
figcure(1)
plot(x,u,x,v,'r');
legend('u','v')
title('Solution of the Fitz-Hugh-Nagumo system')

Akzeptierte Antwort

Stephen23
Stephen23 am 27 Nov. 2018
Bearbeitet: Stephen23 am 27 Nov. 2018
Actually this has nothing to do with division by zero. When you read the NaN help it lists several ways that NaNs can occur:
"These operations produce NaN:"
  • "Any arithmetic operation on a NaN, such as sqrt(NaN)"
  • "Addition or subtraction, such as magnitude subtraction of infinities as (+Inf)+(-Inf)"
  • "Multiplication, such as 0*Inf"
  • "Division, such as 0/0 and Inf/Inf"
  • "Remainder, such as rem(x,y) where y is zero or x is infinity"
In your case your code has increasingly large values in both u and v until on the eighth iteration they include several Inf elements (which critically have the same locations in both arrays):
>> any(any(isinf(u)==isinf(v))) % colocated Inf values.
ans = 1
>> any(isnan(u(:)))
ans = 0
>> any(isnan(v(:)))
ans = 0
Thus on the ninth iteration you calculate this: v - ((Delta * u')'), which is thus a subtraction of infinities and (exactly as the documentation states) produces NaNs:
>> tmp = v - (Delta * u')';
>> any(isnan(tmp(:)))
ans = 1

Weitere Antworten (0)

Kategorien

Mehr zu Loops and Conditional Statements finden Sie in Help Center und File Exchange

Community Treasure Hunt

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

Start Hunting!

Translated by