surface charge columb law
3 visualizzazioni (ultimi 30 giorni)
Mostra commenti meno recenti
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 Commenti
Risposte (0)
Vedere anche
Categorie
Scopri di più su Line Plots in Help Center e File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!