Python help. How do I convert this text to a file that python can execute? When
ID: 3784923 • Letter: P
Question
Python help. How do I convert this text to a file that python can execute? When I save this coding in a blank file, I get errors when I try to run it on python. This is a file a professor gave. I have a similar project and would like to use this code. However, I am not familar at all with python and do not know where indentations, etc. need to go.
import pylab
from mpl_toolkits.mplot3d import axes3d import matplotlib.pyplot as plt
from matplotlib import cm
# Discretization N = 100
nX = 200
nT = 100
# Physical Constants alpha = 1
# T(t)
def T(L,c,t):
return exp(-L*L * c * t)
# A_n for sine series def A(m):
if(m==2):
return 0.5
else:
return -4*sin(0.5*pi*m)/(pi*m**2 - 4*pi)
# X(x)
def X(m,x):
return sin(m*pi*x)
def dX(m,x):
return m*pi*cos(m*pi*x)
def ddX(m,x):
return -(m*pi)**2*sin(m*pi*x)
# u = X(x)T(t)
def U(alpha,x,t,N):
Ys = 0
for n in range(1,N+1):
Ys = Ys + T(n*pi,alpha,t)*A(n)*X(n,x) return Ys
# du/dx = dX(x)/dx T(t)
def dU(alpha,x,t,N):
Ys = 0
for n in range(1,N+1):
Ys = Ys + T(n*pi,alpha,t)*A(n)*dX(n,x) return Ys
# du/dt = ddX(x)/ddx T(t) def ddU(alpha,x,t,N):
Ys = 0
for n in range(1,N+1):
Ys = Ys + T(n*pi,alpha,t)*A(n)*ddX(n,x) return Ys
# Create an x-t grid (appears as a matrix)
[Xm,Tm] = np.mgrid[0:1:nX*1j, 0:0.2:nT*1j]
# Solve U at each x-t point Um = U(alpha,Xm,Tm,N)
def plotXTUsurf():
# Plot the x-t-u surface and a contour above the surface
fig = plt.figure()
ax = fig.gca(projection='3d')
ax.plot_surface(Xm, Tm, Um, rstride=2, cstride=2, alpha=0.9, facecolors=cm.jet(Um)) cset = ax.contour(Xm, Tm, Um, zdir='z', offset=1, cmap=cm.coolwarm)
plt.show()
xlabel("x")
ylabel("time")
# Plot the temperature response at x=0.6 def plotTempX0p6():
t = linspace(0,1,1000)
Ut = U(alpha,0.6,t,50)
fig2 = plt.figure() ax2 = fig2.gca() ax2.plot(t, Ut) xlabel("Time") ylabel("u")
print "Maximum Value at x =0.6 is " + str(max(Ut)) print "Located at " + str(t[argmax(Ut)])
# L squared error versus N
def L2(N):
Error = 0
for n in range(1,N+1): Error += -1.0/2.0*A(n)**2 Error += 0.25
return Error
def PlotL2Error(): n = zeros(N)
L2E = zeros(N) for i in range(1,N): n[i] = i+1
L2E[i] = L2(i+1)
fig3 = plt.figure()
ax3 = fig3.gca() ax3.semilogy(L2E) xlabel("Expansion Terms, N") ylabel("$L_2$ Error")
# Coefficients
def plotCoeffs(): ACoeffs = zeros(N)
for i in range(1,N): ACoeffs[i]=A(i)
fig4 = plt.figure()
ax4 = fig4.gca() ax4.semilogy(ACoeffs,'.') xlabel("n") ylabel("$A_n$")
def FiniteDiffCheck(t,x,h,N):
dt = ( U(alpha,x,t+h,N) - U(alpha,x,t-h,N)) /(2*h)
dxx = ( U(alpha,x+h,t,N) - 2.0* U(alpha,x,t,N) + U(alpha,x-h,t,N)) /(h**2) return dt-alpha*dxx
print(FiniteDiffCheck(0.25,0.25,0.001,100))
# Calculate the energy the hard way
# Create array of du/dx on both boundaries def TotalEnergy():
TimeEffectivelyZero = 4
t = linspace(0,TimeEffectivelyZero,5000) dUtLeft = dU(alpha,0,t,50)
dUtRight = dU(alpha,1.0,t,50)
fig5 = plt.figure()
ax5 = fig5.gca()
ax5.plot(t, dUtLeft)
ax5.plot(t, dUtRight)
xlabel("Time")
ylabel("u")
# Create a rough trapezoidal rule integration routine def TrapIntegrate(array,a,b):
n = size(array)
I=0
for i in range(0,n-1):
I = I + 0.5*(array[i]+array[i+1])
return I*(b-a)/(n-1)
#Energy via flux
print("Energy Transport via Fluxes " + str(TrapIntegrate(dUtLeft-dUtRight,0,TimeEffectivelyZero)))
# Check Initial Condition def CheckIC():
x = linspace(0,1,300) fig6 = plt.figure()
ax6 = fig6.gca()
for i in [1,2,4,8,16,32,64,128,256,512]: U0 = U(alpha,x,0,i)
ax6.plot(x, U0, label="N="+str(i)) xlabel("x")
ylabel("u")
legend(loc='upper right')
# Plot the x-t-u surface and a contour above the surface def FastestChange():
dT = ddU(alpha,Xm,Tm,N)
fig7 = plt.figure()
ax7 = fig7.gca(projection='3d')
ax7.plot_surface(Xm, Tm, dT, rstride=2, cstride=2, alpha=0.9, facecolors=cm.jet(Um)) cset = ax7.contour(Xm, Tm, dT, zdir='z', offset=1, cmap=cm.coolwarm)
plt.show()
xlabel("x")
ylabel("time")
# Plot the x-t-u surface and a contour above the surface def ZeroFlux():
dT = dU(alpha,Xm,Tm,N)
fig7 = plt.figure()
ax7 = fig7.gca()
V = linspace(-1,1,10)
cset = ax7.contourf(Xm, Tm, Um,V, cmap=cm.jet)
V2 = [0]
cset = ax7.contour(Xm, Tm, dT, V2, cmap=cm.binary_r) plt.show()
xlabel("x")
ylabel("time")
Explanation / Answer
"""
AUTHOR :
FILE NAME :
DATE : 30/1/17
"""
from math import pi
import pylab
from mpl_toolkits.mplot3d import axes3d
import matplotlib.pyplot as plt
from matplotlib import cm
import numpy as np
from mpmath import linspace
from numpy.ma import exp, sin, cos, argmax, zeros, size
from scipy.constants import alpha
# Discretization
N = 100
nX = 200
nT = 100
# Physical Constants alpha = 1
# T(t)
def T(L, c, t):
return exp(-L * L * c * t)
# A_n for sine series
def A(m):
if (m == 2):
return 0.5
else:
return -4 * sin(0.5 * pi * m) / (pi * m ** 2 - 4 * pi)
# X(x)
def X(m, x):
return sin(m * pi * x)
def dX(m, x):
return m * pi * cos(m * pi * x)
def ddX(m, x):
return -(m * pi) ** 2 * sin(m * pi * x)
# u = X(x)T(t)
def U(alpha, x, t, N):
Ys = 0
for n in range(1, N + 1):
Ys = Ys + T(n * pi, alpha, t) * A(n) * X(n, x)
return Ys
# du/dx = dX(x)/dx T(t)
def dU(alpha, x, t, N):
Ys = 0
for n in range(1, N + 1):
Ys = Ys + T(n * pi, alpha, t) * A(n) * dX(n, x)
return Ys
# du/dt = ddX(x)/ddx T(t)
def ddU(alpha, x, t, N):
Ys = 0
for n in range(1, N + 1):
Ys = Ys + T(n * pi, alpha, t) * A(n) * ddX(n, x)
return Ys
# Create an x-t grid (appears as a matrix)
[Xm, Tm] = np.mgrid[0:1:nX * 1j, 0:0.2:nT * 1j]
# Solve U at each x-t point
Um = U(alpha, Xm, Tm, N)
def plotXTUsurf():
# Plot the x-t-u surface and a contour above the surface
fig = plt.figure()
ax = fig.gca(projection='3d')
ax.plot_surface(Xm, Tm, Um, rstride=2, cstride=2, alpha=0.9, facecolors=cm.jet(Um))
cset = ax.contour(Xm, Tm, Um, zdir='z', offset=1, cmap=cm.coolwarm)
plt.show()
plt.xlabel("x")
plt.ylabel("time")
# Plot the temperature response at x=0.6
def plotTempX0p6():
t = linspace(0, 1, 1000)
Ut = U(alpha, 0.6, t, 50)
fig2 = plt.figure()
ax2 = fig2.gca()
ax2.plot(t, Ut)
plt.xlabel("Time")
plt.ylabel("u")
print
"Maximum Value at x =0.6 is " + str(max(Ut))
print
"Located at " + str(t[argmax(Ut)])
# L squared error versus N
def L2(N):
Error = 0
for n in range(1, N + 1):
Error += -1.0 / 2.0 * A(n) ** 2
Error += 0.25
return Error
def PlotL2Error():
n = zeros(N)
L2E = zeros(N)
for i in range(1, N): n[i] = i + 1
L2E[i] = L2(i + 1)
fig3 = plt.figure()
ax3 = fig3.gca()
ax3.semilogy(L2E)
plt.xlabel("Expansion Terms, N")
plt.ylabel("$L_2$ Error")
# Coefficients
def plotCoeffs():
ACoeffs = zeros(N)
for i in range(1, N):
ACoeffs[i] = A(i)
fig4 = plt.figure()
ax4 = fig4.gca()
ax4.semilogy(ACoeffs, '.')
plt.xlabel("n")
plt.ylabel("$A_n$")
def FiniteDiffCheck(t, x, h, N):
dt = (U(alpha, x, t + h, N) - U(alpha, x, t - h, N)) / (2 * h)
dxx = (U(alpha, x + h, t, N) - 2.0 * U(alpha, x, t, N) + U(alpha, x - h, t, N)) / (h ** 2)
return dt - alpha * dxx
print(FiniteDiffCheck(0.25, 0.25, 0.001, 100))
# Calculate the energy the hard way
# Create array of du/dx on both boundaries
def TotalEnergy():
TimeEffectivelyZero = 4
t = linspace(0, TimeEffectivelyZero, 5000)
dUtLeft = dU(alpha, 0, t, 50)
dUtRight = dU(alpha, 1.0, t, 50)
fig5 = plt.figure()
ax5 = fig5.gca()
ax5.plot(t, dUtLeft)
ax5.plot(t, dUtRight)
plt.xlabel("Time")
plt.ylabel("u")
# Create a rough trapezoidal rule integration routine
def TrapIntegrate(array, a, b):
n = size(array)
I = 0
for i in range(0, n - 1):
I = I + 0.5 * (array[i] + array[i + 1])
return I * (b - a) / (n - 1)
# Energy via flux
# print("Energy Transport via Fluxes " + str(TrapIntegrate(dUtLeft - dUtRight, 0, TimeEffectivelyZero)))
# Check Initial Condition
def CheckIC():
x = linspace(0, 1, 300)
fig6 = plt.figure()
ax6 = fig6.gca()
for i in [1, 2, 4, 8, 16, 32, 64, 128, 256, 512]:
U0 = U(alpha, x, 0, i)
ax6.plot(x, U0, label="N=" + str(i))
plt.xlabel("x")
plt.ylabel("u")
plt.legend(loc='upper right')
# Plot the x-t-u surface and a contour above the surface
def FastestChange():
dT = ddU(alpha, Xm, Tm, N)
fig7 = plt.figure()
ax7 = fig7.gca(projection='3d')
ax7.plot_surface(Xm, Tm, dT, rstride=2, cstride=2, alpha=0.9, facecolors=cm.jet(Um))
cset = ax7.contour(Xm, Tm, dT, zdir='z', offset=1, cmap=cm.coolwarm)
plt.show()
plt.xlabel("x")
plt.ylabel("time")
# Plot the x-t-u surface and a contour above the surface
def ZeroFlux():
dT = dU(alpha, Xm, Tm, N)
fig7 = plt.figure()
ax7 = fig7.gca()
V = linspace(-1, 1, 10)
cset = ax7.contourf(Xm, Tm, Um, V, cmap=cm.jet)
V2 = [0]
cset = ax7.contour(Xm, Tm, dT, V2, cmap=cm.binary_r)
plt.show()
plt.xlabel("x")
plt.ylabel("time")
Related Questions
Navigate
Integrity-first tutoring: explanations and feedback only — we do not complete graded work. Learn more.