Matlab模拟电荷系电场线和等势面‌

介绍

    本文使用Maltab进行电荷系电场线和等势面的模拟，需要前置知识：电势，电场，梯度以及Matlab的基本语法。通过学习本教程，你能够学会如何绘制任意数量，任意电量，任意位置的电荷系的电场线和电势面。

开始

• 输入点电荷个数clear,number=input('请输点电荷个数:');
• 初始化电荷信息矩阵charge=zeros(number,3);.
• 获取所有电荷信息for i=1:numbers1=['请输入第',num2str(i),'个电荷的电荷量(μC):'];s2=['请输入第',num2str(i),'个电荷的x坐标:'];s3=['请输入第',num2str(i),'个电荷的y坐标:'];charge(i,1)=input(s1)*1e-6;charge(i,2)=input(s2);charge(i,3)=input(s3);end
• 定义真空静电力常数k=9*10^9;
• 定义坐标系范围，并生成网格x=-3:0.01:3; %区间[-3,3]间隔0.01y=-3:0.01:3; %区间[-3,3]间隔0.01[X,Y]=meshgrid(x,y); %X,Y都是601x601的矩阵
• 给电势矩阵分配内存U=zeros(601,601); %存储坐标系中每个位置的电势
• 计算总电势，U与XY的对应关系是U[Y,X]for j=1:numberU=U+k*charge(j,1)./sqrt((X-charge(j,2)).^2+(Y-charge(j,3)).^2);end %这里用到了计算电势的公式，以及Matlab里面的‘.’运算方法
• 对点电荷周边极大的电势进行“削平”（代码太长了，用图片表示^_^）

削平操作

        可能有点难以理解，不过尽量多看看，我就不做解释了，大家都是高材生

• 用Matlab里的等高线画等势面contour(X,Y,U,100,'Color',[1 0.5 0]) %画等势面
• 电荷绝对值求和total=sum(abs(charge(:,1))); %电荷模值求和
• 定义存放电场线数据的容器px=[];py=[]; %存放电场线数据
• 绘制电场线for a=1:number    for b=1:round(abs(charge(a,1))/total*50*number) %电场线数目由归一化电荷决定        [Ex,Ey]=gradient(-U); %电势梯度的相反数是电场强度        j=round(abs(charge(a,1))/total*50*number);        i=1;        %用参数方程分配电场线起点坐标         px(a,b,i)=charge(a,2)+0.05*2^(1/2)*cos(b*2*pi/j);        py(a,b,i)=charge(a,3)+0.05*2^(1/2)*sin(b*2*pi/j);       %绘制电场线，可以试试只写一条的结果，那样只会绘制一半        streamline(X,Y,Ex,Ey,px(a,b,i),py(a,b,i));        streamline(X,Y,-Ex,-Ey,px(a,b,i),py(a,b,i));        i=i+1;    endend
• 美化画图窗口xlabel('X轴'),ylabel('Y轴'),title('电场线和等势面'),text(1.9,-2.20,'电量单位:μC'),text(1.9,-2.80,'长度单位:m')legend('等势面','电场线');
• 给所有点电荷标上正负for a=1:number    if charge(a,1)>0    text(charge(a,2)-0.05,charge(a,3)+0.02,'+','Color','r')    else        text(charge(a,2)-0.03,charge(a,3)+0.02,'-','Color','k')    endend
• 程序编写完毕，大家试试画一张图吧！

5个电荷组成的电荷系

最终程序

clear,number=input('请输点电荷个数:');

charge=zeros(number,3);

%存放电荷信息

for i=1:number

    s1=['请输入第',num2str(i),'个电荷的电荷量(μC):'];

    s2=['请输入第',num2str(i),'个电荷的x坐标:'];

    s3=['请输入第',num2str(i),'个电荷的y坐标:'];

    charge(i,1)=input(s1)*1e-6;

    charge(i,2)=input(s2);

    charge(i,3)=input(s3);

end

k=9*10^9;%真空(静电力常数)

%by-GeekBrohood让我留在你的代码里^_^

x=-3:0.01:3;

y=-3:0.01:3;

[X,Y]=meshgrid(x,y);

U=zeros(601,601);

clf

%计算总电势，U与XY的对应关系是U[Y,X]

for j=1:number

    U=U+k*charge(j,1)./sqrt((X-charge(j,2)).^2+(Y-charge(j,3)).^2);

end

for a=1:number

    for b=-10:10 %横坐标变化

        for c=-10:10 %纵坐标变化

            if abs(U(charge(a,3)*100+301+b,charge(a,2)*100+301+c))>abs(U(charge(a,3)*100+301+5,charge(a,2)*100+301+5))

                   U(charge(a,3)*100+301+b,charge(a,2)*100+301+c)=U(charge(a,3)*100+301+5,charge(a,2)*100+301+5);

            end

        end

    end

end

contour(X,Y,U,100,'Color',[1 0.5 0])%画等势面

hold on

total=sum(abs(charge(:,1)));%电荷模值求和

px=[];

py=[];%存放电场线数据

for a=1:number

    for b=1:round(abs(charge(a,1))/total*50*number)%电场线数目由归一化电荷决定

        [Ex,Ey]=gradient(-U);

        j=round(abs(charge(a,1))/total*50*number);

        i=1;

        px(a,b,i)=charge(a,2)+0.05*2^(1/2)*cos(b*2*pi/j);

        py(a,b,i)=charge(a,3)+0.05*2^(1/2)*sin(b*2*pi/j);%分配电场线起点坐标

        streamline(X,Y,Ex,Ey,px(a,b,i),py(a,b,i));

        streamline(X,Y,-Ex,-Ey,px(a,b,i),py(a,b,i));

        i=i+1;

    end

end

hold on

xlabel('X轴'),

ylabel('Y轴'),

title('电场线和等势面'),

text(1.9,-2.20,'电量单位:μC'),

text(1.9,-2.80,'长度单位:m')

legend('等势面','电场线');

for a=1:number

    if charge(a,1)>0

    text(charge(a,2)-0.05,charge(a,3)+0.02,'+','Color','r')

    else

        text(charge(a,2)-0.03,charge(a,3)+0.02,'-','Color','k')

    end

        

end

总结

    这次我们对电荷系进行了研究，如果大家有新的思路欢迎私戳我。

作者：-42