ラプラス方程式Pythonコード
ラプラス方程式(編集中)
ラプラス方程式とは、2階の線型楕円型偏微分方程式のことになります。領域内においてある境界条件を満たすラプラス方程式を求め、それによりさまざまな解析解を導くことが可能です。
ここでは簡単な例として長方形プレートの平衡温度分布に関して、2次元のラプラス方程式で導き出した解をPythonの3次元描画によって表現します。
ラプラス演算子
ひとまずラプラス方程式に関しての簡単な予備知識を考察していきます。
それぞれの座標とした3次元座標空間において2階の偏微分作用素を
とし、この作用素を次のようにおきます。
![](https://akjp.org/wp-content/uploads/2024/05/laplace_operator_equiv_delta_large_img.png)
ここでナブラと呼ばれる次のもの、
![](https://akjp.org/wp-content/uploads/2024/05/nabla_2024img.png)
これを使って表せば上記式の2階偏微分作用素式の右辺は次のようなものになります。
![](https://akjp.org/wp-content/uploads/2024/05/nabla_square_laplace_operator_img_img2.png)
ラプラス方程式とは、考えうるそれぞれの次元において付与される作用素の2階微分に関して、そのラプラス作用素を用いて表現すると次のようなものになります。
![](https://akjp.org/wp-content/uploads/2024/05/laplace_equation_default_img1.png)
![](https://akjp.org/wp-content/uploads/2024/05/laplaces_equation_default_img.png)
上記式をみてわかるようにラプラス方程式には時間に当たる変数が含まれていません。それにより時間によって変化しない定常状態を表します。時間を反映した変数がないため、このために初期条件はなく、次に示すような境界条件だけが必要となります。
ディリクレ境界条件
ここでは長方形プレートとみなした場合の以下のようなもの、
![](https://akjp.org/wp-content/uploads/2024/05/dirichlet_boundary_condition_img.png)
とした、2次元でのラプラス方程式の考察をしていきます。
境界条件として次のようなものとします。
![](https://akjp.org/wp-content/uploads/2024/05/equilibruim_temperature_distribution_for_a_rectangular_plate_boudary_condition_img1.png)
![](https://akjp.org/wp-content/uploads/2024/05/laplace_equation_boundary_condition_202405img1.jpg)
図示すると次のようなものになります。
このような境界条件から変数分離法などを使って最終的に求まる解が以下のようになります。
![](https://akjp.org/wp-content/uploads/2024/05/laplaces_equation_separation_partial_differential_equation_calculation_result_unxy_202404img-1.png)
ただし上記式のに関しては以下のようになります。
![](https://akjp.org/wp-content/uploads/2024/05/laplace_eq_applied_fourier_sine_series_img2-1.png)
これの導出に関してはこちらを参照ください。
Pythonコード
インポートするのは以下になります。
import numpy
from matplotlib import pyplot
from matplotlib import rcParams
from mpl_toolkits.mplot3d import Axes3D
from matplotlib import cm
次にfigureインスタンスを作成して、add_subplot()メソッドで3D用の軸を追加します。add_subplot()メソッドの第1引数の「111」は、行、列、インデックスを表しています。次にように追記します。
def plot_3D(x, y, p):
fig = pyplot.figure(figsize=(20,10), dpi=80)
ax = fig.add_subplot(111, projection='3d')
X,Y = numpy.meshgrid(x,y)
surf = ax.plot_surface(X,Y,p[:], rstride=1, cstride=1, cmap=cm.gist_heat, linewidth=0, antialiased=False)
matplotlibでグラフの描画範囲の設定
次のようにして範囲指定をします。
ax.set_xlim(0,1)
ax.set_ylim(0,1)
ax.set_xlabel('$x$')
ax.set_ylabel('$y$')
ax.set_zlabel('$z$')
ax.view_init(20,40)