# How to deal with too much recursion

I have put together a very simple climate model based around four equations which define the state at time t based on time t-1, along with a initial state at time t = 0. The problem arises when I wish to evaluate, say, the temperature at year 2000 (t = 2000) or higher. I hit the recursion limit. I could keep increasing the \$RecursionLimit all the time, but I wonder if there is a better way to deal with the problem?

ClearAll[timeStep, waterDepth, gramsPerM3ofWater, joulesPerCalorie,
joulesToHeatWater, solarConstant, albedo, boltzmanConstant,
emissivity, temp, heatContent, incomingFlux, outgoingFlux, flux];

(*Constants*)
timeStep = 31536000;(*One year in seconds*)
waterDepth = 4000.;
gramsPerM3ofWater = 1000000.;
joulesPerCalorie = 4.186;
joulesToHeatWater = waterDepth*gramsPerM3ofWater*joulesPerCalorie;
solarConstant = 1350.;
albedo = 0.3;
incomingFlux = solarConstant (1 - albedo)/4;
boltzmanConstant = 5.67*^-8;
emissivity = 1.;

(*Initial State*)
temp[0] = 0.;
heatContent[0] = temp[0]*joulesToHeatWater;
outgoingFlux[0] = emissivity *boltzmanConstant*temp[0]^4;
flux[0] = (incomingFlux - outgoingFlux[0])*timeStep;

(*Equations*)
heatContent[t_] := heatContent[t - 1] + flux[t - 1]
temp[t_] := heatContent[t]/joulesToHeat
outgoingFlux[t_] := boltzmanConstant*temp[t]^4
flux[t_] := (incomingFlux - outgoingFlux[t])*timeStep

-
Recursion or iteration. That is the question. :) –  belisarius Nov 1 '13 at 19:10
@belisarius Iteration. That is the answer. [Hoping it's okay to answer a rhetorical question.] –  Daniel Lichtblau Nov 1 '13 at 20:29
@DanielLichtblau Recursion or iteration. That is the question. :) –  rm -rf Nov 1 '13 at 20:37
Apparently it's already devolved into reiteration. –  bill s Nov 1 '13 at 21:25
To understand recursion, you must first understand recursion. –  Simon Woods Nov 1 '13 at 22:57

Use Nest, to kill the recursion, as follows:

ClearAll[getNewValues];
getNewValues[{hc_, fl_, out_, temp_}] :=
Module[{newhc, newfl, newout, newtemp},
newhc = hc + fl;
newtemp = newhc/joulesToHeatWater;
newout = boltzmanConstant*newtemp^4;
newfl = (incomingFlux - newout)*timeStep;
{newhc, newfl, newout, newtemp}
];


Then, for example for the 5000 steps:

Nest[getNewValues, {heatContent[0], flux[0], outgoingFlux[0], temp[0]}, 5000] // AbsoluteTiming

(* {0.063185, {4.25409*10^12, 0.000237521, 236.25, 254.066}}  *)

-
This works nicely, thanks! –  Mr Alpha Nov 2 '13 at 17:48
@MrAlpha Was glad to help, and thanks for the accept. This is actually a quite nice application of Nest, I think, and one in which the use of Nest solves what seems to be a rather non-trivial problem in some other approaches. Definitely one of the examples to show the utility of Nest. –  Leonid Shifrin Nov 3 '13 at 1:55