1. 程式人生 > >對Laplace方程的學習——來自流沙公眾號

對Laplace方程的學習——來自流沙公眾號

對Laplace方程的學習
在這裡插入圖片描述
程式碼如下:

import numpy as np
import matplotlib.pyplot as plt
import time

dx = 0.1
dy = 0.1
dx2 = dx*dx
dy2 = dy*dy

def laplace(u):
	nx,ny = u.shape   #檢視陣列的維數
	for i in range(1, nx-1):
		for j in range(1,ny-1):
			u[i,j] = ((u[i+1,j] + u[i-1,j])*dy2 + (u[i,j+1] + u[i,j-1]*dx2)) / (2*(dx2+dy2))

def mat_laplace(u):
	u[1:-1,1:-1] = ((u[2:,1:-1] + u[:-2,1:-1])*dy2 + (u[1:-1,2:] + u[1:-1,:-2])*dx2) / (2*(dx2+dy2))

def calc(N,Niter):
	u = np.zeros([N,N])  # 全域性初始化
	u[0] = 1   #邊界條件
	for i in range(Niter):
		laplace(u)
	return u

def mat_calc(N,Niter):
	u = np.zeros([N,N])  # 全域性初始化
	u[0] = 1   #邊界條件
	for i in range(Niter):
		mat_laplace(u)
	return u

# 執行測試
start = time.clock()
u = calc(100,3000)
end = time.clock()
time_elapsed = end - start   # 普通計算耗時

start = time.clock()
u = mat_calc(100,3000)
end = time.clock()
time_elapsed1 = end - start   # 普通計算耗時

plt.pcolormesh(u)
plt.show()
plt.legend()
print(time_elapsed)
print(time_elapsed1)

疑惑:
執行後編譯器出現了一些警告,暫無法明白。

D:\python\py37\python.exe E:/python/CFD/laplace_equation.py
E:/python/CFD/laplace_equation.py:34: DeprecationWarning: time.clock has been deprecated in Python 3.3 and will be removed from Python 3.8: use time.perf_counter or time.process_time instead
  start = time.clock()
E:/python/CFD/laplace_equation.py:14: RuntimeWarning: overflow encountered in double_scalars
  u[i,j] = ((u[i+1,j] + u[i-1,j])*dy2 + (u[i,j+1] + u[i,j-1]*dx2)) / (2*(dx2+dy2))
E:/python/CFD/laplace_equation.py:36: DeprecationWarning: time.clock has been deprecated in Python 3.3 and will be removed from Python 3.8: use time.perf_counter or time.process_time instead
  end = time.clock()
E:/python/CFD/laplace_equation.py:39: DeprecationWarning: time.clock has been deprecated in Python 3.3 and will be removed from Python 3.8: use time.perf_counter or time.process_time instead
  start = time.clock()
E:/python/CFD/laplace_equation.py:41: DeprecationWarning: time.clock has been deprecated in Python 3.3 and will be removed from Python 3.8: use time.perf_counter or time.process_time instead
  end = time.clock()

結果:
在這裡插入圖片描述
本文僅用作記錄CFD學習過程