Hi,

I've been having a problem with a function definition. Here is all my code.

vm = 5;
r0 = 127.25;
di = r0*2;
de = 273;
v = 0.000000464;
?f = 958;
?t = 7854;
pr = 2.88;
kf = 0.656;
kt = 58;
cf = 4186;
ct = 470;
re = (vm*(2*r0/1000)/v)
e = 0.00015;
f = 1.325 (Log[0.27*(e/(2*r0/1000)) + 5.74*(1/re)^0.9])^(-2)
tinf = 298;
tme = 373;
?f = kf/(?f*cf)
c1 = (-f*?f*vm^2/(4*r0/1000))
g = 9.81;
? = 1/((tme + tinf)/2)
var = 0.000018;
prar = 0.703;
? = 0.000025;
kar = 0.028;

ym[r_] := (((r0/1000) - r)*vm*(f/8)^(0.5))/v;
?m[r_] := v*(0.4*ym[r]/6)*(1 + (r/(r0/1000)))*(1 + 2*(r/(r0/1000))^2);

momento = NDSolve[{(1/r)*D[(v + ?m[r])*r*D[u[r], r], r] == 1/Rho]f)*c1, u'[0.0000001] == 0, u[(r0/1000)] == 0},u[r], {r, -0.0000001, (r0/1000)}, MaxSteps -> 10000];
mom[r_] = u[r] /. First[momento];

hi = (kf/(2*r0/1000))*((re - 1000)*pr*(f/8))/(1 + 12.7*(f/8)^0.5*(pr^(2/3) - 1))
ra = (g*?*((tme + tinf)/2 - tinf)*(de/1000)^3)/(var*?);
he = (kar/(de/1000))*(0.6 + (0.387*ra^(1/6))/(1 + (0.559/prar)^(9/16))^(8/27))^2

Here starts my problem. Mathematica isn't able to solve systems of PDEs with more than one independent variable. But Mathematica is able to solve each one individually. So, I have to solve iteratively both the two equation varying a boundary condition (by the way, it is a partial elliptic equation). But, at the moment this is not the problem.

l = 200;
tempo = 100;
Clear[tmed]; 
Clear[tmedpred]; 
Clear[ttf]; 
Clear[med];
Clear[tpred]; 
tmed[x_, t_] := tinf; (*Initial*)
med[x_, t_] := tinf; (*Initial*)
tmedpred[x_, t_] := tinf; (*Initial*)
tpred[x_, t_] := tinf; (*Initial*)
ttf[x_, r_, t_] := tinf; (*Initial*)
subrel = 0.4;
eps = 10^-4;

(*While[erro[GreaterEqual]eps,*)

enf1 = NDSolve[{(?f*cf*(D[tf[x, r, t], t] + mom[r]*D[tf[x, r, t], x])) == ((1/r)* D[(r*(?f + ?m[r])*D[tf[x, r, t], r]), r]),tf[x, r, 0] == (tme - tinf)*Exp[-1000*x] + tinf,tf[0, r, t] == tme, (D[tf[x, r, t], r] /. r -> 0.0000001) ==  0, -kf*(D[tf[x, r, t], r] /. r -> (r0/1000)) == (1 - Exp[-1000*x])*hi*(tmed[x, t] - ttf[x, r0/1000, t])},tf[x, r, t], {x, 0, l}, {r, 0.0000001, (r0/1000)}, {t, 0, tempo}, Method -> {"MethodOfLines", "SpatialDiscretization" -> {"TensorProductGrid", "MinPoints" -> 13}}];
tff[x_, r_, t_] = tf[x, r, t] /. First[enf1];

med[x_, t_] := (NIntegrate[(mom[r]*tff[x, r, t]*r), {r, 0.0000001, r0/1000}])*2/((2/((r0/1000)^2)*NIntegrate[(mom[r]*r), {r, 0.0000001, r0/1000}])*(r0/   1000)^2));

tpred[x_, t_] := tmedpred[x, t] + subrel*(med[x, t] - tmedpred[x, t]);
tmed[x_?NumericQ, t_?NumericQ] := tpred[x, t];

ent1 = NDSolve[{?t*ct*(D[tt[x, r, t], t]) == (1/r)*(D[(r*kt*D[tt[x, r, t], r]), r]), tt[x, r, 0] == tinf, -kt*(D[tt[x, r, t], r] /. r -> (r0/1000)) == (1 - Exp[-1000*x])*hi*(tmed[x, t] - ttf[x, r0/1000, t]), -kt*(D[tt[x, r, t], r] /. r -> (de/2000)) == (1 - Exp[-1000*t])*he*(ttf[x, de/2000, t] - tinf)}, tt[x, r, t], {x, 0, l}, {r, (r0/1000), (de/2000)}, {t, 0, tempo}];

tmedpred = tmed

Here is the problem. "tmedpred" shall return to the structure of the first NDSolve (enf1) and the code will keep the loop up to an error value. Realize that the While structure is not in the code. I'm doing a manual evaluation.

The error:

$RecursionLimit::reclim: Recursion depth of 1024 exceeded. >> 

The box message:

"A very large output was generated. Here is a sample of it: (<<1>>) + subrel (med[100, 50] - tmedpred[100, 50])

The problem starts with the input tmedpred = tmed

Could some one please help me with this issue? I've been trying to solve this for more than one day.