Commit 5955f958 authored by Lore, Jeremy's avatar Lore, Jeremy
Browse files

Update to preserve input size and unit test

parent b32e92f5

matlab/bfield_library_jdl/Poincare/make_poincare_g_and_rmp.m

+5 −5
Original line number Diff line number Diff line
@@ -2,7 +2,7 @@ clearvars; 
tic;

connect_dots = 0;



gfile_name = 'C:\Work\DIII-D\164723\g164723.03059_410';

gfile_name = 'C:\Work\DIII-D\164273\g164723.03059_410';

g = readg_g3d(gfile_name);



% bfield.type = 'gfile';
@@ -18,7 +18,7 @@ bfield.type = 'gfile+coils'; 
bfield.coil = rmp.coil;

bfield.current = rmp.current;

isaxisym = 0;

nsym = 3;

bfield.nsym = 3;





Rend   = 2.25;    % if [] this is set to max radius of g.bdry
@@ -72,7 +72,7 @@ if connect_dots 
        z(:,i) = z(sort_inds,i);

        phi = phi(sort_inds);

        if ~isaxisym

            phi = mod(phi,2*pi/nsym);

            phi = mod(phi,2*pi/bfield.nsym);

            inds = find(abs(phi - phistart_deg) < 1e-6);

            rfinal(:,i) = r(inds,i);

            zfinal(:,i) = z(inds,i);
@@ -82,8 +82,8 @@ if connect_dots 
    z = zfinal;

else

    if ~isaxisym

        phi = mod(phi,2*pi/nsym);

        phi(abs(phi - 2*pi/nsym) < 1e-6) = 0;

        phi = mod(phi,2*pi/bfield.nsym);

        phi(abs(phi - 2*pi/bfield.nsym) < 1e-6) = 0;

        inds = find(abs(phi - phistart) < 1e-6);

        r = r(inds,:);

        z = z(inds,:);

matlab/bfield_library_jdl/calc_B/bfield_bs_cyl.m

+0 −6
Original line number Diff line number Diff line
@@ -4,12 +4,6 @@ if nargin < 6 
    nowarn = 0;

end

    

if isrow(P_r)

    P_r = P_r.';

    P_phi = P_phi.';

    P_z = P_z.';

end



cp = cos(P_phi);

sp = sin(P_phi);

matlab/bfield_library_jdl/calc_B/bfield_bs_jdl.m

+4 −3
Original line number Diff line number Diff line
@@ -3,9 +3,10 @@ function [Bx,By,Bz,Btot]=bfield_bs_jdl(P_x,P_y,P_z,coil,current) 


I=current.*1e-7; %this is mu0*I/4pi

npts = length(P_x);

Bx=zeros(npts,1);

By=zeros(npts,1);

Bz=zeros(npts,1);

sz = size(P_x);

Bx = zeros(sz);

By = zeros(sz);

Bz = zeros(sz);



coils_x = coil(:,1);

coils_y = coil(:,2);

matlab/bfield_library_jdl/calc_B/bfield_circular_coil_analytic.m

+39 −1
Original line number Diff line number Diff line 
function [Brg,Bzg,Athetag]=bfield_circular_coil_analytic(a,d,r,z)

% Written by R. H. Goulding, 10/1/98

% Calculates the magnetic field components generated by a filamentary current
% loop. It is assumed that the loop is in a plane perpendicular to the z axis
% and is centered on the z axis. All units taken to be MKS (Amperes, meters,
% Tesla). Uses the Matlab function ellipke to calculate the complete elliptic
% integrals of the first and second kind K and E.

