필터 지우기
필터 지우기

surface charge columb law

조회 수: 5 (최근 30일)
abdullah
abdullah 2022년 11월 6일
How Matlab used the columb law to calculate the electric field for points charges , surfaces
does the surface is divided into points charges? What is the approach for this fuction ?
function SurfaceChargeField
% Program to find V and E of surface (rectangular and circular)
% Charge Distribution
%=============================================
% Adjust These Parameters
% Evaluation of electric field at P(x,y,z)
PLoc=[0 0 1] ;
% The discretization element side d
ElementSize=0.001;
Eps0=10e-9/(36*pi);
Eps=Eps0; % Medium is free space
%=====================Surface Parameters
% Charge distribution is assumed to be on a rectangle, Disk or Partial Disk
SurfParam.SurfaceShape= 2; % 1 for rectangle Patch
% 2 for disk, annular disk or partial disk
SurfParam.Rho0=5*10e-9 ; % Unifrom Rho
SurfParam.UseRandChargeDist = 0; % 0 for uninform charge distribution and 1 for random distribution
if SurfParam.SurfaceShape ==1 % Rectangular Patch Shape
SurfParam.X=[0 1]; % X Min and Max
SurfParam.Y=[-1 1]; % Y Min and Max
elseif SurfParam.SurfaceShape ==2
SurfParam.DiskRad =2 ; % Outer Radius of disk
SurfParam.InnerDiskRad=1; % Inner disk rad (Use 0 for full disk)
SurfParam.Phi=[0 2*pi]; % Start and End Angles For Partial Disk
end
%============================================================
%==========================================================
% Use function rho_s_fun below to adjust the charge density
dx=ElementSize;
dy=ElementSize;
ds=dx*dy; % Square Element Areal
format short % prints only 6 decimal places
%====================Check Parameters======================
if SurfParam.SurfaceShape== 1
if (SurfParam.X(2) <=SurfParam.X(1) || SurfParam.Y(2) <=SurfParam.Y(1))
display('Error Check Dimension Values');
return;
end
elseif SurfParam.SurfaceShape== 2
if (SurfParam.DiskRad <=SurfParam.InnerDiskRad)
display('Error DiskRad is Tool Small');
return;
end
%======= Area Around the Disk
SurfParam.X=SurfParam.DiskRad*[-1 1];
SurfParam.Y=SurfParam.DiskRad*[-1 1];
end
tic;
X= SurfParam.X(1) : dx : SurfParam.X(2);
Y= SurfParam.Y(1) : dy : SurfParam.Y(2);
[XLoc,YLoc]=meshgrid(X,Y);
YLoc=flipud(YLoc); % Make First Row For Maximum Y
SurfParam.ArraySize=size(XLoc);
ZLoc=zeros(SurfParam.ArraySize);
% RVec is (r - r')
RVec=cat(3,(PLoc(1)-XLoc),(PLoc(2)-YLoc), (PLoc(3)-ZLoc)); % | r - rprime |
RMag=sqrt(RVec(:,:,1).^2+RVec(:,:,2).^2+RVec(:,:,3).^2);
rho_s=rho_s_fun(XLoc,YLoc,SurfParam);
%=========== Calculate Potential V
VNumerator=ds*rho_s;
VDenominator=4*pi*Eps*RMag;
ElemV=VNumerator./VDenominator;
VField=sum(sum(ElemV));
%=========== Calculate E Field
ENumerator=ds*repmat(rho_s,1,1,3).*RVec;
EDenominator=4*pi*Eps*repmat(RMag.^3,1,1,3); % | r - rprime | ^ 3
ElemField=ENumerator./EDenominator;
EField=[sum(sum(ElemField(:,:,1))),sum(sum(ElemField(:,:,2))),sum(sum(ElemField(:,:,3)))];
% ===============Display results
figure
imagesc(rho_s);
title('Charge Distribution');
colormap(gray)
colorbar
disp('The computed potential in V is');
disp(VField);
disp('The computed electric field in V/m is');
disp(EField);
toc
end
function rho_s=rho_s_fun(XLoc,YLoc,SurfParam)
SurfaceShape=SurfParam.SurfaceShape;
Rho0=SurfParam.Rho0;
if SurfaceShape==1
rho_s=Rho0*ones(SurfParam.ArraySize) ; % uniform charge density
%=================================================
elseif SurfaceShape==2
%For Disk set rho_s to be 0 outsie the disk area
DiskRad=SurfParam.DiskRad;
InnerDiskRad=SurfParam.InnerDiskRad;
ElementRho=sqrt(XLoc.^2+YLoc.^2);
PhiMin=SurfParam.Phi(1);
PhiMax=SurfParam.Phi(2);
ElemPhi=wrapTo2Pi(atan2(YLoc,XLoc)); % Notice atan2 function gives angle from -pi to pi
rho_s=zeros(size(XLoc));
DiskElem=(ElementRho >= InnerDiskRad & ElementRho <= DiskRad & ...
ElemPhi >= PhiMin & ElemPhi <= PhiMax);
% We choose elements inside the disk
rho_s=Rho0*DiskElem ; % Uniform Charge
%======================================================
% Example of Nonuniform charge density
% rho_s=Rho0*ElementRho.^2.*DiskElem;
%===========================================
% Rand Charge Distribution ( TO show elements)
if SurfParam.UseRandChargeDist==1
rho_s=rho_s.*randn(size(rho_s));
end
%============================
end
end

답변 (0개)

카테고리

Help CenterFile Exchange에서 Line Plots에 대해 자세히 알아보기

제품


릴리스

R2022b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by