MATLAB is giving me 0's instead of really small numbers

I don't have code to share, but MATLAB is outputing 0's instead carrying through small numbers. When I change the command window to format long, the calculation show non zero. However, when my code runs it is carrying these zeros through my calculations.

>> Calculation

ans = 

        0.0
        
>> Format long
>> Calculation

ans = 
        
        8.425293911978570e-07

How do I make matlab keep the small numbers?

[Steven Lord]

Can you show us a small segment of the code that appears to return 0's instead of the small numbers you expected?

How small is "really small" in this context? Smaller than realmin?

[Nicholas Morris]

See the comments in the answer below.

The numbers are definitely not smaller than realmin. and I'm not really worried about the fact that my command window is displaying "0" for all values in my matrix, but more the fact that if I use the following:

>>find(array>0)

ans =

  0×1 empty double column vector

so these numbers are actually zero.

I have a similar array in my code with small numbers that you can see display properly

Answer 1

[Nicholas Morris]

Problem solved...Indexing errors.

Answer 2

[Walter Roberson]

MATLAB is keeping the small numbers throughout your calculation.

I recommend that you put

format long g

at the top of your code and run again.

[Nicholas Morris]

I have tried all formats. It still giving me zeros.

When I run my code it gives me an array of all zeros. If I then put format long in the command window and run the line for the i-th step using F9 it will give me a small numberin the right spot.

Why doesn't it do this when I run the code normally?

[Walter Roberson]

Your code might be zeroing the entire array at each step.

But as well, "format long g" does not affect the variable browser. There is a drop-down above the variable name that can select the format for the variable browser; personally though I recommend setting Preferences -> Variables -> Default Array Format to "long g", which controls the default format for the variable browser.

[Nicholas Morris]

I've already tried the "Preferences -> Variables -> Default Array Format to "long g"... : (

No luck.

What do you mean by my code "might be zeroing the entire array at each step"?? Why would it do this?

Answer 3

[vikas singh]

I have written a code and it is giving me H3 and z2 zero. i dont know what is the problem.

clc;
clear all;
close all;
l=10000;
x=0:100:l;
t=[0;1;5;10;15;20;25];

H=30;
alpha=0.00384
Mu=2.8e-1;
D=31;
k=0.01132;
k1=l/H;

for n=1:100
    a=2*H;
    b=(1-(-1)^n*exp(-k*l))/((n*pi+(k^2*l^2)/(n*pi)));
    c=1/(n*pi);
    Cn(n)=a*(b-c);
end
s=0;
for n=1:100
    d= sin(n*pi*x/l);
    q= exp((-n^2*pi^2*t)/(Mu^2*l));
    aa=Cn(n).*d.*q;
    s=s+aa;
    bb=(2*alpha*Mu^2*l^2)/(n*pi)^3;
    s2=bb.*(1-q).*d;
    
end
J=(H+alpha*t).*x/l;
J1=alpha*t;
z1=s-s2-J+J1;


for n=1:100
    syms X T  
    disp(n)
    H0=(Cn(n)/(D-H)).*sin(n*pi*X).*exp(-n^2*pi^2*T/Mu)+(bb/(D-H)).*(1-exp(-n^2*pi^2*T/Mu)).*sin(n*pi*X) +(H/(D-H)+alpha*T*Mu*l/(D-H))*X+(H-alpha*T*Mu*l)/(D-H);
    y=diff(H0*diff(H0,X),X);
    Q=int(y*sin(n*pi*X),X,0,1);
    P=int(Q*exp(n^2*pi^2*T/Mu),T,0,T);
    Hn=-2*P/(exp(n^2*pi^2*T/Mu));
    H1=Hn.*sin(n*pi*X);
    for i=1:size(t,1)
        H2=subs(H1,T,t/(Mu*l));
    end
    format long g
    H3=subs(H2,X,x/l)
    
end

z2=double(((D-H)^2).*H3/l)

[Walter Roberson]

Well, I can get you non-zeros, but it is questionable as to whether the result is what you want.

format long g

l=10000;
x=0:100:l;
t=[0;1;5;10;15;20;25];

H=30;
alpha=0.00384
Mu=2.8e-1;
D=31;
k=0.01132;
k1=l/H;

for n=1:100
    a=2*H;
    b=(1-(-1)^n*exp(-k*l))/((n*pi+(k^2*l^2)/(n*pi)));
    c=1/(n*pi);
    Cn(n)=a*(b-c);
end
s=0;
for n=1:100
    d= sin(n*pi*x/l);
    q= exp((-n^2*pi^2*t)/(Mu^2*l));
    aa=Cn(n).*d.*q;
    s=s+aa;
    bb=(2*alpha*Mu^2*l^2)/(n*pi)^3;
    s2=bb.*(1-q).*d;
    
