fork download
copy 
  1. # your code goes here
Success #stdin #stdout 0.08s 14152KB
comments (?)
stdin
copy 
import numpy as np
import matplotlib.pyplot as plt

# ===== 参数 =====
a = 23e-3          # m
b = 10e-3          # m
f = 36e9           # Hz
c = 3e8            # m/s
m, n = 2, 1        # TM21

# ===== 计算常数 =====
fc = c/2 * np.sqrt((m/a)**2 + (n/b)**2)
lam = c/f
lam_g = lam / np.sqrt(1 - (fc/f)**2)
beta = 2*np.pi / lam_g
kc = np.pi * np.sqrt((m/a)**2 + (n/b)**2)
omega = 2*np.pi*f
eps0 = 8.854e-12
E0 = 1  # 归一化

# ===== 场分量公式（TM21，取实部，t=0） =====
def get_fields(x, y, z):
    # x,y,z 单位均为 m
    sin_m = np.sin(m*np.pi*x/a)
    sin_n = np.sin(n*np.pi*y/b)
    cos_m = np.cos(m*np.pi*x/a)
    cos_n = np.cos(n*np.pi*y/b)
    phase = np.cos(beta*z)  # 取实部

    Ez = E0 * sin_m * sin_n * phase
    Ex = -1j * (beta * m*np.pi/a) / kc**2 * E0 * cos_m * sin_n
    Ex = np.real(Ex * np.exp(-1j*beta*z))
    Ey = -1j * (beta * n*np.pi/b) / kc**2 * E0 * sin_m * cos_n
    Ey = np.real(Ey * np.exp(-1j*beta*z))
    Hx = 1j * (omega*eps0 * n*np.pi/b) / kc**2 * E0 * sin_m * cos_n
    Hx = np.real(Hx * np.exp(-1j*beta*z))
    Hy = -1j * (omega*eps0 * m*np.pi/a) / kc**2 * E0 * cos_m * sin_n
    Hy = np.real(Hy * np.exp(-1j*beta*z))
    Hz = np.zeros_like(x)  # TM模 Hz=0
    return Ex, Ey, Ez, Hx, Hy, Hz

# ===== 网格 =====
N = 100
x = np.linspace(0, a, N)
y = np.linspace(0, b, N)
X, Y = np.meshgrid(x, y)

# ===== 任务2：xy横截面剖面图（z=0, λ/8, λ/4, λ/2） =====
z_list = [0, lam/8, lam/4, lam/2]
titles2 = ['z=0', 'z=λ/8', 'z=λ/4', 'z=λ/2']
comp_names = ['Ex', 'Ey', 'Ez', 'Hx', 'Hy', 'Hz']

for idx, z in enumerate(z_list):
    Ex, Ey, Ez, Hx, Hy, Hz = get_fields(X, Y, z)
    fig, axes = plt.subplots(2, 3, figsize=(12, 8))
    fig.suptitle(f'任务2：xy横截面剖面图  {titles2[idx]}', fontsize=14)
    comps = [Ex, Ey, Ez, Hx, Hy, Hz]
    for i, ax in enumerate(axes.flat):
        im = ax.contourf(X*1000, Y*1000, comps[i], levels=20, cmap='jet')
        ax.set_title(comp_names[i])
        ax.set_xlabel('x (mm)')
        ax.set_ylabel('y (mm)')
        plt.colorbar(im, ax=ax)
    plt.tight_layout()
    plt.show()

# ===== 任务3：大尺寸场分布图（|E| 或 Ez 分布） =====
# 这里画 |E| = sqrt(Ex^2+Ey^2+Ez^2)，更直观展示总场
# 1. xy平面：z=0, λ/8, λ/4, λ/2
z_list = [0, lam/8, lam/4, lam/2]
fig, axes = plt.subplots(2, 2, figsize=(12, 10))
fig.suptitle('任务3：xy平面大尺寸场分布 (|E|)', fontsize=16)
for i, ax in enumerate(axes.flat):
    z = z_list[i]
    Ex, Ey, Ez, Hx, Hy, Hz = get_fields(X, Y, z)
    Emag = np.sqrt(Ex**2 + Ey**2 + Ez**2)
    im = ax.contourf(X*1000, Y*1000, Emag, levels=20, cmap='inferno')
    ax.set_title(f'z = {z_list[i]/lam:.2f}λ')
    ax.set_xlabel('x (mm)')
    ax.set_ylabel('y (mm)')
    plt.colorbar(im, ax=ax)
