ベクトル解析

import numpy as np
import matplotlib.pyplot as plt

# スカラー場の定義
def scalar_field(x, y):
    return np.sin(x) * np.cos(y)

# 勾配の計算
def gradient(x, y):
    dx, dy = np.gradient(scalar_field(x, y), x, y)
    return dx, dy

# グリッドの作成
x = np.linspace(-2*np.pi, 2*np.pi, 100)
y = np.linspace(-2*np.pi, 2*np.pi, 100)
X, Y = np.meshgrid(x, y)

# スカラー場と勾配の計算
Z = scalar_field(X, Y)
dx, dy = gradient(X, Y)

# プロット
plt.figure(figsize=(10, 6))
plt.contourf(X, Y, Z, levels=50, cmap='viridis', alpha=0.7)
plt.colorbar(label='Scalar Field Value')
plt.quiver(X, Y, dx, dy, color='red', scale=50, width=0.002)
plt.title('Gradient of Scalar Field')
plt.xlabel('x')
plt.ylabel('y')
plt.grid()
plt.show()

import numpy as np
import matplotlib.pyplot as plt

# ベクトル場の定義
def vector_field(x, y):
    F_x = -y
    F_y = x
    return F_x, F_y

# 発散の計算
def divergence(x, y):
    F_x, F_y = vector_field(x, y)
    dF_x_dx, dF_y_dy = np.gradient(F_x, x), np.gradient(F_y, y)
    return dF_x_dx + dF_y_dy

# グリッドの作成
x = np.linspace(-2, 2, 100)
y = np.linspace(-2, 2, 100)
X, Y = np.meshgrid(x, y)

# ベクトル場と発散の計算
F_x, F_y = vector_field(X, Y)
div_F = divergence(X, Y)

# プロット
plt.figure(figsize=(10, 6))
plt.contourf(X, Y, div_F, levels=50, cmap='coolwarm', alpha=0.7)
plt.colorbar(label='Divergence')
plt.quiver(X, Y, F_x, F_y, color='black', scale=30, width=0.002)
plt.title('Divergence of Vector Field')
plt.xlabel('x')
plt.ylabel('y')
plt.grid()
plt.show()

import numpy as np
import matplotlib.pyplot as plt

# ベクトル場の定義
def vector_field(x, y):
    F_x = -y
    F_y = x
    return F_x, F_y

# 回転の計算
def curl(x, y):
    F_x, F_y = vector_field(x, y)
    dF_y_dx, dF_x_dy = np.gradient(F_y, x), np.gradient(F_x, y)
    return dF_y_dx - dF_x_dy

# グリッドの作成
x = np.linspace(-2, 2, 100)
y = np.linspace(-2, 2, 100)
X, Y = np.meshgrid(x, y)

# ベクトル場と回転の計算
F_x, F_y = vector_field(X, Y)
curl_F = curl(X, Y)

# プロット
plt.figure(figsize=(10, 6))
plt.contourf(X, Y, curl_F, levels=50, cmap='coolwarm', alpha=0.7)
plt.colorbar(label='Curl')
plt.quiver(X, Y, F_x, F_y, color='black', scale=30, width=0.002)
plt.title('Curl of Vector Field')
plt.xlabel('x')
plt.ylabel('y')
plt.grid()
plt.show()

import numpy as np
import matplotlib.pyplot as plt

# スカラー場とベクトル場の定義
def scalar_field(x, y):
    return np.sin(x) * np.cos(y)

def vector_field(x, y):
    F_x = -y
    F_y = x
    return F_x, F_y

# 勾配の計算
def gradient(x, y):
    dx, dy = np.gradient(scalar_field(x, y), x, y)
    return dx, dy

# グリッドの作成
x = np.linspace(-2*np.pi, 2*np.pi, 100)
y = np.linspace(-2*np.pi, 2*np.pi, 100)
X, Y = np.meshgrid(x, y)

# スカラー場、ベクトル場、勾配の計算
Z = scalar_field(X, Y)
F_x, F_y = vector_field(X, Y)
dx, dy = gradient(X, Y)

# プロット
fig, ax = plt.subplots(1, 2, figsize=(14, 6))

# スカラー場のプロット
contour = ax[0].contourf(X, Y, Z, levels=50, cmap='viridis', alpha=0.7)
ax[0].quiver(X, Y, dx, dy, color='red', scale=50, width=0.002)
fig.colorbar(contour, ax=ax[0], label='Scalar Field Value')
ax[0].set_title('Scalar Field and its Gradient')
ax[0].set_xlabel('x')
ax[0].set_ylabel('y')

# ベクトル場のプロット
div_F = divergence(X, Y)
ax[1].contourf(X, Y, div_F, levels=50, cmap='coolwarm', alpha=0.7)
ax[1].quiver(X, Y, F_x, F_y, color='black', scale=30, width=0.002)
fig.colorbar(contour, ax=ax[1], label='Divergence')
ax[1].set_title('Vector Field and its Divergence')
ax[1].set_xlabel('x')
ax[1].set_ylabel('y')

plt.tight_layout()
plt.show()

import numpy as np
import matplotlib.pyplot as plt

# ベクトル場の定義（流体の速度場）
def velocity_field(x, y):
    U = -y
    V = x
    return U, V

# 発散の計算
def divergence(x, y):
    U, V = velocity_field(x, y)
    dU_dx, dV_dy = np.gradient(U, x), np.gradient(V, y)
    return dU_dx + dV_dy

# グリッドの作成
x = np.linspace(-2, 2, 100)
y = np.linspace(-2, 2, 100)
X, Y = np.meshgrid(x, y)

# ベクトル場と発散の計算
U, V = velocity_field(X, Y)
div = divergence(X, Y)

# プロット
plt.figure(figsize=(12, 6))

# ベクトル場のプロット
plt.subplot(1, 2, 1)
plt.quiver(X, Y, U, V, scale=30, width=0.002)
plt.title('Velocity Field')
plt.xlabel('x')
plt.ylabel('y')
plt.grid()

# 発散のプロット
plt.subplot(1, 2, 2)
plt.contourf(X, Y, div, levels=50, cmap='coolwarm', alpha=0.7)
plt.colorbar(label='Divergence')
plt.title('Divergence of Velocity Field')
plt.xlabel('x')
plt.ylabel('y')
plt.grid()

plt.tight_layout()
plt.show()

import numpy as np
import matplotlib.pyplot as plt

# 層状ベクトル場の定義
def layer_vector_field(x, y):
    U = np.ones_like(x)  # 流れは x 方向に一様
    V = np.zeros_like(y) # y 方向には流れがない
    return U, V

