# Maple code for Newton Cooling maple lab 2, S2009 # #Notes on L2.1 # maple integral table lookup unassign('omega','t','u0','k','k1','k2','k3','x'); integrand:=(79-19*cos(omega*(x-3)))*k3*exp(k*x-k*t); integral:=int(integrand,x=0..t); # Test LHS=RHS for u'+ku=k3 A. unassign('t','C','omega','k','k1','k2','k3'): k:=k1+k2:k3:=k1+0.8*k2: uh:=C*exp(-k*t): up:=(explicit answer from L2.1): myANS:=uh+up: LHS:=diff(myANS,t)+k*myANS: RHS:=k3*(79-19*cos(omega*(t-3))): simplify(expand(LHS-RHS)); #Notes on L2.2 # Extract steady-state solution expr:=(myANS from L2.1): SS:=subs(exp(-k*t)=0,exp(-2*k*t)=0,expr); #Notes on L2.3 unassign('t','u0','k','k1','k2','k3','omega'): AA:=unapply(79-19*cos(omega*(t-3)),(t,omega)); uss:=k3*(79*k^2+79*omega^2-19*k^2*cos(omega*t-3*omega)- 19*k*omega*sin(omega*t-3*omega))/(k*(k^2+omega^2)); uss0:=subs(t=0,uss); U:=unapply((u0-uss0)*exp(-k*t)+uss,(t,u0,k,k3,omega)); k1:=0.32:k2:=0.5: plot({U(t,76,k1+k2,k1+0.8*k2,Pi/12),AA(t,Pi/12)},t=0..48,axes=boxed); #Notes on L2.4 with(plots): unassign('t','u0','k','k1','k2','k3','omega'): uh:=u0*exp(-k*t): up:=(the answer from L2.1): myANS:=uh+up: U:=unapply(myANS,(t,u0,k,k3,omega));k2:=0.5: implicitplot(U(t,76,k1+k2,k1+0.8*k2,Pi/12)=70,t=0..72,k1=0.2..0.48,axes=boxed); plot3d({U(t,76,k1+k2,k1+0.8*k2,Pi/12),70},t=0..72,k1=0.2..0.48);