Error in Untitled (line 47) axis([newx(1) newx(end) -0.1 1.1]) Please help me, appreciate that!

% This program implements an explicit scheme to solve the Burgers−Huxley
% eqution, with the addition of shifting the profile. The code ... outputs the
% travelling wave profile, with its corresponding speeds with respect to
% variation in time.
clear all; clc;
% Spacial increment.
dx = 0.1;
% Time increment. 
dt = (0.95)*(1/2)*dx^2;
% Initial spacial domain.
x = [-20:dx:20];
% Time frame for simulation.
T = 600;
% Initial conditions.
u = (1-tanh(x))/2;
% Parameters k,m and n.
pk = 3; pm = 3; pn = 3;
% Defining initail states of paramters that will be used in the program.
iter = 1;
cumshift = 0;
xoffset = 0;
xpos = 0;
tm = 0;
% Implementing a shift every 300th iteration.
for j = 0:dt:T
% Applying a numerical shift every .
if (mod(iter,300)==0)
figure(1)
% Returns indices of the vector "u" that agree with the condition.
ind = find(u>0.1 & u<0.9);
% Returns the x value that corresponds to the u value at 0.5.
xc = interp1(u(ind),x(ind),0.5);
% Determining the number of steps xc is away from 0.
N = floor(xc/dx);
% Creating shift variable.
xshift = N*dx;
% Cumulated shift.
cumshift = cumshift + xshift;
% Rounding error.
xoffset = xc - xshift;
% Creating shifted u.
%u = [u(N+1:end) 1,N];
% Shifting x domain.
newx = x+cumshift;
% Plotting shifting wave solution.
plot(newx,u,'b');
axis([newx(1) newx(end) -0.1 1.1])
pause(0.01)
% Adding the total shift to a vector.
xpos = [xpos cumshift+xoffset];
% Adding the accumulated time frames between shits.
tm = [tm j];
end
% Explicit scheme.
u(2:end-1) = u(2:end-1) + ... (dt/(dx)ˆ2).*(u(3:end)−2.*u(2:end−1)+u(1:end−2)) ...
- (dt/(2*dx)).*u(2:end-1).^pk.*(u(3:end) - u(1:end-2)) ... 
+ dt.*u(2:end-1).^pm.*(1-u(2:end-1).^pn);
% Boundary conditions.
u(1) = 1;
u(length(x)) = 0;
iter = iter+1;
end
% Defining parameter values.
sum = 0;
k = 0;
speed(1) = 0;
% Computing the speeds with respect to time.
for j = 2:length(xpos)
    g = (xpos(j)-xpos(j-1))/(tm(j)-tm(j-1));
    sum = sum + g;
k = k+1;
speed(j) = sum/k;
end
figure (4)
% Plotting a Speed vs Time graph.
plot(tm,speed)
fprintf('final speed = %1.9f\n', speed(end));

Answer 1

[Shubhankar Poundrik]

Hi Jin,

I understand that you are trying to set axis limits while plotting some values, and are getting an error in the process.

[newx(1) newx(end) -0.1 1.1]

By adding the above line of code just above the line on which the error occurs(line 47) and observing the output it is clear that at one point in the loop the values of news(1) and news(end) become NaN. This leads to an error as they must be numeris and the second value must be greater than the first.

The code has to be modified to make sure that the values of news(1) and news(end) do not become NaN and are numeric.

Regards,

Shubhankar

[Jin Hao Yen]

Thanks for your time. Any good suggestions on modifying the values of newx(1) and newx(end) ? Or I just modify it with trial and error?

Answer 2

[Walter Roberson]

% Returns indices of the vector "u" that agree with the condition.
ind = find(u>0.1 & u<0.9);

3 of your points satisfy that condition.

% Returns the x value that corresponds to the u value at 0.5.
xc = interp1(u(ind),x(ind),0.5);

In order to be able to interp1() at position 0.5, then u(ind) must include at least one number >= 0.5 and one number <= 0.5 . But it doesn't -- the three numbers found are all > 0.5. So interp1() cannot interpolate at 0.5 and so returns nan. That makes lots of your other variables nan; in particular it makes all your newx nan, and then you try to axes([nan nan -0.1 1.1]) which is not permitted.

The problem does not (immediately?) occur if you change your dx from 0.1 to 0.01

[Jin Hao Yen]

Hi Walter, may I know which three numbers are found >0.5? Sorry, I am quite new to this, haven't learn it before.