# グリッドの作成
x = np.linspace(-2, 2, 100)
y = np.linspace(-2, 2, 100)
X, Y = np.meshgrid(x, y)

# 層状ベクトル場の計算
U, V = layer_vector_field(X, Y)

# プロット
plt.figure(figsize=(6, 6))
plt.quiver(X, Y, U, V, scale=10, width=0.002)
plt.title('Layer Vector Field')
plt.xlabel('x')
plt.ylabel('y')
plt.grid()
plt.show()

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

# 管状ベクトル場の定義
def tube_vector_field(x, y, z):
    U = -y
    V = x
    W = np.zeros_like(z)
    return U, V, W

# グリッドの作成
x = np.linspace(-2, 2, 10)
y = np.linspace(-2, 2, 10)
z = np.linspace(-1, 1, 5)
X, Y, Z = np.meshgrid(x, y, z)

# 管状ベクトル場の計算
U, V, W = tube_vector_field(X, Y, Z)

# プロット
fig = plt.figure(figsize=(10, 8))
ax = fig.add_subplot(111, projection='3d')
ax.quiver(X, Y, Z, U, V, W, length=0.2, normalize=True)
ax.set_title('Tube Vector Field')
ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_zlabel('z')
plt.show()

import numpy as np
import matplotlib.pyplot as plt

# 層状ベクトル場の定義（例：ポテンシャル関数が存在する流れ場）
def layer_potential_field(x, y):
    # ベクトル場が一定である場合、ポテンシャル関数の勾配はベクトル場になります
    phi = -y
    return phi

# グリッドの作成
x = np.linspace(-2, 2, 100)
y = np.linspace(-2, 2, 100)
X, Y = np.meshgrid(x, y)

# ポテンシャル関数の計算
phi = layer_potential_field(X, Y)

# プロット
plt.figure(figsize=(6, 6))
plt.contourf(X, Y, phi, levels=50, cmap='plasma', alpha=0.7)
plt.colorbar(label='Potential Function Value')
plt.title('Layer Potential Field')
plt.xlabel('x')
plt.ylabel('y')
plt.grid()
plt.show()

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

# ベクトル場の定義
def vector_field(x, y, z):
    U = np.sin(y) * np.cos(z)
    V = np.sin(z) * np.cos(x)
    W = np.sin(x) * np.cos(y)
    return U, V, W

# グリッドの作成
x = np.linspace(-2, 2, 10)
y = np.linspace(-2, 2, 10)
z = np.linspace(-2, 2, 10)
X, Y, Z = np.meshgrid(x, y, z)

# ベクトル場の計算
U, V, W = vector_field(X, Y, Z)

# 3Dプロット
fig = plt.figure(figsize=(12, 8))
ax = fig.add_subplot(111, projection='3d')
ax.quiver(X, Y, Z, U, V, W, length=0.5, normalize=True, color='b')
ax.set_title('Vector Field and Divergence')
ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_zlabel('Z')
plt.show()

import numpy as np
import matplotlib.pyplot as plt

# スカラー場ポテンシャルの定義
def potential(x, y):
    return np.sin(x) * np.cos(y)

# グリッドの作成
x = np.linspace(-2*np.pi, 2*np.pi, 100)
y = np.linspace(-2*np.pi, 2*np.pi, 100)
X, Y = np.meshgrid(x, y)

# ポテンシャルの計算
Z = potential(X, Y)

# プロット
plt.figure(figsize=(8, 6))
plt.contourf(X, Y, Z, levels=50, cmap='viridis', alpha=0.7)
plt.colorbar(label='Potential')
plt.title('Potential Function')
plt.xlabel('x')
plt.ylabel('y')
plt.grid()
plt.show()

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

# ベクトル・ポテンシャルの定義
def vector_potential(x, y, z):
    A_x = -y
    A_y = x
    A_z = np.zeros_like(z)
    return A_x, A_y, A_z

# グリッドの作成
x = np.linspace(-2, 2, 10)
y = np.linspace(-2, 2, 10)
z = np.linspace(-2, 2, 10)
X, Y, Z = np.meshgrid(x, y, z)

# ベクトル・ポテンシャルの計算
A_x, A_y, A_z = vector_potential(X, Y, Z)

# 3Dプロット
fig = plt.figure(figsize=(12, 8))
ax = fig.add_subplot(111, projection='3d')
ax.quiver(X, Y, Z, A_x, A_y, A_z, length=0.5, normalize=True, color='r')
ax.set_title('Vector Potential')
ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_zlabel('Z')
plt.show()

import numpy as np
import matplotlib.pyplot as plt

# スカラー場の定義
def scalar_field(x, y):
    return np.sin(x) * np.cos(y)

# 直交座標系でのスカラー場
x = np.linspace(-2*np.pi, 2*np.pi, 50)
y = np.linspace(-2*np.pi, 2*np.pi, 50)
X, Y = np.meshgrid(x, y)
Z = scalar_field(X, Y)

# 円筒座標系に変換
R = np.sqrt(X**2 + Y**2)
Theta = np.arctan2(Y, X)
Z_cylindrical = Z  # z座標のスカラー場

# スカラー場の微分 (直交座標系)
dphi_dx, dphi_dy = np.gradient(Z, x, y)

# 円筒座標系での微分
dphi_dr = np.gradient(Z, R, axis=0)
dphi_dtheta = np.gradient(Z, Theta, axis=1)

# プロット
fig, axs = plt.subplots(2, 2, figsize=(12, 10))

# 直交座標系のスカラー場
c1 = axs[0, 0].contourf(X, Y, Z, levels=50, cmap='viridis')
fig.colorbar(c1, ax=axs[0, 0])
axs[0, 0].set_title('Scalar Field in Cartesian Coordinates')

# スカラー場の微分
axs[0, 1].quiver(X, Y, dphi_dx, dphi_dy, color='r')
axs[0, 1].set_title('Gradient of Scalar Field (Cartesian)')

# 円筒座標系のスカラー場
axs[1, 0].contourf(R, Theta, Z_cylindrical, levels=50, cmap='viridis')
axs[1, 0].set_title('Scalar Field in Cylindrical Coordinates')

# 円筒座標系での微分
axs[1, 1].quiver(R, Theta, dphi_dr, dphi_dtheta, color='b')
axs[1, 1].set_title('Gradient of Scalar Field (Cylindrical)')

plt.tight_layout()
plt.show()

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.colors import Normalize

# 円柱座標系のベクトル場を定義
def velocity_field(r, theta):
    u_r = np.zeros_like(r)  # 半径方向の速度
    u_theta = np.ones_like(r)  # 角度方向の速度（単位円）
    return u_r, u_theta

