Iterating with NDSolveValue in a Do-loop [duplicate]

I am trying to get a differential solution to converge. There are two separate time periods that I am solving for and the entire solution is periodic. I want to use the output from one as the input for the next and so on and so forth. I have attempted to do this with a Do loop but it does not seem to be working. I don't know if the Do loop is set up to work this way or if I should use another approach. Please let me know if you have any suggestions on how to get this to work.

Here is the code:

<< "NDSolveFEM"

ts = .000250; tsl = .000250; (* m *)
ρs = 3980;  ρsl = 958; (* kg/m3 *)
ks = .035;  ksl = .00067; (* kW/m/K *)
cs = .75; csl = 4.22; (* kJ/kg/K *)

ρ = If[0 <= z < ts, ρs, ρsl];
k = If[0 <= z < ts, ks, ksl];
c =  If[0 <= z < ts, cs, csl];
q = 28;

td = 0.060; tg = 0.01;

eqn1 = k*\!$$\*SubscriptBox[\(∂$$, $$z$$]$$T1[z]$$\) + q;
Tbl = 100;
Subscript[Γ1, D] =
DirichletCondition[T1[z] == Tbl, z == ts + tsl];

BCr = NDSolveValue[{eqn1 == 0, Subscript[Γ1, D]},
T1, {z, 0, ts + tsl}];

Ti[z_] := \[Piecewise] {
{BCr[z], 0 <= z <= ts},
{100, True}
};

eqn2 = ρ*c*\!$$\*SubscriptBox[\(∂$$, $$t$$]$$T[t, z]$$\) - k*\!$$\*SubscriptBox[\(∂$$, $$z, z$$]$$T[t, z]$$\);
Subscript[Γ, D] =
DirichletCondition[T[t, z] == Tbl, z == ts + tsl];
Subscript[Γ, N] = NeumannValue[-q, z == 0];

solnd = NDSolveValue[{eqn2 == Subscript[Γ, N],
Subscript[Γ, D], T[0, z] == Ti[z]},
T, {t, 0, td}, {z, 0, ts + tsl},
Method -> {"PDEDiscretization" -> {"MethodOfLines",
"SpatialDiscretization" -> {"FiniteElement",
"MeshOptions" -> {"MaxCellMeasure" -> (ts + tsl)/1011}}}}];

eqn3 = ρ*c*\!$$\*SubscriptBox[\(∂$$, $$t$$]$$Tg[t, z]$$\) - k*\!$$\*SubscriptBox[\(∂$$, $$z, z$$]$$Tg[t, z]$$\);
TiG[z_?NumericQ] := solnd[td, z];
Tits := solnd[td, ts]
Subscript[Γd, N] = NeumannValue[-q, z == 0];

solng = NDSolveValue[{eqn3 == Subscript[Γd,
N], (Tg[t, z] /. z -> ts) == (-300*t) + Tits,
Tg[0, z] == TiG[z]}, Tg, {t, 0, tg}, {z, 0, ts},
Method -> {"PDEDiscretization" -> {"MethodOfLines",
"SpatialDiscretization" -> {"FiniteElement",
"MeshOptions" -> {"MaxCellMeasure" -> (ts)/1011}}}}];

Do[
Ti[z_?NumericQ] := \[Piecewise] {
{solng[tg, z], 0 <= z <= ts},
{100, True}
};

solnd = NDSolveValue[{eqn2 == Subscript[Γ, N],
Subscript[Γ, D], T[0, z] == Ti[z]},
T, {t, 0, td}, {z, 0, ts + tsl},
Method -> {"PDEDiscretization" -> {"MethodOfLines",
"SpatialDiscretization" -> {"FiniteElement",
"MeshOptions" -> {"MaxCellMeasure" -> (ts + tsl)/1011}}}}];

TiG[z_?NumericQ] := solnd[td, z];
Tits := solnd[td, ts];

solng = NDSolveValue[{eqn3 == Subscript[Γd,
N], (Tg[t, z] /. z -> ts) == (-300*t) + Tits,
Tg[0, z] == TiG[z]}, Tg, {t, 0, tg}, {z, 0, ts},
Method -> {"PDEDiscretization" -> {"MethodOfLines",
"SpatialDiscretization" -> {"FiniteElement",
"MeshOptions" -> {"MaxCellMeasure" -> (ts)/1011}}}}];
, {i, 10}]

Animate[Plot[solnd[t, z], {z, 0, ts + tsl},
GridLines -> {{ts, ts + .00005}, {0}}, PlotRange -> {99, 111}], {t,
0, td}]
Plot[solnd[t, ts], {t, 0, td}, GridLines -> {{0}, {107}},
PlotRange -> Full]

Animate[Plot[solng[t, z], {z, 0, ts}, PlotRange -> {100, 111}], {t, 0,
tg}]
Plot[solng[t, ts], {t, 0, tg}, PlotRange -> Full]
• How does this question differ from your 101880? Dec 14 '15 at 3:00
• Your definitions of ρ, k, c are suspect. They are fixed in value for the duration of your computation at the first time that global z gets bound to any numerical value. If you want these variables to bound differently depending on the current value of z, you need to use SetDelayed ( := ) not Set ( = ). Dec 14 '15 at 6:03