plt.tight_layout()
plt.show()

# 2. xz平面（y=0, b/8, b/4, b/2）
y_list = [0, b/8, b/4, b/2]
x_z = np.linspace(0, a, N)
z_xz = np.linspace(0, lam/2, N)  # 取半个波长长度
X_xz, Z_xz = np.meshgrid(x_z, z_xz)

fig, axes = plt.subplots(2, 2, figsize=(12, 10))
fig.suptitle('任务3：xz平面大尺寸场分布 (|E|)', fontsize=16)
for i, ax in enumerate(axes.flat):
    yy = y_list[i]
    # 重新计算（yy固定）
    sin_m = np.sin(m*np.pi*X_xz/a)
    sin_n = np.sin(n*np.pi*yy/b)
    cos_m = np.cos(m*np.pi*X_xz/a)
    cos_n = np.cos(n*np.pi*yy/b)
    phase = np.cos(beta*Z_xz)
    Ez = E0 * sin_m * sin_n * phase
    Ex = np.real(-1j * (beta * m*np.pi/a) / kc**2 * E0 * cos_m * sin_n * np.exp(-1j*beta*Z_xz))
    Ey = np.real(-1j * (beta * n*np.pi/b) / kc**2 * E0 * sin_m * cos_n * np.exp(-1j*beta*Z_xz))
    Emag = np.sqrt(Ex**2 + Ey**2 + Ez**2)
    im = ax.contourf(X_xz*1000, Z_xz*1000, Emag, levels=20, cmap='inferno')
    ax.set_title(f'y = {yy/b:.2f}b')
    ax.set_xlabel('x (mm)')
    ax.set_ylabel('z (mm)')
    plt.colorbar(im, ax=ax)
plt.tight_layout()
plt.show()

# 3. zy平面（x=0, a/8, a/4, a/2）
x_list = [0, a/8, a/4, a/2]
y_zy = np.linspace(0, b, N)
z_zy = np.linspace(0, lam/2, N)
Y_zy, Z_zy = np.meshgrid(y_zy, z_zy)

fig, axes = plt.subplots(2, 2, figsize=(12, 10))
fig.suptitle('任务3：zy平面大尺寸场分布 (|E|)', fontsize=16)
for i, ax in enumerate(axes.flat):
    xx = x_list[i]
    sin_m = np.sin(m*np.pi*xx/a)
    sin_n = np.sin(n*np.pi*Y_zy/b)
    cos_m = np.cos(m*np.pi*xx/a)
    cos_n = np.cos(n*np.pi*Y_zy/b)
    phase = np.cos(beta*Z_zy)
    Ez = E0 * sin_m * sin_n * phase
    Ex = np.real(-1j * (beta * m*np.pi/a) / kc**2 * E0 * cos_m * sin_n * np.exp(-1j*beta*Z_zy))
    Ey = np.real(-1j * (beta * n*np.pi/b) / kc**2 * E0 * sin_m * cos_n * np.exp(-1j*beta*Z_zy))
    Emag = np.sqrt(Ex**2 + Ey**2 + Ez**2)
    im = ax.contourf(Y_zy*1000, Z_zy*1000, Emag, levels=20, cmap='inferno')
    ax.set_title(f'x = {xx/a:.2f}a')
    ax.set_xlabel('y (mm)')
    ax.set_ylabel('z (mm)')
    plt.colorbar(im, ax=ax)
plt.tight_layout()
plt.show()

print(f"截止频率 fc = {fc/1e9:.2f} GHz")
print(f"自由空间波长 λ = {lam*1000:.2f} mm")
print(f"波导波长 λ_g = {lam_g*1000:.2f} mm")
stdout
copy 
Standard output is empty
https://ideone.com/JGpuoA
language:
Python 3 (python  3.12)
created:
3 days ago
visibility:
public

Share or Embed source code

Discover > Sphere Engine API

The brand new service which powers Ideone!

Discover > IDE Widget

Widget for compiling and running the source code in a web browser!