Tell me more ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.
up vote 0 down vote favorite
1

I've used Wolfram Mathematica 8 Trial version for two days. Any pointers are appreciated. Thanks in advance!

EDIT: There are supposed to be 4 plots, but only 3 are showing. Everything works except last Plot. The input is as follows:

wallThickness = 0.2
wallArea = 2.2*2.2

c1 = c2 = 1.1
v1 = v2 = wallArea*4.4
p1 = p2 = 1.15

emissivity = 0.9
stefan = 5.6704*10^(-8)

hAir = 100
kWall = 10^(-1.4)

R1 = R3 = 1/(hAir*wallArea)
R2 = wallThickness/(kWall*wallArea)
R = 6 (R1 + R2 + R3)

a = p1*v1*c1
b = p2*v2*c2

fullyTime = 358.4
decayTime = 636
extTime = 1100

growCoef = 0.0029
decayCoef = 0.00003

grow[t_] = growCoef*t^2
fully[t_] = grow[fullyTime]
decay[t_] = grow[fullyTime]*Exp[-((t - decayTime)^2)*decayCoef]

heatRate[t_] = 
 2*1000*Piecewise[{{grow[t], 0 <= t < fullyTime}, {fully[t], 
     fullyTime <= t < decayTime}, {decay[t], 
     decayTime <= t < extTime}}]

Plot[heatRate[t], {t, 0, 1100}]

heat[t_] = Integrate[heatRate[t], t]
Plot[heat[t], {t, 0, 1100}]

sol1 = DSolve[{Derivative[1][q][t] + 
     q[t]/(R*b) == (heat[t] - 5*q[t])/(R*a), q[0] == 0}, q, t]
Plot[q[t]/1000 /. sol1, {t, 0, 1100}]

j = q[t] /. sol1[[1]]

sol2 = DSolve[{Derivative[1][w][t] == j - w[t]/(6*R*b), w[0] == 0}, w,
   t]
Plot[Evaluate[w[o] /. sol2], {o, 0, 1100}]

Quit[]
share|improve this question
In the Dsolve for sol2 there should be == instead of =. – Andrew May 16 '12 at 17:48
Also, sol2 is the solution for w so your last Plot should be Plot[Evaluate[w[o] /. sol2] ...]. To get rid of the warnings produced by DSolve you could Rationalize your differential equations first, i.e. do sol1 = DSolve[Rationalize[{...}], q, t] etc. – Heike May 16 '12 at 18:12
Thank you @Andrew and (a)Heike. I have corrected mentioned errors, but new replacement errors have occured. – loland May 16 '12 at 19:26
1  
This gives me some Reduce::ratnz, but otherwise works fine. @loland Quit the kernel using Quit[] to be sure you have no stale definitions and try again – Ajasja May 16 '12 at 19:50
@Ajasja after restarting and adding Quit[] it works, except the last Plot[] is not working. – loland May 16 '12 at 20:48
show 8 more comments

closed as too localized by rm -rf♦, rcollyer, acl, Szabolcs, Oleksandr R. May 25 '12 at 3:01

This question is unlikely to help any future visitors; it is only relevant to a small geographic area, a specific moment in time, or an extraordinarily narrow situation that is not generally applicable to the worldwide audience of the internet. For help making this question more broadly applicable, see the FAQ.

1 Answer

up vote 1 down vote 
wallThickness = 0.2;
wallArea = 2.2*2.2;

c1 = c2 = 1.1;
v1 = v2 = wallArea*4.4;
p1 = p2 = 1.15;

emissivity = 0.9;
stefan = 5.6704*10^(-8);

hAir = 100;
kWall = 10^(-1.4);

R1 = R3 = 1/(hAir*wallArea);
R2 = wallThickness/(kWall*wallArea);
R = 6 (R1 + R2 + R3);

a = p1*v1*c1;
b = p2*v2*c2;

fullyTime = 358.4;
decayTime = 636;
extTime = 1100;

growCoef = 0.0029;
decayCoef = 0.00003;

grow[t_] := growCoef*t^2
fully[t_] := grow[fullyTime]
decay[t_] := grow[fullyTime]*Exp[-((t - decayTime)^2)*decayCoef]

heatRate[t_] := 
  2*1000*Piecewise[{{grow[t], 0 <= t < fullyTime}, {fully[t], 
      fullyTime <= t < decayTime}, {decay[t], 
      decayTime <= t < extTime}}];

Plot[heatRate[t], {t, 0, 1100}]

heat[t_] := Integrate[heatRate[t], t];
Plot[Evaluate[heat[t]], {t, 0, 1100}]

sol1 = NDSolve[{Derivative[1][q][t] + 
      q[t]/(R*b) == (heat[t] - 5*q[t])/(R*a), q[0] == 0}, 
   q, {t, 0, 1100}];
Plot[Evaluate[q[t]/1000 /. sol1], {t, 0, 1100}]

tj = q[t] /. sol1[[1]];

sol2 = NDSolve[{Derivative[1][w][t] == tj - w[t]/(6*R*b), w[0] == 0}, 
   w, {t, 0, 1100}];
Plot[Evaluate[w[o] /. sol2], {o, 0, 1100}]

try above code.

share|improve this answer
3  
Might I ask, what did you change from the OPs code? – rcollyer May 17 '12 at 13:47

Not the answer you're looking for? Browse other questions tagged or ask your own question.