利用MATLAB模拟点电荷的电场分布
一、实验目的：
1．熟悉点电荷的电场分布情况；
2．学会使用Matlab绘图；
二、实验内容：
1．根据库伦定律，利用Matlab强大的绘图功能画出单个点电荷的电场分布情况，包括电力线和等势面。
2．根据库伦定律,利用Matlab强大的绘图功能画出一对点电荷的电场分布情况,包括电力线的分布和等势面。
3．实验内容1中，可以在正电荷和负电荷中任选一组画出其电场分布，实验内容2中，可以在一对正电荷，一对负电荷和一正一负一对电荷中选择一组画出其电场分布情况。
三、实验步骤：
一．对于单个点荷的电力线和等势线：
真空中点电荷的场强大小是：
（式1）
其中k=
为静电力恒量，q为点电荷的电量，r为点电荷到场点
(x,y)的距离。电场呈球对称分布，本实验中，取点电荷为正电荷，电力线是以电荷为起点的射线簇。以无穷远处为零势点，点电荷的电势为：
（式2）
当U取常数时，此式就是等势面方程。等势面是以电荷中心，以r为半径的球面。
（1） 平面电力线的画法：
在平面上，电力线是等角平分布的射线簇，取射线的半径为
=0.12。其程序如下：
r0=0.12; % 射线的半径
th=linspace(0,2*pi,13); % 电力线的角度
[x,y]=pol2cart(th,r0); % 将极坐标转化为直角坐标
x=[x;0.1*x]; % 插入x的起始坐标
y=[y;0.1*y]; % 插入y的起始坐标
plot(x,y,‘b’) % 用蓝色画出所有电力线
grid on % 加网格
Hold on % 保持图像
plot(0,0,‘o’,‘MarkerSize’,12) % 画电荷
xlabel(‘x’,‘fontsize’,16) % 用16号字体标出X轴
ylabel(‘y’,‘fontsize’,16) % 用16号字体标出Y轴
title(‘正电荷的电力线’,‘fontsize’,20) % 添加标题
图1 正电荷的电力线
（2）平面等势面的画法
在过电荷的截面上，等势线就是以电荷为中心的圆簇。此实验中，由于r0=0.12，k=9x109，考虑到电势的大小，取q=1x10-9C，且最大的等势线的半径应该比射线的半径小一点，取r0=0.1，其电势为
。等势线共取7条，且最大的电势为最小电势的3倍。在电场线的基础上画出点电荷的等势线图，可以省略一些基本参数的设置，其图如图2所示，其程序如下：
k=9e9; % 设定k值
q=1e-9; % 设定电荷电量
r0=0.1; % 设定最大等势线的半径
u0=k*q/r0; % 算出最小的电势
u=linspace(1,3,7)*u0; % 求出各条等势线的电势大小
x=linspace(-r0,r0,100); % 将X坐标分成100等份
[X,Y]=meshgrid(x); % 在直角坐标中形成网格坐标
r=sqrt(X.2+Y.2); % 各个网格点到电荷点的距离
U=k*q./r; % 各点的电势
contour(X,Y,U,u) % 画出点电荷的电势面
title(‘正电荷的电场线和等势线’,‘fontsize’,20) %显示标题
图2 正电荷的电场线和等势线
(3) 点电荷的立体电力线
点电荷的立体等势线呈球形发射状的射线簇，因此要先形成三维单位球面坐标，参数还是用前面画平面图的参数。因此其程序如下：
r0=0.12 % 重新设定电力线的半径
[X,Y,Z]=sphere(8); % 形成三维单位球面坐标，绕Z轴一周有8条电力线
x=r0*X(😃'; % 将X化成行向量
y=r0*Y(😃'; % 将Y化成行向量
z=r0*Z(😃'; % 将Z化成行向量
x=[x;zeros(size(x))]; % 对x坐标插入原点
y=[y;zeros(size(y))]; % 对y坐标插入原点
z=[z;zeros(size(z))]; % 对z坐标插入原点
plot3(x,y,z，‘b’) % 画出所有电力线
Hold on % 保持图像
xlabel(‘x’,‘fontsize’,16) % 用16号字体标出X轴
ylabel(‘y’,‘fontsize’,16) % 用16号字体标出Y轴
zlabel(‘z’,‘fontsize’,16) % 用16号字体标出Z轴
title(‘正电荷电场线的三维图形’,‘fontsize’,20) % 添加标题
其图形如下：
图3 正电荷电场线的三维图形
(4) 点电荷的等势面
画等势面时同样要先形成球面，不同的等势面对应不同的半径，而坐标所形成的一个一维的行向量，而三维单位球面的每一维都是21*21的网格矩阵，矩阵的维度不一样，不能直接相乘。因此为减少计算量，只画5条等势面。其程序如下：
u=linspace(1,3,5)*u0; % 计算各面的电势
r=k*q./u; % 计算各等势面的半径
[X,Y,Z]=sphere; % 形成三维的单位球
Z(X<0&Y<0)=nan; % 把球面的四分之一设为非数，便于观察
surf(r(1)*X,r(1)*Y,r(1)*Z); % 画最外面的等势面
hold on; % 保持图形
surf(r(2)*X,r(2)*Y,r(2)*Z);
hold on;
surf(r(3)*X,r(3)*Y,r(3)*Z);
hold on;
surf(r(4)*X,r(4)*Y,r(4)*Z);
hold on;
surf(r(5)*X,r(5)*Y,r(5)*Z);
shading interp % 将各球面的颜色设置成浓淡变化的
xlabel(‘x’,‘fontsize’,16) % 标记X坐标轴
ylabel(‘y’,‘fontsize’,16) % 标记X坐标轴
zlabel(‘z’,‘fontsize’,16) % 标记X坐标轴
title(‘正电荷等势面的三维图形’,‘fontsize’,20) % 添加标题
图四 正电荷等势面的三维图形
二．对于一对点电荷的电力线与等势线
到于两个点电荷的电场分布，比一个点电荷的电场分布要复杂得多，电场线的切线为该点电场强度E的方向。因此画电场线需要先计算出当前点的电场强度E方向，而E又是一个矢量，没有像电势U那样可以直接进行标量计算。因此对于多个点电荷的电场来说，先画出其等势线会更方便一些。
(1) 一对点电荷的平面等势线
对于两个点电荷，不妨取，正
电荷在x轴的正方向，负电荷在x轴的负方向，它们到原点的距离定为a=0.02；假设平面的范围为
=0.05，
=0.04。则其程序如下：
q1=1e-9; % 设置正电荷电量
q2=-1e-9; % 设置负电荷电量
a=0.02; % 设置电荷到原点的距离
xx0=0.05; % 设置X轴的范围
yy0=0.04; % 设置Y轴的范围
x=linspace(-xx0,xx0,20); % 将X轴进行20等分
y=linspace(-yy0,yy0,50); % 将Y轴进行50等分
[X,Y]=meshgrid(x); % 形成网格坐标
r1=sqrt((X-a).2+Y.2); % 各点到正电荷的距离
r2=sqrt((X+a).2+Y.2); % 各点到负电荷的距离
U=kq1./r1+kq2./r2; % 各点的电势
u0=500; % 设定最大电势的大小
u=linspace(u0,-u0,11); % 计算各等势线的电势
contour(X,Y,U,u,‘k-’); % 画出所有的等势线
Grid on % 形成网格
Hold on % 保持图形
plot(0,0,‘o’,‘MarkerSize’,12) % 画电荷
xlabel(‘x’,‘fontsize’,16) % 用16号字体标出X轴
ylabel(‘y’,‘fontsize’,16) % 用16号字体标出Y轴
title(‘一对相异电荷的等势线图’,‘fontsize’,20) % 添加标题
图五 一对相异电荷的等势线图
(2) 一对点电荷的平面电场线
各点的电场强度方向代替电力线。根据电势的梯度可以求出各点的场强的两个分量再在此方向上标上箭头。其程序如下所示：
[Ex,Ey]=gradient(-U); % 各点的场强的两个分量
E=sqrt(Ex.2+Ey.2); % 各点的合场强
Ex=Ex./E; % 为使箭头等长，将场强归一化
Ey=Ey./E;
quiver(X,Y,Ex,Ey)； % 标出各网点的电场强度方向
title(‘一对相异电荷的等势线图和电场线图’,‘fontsize’,20)； 标出标题
其图如图六所示：
图六 一对相异电荷的等势线图和电场线图
四、 实验总结：
1、以上方法同样适合不同极性和大小的电荷的组合，只需要按要求设置不同的电荷电量即可。
2、以Matlab强大的数值计算和绘图功能，将抽象的电场形象化，大大增加了我们对电场的感性认识，为今后的学习打了牢固的基础。
3、由于对Matlab 的不熟悉，对多点电荷的三维电场分布情况，无法将具体的思路转化成实际的程序语言，所以无法绘出它的三维图形。今后一定要多学习Matlab这方面的知识。