end
J=(H+alpha*t).*x/l;
J1=alpha*t;
z1=s-s2-J+J1;

syms X T 
for n=1:100
    disp(n)
    H0=(Cn(n)/(D-H)).*sin(n*pi*X).*exp(-n^2*pi^2*T/Mu)+(bb/(D-H)).*(1-exp(-n^2*pi^2*T/Mu)).*sin(n*pi*X) +(H/(D-H)+alpha*T*Mu*l/(D-H))*X+(H-alpha*T*Mu*l)/(D-H);
    y=diff(H0*diff(H0,X),X);
    Q=int(y*sin(n*pi*X),X,0,1);
    P=int(Q*exp(n^2*pi^2*T/Mu),T,0,T);
    Hn=-2*P/(exp(n^2*pi^2*T/Mu));
    H1=Hn.*sin(n*pi*X);
    for i=1:size(t,1)
        H2(i,1)=subs(H1,T,t(i)/(Mu*l));
    end
    H3(:,:,n)=subs(H2,X,x/l);
end

z2=double(((D-H)^2).*H3/l);

z2 is now length(t) (=7) by length(x) (=101) by 100

I suspect that the "for n" loop could be optimized somewhat, but I would need to test that later.

[vikas singh]

I want to get z=(z1-z2) and then plot between (x,z).

[vikas singh]

and z1 is coming length(t)(=7)by length(x)(=101). so we have to change the same size in z2 to get by adding all those 100 matrix. please tell me here how to add those 100 matrix.

.

[Walter Roberson]

sum(z2, 3)

would add along the third dimension.

[vikas singh]

Thanks a lot walter sir

[vikas singh]

I'm getting H1 value in the following equation 0. I'm not getting where it is wrong.

clc;
close all;
clear all;
format long g
t=365;
k1=0.025;
l=1000;
k=0.0864;
x=0:0.1:l;
epc=0.3;
H=30;
a1=0.0138*l/(k*k1*(1+k1));
mu=epc/(k*k1*H)
for n=1:100
    a2(n)=(1-cos(1000*n*pi/l))/(n^2*pi^2);
end
for n=1:100
    a3=(1-exp(-n^2*pi^2*t/(l^2*mu))).*sin(n*pi*x/l);
end
H1=a1.*a2(n).*a3
plot(x,H1)

[Walter Roberson]

you are overwriting all of a3 in the last loop instead of doing indexed assignments like in the other loop

[vikas singh]

how to correct it?

[Walter Roberson]

a3(n)=(1-exp(-n^2*pi^2*t/(l^2*mu))).*sin(n*pi*x/l);

Same technique as you used in the loop calculating a2.

However:

H1=a1.*a2(n).*a3

That line is after all of the loops, so n is going to have the last value it was assigned, which in this case would be the 100 from the for n=1:100 loop. So that line of code is equivalent to

H1=a1.*a2(100).*a3

Are you sure that is what you want to do? If so what is the point of calculating the other entries in a2 ?

[vikas singh]

I need all entries for a2 to multiply one by one and get 100 values of H1.

[vikas singh]

I'm using the following code to get but H2 coming zero

clc;
close all;
clear all;
format long g
t=365;
k1=0.025;
l=1000;
k=0.0864;
x=0:0.1:l;
epc=0.3;
H=30;
a1=0.0138*l/(k*k1*(1+k1));
mu=epc/(k*k1*H)

mu = 
          4.62962962962963

H1=0

H1 = 
     0

for n=1:100
    a2(n)=(1-cos(1000*n*pi/l))/(n^2*pi^2);
    a3=(1-exp(-n^2*pi^2*t/(l^2*mu))).*sin(n*pi*x/l);
    H2=H1+a1.*a2(n).*a3;
end

plot(x,H2)

[Walter Roberson]

You initialize H1 to 0, and you do not change it after that. That makes it pointless to calculate H1 + something -- it would be adding 0 each time, which you might as well skip doing.

You are overwriting all of H2 every time through your for n=1:100 loop.

I wonder if what you mean is

clc;
close all;
clear all;
format long g
t=365;
k1=0.025;
l=1000;
k=0.0864;
x=0:0.1:l;
epc=0.3;
H=30;
a1=0.0138*l/(k*k1*(1+k1));
mu=epc/(k*k1*H)

mu = 
          4.62962962962963

