Ch10mini_project.sagews

# Mini project: IVP for heat flow in a thin circular ring
#1. Solve the IVP for the ODE u'_t=u''_xx, -pi<=x<pi,t>=0, u(x,0)=pi-|x|.
# Use the general solution to one-dimentional heat problem in a thin ring derived in the notes.
# Plot three snapshots of the solution to the problem in one figure.
###################################
#2. Write a CAS function heat_in_ring(f) that takes a continuous 2pi-periodic function f and returns the solution to the IVP for the heat equation with initial condition u(x,0)=f(x). Test your function on the problem in part 1.

reset

<function reset at 0x7f6d788e30c8>

#1. The formula for solution involves Fourier series for the initial function. The function pi-|x| is even, only "a" coefficients are needed. Use your acoef function to find some number of a's.
f(x)=pi-abs(x)
def acoef(f,n,x0,X):
    L=[integral(f(x)*cos(k*x),x,-pi,pi) for k in range(n) ]
    return [N(c/pi) for c in L]
# Replace t variable with x in your original code

A=acoef(f,5,-pi,pi);A

[3.14159265358979, 1.27323954473516, 0.000000000000000, 0.141471060526129, 0.000000000000000]

#Encode the formula for solution to the IVP.
u(x,t)=A[0]/2+sum(A[k]*cos(k*x)*exp(-k^2*t) for k in range(1,5))

N(u(1,1.5))

1.72429501682374

# Create a graphics object, like G=Graphics()
# Make three plots of u(x,t_j), add them to G
#Use the command G.show() to plot all three curves in one figure. Add legend (like in Chapter 10, Fig. 10.9)
G=Graphics()
p0=plot(u(x,0),(x,-pi,pi),color='red',figsize=[5,5],legend_label='u(x,0)');p0

f_plot=plot(f,(x,-pi,pi),figsize=[4,3],color='red',legend_label=('f(x)=pi-|x|'))
u_plot1=plot(u(x,0),(x,-pi,pi),figsize=[5,4],color='green',legend_label=('u(x,0)'));f_plot+u_plot1
u_plot2=plot(u(x,1),(x,-pi,pi),figsize=[5,4],color='blue',legend_label=('u(x,1)'));f_plot+u_plot2
u_plot3=plot(u(x,3),(x,-pi,pi),figsize=[5,4],color='green',legend_label=('u(x,3)'));f_plot+u_plot3

fmh

(u_plot1+u_plot2+u_plot3).show()

#2.Put your interactive solution together as a function.
#Assumptions: f is a 2pi-periodic function, n is the number of terms in Fourier series for f; t2, t3 are time moment after initial time t0=0
f(x)=pi-abs(x)
def heat_in_ring(f,n,t1,t2):
    A=acoef(f,2,0,2*pi)
    u(x,t)=A[0]/2+sum(A[k]*cos(k*x)*exp(-k^2*t) for k in range(1,n+1))
    G=Graphics()
    plot0=plot(u(x,0),(x,0,2*pi),thickness=3,color='red',legend_label=('u(x,0)'))
    G+=plot0
    plot1=plot(u(x,1),(x,0,2*pi),thickness=3,color='green',legend_label=('u(x,1)'))
    G+=plot1
    plot2=plot(u(x,3),(x,0,2*pi),thickness=3,color='blue',legend_label=('u(x,3)'))
    G+=plot2
    return G

heat_in_ring(f,0.5,0.5,2)

#Score: 5/10. No Part 2.