| Author |
Comment/Response |
angel nyirenda
|
08/27/12 04:56am
i cant figure out why my code does output empty graphs when executed. the code models a TL irradiation process by solving a set of kinetic equations.
Remove["Global`*"];
Off[General::stop];
programMain := (A = 10^-17; Ar = 10^-13; Ah = 10^-15; N1 = 10^15;
Nh = 3*10^14;
solveDiff[n10_, nh0_, nv0_, nc0_, R_, tfinal_] := Module[{t},
sol = NDSolve[{n1'[t] == A*(N1 - n1[t])*nc[t],
nc'[t] == nh'[t] + nv'[t] - n1'[t],
nh'[t] == -Ar*nh[t]*nc[t] + Ah*nv[t]*(Nh - nh[t]),
nv'[t] == R - Ah*nv[t]*(Nh - nh[t]),
n1[0] == n10, nh[0] == nh0, nv[0] == nv0, nc[0] == nc0},
{n1[t], nh[t], nv[t], nc[t]},
{t, 0, tfinal}, MaxSteps -> 50000]]
initValues[b_, d_] := Module[{}, n10 == Last[n1[d] /. b];
nh0 == Last[nh[d] /. b];
nv0 == Last[nv[d] /. b];
nc0 == Last[nc[d] /. b]];(*irradiation*)
R = 10^14;
tfinal9 = irrTime;
solveDiff[n10, nh0, nv0, nc0, R, tfinal9];
sol9 = sol;
initValues[sol9, tfinal9];(*Relaxation stage*)
R = 0;
tfinal10 = 60;
solveDiff[n10, nh0, nv0, nc0, R, tfinal10];
sol10 = sol;)
tstart = 0; tend = 20; tstep = 1; n1List1 = {}; n1List2 = {};
For[tloop = tstart, tloop <= tend, tloop += tstep,
n10 = 0; nh0 = 0; nv0 = 0; nc0 = 0;
irrTime = tloop;
programMain;
initValues[sol9,
tfinal9];(*find n(t) at the end of irradiation stage*)
AppendTo[n1List1, {tloop, n10}];
initValues[sol10, tfinal10];
AppendTo[n1List2, {tloop, n10}];]
gr1 = ListPlot[n1List1, PlotRange -> All, PlotJoined -> True,
PlotLabel -> "n(t) vs t"]
gr2 = ListPlot[n1List2, PlotRange -> All, PlotJoined -> True,
PlotLabel -> "n(t+T) vs t"]
URL: www.ru.ac.za, |
|