# 円柱座標系のグリッド作成
r = np.linspace(0.1, 2, 20)  # 半径方向（0で割り算を防ぐため0.1から開始）
theta = np.linspace(0, 2*np.pi, 20)  # 角度方向
R, Theta = np.meshgrid(r, theta)

# 速度場の計算
U_r, U_theta = velocity_field(R, Theta)

# カルティナル座標系のプロット
X = R * np.cos(Theta)  # x座標
Y = R * np.sin(Theta)  # y座標

# プロット
fig, ax = plt.subplots(figsize=(8, 8))
quiver = ax.quiver(X, Y, U_r * np.cos(Theta) - U_theta * np.sin(Theta), 
                    U_r * np.sin(Theta) + U_theta * np.cos(Theta), 
                    angles='xy', scale_units='xy', scale=1, color='b')

# スケールとタイトル
ax.set_aspect('equal')
ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_title('Velocity Field in Cylindrical Coordinates')

plt.grid(True)
plt.show()

import numpy as np
import matplotlib.pyplot as plt

# ベクトル場の定義
def vector_field(x, y):
    return -y, x

# 曲線の定義
t = np.linspace(0, 2 * np.pi, 100)
x_curve = np.cos(t)
y_curve = np.sin(t)

# ベクトル場の計算
X, Y = np.meshgrid(np.linspace(-2, 2, 20), np.linspace(-2, 2, 20))
U, V = vector_field(X, Y)

# プロット設定
fig, axs = plt.subplots(1, 2, figsize=(14, 7))

# ベクトル場と曲線
axs[0].quiver(X, Y, U, V, color='b')
axs[0].plot(x_curve, y_curve, color='r', label='Curve C')
axs[0].set_title('Vector Field with Curve C')
axs[0].set_xlabel('X')
axs[0].set_ylabel('Y')
axs[0].legend()

# 面積分のプロット（単位面積に沿ったベクトル場のフラックス）
X, Y = np.meshgrid(np.linspace(-1, 1, 20), np.linspace(-1, 1, 20))
U, V = vector_field(X, Y)
Z = np.zeros_like(X)

# 面の定義
X_surface, Y_surface = np.meshgrid(np.linspace(-1, 1, 10), np.linspace(-1, 1, 10))

# プロット設定
contour = axs[1].quiver(X_surface, Y_surface, U, V, scale=20, color='g')
axs[1].set_title('Surface Integral of Vector Field')
axs[1].set_xlabel('X')
axs[1].set_ylabel('Y')

plt.tight_layout()
plt.show()

import numpy as np
import matplotlib.pyplot as plt
from scipy.linalg import lapack

# Define a 2D grid
x = np.linspace(-5, 5, 100)
y = np.linspace(-5, 5, 100)
X, Y = np.meshgrid(x, y)

# Define a function f(x, y)
Z = np.sin(np.sqrt(X**2 + Y**2))

# Calculate Laplacian
def laplacian(Z, dx):
    Z = np.pad(Z, pad_width=1, mode='constant', constant_values=0)
    Z_laplacian = np.zeros_like(Z)
    Z_laplacian[1:-1, 1:-1] = (Z[:-2, 1:-1] + Z[2:, 1:-1] + Z[1:-1, :-2] + Z[1:-1, 2:] - 4 * Z[1:-1, 1:-1]) / dx**2
    return Z_laplacian[1:-1, 1:-1]

dx = x[1] - x[0]
Z_laplacian = laplacian(Z, dx)

# Plot the result
plt.figure(figsize=(12, 6))
plt.subplot(1, 2, 1)
plt.title('Original Function')
plt.contourf(X, Y, Z, cmap='viridis')
plt.colorbar()

plt.subplot(1, 2, 2)
plt.title('Laplacian')
plt.contourf(X, Y, Z_laplacian, cmap='viridis')
plt.colorbar()

plt.show()

import numpy as np
import matplotlib.pyplot as plt

# Define a vector field
x = np.linspace(-2, 2, 20)
y = np.linspace(-2, 2, 20)
X, Y = np.meshgrid(x, y)
U = -Y
V = X

# Define a scalar potential field
phi = np.sin(np.sqrt(X**2 + Y**2))

# Calculate gradient of scalar potential
U_phi, V_phi = np.gradient(phi)

plt.figure(figsize=(10, 10))
plt.quiver(X, Y, U, V, color='blue', alpha=0.5, label='Vector Field')
plt.quiver(X, Y, U_phi, V_phi, color='red', alpha=0.5, label='Gradient of Scalar Field')
plt.title('Vector Field and Its Gradient')
plt.xlabel('x')
plt.ylabel('y')
plt.legend()

plt.show()

import numpy as np
import matplotlib.pyplot as plt

# Define a vector field
x = np.linspace(-2, 2, 20)
y = np.linspace(-2, 2, 20)
X, Y = np.meshgrid(x, y)
U = -Y
V = X

# Define a contour (e.g., circle)
theta = np.linspace(0, 2*np.pi, 100)
R = 1
Cx = R * np.cos(theta)
Cy = R * np.sin(theta)

plt.figure(figsize=(10, 10))
plt.quiver(X, Y, U, V, color='blue', alpha=0.5)
plt.plot(Cx, Cy, color='red')
plt.title('Vector Field with Contour (Circle)')
plt.xlabel('x')
plt.ylabel('y')

plt.show()

import numpy as np
import matplotlib.pyplot as plt

# Define a vector field
x = np.linspace(-5, 5, 20)
y = np.linspace(-5, 5, 20)
X, Y = np.meshgrid(x, y)
U = -Y
V = X

# Calculate curl (z-component)
curl = np.gradient(V, axis=0) - np.gradient(U, axis=1)

# Plot the vector field
plt.figure(figsize=(8, 8))
plt.quiver(X, Y, U, V)
plt.title('Vector Field')

plt.figure(figsize=(8, 8))
plt.contourf(X, Y, curl, cmap='viridis')
plt.colorbar()
plt.title('Curl of the Vector Field')

plt.show()

import numpy as np
import matplotlib.pyplot as plt

# Define a grid in Cartesian coordinates
x = np.linspace(-5, 5, 20)
y = np.linspace(-5, 5, 20)
X, Y = np.meshgrid(x, y)
U = -Y
V = X

# Plot the vector field
plt.figure(figsize=(8, 8))
plt.quiver(X, Y, U, V, scale=5)
plt.title('Vector Field in Cartesian Coordinates')
plt.xlabel('x')
plt.ylabel('y')
plt.grid()
plt.show()

import numpy as np
import matplotlib.pyplot as plt

