POLYMATH Report          DEQ 
 Ordinary Differential Equations       2016-May-31 

Calculated values of DEQ variables
  Variable Initial value Final value Minimal value Maximal value
1 dTdt 0 -0.004826 -2.26445 0
2 errsum 0 720.924 -0.016445 720.924
3 Kc 0 0 0 0
4 q 10000 10000 10000 10000
5 rhoVCp 4000 4000 4000 4000
6 step 0 1 0 1
7 t 0 60 0 60
8 T 80 60.0386 60.0386 80
9 T0 80 60.0438 60.0438 80
10 taud 1 1 1 1
11 tauI 2 2 2 2
12 taum 5 5 5 5
13 Ti 60 40 40 60
14 Tm 80 60.1148 60.1148 80.0349
15 Tr 80 80 80 80
16 WC 500 500 500 500

Differential equations
1 d(T)/d(t) = (WC*(Ti-T)+q)/rhoVCp
2 d(T0)/d(t) = (T-T0-(taud/2)*dTdt)*2/taud
3 d(Tm)/d(t) = (T0-Tm)/taum
4 d(errsum)/d(t) = Tr-Tm

Explicit equations
1 WC = 500
2 rhoVCp = 4000
3 taud = 1
4 taum = 5
5 Tr = 80
6 Kc = 0
7 tauI = 2
8 step = if (t<10) then (0) else (1)
9 Ti = 60+step*(-20)
10 q = 10000+Kc*(Tr-Tm)+Kc/tauI*errsum
11 dTdt = (WC*(Ti-T)+q)/rhoVCp

Problem source text
# S. 23(a)* - ODE System
# Dynamics of Heated Tank-Open loop
# Verified Final Values: T = 60.0386,T0 = 60.0438, Tm = 60.1148
# Ref.: Comput. Appl. Eng. Educ. 6: 178-179, 1998
d(T)/d(t)=(WC*(Ti-T)+q)/rhoVCp
d(T0)/d(t)=(T-T0-(taud/2)*dTdt)*2/taud
d(Tm)/d(t)=(T0-Tm)/taum
d(errsum)/d(t)=Tr-Tm
WC=500
rhoVCp=4000
taud=1
taum=5
Tr=80
Kc=0
tauI=2
step=if (t<10) then (0) else (1)
Ti=60+step*(-20)
q=10000+Kc*(Tr-Tm)+Kc/tauI*errsum
dTdt=(WC*(Ti-T)+q)/rhoVCp
t(0)=0
T(0)=80
T0(0)=80
Tm(0)=80
errsum(0)=0
t(f)=60

Matlab formatted problem
Create m file called PolyOde.m and paste the following text into it.
% S. 23(a)* - ODE System
% Dynamics of Heated Tank-Open loop
% Verified Final Values: T = 60.0386,T0 = 60.0438, Tm = 60.1148
% Ref.: Comput. Appl. Eng. Educ. 6: 178-179, 1998
function PolyOde
   tspan = [0 60]; % Range for the independent variable
   y0 = [80; 80; 80; 0]; % Initial values for the dependent variables
   [t,y]=ode45(@ODEfun,tspan, y0);
   plot (t,y);
   xlabel('t');
   legend('T','T0','Tm','errsum');
   fprintf('T = %16.6f \n',y(length(y),1));
   fprintf('T0 = %16.6f \n',y(length(y),2));
   fprintf('Tm = %16.6f \n',y(length(y),3));
   fprintf('errsum = %16.6f \n',y(length(y),4));
end

function dYfuncvecdt = ODEfun(t,Yfuncvec)
   T = Yfuncvec(1);
   T0 = Yfuncvec(2);
   Tm = Yfuncvec(3);
   errsum = Yfuncvec(4);
   WC = 500;
   rhoVCp = 4000;
   taud = 1;
   taum = 5;
   Tr = 80;
   Kc = 0;
   tauI = 2;
   if (t < 10)
       step = 0;
   else
       step = 1;
   end

   Ti = 60 + step * -20;
   q = 10000 + Kc * (Tr - Tm) + Kc / tauI * errsum;
   dTdt = (WC * (Ti - T) + q) / rhoVCp;
   dTdt = (WC * (Ti - T) + q) / rhoVCp;
   dT0dt = (T - T0 - (taud / 2 * dTdt)) * 2 / taud;
   dTmdt = (T0 - Tm) / taum;
   derrsumdt = Tr - Tm;
   dYfuncvecdt = [dTdt; dT0dt; dTmdt; derrsumdt];
end

General Settings
Total number of equations 15
Number of differential equations 4
Number of explicit equations 11
Reporting digits 8
Elapsed time 0.27 sec
Solution method RKF_45
Step size guess. h 1E-06
Truncation error tolerance. eps 1E-06
Calculated Intermediate data points 50

