admin 管理员组

文章数量: 1184232

 
# Optimized Calculation and Visualization of Electric Potential
# and Field with Numba and SymPy
import numpy as np
import matplotlib.pyplot as plt
from numba import jit
from sympy import symbols, sqrt, lambdify

# Constants
k = 8.9875517873681764e9  # Coulomb's constant in N·m²/C²

@jit(nopython=True)
def compute_potential(x_q, y_q, z_q, q, x_p, y_p, z_p):
    r = np.sqrt((x_p - x_q)**2 + (y_p - y_q)**2 + (z_p - z_q)**2)
    return k * q / r  # 注意：若观测点与电荷重合会除零

def plot_electric_field(point_charge_pos, point_charge_value):
    """
    Render the electric field vector field with optimized calculations.
    :param point_charge_pos: (x0, y0, z0)
    :param point_charge_value: q
    """
    # ---- Symbolic scalar formulation (avoid SymPy vector module) ----
    x, y, z = symbols('x y z', real=True)
    x0, y0, z0 = map(float, point_charge_pos)

    r = sqrt((x - x0)**2 + (y - y0)**2 + (z - z0)**2)
    V_sym = k * point_charge_value / r

    # Electric field components (2D slice at z=0)
    Ex_sym = (-V_sym.diff(x)).subs(z, 0)
    Ey_sym = (-V_sym.diff(y)).subs(z, 0)

    # Vectorized numeric functions
    Ex_fn = lambdify((x, y), Ex_sym, 'numpy')
    Ey_fn = lambdify((x, y), Ey_sym, 'numpy')

    # ---- Grid & evaluate ----
    lin = np.linspace(-10.0, 10.0, 25)
    X, Y = np.meshgrid(lin, lin)

    # 避免在电荷所在点除零：对该点做个很小位移或置 0
    eps = 1e-12
    X_safe = X.copy()
    Y_safe = Y.copy()
    mask = (np.isclose(X, x0) & np.isclose(Y, y0))
    X_safe[mask] += eps

    U = Ex_fn(X_safe, Y_safe).astype(float)
    W = Ey_fn(X_safe, Y_safe).astype(float)

    # 可选：对箭头长度做裁剪，避免中心过大
    mag = np.hypot(U, W)
    cap = np.percentile(mag[~np.isnan(mag)], 95)
    cap = cap if cap > 0 else 1.0
    U = np.clip(U, -cap, cap)
    W = np.clip(W, -cap, cap)

    # ---- Plot ----
    plt.figure(figsize=(10, 8))
    plt.quiver(X, Y, U, W, scale_units='xy', angles='xy')
    plt.scatter([x0], [y0], color='red', s=100, label='Point Charge')
    plt.title('Optimized Electric Field Visualization (z=0 plane)')
    plt.xlabel('X axis')
    plt.ylabel('Y axis')
    plt.legend()
    plt.grid(True, alpha=0.3)
    plt.axis('equal')
    plt.tight_layout()
    plt.show()

# Example
point_charge_position = (0.0, 0.0, 0.0)
point_charge_magnitude = 1e-9
observation_position = (3.0, 4.0, 5.0)

potential_at_obs = compute_potential(*point_charge_position,
                                     point_charge_magnitude, *observation_position)
print(f"Electric Potential at {observation_position}: {potential_at_obs:.5f} V")

plot_electric_field(point_charge_position, point_charge_magnitude)

The computed result should be 

Electric Potential at (3.0, 4.0, 5.0): 1.27103 V

Here is the generated plot:

本文标签： Numpy python Compute potential Electric

版权声明：本文标题：Using Python and Numpy to Compute Electric potential 内容由网友自发贡献，该文观点仅代表作者本人， 转载请联系作者并注明出处：http://www.roclinux.cn/b/1765518981a3388738.html， 本站仅提供信息存储空间服务，不拥有所有权，不承担相关法律责任。如发现本站有涉嫌抄袭侵权/违法违规的内容，一经查实，本站将立刻删除。