plotedfield.m

解决麦克斯韦的matlab源文件

function plotedfield(mesh,ev,titstr)
% Quiver plot of edge element vector field on midpoints of triangles
%
% mesh -> data structure for 2D unstructured grid
% ev -> Column vector of length mesh.Ne containing the edge dofs
%

if (size(ev,1) ~= mesh.Ne), error('Size mismatch for argument vector'); end
if (size(ev,2) ~= 1), error('ev must be a column vector'); end
if (nargin < 3), titstr = 'Center values of edge element v.f.'; end

cntr = zeros(mesh.Nt,2);
vals = zeros(mesh.Nt,2);

for i=1:mesh.Nt
% Global indices of vertices and edges of triangle
  vidx = mesh.trv(i,:);
  eidx = mesh.tre(i,:);
% Compute information about geometry of current triangle
  tgeo = mesh.vt(vidx,:)';
  x = tgeo(1,:);
  y = tgeo(2,:);
  c = sum(tgeo')/3.0;
  vofs = tgeo-[c' c' c'];
  cntr(i,:) = c;
  d = (x(2)-x(1))*(y(3)-y(1))-(x(3)-x(1))*(y(2)-y(1));
  area = abs(d)/2;
% Get value of edge element function at center of gravity
% (Uses gradients of barycentric coordinate functions)
  G = 1/(2*area)*[y(2)-y(3) , y(3)-y(1) , y(1)-y(2);...
		  x(3)-x(2) , x(1)-x(3) , x(2)-x(1)];
  vals(i,:) = ([G(:,3)-G(:,2) , G(:,1)-G(:,3) , G(:,2)-G(:,1)]/3*...
	       (ev(eidx,1).*mesh.treo(i,:)'))';
end

hold on;
title(titstr);
% Plot boundary of domain
ebd_idx  = find(mesh.ebfl(:,1) ~= 0);
for i=ebd_idx'
  ev = [mesh.vt(mesh.ep(i,1),:)' mesh.vt(mesh.ep(i,2),:)']; 
  plot(ev(1,:),ev(2,:),'k-');
end

quiver(cntr(:,1),cntr(:,2),vals(:,1),vals(:,2),0.75,'b-');
hold off;