import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline
T = np.zeros([100,10])
for i in range(len(T)):
T[i][-1] = 10
i = 0
while i < len(T)-1:
for j in range(1,len(T[0])-1):
T[i+1][j] = T[i][j] +(1./10)*(T[i][j-1] +T[i][j+1] -2*T[i][j])
i += 1
plt.hold(True)
for i in range(len(T)):
plt.plot(T[i])
np.linspace(0,np.pi,10)
array([ 0. , 0.34906585, 0.6981317 , 1.04719755, 1.3962634 , 1.74532925, 2.0943951 , 2.44346095, 2.7925268 , 3.14159265])
T = np.zeros([100,20])
T[0] = np.sin(np.linspace(0,np.pi,20))**2
i = 0
while i < len(T)-1:
for j in range(1,len(T[0])-1):
T[i+1][j] = T[i][j] +(1./10)*(T[i][j-1] +T[i][j+1] -2*T[i][j])
i += 1
plt.hold(True)
for i in range(len(T)):
plt.plot(T[i])
T = np.zeros([100,20])
T[0][1:-1] = 100.
T
array([[ 0., 100., 100., ..., 100., 100., 0.], [ 0., 0., 0., ..., 0., 0., 0.], [ 0., 0., 0., ..., 0., 0., 0.], ..., [ 0., 0., 0., ..., 0., 0., 0.], [ 0., 0., 0., ..., 0., 0., 0.], [ 0., 0., 0., ..., 0., 0., 0.]])
i = 0
while i < len(T)-1:
for j in range(1,len(T[0])-1):
T[i+1][j] = T[i][j] +(1./10)*(T[i][j-1] +T[i][j+1] -2*T[i][j])
i += 1
plt.hold(True)
for i in range(len(T)):
plt.plot(T[i])
from ipywidgets import StaticInteract, RangeWidget, RadioWidget
def update(T,t):
for i in range(0,t):
for j in range(1,len(T[0])-1):
T[i+1][j] = T[i][j] +(1./10)*(T[i][j-1] +T[i][j+1] -2*T[i][j])
return T[t]
def plot(t):
T = np.zeros([10,10])
for i in range(len(T)):
T[i][-1] = 10
fig, ax = plt.subplots(figsize=(4, 3),
subplot_kw={'axisbg':'#EEEEEE',
'axisbelow':True})
ax.grid(color='w', linewidth=2, linestyle='solid')
array = update(T,t)
ax.plot(array, lw=5, alpha=0.4, label = t)
ax.legend(loc='upper right')
return fig
StaticInteract(plot, t=RangeWidget(1., 4., 1.))
<ipywidgets.interact.StaticInteract at 0xb104514c>
for t in range(len(T)):
plt.hold(True)
array = update(T,t)
plt.plot(array)
T = np.zeros([10,10])
T[0][1:-1] = 100
op = np.zeros([8,8])
a = 1./2
for i in range(len(op)):
for j in range(len(op[0])):
if i==j:
op[i][j] = 1+2*a
elif abs(i-j) == 1:
op[i][j] = -a
op
array([[ 2. , -0.5, 0. , 0. , 0. , 0. , 0. , 0. ], [-0.5, 2. , -0.5, 0. , 0. , 0. , 0. , 0. ], [ 0. , -0.5, 2. , -0.5, 0. , 0. , 0. , 0. ], [ 0. , 0. , -0.5, 2. , -0.5, 0. , 0. , 0. ], [ 0. , 0. , 0. , -0.5, 2. , -0.5, 0. , 0. ], [ 0. , 0. , 0. , 0. , -0.5, 2. , -0.5, 0. ], [ 0. , 0. , 0. , 0. , 0. , -0.5, 2. , -0.5], [ 0. , 0. , 0. , 0. , 0. , 0. , -0.5, 2. ]])
opinv = np.linalg.inv(op)
len(opinv)
8
for i in range(1,len(T)):
T[i][1:-1] = np.dot(opinv,T[i-1][1:-1])
plt.hold(True)
for i in range(len(T)):
plt.plot(T[i])
len(T[0][1:-1])
8
I = np.identity(8)
B = np.zeros([8,8])
a = 1./2
for i in range(len(B)):
for j in range(len(B[0])):
if i==j:
B[i][j] = 2
elif abs(i-j) == 1:
B[i][j] = -1
op = np.dot(np.linalg.inv(2*I+a*B),2*I-a*B)
for i in range(1,len(T)):
T[i][1:-1] = np.dot(opinv,T[i-1][1:-1])
plt.hold(True)
for i in range(len(T)):
plt.plot(T[i])
from mpl_toolkits.mplot3d import Axes3D
fig = plt.figure()
ax = plt.axes(projection='3d')
#ax.scatter(x,y,T)
ax.plot_surface(x,y,T)
<mpl_toolkits.mplot3d.art3d.Poly3DCollection at 0xae1bcbec>
x = np.arange(len(T))
y = np.arange(len(T))
x,y = np.meshgrid(x,y)