0

我正在尝试使用有限差分法在 python 中进行迭代计算。我从中找到了有限差分法：

http://depa.fquim.unam.mx/amyd/archivero/DiferenciasFinitas3_25332.pdf

代码计算的值是正确的。问题是它只显示最终值。我想要的是提取 x 方向上任何点的值，以便我可以绘制它们，以及提取任何时间点的值，例如计算中途点的值。这是进行迭代计算的正确方法吗？代码如下所示：

import numpy as np
import scipy as sp
import time
import matplotlib as p
L=0.005
Nx=3
T=5
N1=5
k=0.5
rho=1200
c=1000
a=(k/(rho*c))
x = np.linspace(0, L, Nx+1)   # mesh points in space
dx = x[1] - x[0]
t = np.linspace(0, T, N1)    # time
dt = t[1] - t[0]
toutside=5
Coefficient = a*dt/dx**2
bi=0.5
ui = sp.zeros(Nx+1)
u = sp.zeros(Nx+1)
for i in range(Nx+1):
  ui[i] = 50 # initial values

for n in range(0, N1):

   for i in range(0,1):
      u[i] = 2*Coefficient*(ui[i+1]+bi*toutside)+(1-2*Coefficient-2*bi*Coefficient)*ui[i]
   for i in range(1,Nx):
      u[i] = Coefficient*(ui[i+1]+ui[i-1])+(1-2*Coefficient)*ui[i]
   for i in range(Nx,Nx+1):
      u[i] = 2*Coefficient*(ui[i-1])+(1-2*Coefficient)*ui[i]
      ui[:]= u #updates matrix for next loop
 print ui

我已经根据 danodonovan 对以下问题的回答修改了我的代码：

for n in range(0, N1):

for i in range(0,1):
    u[i] = 2*Coefficient*(ui[i+1]+bi*toutside)+(1-2*Coefficient-2*bi*Coefficient)*ui[i]
for i in range(1,Nx):
    u[i] = Coefficient*(ui[i+1]+ui[i-1])+(1-2*Coefficient)*ui[i]
for i in range(Nx,Nx+1):
    u[i] = 2*Coefficient*(ui[i-1])+(1-2*Coefficient)*ui[i]
    ui=u
    a=list(ui)

print a

当我尝试将整个列表从循环中取出时，只会产生最终值。如何提取整个列表？这是使用前一行的值来计算新行的值进行迭代计算的正确方法吗？

4

2 回答 2

0

首先缩进print ui以查看uifor all的值N1。然后将结果附加到列表中：

res = []
for n in range(0, N1):

   for i in range(0,1):
      u[i] = 2*Coefficient*(ui[i+1]+bi*toutside)+(1-2*Coefficient-2*bi*Coefficient)*ui[i]
   for i in range(1,Nx):
      u[i] = Coefficient*(ui[i+1]+ui[i-1])+(1-2*Coefficient)*ui[i]
   for i in range(Nx,Nx+1):
      u[i] = 2*Coefficient*(ui[i-1])+(1-2*Coefficient)*ui[i]
      ui[:]= u #updates matrix for next loop
   print ui
   res.append(ui.copy())

print res

产生这个结果：

[array([ 41.5625,  50.    ,  50.    ,  50.    ]),
 array([ 37.87109375,  48.41796875,  50.        ,  50.        ]),
 array([ 35.6628418 ,  46.73706055,  49.70336914,  50.        ]),
 array([ 34.06639099,  45.21682739,  49.20280457,  49.88876343]),
 array([ 32.79785633,  43.87349129,  48.58405113,  49.63152885])]
于 2013-05-30T14:14:04.207 回答
0

（如果我正确理解你的问题）使用matplotlib你可以做

 import matplotlib as p

 # and after your loop has completed

 p.pyplot.plot(range(0, N1), ui, 'o-')
 p.pyplot.show()

得到一个简单的数据u图range(0, N1)。

我不知道您期望ui成为什么ui[:]= u是一件奇怪的事情-它设置了 to 的副本，ui 但您u不保留正在生成的副本。uiui[:]

提示：认为ui[:]是一样的list(ui)

于 2013-05-30T13:52:52.743 回答