# Define a grid in cylindrical coordinates
r = np.linspace(0, 5, 20)
theta = np.linspace(0, 2*np.pi, 20)
R, Theta = np.meshgrid(r, theta)
X = R * np.cos(Theta)
Y = R * np.sin(Theta)

# Define a vector field in cylindrical coordinates
U = -R * np.sin(Theta)
V = R * np.cos(Theta)

# Plot the vector field
plt.figure(figsize=(8, 8))
plt.quiver(X, Y, U, V, scale=5)
plt.title('Vector Field in Cylindrical Coordinates')
plt.xlabel('x')
plt.ylabel('y')
plt.grid()
plt.show()

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

# Define a grid in spherical coordinates
r = np.linspace(0, 5, 10)
theta = np.linspace(0, 2*np.pi, 20)
phi = np.linspace(0, np.pi, 10)
R, Theta, Phi = np.meshgrid(r, theta, phi)

X = R * np.sin(Phi) * np.cos(Theta)
Y = R * np.sin(Phi) * np.sin(Theta)
Z = R * np.cos(Phi)

# Define a vector field in spherical coordinates
U = np.sin(Phi) * np.cos(Theta)
V = np.sin(Phi) * np.sin(Theta)
W = np.cos(Phi)

# Plot the vector field
fig = plt.figure(figsize=(12, 8))
ax = fig.add_subplot(111, projection='3d')
ax.quiver(X, Y, Z, U, V, W, length=0.5)
ax.set_title('Vector Field in Spherical Coordinates')
ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_zlabel('z')
plt.show()

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

# Define a grid
x = np.linspace(-2, 2, 20)
y = np.linspace(-2, 2, 20)
X, Y = np.meshgrid(x, y)
Z = np.sin(X) * np.cos(Y)

# Define a vector field
U = -Y
V = X

# Calculate curl (z-component)
curl = np.gradient(V, axis=0) - np.gradient(U, axis=1)

# Plot the vector field and surface
fig = plt.figure(figsize=(14, 7))

# Surface plot
ax1 = fig.add_subplot(121, projection='3d')
ax1.plot_surface(X, Y, Z, cmap='viridis', alpha=0.7)
ax1.set_title('Surface Plot')
ax1.set_xlabel('x')
ax1.set_ylabel('y')
ax1.set_zlabel('z')

# Vector field plot
ax2 = fig.add_subplot(122)
ax2.quiver(X, Y, U, V, curl, cmap='viridis')
ax2.set_title('Curl of the Vector Field')
ax2.set_xlabel('x')
ax2.set_ylabel('y')
plt.colorbar(ax2.quiver(X, Y, U, V, curl), ax=ax2)

plt.show()

import numpy as np
import matplotlib.pyplot as plt

# Define the grid
x = np.linspace(-5, 5, 100)
y = np.linspace(-5, 5, 100)
X, Y = np.meshgrid(x, y)

# Define the potential energy function
U = 0.5 * (X**2 + Y**2)

# Compute the gradient (force)
dU_dx, dU_dy = np.gradient(U, x, y)

# Plot the potential energy
plt.figure(figsize=(12, 6))

plt.subplot(1, 2, 1)
plt.contourf(X, Y, U, cmap='viridis')
plt.colorbar(label='Potential Energy U')
plt.title('Potential Energy U(x, y)')

# Plot the gradient (force)
plt.subplot(1, 2, 2)
plt.quiver(X, Y, -dU_dx, -dU_dy, color='r')
plt.title('Gradient of Potential Energy (Force)')
plt.xlabel('x')
plt.ylabel('y')
plt.grid()

plt.tight_layout()
plt.show()

import numpy as np
import matplotlib.pyplot as plt

# Define the grid
x = np.linspace(-5, 5, 100)
y = np.linspace(-5, 5, 100)
X, Y = np.meshgrid(x, y)

# Define the velocity field
u = -Y
v = X

# Compute the divergence
div_v = np.gradient(u, axis=0) + np.gradient(v, axis=1)

# Plot the velocity field
plt.figure(figsize=(12, 6))

plt.subplot(1, 2, 1)
plt.quiver(X, Y, u, v)
plt.title('Velocity Field')
plt.xlabel('x')
plt.ylabel('y')
plt.grid()

# Plot the divergence
plt.subplot(1, 2, 2)
plt.contourf(X, Y, div_v, cmap='viridis')
plt.colorbar(label='Divergence of Velocity Field')
plt.title('Divergence (div)')
plt.xlabel('x')
plt.ylabel('y')

plt.tight_layout()
plt.show()

import numpy as np
import matplotlib.pyplot as plt
from matplotlib import cm

# 空間の範囲を定義
x = np.linspace(-5, 5, 100)
y = np.linspace(-5, 5, 100)
X, Y = np.meshgrid(x, y)

# 電場と磁場を定義 (例: 簡単な静的場)
E_x = -Y
E_y = X
B_z = np.ones_like(X)  # 簡単化のため磁場を一定

# 電場のプロット
plt.figure(figsize=(12, 6))

plt.subplot(1, 2, 1)
plt.quiver(X, Y, E_x, E_y, color='r')
plt.title('Electric Field')
plt.xlabel('x')
plt.ylabel('y')
plt.grid(True)

# 磁場のプロット (例: 簡単な定常磁場)
plt.subplot(1, 2, 2)
plt.imshow(B_z, extent=(-5, 5, -5, 5), origin='lower', cmap=cm.viridis)
plt.colorbar(label='B_z')
plt.title('Magnetic Field')
plt.xlabel('x')
plt.ylabel('y')
plt.grid(True)

plt.tight_layout()
plt.show()

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import quad

# 空間曲線のベクトル関数
def r(t):
    x = np.sin(t)
    y = np.cos(t)
    z = t
    return np.array([x, y, z])

# 曲線の長さを計算するための微分
def dr_dt(t):
    return np.array([np.cos(t), -np.sin(t), 1])

# 曲線の長さを計算する関数
def curve_length(a, b):
    integral, _ = quad(lambda t: np.linalg.norm(dr_dt(t)), a, b)
    return integral

# ベクトル関数の線積分（例: 定数ベクトル場 F(x, y, z) = (1, 1, 1)）
def F(x, y, z):
    return np.array([1, 1, 1])

# 線積分を計算する関数
def line_integral(a, b):
    integral, _ = quad(lambda t: np.dot(F(*r(t)), np.linalg.norm(dr_dt(t))), a, b)
    return integral

# パラメータ範囲
t_values = np.linspace(0, 2 * np.pi, 100)
r_values = np.array([r(t) for t in t_values])