% Input parameters: 
% 
%  a   radius of current loop 
%  d   axial distance of loop from z=0 point ( can be >0 or <0 ) %  I   current in loop  %  r   radius at which B components are calculated %  z   axial location at which B components are calculated
% Output parameters:% Br   radial B component% Bz   axial B component% Atheta  azimuthal vector potential% J.D. Lore modified from Gouldingcnst=1.e-7;   % mu0/(4*pi)
u2=4*r.*a./((r+a).^2+(z-d).^2);  % argument of elliptical integral functionsu=sqrt(u2);[K,E]=ellipke(u2);
% Brg=-2*cnst./r.^(1.5).*(z-d)./sqrt(a).*(u/2.*K-u/4.*(u2-2)./(u2-1).*E);% Bzg=2*cnst./r.*(sqrt(r./a).*u/2.*K-u./(4*(u2-1)).*(sqrt(r./a).*(u2-2)+sqrt(a./r).*u2).*E);% Athetag=2*cnst./sqrt(r).*sqrt(a).*((2./u-u).*K-2./u.*E);

sa = sqrt(a);sroa = sqrt(r./a);Brg=-2*cnst./r.^(1.5).*(z-d)./sa.*(u/2.*K-u/4.*(u2-2)./(u2-1).*E);Bzg=2*cnst./r.*(sroa.*u/2.*K-u./(4*(u2-1)).*(sroa.*(u2-2)+u2./sroa).*E);if nargout > 2
    Athetag=2*cnst./sroa.*((2./u-u).*K-2./u.*E);
end
 
 No newline at end of file
 
function [Brg,Bzg,Athetag]=bfield_circular_coil_analytic(a,d,r,z)



% Calculates the magnetic field components generated by a filamentary current

% loop. It is assumed that the loop is in a plane perpendicular to the z axis

% and is centered on the z axis. All units taken to be MKS (Amperes, meters,

% Tesla). Uses the Matlab function ellipke to calculate the complete elliptic

% integrals of the first and second kind K and E.



% Input parameters: 

% 

%  a   radius of current loop 

%  d   axial distance of loop from z=0 point ( can be >0 or <0 ) 

%  r   radius at which B components are calculated 

%  z   axial location at which B components are calculated



% Output parameters:

% Br   radial B component

% Bz   axial B component

% Atheta  azimuthal vector potential

% Written by R. H. Goulding, 10/1/98

% J.D. Lore modified from Goulding

cnst=1.e-7;   % mu0/(4*pi)



u2=4*r.*a./((r+a).^2+(z-d).^2);  % argument of elliptical integral functions

u=sqrt(u2);

[K,E]=ellipke(u2);



% Brg=-2*cnst./r.^(1.5).*(z-d)./sqrt(a).*(u/2.*K-u/4.*(u2-2)./(u2-1).*E);

% Bzg=2*cnst./r.*(sqrt(r./a).*u/2.*K-u./(4*(u2-1)).*(sqrt(r./a).*(u2-2)+sqrt(a./r).*u2).*E);

% Athetag=2*cnst./sqrt(r).*sqrt(a).*((2./u-u).*K-2./u.*E);





sa = sqrt(a);

sroa = sqrt(r./a);

Brg=-2*cnst./r.^(1.5).*(z-d)./sa.*(u/2.*K-u/4.*(u2-2)./(u2-1).*E);

Bzg=2*cnst./r.*(sroa.*u/2.*K-u./(4*(u2-1)).*(sroa.*(u2-2)+u2./sroa).*E);

if nargout > 2

    Athetag=2*cnst./sroa.*((2./u-u).*K-2./u.*E);

end

matlab/bfield_library_jdl/fieldline_following/dphi/fl_derivs_dphi_gfile_coils.m

+2 −2
Original line number Diff line number Diff line
@@ -20,6 +20,6 @@ Bout.br = Bout.br + Br; 
Bout.bphi = Bout.bphi + Bphi;

Bout.bz = Bout.bz + Bz;

   

df(1:2:N-1) = RZ(1:2:N-1).'.*Bout.br./Bout.bphi;

df(2:2:N)   = RZ(1:2:N-1).'.*Bout.bz./Bout.bphi;

df(1:2:N-1) = RZ(1:2:N-1).*Bout.br./Bout.bphi;

df(2:2:N)   = RZ(1:2:N-1).*Bout.bz./Bout.bphi;

ierr = 0;
 
 No newline at end of file