import numpy as np
import matplotlib.pyplot as plt
from matplotlib import cm
def solve_heat_distribution(known_temperatures, max_iter=1000, tolerance=1e-6,
grid_size=7):
"""
5×5メッシュ上の熱分布を解く関数
Parameters:
-----------
known_temperatures : dict
既知の温度を持つセルの位置と温度値 {(i, j): temperature}
max_iter : int
最大反復回数
tolerance : float
収束判定の閾値
Returns:
--------
temperature_grid : numpy.ndarray
計算された温度分布
iterations : int
収束までの反復回数
"""
# 5×5メッシュの初期化
temperature_grid = np.zeros((grid_size, grid_size))
# 既知の温度値を設定
for (i, j), temp in known_temperatures.items():
temperature_grid[i, j] = temp
# 既知の位置のマスク(True=既知、False=未知)
mask = np.zeros((grid_size, grid_size), dtype=bool)
for i, j in known_temperatures.keys():
mask[i, j] = True
# 反復計算(ガウス・ザイデル法)
iterations = 0
for iter_count in range(max_iter):
old_temperatures = temperature_grid.copy()
# 各グリッドポイントを更新
for i in range(grid_size):
for j in range(grid_size):
if not mask[i, j]: # 未知の温度点のみ更新
# 隣接セルの温度を取得(境界を考慮)
neighbors = []
if i > 0: # 上
neighbors.append(temperature_grid[i-1, j])
if i < grid_size-1: # 下
neighbors.append(temperature_grid[i+1, j])
if j > 0: # 左
neighbors.append(temperature_grid[i, j-1])
if j < grid_size-1: # 右
neighbors.append(temperature_grid[i, j+1])
# 隣接セルの平均を計算
temperature_grid[i, j] = sum(neighbors) / len(neighbors)
# 収束判定
diff = np.max(np.abs(temperature_grid - old_temperatures))
iterations = iter_count + 1
if diff < tolerance:
print(f"収束しました(反復回数: {iterations})")
break
return temperature_grid, iterations
def visualize_temperature(temperature_grid, known_positions=None, grid_size=7):
"""
温度分布を可視化する関数
Parameters:
-----------
temperature_grid : numpy.ndarray
計算された温度分布
known_positions : list
既知の温度を持つセルの位置 [(i, j), ...]
"""
fig, ax = plt.subplots(figsize=(10, 8))
# ヒートマップの描画
im = ax.imshow(temperature_grid, cmap=cm.jet)
# カラーバーの追加
cbar = fig.colorbar(im, ax=ax)
cbar.set_label('温度')
# グリッド線の表示
ax.set_xticks(np.arange(-.5, grid_size, 1), minor=True)
ax.set_yticks(np.arange(-.5, grid_size, 1), minor=True)
ax.grid(which="minor", color="black", linestyle='-', linewidth=2)
# 座標表示
ax.set_xticks(np.arange(0, grid_size, 1))
ax.set_yticks(np.arange(0, grid_size, 1))
# 既知の温度位置をマーク
if known_positions:
for i, j in known_positions:
ax.text(j, i, f"{temperature_grid[i, j]:.1f}", ha="center", va="center",
color="black", fontweight="bold", bbox=dict(facecolor='white', alpha=0.7))
# 値の表示
for i in range(grid_size):
for j in range(grid_size):
if known_positions and (i, j) not in known_positions:
ax.text(j, i, f"{temperature_grid[i, j]:.1f}", ha="center", va="center")
plt.title('5×5メッシュの熱分布')
plt.show()
# 使用例
if __name__ == "__main__":
grid_size = 7
# 既知の温度値を設定(例)
known_temperatures = {
(1, 1): 100, # 左上角100度
(2, 1): 80, # 右上角: 100度
(3, 1): 50, # 右上角: 100度
(4, 1): 80, # 右上角: 100度
(5, 1): 100, # 右上角: 100度
(1, 3): 30, # 左下角: 0度
(2, 3): 40, # 右下角: 0度
# (3, 2): 30, # 中央: 50度
}
# 熱分布を計算
temperature_grid, iterations = solve_heat_distribution(known_temperatures)
# 結果を表示
print("計算された温度分布:")
print(temperature_grid)
print(f"反復回数: {iterations}")
# 結果を可視化
visualize_temperature(temperature_grid, list(known_temperatures.keys()),
grid_size)
# # 直接法(行列解法)による解法
# print("\n直接法による解法:")
# def solve_direct_method(known_temperatures, grid_size=5):
# """行列方程式を解いて熱分布を計算する関数"""
# # 未知の温度ポイントのインデックスを取得
# unknown_indices = []
# for i in range(grid_size):
# for j in range(grid_size):
# if (i, j) not in known_temperatures:
# unknown_indices.append((i, j))
# n_unknowns = len(unknown_indices)
# # 係数行列と右辺ベクトルを構築
# A = np.zeros((n_unknowns, n_unknowns))
# b = np.zeros(n_unknowns)
# # インデックスマッピングを作成
# idx_map = {idx: k for k, idx in enumerate(unknown_indices)}
# # 方程式を構築
# for k, (i, j) in enumerate(unknown_indices):
# # 当該セルに対応する行に1を設定
# A[k, k] = 1
# # 隣接セルの貢献を計算
# neighbors = []
# neighbor_indices = []
# if i > 0: # 上
# neighbors.append((i-1, j))
# if i < grid_size-1: # 下
# neighbors.append((i+1, j))
# if j > 0: # 左
# neighbors.append((i, j-1))
# if j < grid_size-1: # 右
# neighbors.append((i, j+1))
# for ni, nj in neighbors:
# if (ni, nj) in known_temperatures:
# # 既知の境界点からの貢献
# b[k] += known_temperatures[(ni, nj)] / len(neighbors)
# else:
# # 未知点からの貢献
# nk = idx_map[(ni, nj)]
# A[k, nk] = -1.0 / len(neighbors)
# # 連立方程式を解く
# x = np.linalg.solve(A, b)
# # 結果を元のグリッドに戻す
# result = np.zeros((grid_size, grid_size))
# for (i, j), temp in known_temperatures.items():
# result[i, j] = temp
# for k, (i, j) in enumerate(unknown_indices):
# result[i, j] = x[k]
# return result
# # 直接法で解く
# direct_solution = solve_direct_method(known_temperatures)
# print(direct_solution)
# # 直接法の結果を可視化
# plt.figure(figsize=(10, 8))
# plt.title('直接法による5×5メッシュの熱分布')
# visualize_temperature(direct_solution, list(known_temperatures.keys()))