# 空間曲線のプロット
fig = plt.figure(figsize=(12, 6))
ax = fig.add_subplot(111, projection='3d')
ax.plot(r_values[:, 0], r_values[:, 1], r_values[:, 2], label='Curve r(t)')
ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_zlabel('Z')
ax.set_title('Space Curve r(t)')
ax.legend()

# 曲線の長さと線積分を計算
length = curve_length(0, 2 * np.pi)
integral = line_integral(0, 2 * np.pi)

print(f"Curve Length: {length:.2f}")
print(f"Line Integral: {integral:.2f}")

plt.show()

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

# ベクトル場の定義
def F(x, y, z):
    return np.array([x, y, z])

# 発散を計算
def divergence(x, y, z):
    return 1  # 発散が一定のベクトル場（例: F = (x, y, z)）

# 体積を定義
x = np.linspace(-1, 1, 10)
y = np.linspace(-1, 1, 10)
z = np.linspace(-1, 1, 10)
X, Y, Z = np.meshgrid(x, y, z)

# 発散のプロット
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')

# メッシュデータのプロット
sc = ax.scatter(X, Y, Z, c=np.ones_like(X), cmap='viridis')
plt.colorbar(sc)

ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_zlabel('Z')
ax.set_title('Divergence of the Vector Field')

plt.show()

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

# 二変数関数の定義
def f(x, y):
    return x**2 + y**2

# 勾配（完全微分）の計算
def grad_f(x, y):
    df_dx = 2*x
    df_dy = 2*y
    return np.array([df_dx, df_dy])

# 空間の範囲を定義
x = np.linspace(-2, 2, 20)
y = np.linspace(-2, 2, 20)
X, Y = np.meshgrid(x, y)
Z = f(X, Y)

# 勾配の計算
U, V = grad_f(X, Y)

# 3Dプロット
fig = plt.figure(figsize=(14, 7))

# 関数のプロット
ax = fig.add_subplot(121, projection='3d')
ax.plot_surface(X, Y, Z, cmap='viridis', alpha=0.7)
ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_zlabel('f(x, y)')
ax.set_title('Function f(x, y) = x^2 + y^2')

# 勾配ベクトルのプロット
ax2 = fig.add_subplot(122)
ax2.quiver(X, Y, U, V, color='r', scale=10)
ax2.set_xlabel('X')
ax2.set_ylabel('Y')
ax2.set_title('Gradient of f(x, y)')
ax2.set_aspect('equal')

plt.tight_layout()
plt.show()

# スカラー場（例: 高さ関数 z = x^2 + y^2）
def scalar_field(x, y):
    return x**2 + y**2

# 勾配ベクトルの計算
def gradient_field(x, y):
    df_dx = 2 * x
    df_dy = 2 * y
    return np.array([df_dx, df_dy])

# スカラー場と勾配の計算
Z_scalar = scalar_field(X, Y)
grad = np.array([gradient_field(x, y) for x, y in zip(X.flatten(), Y.flatten())])
Gx, Gy = grad[:, 0].reshape(X.shape), grad[:, 1].reshape(Y.shape)

# スカラー場のプロットと勾配ベクトル
fig = plt.figure(figsize=(12, 6))
ax = fig.add_subplot(121, projection='3d')
ax.plot_surface(X, Y, Z_scalar, cmap='viridis', alpha=0.6)
ax.quiver(X, Y, Z_scalar, Gx, Gy, 0, length=0.1, color='r', normalize=True)
ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_zlabel('Z')
ax.set_title('Scalar Field with Gradient Vectors')

plt.show()

# 曲線のパラメータ表示 (例: 緯度線)
u = np.linspace(0, 2 * np.pi, 100)
v = np.pi / 4  # 固定された緯度

# 曲線のパラメータ表示
X_curve = np.sin(v) * np.cos(u)
Y_curve = np.sin(v) * np.sin(u)
Z_curve = np.cos(v)

# 曲面と曲線のプロット
fig = plt.figure(figsize=(12, 7))
ax = fig.add_subplot(111, projection='3d')
ax.plot_surface(X, Y, Z, alpha=0.6, cmap='viridis')
ax.plot(X_curve, Y_curve, Z_curve, color='r', lw=2, label='Curve on Surface')
ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_zlabel('Z')
ax.set_title('Curve on the Surface')
ax.legend()

plt.show()

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

# パラメータ曲面の定義 (例: 球面)
u = np.linspace(0, 2 * np.pi, 100)
v = np.linspace(0, np.pi, 100)
U, V = np.meshgrid(u, v)
X = np.sin(V) * np.cos(U)
Y = np.sin(V) * np.sin(U)
Z = np.cos(V)

# 曲面のプロット
fig = plt.figure(figsize=(10, 7))
ax = fig.add_subplot(111, projection='3d')
ax.plot_surface(X, Y, Z, cmap='viridis')
ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_zlabel('Z')
ax.set_title('Parametric Surface (Sphere)')
plt.show()

import numpy as np
import matplotlib.pyplot as plt

# 2次元スカラー場の定義
def scalar_field(x, y):
    return np.sin(x) * np.cos(y)

# 2次元ベクトル場の定義
def vector_field(x, y):
    u = -np.sin(x) * np.cos(y)
    v = np.cos(x) * np.sin(y)
    return u, v

# グリッドの設定
x = np.linspace(-2*np.pi, 2*np.pi, 100)
y = np.linspace(-2*np.pi, 2*np.pi, 100)
X, Y = np.meshgrid(x, y)

# スカラー場の計算
Z = scalar_field(X, Y)

# ベクトル場の計算
U, V = vector_field(X, Y)

# スカラー場のプロット
plt.figure(figsize=(12, 6))

# スカラー場の等高線プロット
plt.subplot(1, 2, 1)
contour = plt.contourf(X, Y, Z, cmap='viridis')
plt.colorbar(contour)
plt.title('Scalar Field')
plt.xlabel('x')
plt.ylabel('y')

# ベクトル場のプロット
plt.subplot(1, 2, 2)
plt.quiver(X, Y, U, V, scale=20, color='r')
plt.title('Vector Field')
plt.xlabel('x')
plt.ylabel('y')

plt.tight_layout()
plt.show()

import numpy as np
import matplotlib.pyplot as plt

# 2次元スカラー場の定義
def scalar_field(x, y):
    return np.sin(np.sqrt(x**2 + y**2))

# グリッドの設定
x = np.linspace(-5, 5, 100)
y = np.linspace(-5, 5, 100)
X, Y = np.meshgrid(x, y)

# スカラー場の計算
Z = scalar_field(X, Y)

