simple for loop returns NaN somehow

Hello,

I'm trying to make a simple nested for loop where 2 for loops are run 100-10000 times. The first 4 or so iterations it runs just fine, but after that it returns a NaN on one of the variables and that of course ruins the whole thing and turns everything else into NaN.

Relevant code:

for k = 1:STEPS

for i = 1:FE

      REACTAB(i)=STEP*k01*exp(-Ea1/(R*TEMP1(i)))*COMPA(i);       % reaction per step
      REACTBC(i)=STEP*k02*exp(-Ea2/(R*TEMP1(i)))*COMPB(i);
      COMPA(i)=COMPA(i)-REACTAB(i);                               % composition change BEFORE moving
      COMPB(i)=COMPB(i)+REACTAB(i)-REACTBC(i);
      COMPC(i)=COMPC(i)+REACTBC(i);

      QFLOW23(i)=(TEMP2(i)-TEMP3(i))*STEP*SUR23*(HCtemporary);
      TEMP1(i)=TEMP1(i)+(REACTAB(i)*DH1+REACTBC(i)*DH2)/CAP1;
      QFLOW12(i)=(TEMP1(i)-TEMP2(i))*STEP*SUR12*(lambdaw/w);
      TEMP2(i)=TEMP2(i)+(QFLOW12(i)/CAP2)-(QFLOW23(i)/CAP2);
      TEMP1(i)=TEMP1(i)+(QFLOW12(i)/CAP1);
      TEMP3(i)=TEMP3(i)+(QFLOW23(i)/CAP3);

end

for i = FE:-1:1
    COMPA(i+1) = COMPA(i);  % move actual mass/composition one step
    COMPB(i+1) = COMPB(i);
    COMPC(i+1) = COMPC(i);
    TEMP1(i+1) = TEMP1(i);  % move temperature one step
    % TEMP2 is static
    TEMP3(i+1) = TEMP3(i);  % move temperature one step
end
    % new flow input
    COMPA(1)=FEED;               %
    COMPB(1)=0;
    COMPC(1)=0;
    TEMP1(1)=Treactinlet;        %
    TEMP3(FE:FE+10)=Tcoolinlet;        %

end

Using "dbstop if naninf" it stopped at the first TEMP1(i) iteration, so specifically:

      TEMP1(i)=TEMP1(i)+(REACTAB(i)*DH1+REACTBC(i)*DH2)/CAP1;

returns NaN. What I understand from NaN in Matlab it can only come from divide or multiply by zero or infinity and some more rarer cases, yet in this code, DH1 is a constant, as is DH2 and CAP1. The other arrays also start at nonzero values and I don't see how other parts of the code could introduce a zero or inf anywhere.

ALL of the variables in the code above are non-zero, real numbers. I've got no idea how a NaN could be produced.

Any help would be greatly appreciated!

Cheers

I figured it out. I had been uncarefull with some subtraction and adding which ,combined with the exponential, resulted in some -inf and +inf variables after enough iterations. Thanks anyway!

[Walter Roberson]

Perhaps you could save() the variables into a .mat file to allow us to test ?

[Matt J]

Using "dbstop if naninf" it stopped at the first TEMP1(i) iteration, so specifically:

When it stops, what are the values on the right hand side of

TEMP1(i)=TEMP1(i)+(REACTAB(i)*DH1+REACTBC(i)*DH2)/CAP1

[Simon Schouten]

I figured it out; turned out it was just a casual mistake, but thanks anyway!

[the cyclist]

I am glad it was not a formal mistake.