<Text-field layout="Heading 1" style="Heading 1"><Font bold="false" executable="false">Excercise 1</Font></Text-field><Text-field layout="Heading 2" style="Heading 2"><Font bold="false" italic="true" size="12">Critical parameters</Font></Text-field>restart;
p:=R*T/(v-b)-a/v^2;
diff(p,v);
diff(p,v$2);
eqn:={subs({pa=1,T=1,v=1},pa-p),subs({pa=1,T=1,v=1},diff(p,v)),subs({pa=1,T=1,v=1},diff(p,v$2))};solve(eqn,{a,b,R});<Text-field layout="Heading 1" style="Heading 1"><Font bold="false" executable="false" italic="true" size="12">Figures</Font></Text-field>restart;p:=R*T/(v-b)-a/v^2;
pT:=subs({a=3,b=1/3,R=8/3},p);plot({subs(T=1,pT),subs(T=.9,pT),subs(T=.8,pT)},v=.5..4);plot3d(pT,T=.5..1.2,v=.5..4);<Text-field layout="Heading 1" style="Heading 1"><Font bold="false" size="14">CompHW1</Font></Text-field>restart;
p:=R*T/(v-b)-a/v^2;
pT:=subs({T=Tf},p);
T:=e/c+a/v/c;
s:=int(p/T,v);
solve({subs({a=3,b=1/3,R=8/3,c=1},p)=1,subs({a=3,b=1/3,R=8/3,c=1},T)=1},{e,v});subs({a=3,b=1/3,R=8/3,c=1},s);
plot3d(%,e=-4.5..2,v=.35..4); mu:=c*T-2*a/v+R*T*(v/(v-b)-log(v-b))-c*T*log(T);ef:=solve(T=Tf,e);
mu:=subs({e=ef},ef+p*v-s*Tf);
muf:=subs({a=3,b=1/3,R=8/3,c=1},mu);
puf:=subs({a=3,b=1/3,R=8/3,c=1},pT);
plot(subs(Tf=.9,muf),v=.57..2.9);plot({[subs(Tf=.9,puf),.2*v,v=.57..2.9],[subs(Tf=.9,puf),subs(Tf=.9,muf),v=.57..2.9]});<Text-field layout="Heading 1" style="Heading 1"><Font bold="false" italic="false">Excercise 2</Font></Text-field><Text-field layout="Heading 2" style="Heading 2"><Font bold="false" italic="true" size="14">Ordinary reversible simple system</Font></Text-field>R:=1:c:=3/2:
T0:=0.9:p0:=0.6:
alpha:=100:beta:=100:
p:=R*T(t)/v(t):
de1:=c*diff(T(t),t)+alpha*(T(t)-T0)+p*beta*(p-p0);
de2:=diff(v(t),t)-beta*(p-p0);solth:=dsolve({de1,de2,T(0)=0.48,v(0)=0.4},{T(t),v(t)},type=numeric);
solth1:=dsolve({de1,de2,T(0)=0.96,v(0)=0.8},{T(t),v(t)},type=numeric);
solth2:=dsolve({de1,de2,T(0)=1.44,v(0)=1.2},{T(t),v(t)},type=numeric);print(solth(1));
print(rhs(solth(1)[3]));#plot(R*T0/v,v=0.3..3):
#plot('rhs(solth(t)[3])',t=0..1.30);
#plot('rhs(solth(t)[2])',t=0..0.30);
plot({['rhs(solth(t)[3])','rhs(solth(t)[2])'/'rhs(solth(t)[3])',
t=0..5],['rhs(solth1(t)[3])','rhs(solth1(t)[2])'/'rhs(solth1(t)[3])',
t=0..5],['rhs(solth2(t)[3])','rhs(solth2(t)[2])'/'rhs(solth2(t)[3])',
t=0..5],[v,R*T0/v,v=0.7..2]}, color=[red,blue,green,black]);<Text-field layout="Heading 1" style="Heading 1"><Font bold="false" size="14">CompHW2</Font></Text-field><Text-field layout="Heading 1" style="Heading 1"><Font bold="false" italic="true" size="14">Reversible</Font></Text-field>restart;
p:=R*T/(v-b)-a/v^2:
e:=c*T-a/v:
dTe:=subs({a=3,b=1/3,c=4,R=8/3},diff(e,T)):dve:=subs({a=3,b=1/3,R=8/3,v=v(t)},diff(e,v));
subs(T=solve(diff(p,v)=0,T),p):
binc:=subs({a=3,b=1/3,R=8/3},%);
pt:=subs({a=3,b=1/3,R=8/3,T=T(t),v=v(t)},p);
T0:=0.9:p0:=0.6:
vi:=[0.4,0.8,1.2]:
for i from 1 to 3 do
Ti[i]:=solve(subs({a=3,b=1/3,R=8/3,v=vi[i]},p)=1.2,T):
od:
alpha:=10:beta:=10:
de1:=dTe*diff(T(t),t)+dve*diff(v(t),t)+alpha*(T(t)-T0)+pt*beta*(pt-p0);
de2:=diff(v(t),t)-beta*(pt-p0);
solth1:=dsolve({de1,de2,T(0)=Ti[1],v(0)=vi[1]},{T(t),v(t)},type=numeric);
solth2:=dsolve({de1,de2,T(0)=Ti[2],v(0)=vi[2]},{T(t),v(t)},type=numeric);
solth3:=dsolve({de1,de2,T(0)=Ti[3],v(0)=vi[3]},{T(t),v(t)},type=numeric);print(solth1(1));tend:=5;
plot({['rhs(solth1(t)[3])',8/3*'rhs(solth1(t)[2])'/('rhs(solth1(t)[3])'-1/3)-3/'rhs(solth1(t)[3])'^2,t=0..tend],
['rhs(solth2(t)[3])',8/3*'rhs(solth2(t)[2])'/('rhs(solth2(t)[3])'-1/3)-3/'rhs(solth2(t)[3])'^2,t=0..tend],
['rhs(solth3(t)[3])',8/3*'rhs(solth3(t)[2])'/('rhs(solth3(t)[3])'-1/3)-3/'rhs(solth3(t)[3])'^2,t=0..tend],
[v,binc,v=0.7..5]}, color=[blue,blue,blue,red],labels=[v,"p"],thickness=2);<Text-field layout="Heading 1" style="Heading 1"><Font bold="false" italic="true" size="14">Extended</Font></Text-field>restart;
p:=R*T/(v-b)-a/v^2:
e:=c*T-a/v:
dTe:=subs({a=3,b=1/3,c=4,R=8/3},diff(e,T)):dve:=subs({a=3,b=1/3,R=8/3,v=v(t)},diff(e,v));
subs(T=solve(diff(p,v)=0,T),p):
binc:=subs({a=3,b=1/3,R=8/3},%);
pt:=subs({a=3,b=1/3,R=8/3,T=T(t),v=v(t)},p)-k*u(t);
T0:=0.9:p0:=0.6:
vi:=[0.4,0.8,1.2]:
for i from 1 to 3 do
Ti[i]:=solve(subs({a=3,b=1/3,R=8/3,v=vi[i]},p)=1.2,T):
od:
alpha:=10:beta:=10:k:=0;
de1:=dTe*diff(T(t),t)+dve*diff(v(t),t)+alpha*(T(t)-T0)+pt*diff(v(t),t);
de2:=diff(v(t),t)-u(t);
de3:=diff(u(t),t)-beta*(pt-p0);
solth1:=dsolve({de1,de2,de3,T(0)=Ti[1],v(0)=vi[1],u(0)=0},{T(t),v(t),u(t)},type=numeric);
solth2:=dsolve({de1,de2,de3,T(0)=Ti[2],v(0)=vi[2],u(0)=0},{T(t),v(t),u(t)},type=numeric);
solth3:=dsolve({de1,de2,de3,T(0)=Ti[3],v(0)=vi[3],u(0)=0},{T(t),v(t),u(t)},type=numeric);tend:=15;
plot({['rhs(solth1(t)[4])',8/3*'rhs(solth1(t)[2])'/('rhs(solth1(t)[4])'-1/3)-3/'rhs(solth1(t)[4])'^2,t=0..tend],
['rhs(solth2(t)[4])',8/3*'rhs(solth2(t)[2])'/('rhs(solth2(t)[4])'-1/3)-3/'rhs(solth2(t)[4])'^2,t=0..tend],
['rhs(solth3(t)[4])',8/3*'rhs(solth3(t)[2])'/('rhs(solth3(t)[4])'-1/3)-3/'rhs(solth3(t)[4])'^2,t=0..tend],
[v,binc,v=0.67..5]}, color=[blue,blue,blue,red],labels=[v,"p"],thickness=2);<Text-field layout="Heading 1" style="Heading 1"><Font bold="false" executable="false" italic="false">Excercise 3</Font></Text-field>a:=.1;Ta:=1:tp:=0.1:
heatcond := diff(T(x,t),t)=a*diff(T(x,t),x,x);
ibc := {T(x,0)=1,T(0,t)=0, D[1](T)(1,t)=0}:
ibc := {T(x,0)=1,T(0,t)=1, D[1](T)(1,t)=0}:
ibc := {T(x,0)=1,T(0,t)=1, T(1,t)=0}: ibc := {T(x,0)=0,T(0,t)=Ta*(Heaviside(t)-Heaviside(t-tp)), T(1,t)=0};pdsol := pdsolve(heatcond,ibc,numeric);pdsol:-value(t=1);
print(%(0.5));p1:= pdsol:-plot(t=0.001):
p2:= pdsol:-plot(t=0.1):
p3:= pdsol:-plot(t=.5):
plots[display]({p1,p2,p3},title=`Temperature distribution`);pdsol:-plot3d(t=0..3,x=0..1,axes=boxed);pdsol:-animate(t=1,frames=50);