# ラプラシアンの計算
def laplacian_2d(Z, dx):
    d2Z_dx2 = np.gradient(np.gradient(Z, axis=0), axis=0)
    d2Z_dy2 = np.gradient(np.gradient(Z, axis=1), axis=1)
    return d2Z_dx2 + d2Z_dy2

dx = x[1] - x[0]  # グリッドの間隔
Laplacian_Z = laplacian_2d(Z, dx)

# グラフの描画
plt.figure(figsize=(12, 6))

# スカラー場のプロット
plt.subplot(1, 2, 1)
plt.contourf(X, Y, Z, cmap='viridis')
plt.colorbar()
plt.title('Scalar Field')
plt.xlabel('x')
plt.ylabel('y')

# ラプラシアンのプロット
plt.subplot(1, 2, 2)
plt.contourf(X, Y, Laplacian_Z, cmap='viridis')
plt.colorbar()
plt.title('Laplacian of Scalar Field')
plt.xlabel('x')
plt.ylabel('y')

plt.tight_layout()
plt.show()

from mpl_toolkits.mplot3d import Axes3D

# 3次元スカラー場の定義
def scalar_field_3d(x, y, z):
    return np.sin(np.sqrt(x**2 + y**2 + z**2))

# グリッドの設定
x = np.linspace(-5, 5, 30)
y = np.linspace(-5, 5, 30)
z = np.linspace(-5, 5, 30)
X, Y, Z = np.meshgrid(x, y, z)

# スカラー場の計算
F = scalar_field_3d(X, Y, Z)

# ラプラシアンの計算
def laplacian_3d(F, dx):
    d2F_dx2 = np.gradient(np.gradient(F, axis=0), axis=0)
    d2F_dy2 = np.gradient(np.gradient(F, axis=1), axis=1)
    d2F_dz2 = np.gradient(np.gradient(F, axis=2), axis=2)
    return d2F_dx2 + d2F_dy2 + d2F_dz2

dx = x[1] - x[0]  # グリッドの間隔
Laplacian_F = laplacian_3d(F, dx)

# グラフの描画
fig = plt.figure(figsize=(12, 6))