H1=0;
for n=1:100
    a2(n)=(1-cos(1000*n*pi/l))/(n^2*pi^2);
    a3=(1-exp(-n^2*pi^2*t/(l^2*mu))).*sin(n*pi*x/l);
    H1=H1+a1.*a2(n).*a3;
end

plot(x,H1)

[vikas singh]

Sir the graph should be exactly this only.

but why it is not coming when we take

H2=H1+a1.*a2(n).*a3;

[vikas singh]

thank you sir for helping. I'm first time using matlab, so sometimes i miss the things.

[Walter Roberson]

Your H1 was initialized to 0. You had

H2=H1+a1.*a2(n).*a3;

with H1 being 0, that becomes equivalent to

H2=a1.*a2(n).*a3;

which overwrites all of H2 in each iteration.

Now suppose in the last iteration, H2 comes out as 0. The values in the other iterations would not matter, because you do not save the H2 values for the other values of n, and if the last entry came out 0, then after the loop, H2 would be only 0.

[vikas singh]

Sir I have a following code. Now I want the integration of Q*sin(n*pi*x/l) over 0 to l. how to write it in another for loop.

clc;
close all;
clear all;
format long g
H=50;
k=10^-4;           %in m/sec
k1=0.025;          %delta value
epc=0.3            %epcilon value
l=600;
x=0:1:l;
nu=k*k1*H/epc      % nu bar value
omg1=1*10^-4;      % omega1 recharge value for x<=250 and x>=350
omg2=2*10^-4;      %omega2 recharge value for 250<x<350
for i=1:length(x)
    if i<=251
        omg(i)=omg1;        
    elseif i<=351
        omg(i)=omg2;
    elseif i<=601
        omg(i)=omg1;         
    end
    Q(i)=2*omg(i)*nu/(k1*(k1+1)*k);
end

[vikas singh]

I have added the following code to get it but it is giving me zero.

syms x
for n=1:100
    q1=2/l;
    q2=sin(n*pi*x);
    q3(i)=Q(i).*q2;
    q4=int(q3,x,0,l)
    
end

[Walter Roberson]

q2=sin(n*pi*x);

n is an integer. So you have sin() of an integer multiple of pi times x. You are integrating that over 0 to 600. But that is an integer multiple of 2 pi, and the integral of sine over an integer multiple of 2 pi is 0.

[vikas singh]

sorry

q2=sin(n*pi*x/l). now also it is giving zero. the entire code is given following.

clc;
close all;
clear all;
format long g
H=50;
k=10^-4;           %in m/sec
k1=0.025;          %delta value
epc=0.3            %epcilon value
l=600;
x=0:1:l;
nu=k*k1*H/epc      %nu bar value
omg1=1*10^-4;      %omega1 recharge value for x<=250 and x>=750
omg2=2*10^-4;      %omega2 recharge value for 250<x<750
for i=1:length(x)
    if i<=251
        omg(i)=omg1;        
    elseif i<=351
        omg(i)=omg2;
    elseif i<=601
        omg(i)=omg1;         
    end
    Q(i)=2*omg(i)*nu/(k1*(k1+1)*k);
end
syms x
for n=1:100
    q1=2/l;
    q2=sin(n*pi*x/l);
    q3(i)=Q(i).*q2;
    q4=int(q3,x,0,l)
    
end

[Walter Roberson]

format long g

S = @(v) sym(v);

H = S(50);
k = S(10)^-4;           %in m/sec
k1 = S(25)/1000;        %delta value
epc = S(3)/10;          %epcilon value
l = S(600);
x=0:1:l;
nu=k*k1*H/epc      %nu bar value

nu = 

omg1 = S(1)/10^4;      %omega1 recharge value for x<=250 and x>=750
omg2 = S(2)/10^4;      %omega2 recharge value for 250<x<750

syms x
syms n integer
omg = piecewise(250 <x & x < 750, omg2, omg1);
Q = 2 * omg * nu/(k1*(k1+1)*k);

q1 = 2/l;
q2 = sin(n*sym(pi)*x/l);
q3 = Q*q2;
q4 = int(q3, x, 0, l);
disp(q4)

q4s = simplify(q4, 'steps', 50)

q4s = 

This would not always be 0, but when n is a multiple of 24 then you would have cos() of a multiple of 5 pi which would give you +/- 1, and square of that would be 1, and (-1) to a multiple of 24 is going to be 1, so you are going to have 1600 - 1600 which would be 0.

Notice that you do not save the individual q4 results in your code, so you are only going to see the result for n = 600 and 600 happens to be a multiple of 24...

[vikas singh]

Sir I need the value of q4 for each value of n=1:100.

[Walter Roberson]

Take the output of q4s in the code I show and do

output = subs(q4s, n, 1:100);