Hello everyone,

I am solving a simple heat equation of variable u with the source term f.

At first, function f only depends on u so I can define f outside the time loops and it gives a good result. But then I want f to depend on time t and I try to define it inside the time loops. I found that even though I haven’t changed the function f, the result I obtained was different from the previous case.

Does anyone have an idea why I got this difference and how I can avoid it?

Here comes my code:

```
real r=0.5, D = 1., nuE = 0.08, muE = 0.05 ,muF = 0.1, b = 10 , K = 200, dt=0.01, Tf=5. ;
real r1 = 1, r2 = 5, c= 0.05, lbda = 500;
border C(z=0, 2*pi){x=10*cos(z); y=10*sin(z);}
// The triangulated domain Th is on the left side of its boundary
mesh Th = buildmesh(C(100));
fespace Vh(Th,P1);
Vh uh,vh,uh0,f;
//Function f outside the time loop
f = r*nuE*b*uh/(b*uh/K - muE - nuE)*(1-exp(-uh))- muF*uh;
func init = 200*(sqrt(x*x+y*y) >= 7);
macro Grad(u)[dx(u),dy(u)]//
problem chaleur(uh,vh) = int2d(Th)(uh*vh/dt + Grad(uh)'*Grad(vh)*D) - int2d(Th, optimize = 0)(uh0*vh/dt + f*vh);
// Time loop
real t = 0;
uh0 = init;
for (int m = 0; m <= Tf/dt; m++){
// Update
t = t+dt;
// Function f inside the time loop
//f = r*nuE*b*uh/(b*uh/K - muE - nuE)*(1-exp(-uh))- muF*uh;
// Solve
chaleur;
uh0 = uh;
// Plot
if (!(m % 10)) {
plot(uh, cmm="t="+(t-dt)+"[days]", value = true, fill = true, wait=1);
}
}
```

Thank you very much!