# スカラー場のスライスプロット
ax1 = fig.add_subplot(1, 2, 1, projection='3d')
X2, Y2 = np.meshgrid(x, y)
ax1.plot_surface(X2, Y2, F[:, :, F.shape[2] // 2], cmap='viridis')
ax1.set_title('Scalar Field (z = 0)')

# ラプラシアンのスライスプロット
ax2 = fig.add_subplot(1, 2, 2, projection='3d')
ax2.plot_surface(X2, Y2, Laplacian_F[:, :, F.shape[2] // 2], cmap='viridis')
ax2.set_title('Laplacian of Scalar Field (z = 0)')

plt.tight_layout()
plt.show()

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

# Define vectors A and B
A = np.array([2, 3, 4])
B = np.array([5, 6, 7])

# Dot product
dot_product = np.dot(A, B)
print(f"Dot Product of A and B: {dot_product}")

# Cross product
cross_product = np.cross(A, B)
print(f"Cross Product of A and B: {cross_product}")

# Visualization
fig = plt.figure(figsize=(10, 8))

# 3D plot
ax = fig.add_subplot(111, projection='3d')

# Plot vector A
ax.quiver(0, 0, 0, A[0], A[1], A[2], color='r', label='A', arrow_length_ratio=0.1)

# Plot vector B
ax.quiver(0, 0, 0, B[0], B[1], B[2], color='b', label='B', arrow_length_ratio=0.1)

# Plot cross product (vector C)
ax.quiver(0, 0, 0, cross_product[0], cross_product[1], cross_product[2], color='g', label='A × B', arrow_length_ratio=0.1)

# Setting the axes properties
ax.set_xlim([0, 7])
ax.set_ylim([0, 7])
ax.set_zlim([0, 7])
ax.set_xlabel('X axis')
ax.set_ylabel('Y axis')
ax.set_zlabel('Z axis')

# Display the legend
ax.legend()

# Show plot
plt.show()

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

# Define a vector field F(x, y, z) = [x, y, z]
def vector_field(x, y, z):
    return np.array([x, y, z])

# Calculate the divergence of the vector field
def divergence(F, x, y, z):
    dFx_dx = np.gradient(F[0], x, axis=0)
    dFy_dy = np.gradient(F[1], y, axis=1)
    dFz_dz = np.gradient(F[2], z, axis=2)
    return dFx_dx + dFy_dy + dFz_dz

# Create a 3D grid
x = y = z = np.linspace(-2, 2, 10)
X, Y, Z = np.meshgrid(x, y, z)

# Compute the vector field and its divergence
F = vector_field(X, Y, Z)
div_F = divergence(F, X, Y, Z)

# Plot the vector field and its divergence
fig = plt.figure(figsize=(12, 8))
ax = fig.add_subplot(111, projection='3d')

# Quiver plot for vector field
ax.quiver(X, Y, Z, F[0], F[1], F[2], length=0.2, color='b', label='Vector Field')

# Scatter plot for divergence
ax.scatter(X, Y, Z, c=div_F, cmap='coolwarm', label='Divergence')

ax.set_xlabel('X axis')
ax.set_ylabel('Y axis')
ax.set_zlabel('Z axis')
ax.set_title("Gauss's Theorem - Vector Field and Divergence")
ax.legend()
plt.show()

# Define a 2D vector field F(x, y) = [-y, x]
def vector_field_2d(x, y):
    return np.array([-y, x])

# Calculate the curl of the vector field (in 2D, this reduces to a scalar)
def curl_2d(F, x, y):
    dFx_dy = np.gradient(F[0], y, axis=0)
    dFy_dx = np.gradient(F[1], x, axis=1)
    return dFy_dx - dFx_dy

# Create a 2D grid
x = y = np.linspace(-2, 2, 20)
X, Y = np.meshgrid(x, y)

# Compute the vector field and its curl
F = vector_field_2d(X, Y)
curl_F = curl_2d(F, X, Y)

# Plot the vector field and its curl
fig, ax = plt.subplots(figsize=(8, 6))

# Quiver plot for vector field
ax.quiver(X, Y, F[0], F[1], color='b', label='Vector Field')

# Contour plot for curl
contour = ax.contourf(X, Y, curl_F, cmap='coolwarm', alpha=0.7)
plt.colorbar(contour)

ax.set_xlabel('X axis')
ax.set_ylabel('Y axis')
ax.set_title("Stokes' Theorem - Vector Field and Curl")
ax.legend()
plt.show()

# Define a simple vector field F(x, y) = [x, y]
def vector_field_2d_green(x, y):
    return np.array([x, y])

# Define the region and curve
theta = np.linspace(0, 2*np.pi, 100)
x_circle = np.cos(theta)
y_circle = np.sin(theta)

# Compute vector field on the curve
F_circle = vector_field_2d_green(x_circle, y_circle)

# Line integral (approximation by summing vectors along the circle)
line_integral = np.sum(F_circle[0] * np.gradient(x_circle) + F_circle[1] * np.gradient(y_circle))

# Compute the curl of the vector field over the region
X, Y = np.meshgrid(np.linspace(-1, 1, 50), np.linspace(-1, 1, 50))
F = vector_field_2d_green(X, Y)
curl_F = curl_2d(F, X, Y)

# Area integral (approximation by summing the curl over the region)
area_integral = np.sum(curl_F) * (4 / 50**2)

# Plot the vector field and the curve
fig, ax = plt.subplots(figsize=(8, 6))

# Vector field plot
ax.quiver(X, Y, F[0], F[1], color='b', alpha=0.5, label='Vector Field')

# Curve plot
ax.plot(x_circle, y_circle, 'r', label='Curve (C)')

# Curl plot
contour = ax.contourf(X, Y, curl_F, cmap='coolwarm', alpha=0.7)
plt.colorbar(contour)

ax.set_xlabel('X axis')
ax.set_ylabel('Y axis')
ax.set_title("Green's Theorem - Vector Field and Curl")
ax.legend()
plt.show()

print(f"Line Integral around C: {line_integral}")
print(f"Area Integral of Curl: {area_integral}")

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.colors import Normalize

# グリッドの設定
x = np.linspace(-5, 5, 100)
y = np.linspace(-5, 5, 100)
X, Y = np.meshgrid(x, y)

# 電場の定義（静電場）
def electric_field(X, Y):
    k = 1  # 定数（例: クーロン定数）
    r = np.sqrt(X**2 + Y**2)
    E_x = k * X / (r**3 + 1e-9)  # x成分
    E_y = k * Y / (r**3 + 1e-9)  # y成分
    return E_x, E_y

# 磁場の定義（定常磁場）
def magnetic_field(X, Y):
    B_x = np.zeros_like(X)  # x成分はゼロ
    B_y = np.zeros_like(Y)  # y成分はゼロ
    return B_x, B_y

# 電場と磁場の計算
E_x, E_y = electric_field(X, Y)
B_x, B_y = magnetic_field(X, Y)

# プロット設定
fig, axs = plt.subplots(1, 2, figsize=(12, 6))

# 電場のプロット
strm = axs[0].streamplot(X, Y, E_x, E_y, color='blue', linewidth=1.5, density=2)
axs[0].set_title('Electric Field')
axs[0].set_xlabel('x')
axs[0].set_ylabel('y')
axs[0].set_aspect('equal')
fig.colorbar(strm.lines, ax=axs[0])

# 磁場のプロット
strm = axs[1].streamplot(X, Y, B_x, B_y, color='red', linewidth=1.5, density=2)
axs[1].set_title('Magnetic Field')
axs[1].set_xlabel('x')
axs[1].set_ylabel('y')
axs[1].set_aspect('equal')
fig.colorbar(strm.lines, ax=axs[1])

plt.tight_layout()
plt.show()

import numpy as np
import matplotlib.pyplot as plt
from scipy import ndimage

# グリッドの設定
n = 100
x = np.linspace(-5, 5, n)
y = np.linspace(-5, 5, n)
X, Y = np.meshgrid(x, y)

# ポアソン方程式の右辺（体積電荷密度の例）
rho = np.zeros_like(X)
rho[np.sqrt(X**2 + Y**2) < 1] = 1  # 点電荷を模した体積電荷密度

# ポテンシャルの初期化
phi = np.zeros_like(X)

# 数値解法（ガウス・ザイデル法）
def solve_poisson(phi, rho, iterations=1000):
    for _ in range(iterations):
        phi_new = np.copy(phi)
        phi_new[1:-1, 1:-1] = 0.25 * (phi[1:-1, :-2] + phi[1:-1, 2:] + phi[:-2, 1:-1] + phi[2:, 1:-1] - rho[1:-1, 1:-1])
        phi = np.copy(phi_new)
    return phi

# ポテンシャルの計算
phi = solve_poisson(phi, rho)

# プロット
plt.figure(figsize=(10, 6))
plt.contourf(X, Y, phi, levels=50, cmap='inferno')
plt.colorbar(label='Potential $\phi$')
plt.title('Potential Field from Point Charge (Poisson\'s Equation)')
plt.xlabel('x')
plt.ylabel('y')
plt.show()

import numpy as np
import matplotlib.pyplot as plt

# グリッドの設定
n = 100
x = np.linspace(0, 1, n)
y = np.linspace(0, 1, n)
X, Y = np.meshgrid(x, y)

# ポテンシャルの初期化
phi = np.zeros_like(X)

# 境界条件
phi[:, 0] = 0   # 左端
phi[:, -1] = 0  # 右端
phi[0, :] = 0   # 下端
phi[-1, :] = 1  # 上端

# 数値解法（ガウス・ザイデル法）
def solve_laplace(phi, iterations=1000):
    for _ in range(iterations):
        phi_new = np.copy(phi)
        phi_new[1:-1, 1:-1] = 0.25 * (phi[1:-1, :-2] + phi[1:-1, 2:] + phi[:-2, 1:-1] + phi[2:, 1:-1])
        phi = np.copy(phi_new)
    return phi

# ポテンシャルの計算
phi = solve_laplace(phi)

# プロット
plt.figure(figsize=(10, 6))
plt.contourf(X, Y, phi, levels=50, cmap='viridis')
plt.colorbar(label='Potential $\phi$')
plt.title('Potential Field (Laplace\'s Equation)')
plt.xlabel('x')
plt.ylabel('y')
plt.show()

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.colors import Normalize

# グリッドサイズとパラメータ設定
n = 100
x = np.linspace(0, 1, n)
y = np.linspace(0, 1, n)
X, Y = np.meshgrid(x, y)

# 空気の速度場の初期条件（例: ダムブレイク問題のシンプルなモデル）
U = np.sin(np.pi * X) * np.cos(np.pi * Y)
V = -np.cos(np.pi * X) * np.sin(np.pi * Y)

# 圧力場の初期条件（速度場に基づくシンプルなモデル）
P = np.sin(np.pi * X) * np.sin(np.pi * Y)

# プロット設定
fig, axs = plt.subplots(1, 2, figsize=(14, 6))

# 速度場のプロット
strm = axs[0].streamplot(X, Y, U, V, color='blue', linewidth=1.5, density=2)
axs[0].set_title('Velocity Field')
axs[0].set_xlabel('x')
axs[0].set_ylabel('y')
axs[0].set_aspect('equal')
fig.colorbar(strm.lines, ax=axs[0])

# 圧力場のプロット
contour = axs[1].contourf(X, Y, P, levels=50, cmap='viridis')
fig.colorbar(contour, ax=axs[1], label='Pressure')
axs[1].set_title('Pressure Field')
axs[1].set_xlabel('x')
axs[1].set_ylabel('y')
axs[1].set_aspect('equal')

plt.tight_layout()
plt.show()

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.colors import Normalize

# グリッドサイズとパラメータ設定
n = 100
x = np.linspace(0, 1, n)
y = np.linspace(0, 1, n)
X, Y = np.meshgrid(x, y)

# 初期条件設定
U = np.zeros((n, n))
V = np.zeros((n, n))
P = np.zeros((n, n))

# 初期の速度場（例: 初期速度がある領域）
U[40:60, 40:60] = 1.0
V[40:60, 40:60] = 0.5

# 圧力場の初期条件
# ここではシンプルに圧力がx座標に依存する場合を考えます
P = 1.0 - X

# プロット設定
fig, axs = plt.subplots(1, 2, figsize=(14, 6))

# 速度場のプロット
strm = axs[0].streamplot(X, Y, U, V, color='blue', linewidth=1.5, density=2)
axs[0].set_title('Velocity Field')
axs[0].set_xlabel('x')
axs[0].set_ylabel('y')
axs[0].set_aspect('equal')
fig.colorbar(strm.lines, ax=axs[0])

# 圧力場のプロット
contour = axs[1].contourf(X, Y, P, levels=50, cmap='viridis')
fig.colorbar(contour, ax=axs[1], label='Pressure')
axs[1].set_title('Pressure Field')
axs[1].set_xlabel('x')
axs[1].set_ylabel('y')
axs[1].set_aspect('equal')

plt.tight_layout()
plt.show()

import numpy as np
import matplotlib.pyplot as plt

# グリッドサイズとパラメータ設定
n = 100
x = np.linspace(-2, 2, n)
y = np.linspace(-2, 2, n)
X, Y = np.meshgrid(x, y)

# スカラー場の定義
phi = np.sin(np.pi * X) * np.cos(np.pi * Y)

# 勾配の計算
dphi_dx, dphi_dy = np.gradient(phi, x, y)
grad_phi_x, grad_phi_y = dphi_dx, dphi_dy

# ベクトル場の定義（例: 流体場の簡単なモデル）
U = -Y
V = X

# 発散の計算
div_V = np.gradient(U, x, axis=0) + np.gradient(V, y, axis=1)

# 回転の計算
curl_V_x = np.gradient(V, x, axis=0) - np.gradient(U, y, axis=1)
curl_V_y = np.gradient(U, x, axis=1) - np.gradient(V, y, axis=0)

# プロット設定
fig, axs = plt.subplots(2, 3, figsize=(18, 12))

# 勾配のプロット
axs[0, 0].quiver(X, Y, grad_phi_x, grad_phi_y, color='blue')
axs[0, 0].set_title('Gradient of φ')
axs[0, 0].set_xlabel('x')
axs[0, 0].set_ylabel('y')
axs[0, 0].set_aspect('equal')

# スカラー場のプロット
contour_phi = axs[0, 1].contourf(X, Y, phi, levels=50, cmap='viridis')
fig.colorbar(contour_phi, ax=axs[0, 1])
axs[0, 1].set_title('Scalar Field φ')
axs[0, 1].set_xlabel('x')
axs[0, 1].set_ylabel('y')
axs[0, 1].set_aspect('equal')

# 発散のプロット
contour_div_V = axs[0, 2].contourf(X, Y, div_V, levels=50, cmap='plasma')
fig.colorbar(contour_div_V, ax=axs[0, 2])
axs[0, 2].set_title('Divergence of V')
axs[0, 2].set_xlabel('x')
axs[0, 2].set_ylabel('y')
axs[0, 2].set_aspect('equal')

# 回転のプロット
axs[1, 0].quiver(X, Y, curl_V_x, curl_V_y, color='green')
axs[1, 0].set_title('Curl of V')
axs[1, 0].set_xlabel('x')
axs[1, 0].set_ylabel('y')
axs[1, 0].set_aspect('equal')

# ベクトル場のプロット
strm = axs[1, 1].streamplot(X, Y, U, V, color='red', linewidth=1.5, density=2)
axs[1, 1].set_title('Vector Field V')
axs[1, 1].set_xlabel('x')
axs[1, 1].set_ylabel('y')
axs[1, 1].set_aspect('equal')
fig.colorbar(strm.lines, ax=axs[1, 1])

# カラーバーのプロット
contour_curl_V = axs[1, 2].contourf(X, Y, np.sqrt(curl_V_x**2 + curl_V_y**2), levels=50, cmap='cividis')
fig.colorbar(contour_curl_V, ax=axs[1, 2])
axs[1, 2].set_title('Magnitude of Curl of V')
axs[1, 2].set_xlabel('x')
axs[1, 2].set_ylabel('y')
axs[1, 2].set_aspect('equal')

plt.tight_layout()
plt.show()

import numpy as np
import matplotlib.pyplot as plt
from scipy.sparse import diags
from scipy.sparse.linalg import spsolve

# グリッド設定
n = 20
h = 1.0 / (n - 1)
x = np.linspace(0, 1, n)
y = np.linspace(0, 1, n)
X, Y = np.meshgrid(x, y)

# メッシュ作成
N = (n - 2) * (n - 2)
A = np.zeros((N, N))
b = np.zeros(N)

def index(i, j):
    return (i - 1) * (n - 2) + (j - 1)

# スタディアム
for i in range(1, n - 1):
    for j in range(1, n - 1):
        k = index(i, j)
        A[k, k] = -4
        b[k] = -1  # ソース項 f(x, y) = -1
        if i > 1:
            A[k, index(i-1, j)] = 1
        if i < n - 2:
            A[k, index(i+1, j)] = 1
        if j > 1:
            A[k, index(i, j-1)] = 1
        if j < n - 2:
            A[k, index(i, j+1)] = 1

# 線形システムの解
u = spsolve(A, b)

# 結果のプロット
U = np.zeros((n, n))
U[1:-1, 1:-1] = u.reshape((n - 2, n - 2))
plt.contourf(X, Y, U, cmap='viridis')
plt.colorbar()
plt.title('2D Poisson Equation Solution')
plt.xlabel('x')
plt.ylabel('y')
plt.show()