Calculated data points
    t T T0 Tm errsum step Ti q dTdt
1 0 80 80 80 0 0 60 10000 0
2 3.00046 80 80 80 0 0 60 10000 0
3 3.96046 80 80 80 0 0 60 10000 0
4 4.92046 80 80 80 0 0 60 10000 0
5 6.84046 80 80 80 0 0 60 10000 0
6 7.80046 80 80 80 0 0 60 10000 0
7 8.76046 80 80 80 0 0 60 10000 0
8 9.72046 80 80 80 0 0 60 10000 0
9 10.8191 78.0536 79.9426 80.0349 -0.016445 1 40 10000 -2.26445
10 12.0603 75.4591 77.477 79.7436 0.064128 1 40 10000 -1.94579
11 13.3841 73.1015 74.8453 78.8901 0.916614 1 40 10000 -1.67542
12 14.679 71.1436 72.6292 77.6968 3.09797 1 40 10000 -1.42329
13 15.8502 69.6259 70.9093 76.4475 6.51647 1 40 10000 -1.23489
14 17.2694 68.0612 69.136 74.8431 12.6908 1 40 10000 -1.03991
15 18.3817 67.0148 67.9501 73.5794 19.1315 1 40 10000 -0.908924
16 19.6524 65.9845 66.7825 72.176 28.1884 1 40 10000 -0.779443
17 21.1112 64.987 65.652 70.6569 40.7232 1 40 10000 -0.653539
18 21.9233 64.5056 65.1063 69.8638 48.6363 1 40 10000 -0.592564
19 23.7481 63.5867 64.0649 68.2342 68.6516 1 40 10000 -0.47577
20 24.7081 63.1811 63.6052 67.4636 80.3214 1 40 10000 -0.422272
21 25.6681 62.8214 63.1976 66.7521 92.7025 1 40 10000 -0.374522
22 26.6281 62.5023 62.836 66.098 105.739 1 40 10000 -0.332171
23 28.5481 61.9684 62.2309 64.9513 133.565 1 40 10000 -0.261295
24 29.5081 61.7458 61.9786 64.4529 148.255 1 40 10000 -0.231748
25 30.4681 61.5484 61.7549 64.0001 163.401 1 40 10000 -0.205542
26 31.4281 61.3733 61.5564 63.5896 178.961 1 40 10000 -0.182299
27 33.3481 61.0803 61.2243 62.8829 211.17 1 40 10000 -0.143402
28 34.3081 60.9581 61.0859 62.5804 227.75 1 40 10000 -0.127186
29 35.2681 60.8498 60.9631 62.308 244.606 1 40 10000 -0.112804
30 36.2281 60.7537 60.8542 62.0631 261.71 1 40 10000 -0.100048
31 38.1481 60.5929 60.6719 61.6455 296.564 1 40 10000 -0.078701
32 39.1081 60.5258 60.5959 61.4685 314.271 1 40 10000 -0.069801
33 40.0681 60.4664 60.5286 61.3099 332.139 1 40 10000 -0.061908
34 41.0281 60.4136 60.4688 61.1679 350.151 1 40 10000 -0.054907
35 42.9481 60.3254 60.3688 60.9274 386.548 1 40 10000 -0.043192
36 43.9081 60.2886 60.3271 60.8259 404.907 1 40 10000 -0.038308
37 44.8681 60.2559 60.2901 60.7354 423.359 1 40 10000 -0.033976
38 45.8281 60.227 60.2573 60.6546 441.892 1 40 10000 -0.030134
39 47.7481 60.1786 60.2024 60.5182 479.172 1 40 10000 -0.023704
40 48.7081 60.1584 60.1795 60.461 497.902 1 40 10000 -0.021024
41 49.6681 60.1405 60.1592 60.4099 516.684 1 40 10000 -0.018646
42 50.6281 60.1246 60.1412 60.3645 535.513 1 40 10000 -0.016538
43 52.5481 60.098 60.1111 60.288 573.29 1 40 10000 -0.013009
44 53.5081 60.0869 60.0985 60.2559 592.229 1 40 10000 -0.011538
45 54.4681 60.0771 60.0874 60.2274 611.197 1 40 10000 -0.010233
46 55.4281 60.0684 60.0775 60.202 630.191 1 40 10000 -0.009076
47 57.3481 60.0538 60.061 60.1594 668.246 1 40 10000 -0.00714
48 58.3081 60.0477 60.0541 60.1416 687.302 1 40 10000 -0.006332
49 59.2681 60.0423 60.0479 60.1257 706.374 1 40 10000 -0.005616
50 60 60.0386 60.0438 60.1148 720.924 1 40 10000 -0.004826