你的程序呢亲？

很抱歉，我是一名语言模型AI，无法在本地运行程序。但是我可以给你提供一个大致的程序框架，供你参考：

import numpy as np
import matplotlib.pyplot as plt

# 设置计算参数
L = 1.0 # 平行板电容器的长度
W = 1.0 # 平行板电容器的宽度
h = 0.01 # 网格步长
N = int(L/h) # 网格数
d = 0.01 # 计算步长
epsilon_r = 1 # 相对介电常数
epsilon_0 = 8.854e-12 # 真空介电常数

# 初始化电势和电荷密度
U = np.zeros((N,N))
rho = np.zeros((N,N))

# 设置边界条件，上下两个平行板的电势为100V和-100V
U[0,:] = 100
U[-1,:] = -100

# 迭代求解拉普拉斯方程
for k in range(1000):
    U_old = U.copy() # 保存上一步的电势值
    for i in range(1,N-1):
        for j in range(1,N-1):
            U[i,j] = (U_old[i-1,j] + U_old[i+1,j] + U_old[i,j-1] + U_old[i,j+1]) / 4
    # 判断是否收敛
    diff = np.abs(U - U_old)
    if np.max(diff) < 1e-6:
        break

# 计算电荷密度
for i in range(1,N-1):
    for j in range(1,N-1):
        rho[i,j] = -epsilon_0/epsilon_r * (U[i-1,j] + U[i+1,j] + U[i,j-1] + U[i,j+1] - 4*U[i,j]) / d**2

# 绘制电荷密度分布图
plt.imshow(rho, cmap='jet', extent=[0, L, 0, W])
plt.colorbar()
plt.xlabel('x')
plt.ylabel('y')
plt.show()

需要注意的是，由于本程序是一种数值计算方法，因此求解结果可能会受到计算精度、迭代次数等因素的影响，因此需要根据具体情况进行调整，以达到较为准确的计算结果。