1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 |
import os import numpy as np import matplotlib.pyplot as plt import scipy.sparse as sp import matplotlib as mpl import scipy.sparse.linalg os.chdir(r'D:\projects\wordpress\ex42') os.getcwd() N = 100 u0_t, u0_b = 5, -5 u0_l, u0_r = 3, -1 dx = 1. / (N+1) A_1d = (sp.eye(N, k=-1) + sp.eye(N, k=1) - 4 * sp.eye(N))/dx**2 A = sp.kron(sp.eye(N), A_1d) + (sp.eye(N**2, k=-N) + sp.eye(N**2, k=N))/dx**2 print(A) print(A.nnz * 1.0/ np.prod(A.shape) * 2000) d = np.zeros((N, N)) d[0, :] += -u0_b d[-1, :] += -u0_t d[:, 0] += -u0_l d[:, -1] += -u0_r d = d.reshape(N**2) / dx**2 u = sp.linalg.spsolve(A, d).reshape(N, N) U = np.vstack([np.ones((1, N+2)) * u0_b, np.hstack([np.ones((N, 1)) * u0_l, u, np.ones((N, 1)) * u0_r]), np.ones((1, N+2)) * u0_t]) fig, ax = plt.subplots(1, 1, figsize=(8, 6)) x = np.linspace(0, 1, N+2) X, Y = np.meshgrid(x, x) c = ax.pcolor(X, Y, U, vmin=-5, vmax=5, cmap=mpl.cm.get_cmap('RdBu_r')) cb = plt.colorbar(c, ax=ax) ax.set_xlabel(r"$x_1$", fontsize=18) ax.set_ylabel(r"$x_2$", fontsize=18) cb.set_label(r"$u(x_1, x_2)$", fontsize=18) fig.savefig("ch11-fdm-2d.pdf") fig.tight_layout() plt.savefig("example42.png", dpi=100) plt.show() plt.close() <img class=" wp-image-204 aligncenter" src="https://gantovnik.com/bio-tips/wp-content/uploads/2019/01/example42.png" alt="example42" width="618" height="463" /> x = np.linspace(0, 1, N+2) X, Y = np.meshgrid(x, x) fig = plt.figure(figsize=(12, 5.5)) cmap = mpl.cm.get_cmap('RdBu_r') ax = fig.add_subplot(1, 2, 1) p = ax.pcolor(X, Y, U, vmin=-5, vmax=5, cmap=cmap) ax.set_xlabel(r"$x_1$", fontsize=18) ax.set_ylabel(r"$x_2$", fontsize=18) ax = fig.add_subplot(1, 2, 2, projection='3d') p = ax.plot_surface(X, Y, U, vmin=-5, vmax=5, rstride=3, cstride=3, linewidth=0, cmap=cmap) ax.set_xlabel(r"$x_1$", fontsize=16) ax.set_ylabel(r"$x_2$", fontsize=16) cb = plt.colorbar(p, ax=ax, shrink=0.75) cb.set_label(r"$u(x_1, x_2)$", fontsize=18) plt.savefig("example42_2.png", dpi=100) plt.show() plt.close() <img class=" wp-image-208 aligncenter" src="https://gantovnik.com/bio-tips/wp-content/uploads/2019/01/example42_2.png" alt="example42_2" width="636" height="291" /> |
Recent Comments