下面的代码求解一维热方程,该方程代表了一根棒,其端部在初始条件10 * np.sin(np.pi * x)下保持零温度。
Dirichlet边界条件(两端为零温度)如何用于此计算?有人告诉我矩阵A的上,下两行包含两个非零元素,而缺少的第三个元素is Dirichlet条件。但是我不知道这种情况会通过哪种机制影响计算。在A中缺少元素的情况下,u_ {0}或u_ {n}如何为零?
下面的有限差分法使用Crank-Nicholson。
import numpy as np
import scipy.linalg
# Number of internal points
N = 200
# Calculate Spatial Step-Size
h = 1/(N+1.0)
k = h/2
x = np.linspace(0,1,N+2)
x = x[1:-1] # get rid of the '0' and '1' at each end
# Initial Conditions
u = np.transpose(np.mat(10*np.sin(np.pi*x)))
# second derivative matrix
I2 = -2*np.eye(N)
E = np.diag(np.ones((N-1)), k=1)
D2 = (I2 + E + E.T)/(h**2)
I = np.eye(N)
TFinal = 1
NumOfTimeSteps = int(TFinal/k)
for i in range(NumOfTimeSteps):
# Solve the System: (I - k/2*D2) u_new = (I + k/2*D2)*u_old
A = (I - k/2*D2)
b = np.dot((I + k/2*D2), u)
u = scipy.linalg.solve(A, b)
让我们看一个简单的例子。我们假设N = 3
,即三个内部点,但我们首先还将在矩阵D2
中包括描述近似二阶导数的边界点:
1 / 1 -2 1 0 0 \
D2 = --- | 0 1 -2 1 0 |
h^2 \ 0 0 1 -2 1 /
[第一行表示x_1
处的近似二阶导数为1/h^2 * (u_0 - 2*u_1 + u_2)
。我们知道u_0 = 0
,但是由于Dirichlet边界条件是同质的,因此我们可以简单地从方程中省略它,并且对于矩阵e可以得到相同的结果
1 / 0 -2 1 0 0 \
D2 = --- | 0 1 -2 1 0 |
h^2 \ 0 0 1 -2 0 /
由于u_0
和u_{n+1}
不是真正的未知数-已知它们为零-我们可以将它们从矩阵中完全删除,然后得到
1 / 2 1 0 \
D2 = --- | 1 -2 1 |
h^2 \ 0 1 -2 /
矩阵中缺少的条目确实对应于边界条件为零的事实。
有人告诉我,矩阵A包含两个非零元素,以及缺少的第三部分元素(即零)是Dirichlet条件。
我不确定您被告知了什么,但这是我对这个问题的了解。
狄利克雷边界条件确定电位值(在这种情况下为温度)。诺伊曼边界条件将指定点处的通量或一阶导数。如果表面上有对流边界条件,则需要此。
关于处理Dirichlet边界条件,您在不首先考虑边界条件的情况下制定系统矩阵。如果给定节点的温度固定,则可以通过以下方式进行处理:
我这个答案也没有任何意义,我们要求的是在计算过程中固定值不变,这会影响整个解决方案和整个结果吗?那么如何使值不变?在此过程中,开头和结尾的零不应该改变?怎么样?