277
u(i, j, 0) = mu(x, y): Next j: Next i
'1.k>0
For k = 0 To M - 1
'1
tk = tau * k: tk1 = tk + tau: tk12 = tk + tau / 2
For j = 1 To Ny - 1: y = hy * j
'1.1
For i = 1 To Nx - 1: x = hx * i
fi(i) = hxy * u(i, j, k) - hx2 * ((u(i, j - 1, k) + u(i, j + 1, k)) / hy2 + f(x, y, tk12))
Next i
'1.2
d = 2 + 2 * hx2 / tau
m10 = (mu_1(y, tk1) + mu_1(y, tk)) / 2
m10 = m10 - (tau / 4) * (mu_1(y - hy, tk1) - 2 * mu_1(y, tk1) + mu_1(y + hy, tk1)) / hy2
m10 = m10 + (tau / 4) * (mu_1(y - hy, tk) - 2 * mu_1(y, tk) + mu_1(y + hy, tk)) / hy2
m2N = (mu_2(y, tk1) + mu_2(y, tk)) / 2
m2N = m2N - (tau / 4) * (mu_2(y - hy, tk1) - 2 * mu_2(y, tk1) + mu_2(y + hy, tk1)) / hy2
m2N = m2N + (tau / 4) * (mu_2(y - hy, tk) - 2 * mu_2(y, tk) + mu_2(y + hy, tk)) / hy2
a(1) = 1 / d: b(1) = (fi(1) - m10) / (-d)
For i = 2 To Nx - 2
a(i) = -1 / (a(i - 1) - d):
b(i) = (fi(i) - b(i - 1)) / (a(i - 1) - d)
Next i
'1.3
uk12(Nx - 1, j) = (fi(Nx - 1) - m2N - b(Nx - 2)) / (a(Nx - 2) - d)
For i = Nx - 2 To 1 Step -1: uk12(i, j) = a(i) * uk12(i + 1, j) + b(i):
Next i
uk12(0, j) = m10: uk12(Nx, j) = m2N
Next j
'2
For i = 1 To Nx - 1: x = hx * i
'2.1
For j = 1 To Ny - 1: y = j * hy
fi(j) = hyx * uk12(i, j) - hy2 * ((uk12(i - 1, j) + uk12(i + 1, j)) / hx2 + f(x, y, tk12))
Next j
'2.2
d = 2 + 2 * hy2 / tau
a(1) = 1 / d: b(1) = (fi(1) - mu_3(x, tk1)) / (-d)
For j = 2 To Ny - 2
a(j) = -1 / (a(j - 1) - d): b(j) = (fi(j) - b(j - 1)) / (a(j - 1) - d)
Next j
'2.3
u(i, Ny - 1, k + 1) = (fi(Ny - 1) - mu_4(x, tk1) - b(Ny - 2)) / (a(Ny - 2) - d)
For j = Ny - 2 To 1 Step -1: u(i, j, k + 1) = a(j) * u(i, j + 1, k + 1) + b(j): Next j
Next i
For j = 0 To Ny: y = hy * j
u(0, j, k + 1) = mu_1(y, tk1): u(Nx, j, k + 1) = mu_2(y, tk1)
Next j
For i = 0 To Nx: x = hx * i
u(i, 0, k + 1) = mu_3(x, tk1): u(i, Ny, k + 1) = mu_4(x, tk1)
Next i
