diff --git a/matlab/Magnetic_Field_GLS_2020/Cryogenic/Draw_B_Ellip_170_geom_proto.m b/matlab/Magnetic_Field_GLS_2020/Cryogenic/Draw_B_Ellip_170_geom_proto.m index e919daa..107ced4 100644 --- a/matlab/Magnetic_Field_GLS_2020/Cryogenic/Draw_B_Ellip_170_geom_proto.m +++ b/matlab/Magnetic_Field_GLS_2020/Cryogenic/Draw_B_Ellip_170_geom_proto.m @@ -1,135 +1,169 @@ % Display the geometry of the gt170 gyrotron gun for the first prototype. % The magnetic field is the one of the asg magnet % This also creates the h5 files of the geometry and of the magnetic field % to be used in FENNECS +%% Define and calculate the magnetic field z=linspace(-0.05,0.2,100); r=linspace(0.02,0.09,100); [Z,R] = meshgrid(z,r); % I = [ 9.7 7.7 96.2499 91.0124 87.4723]; % I = [0 0 6.79541 1.967 88.2799]; % I = [0 0 5 1.967 88.2799]; -I = [61.56 66.45 113.43 108.09 ]; +% Current in coils +I = 0.6*[61.56 66.45 113.43 108.09 ]; res=zeros([size(Z),3]); -%Bz=zeros(size(Z)); -%Br=zeros(size(Z)); + for i=1:size(R,1) res(i,:,:) = B_Ellip_Cryogenic_170({'aphi','bz','br'},'cryogenic',I,R(i,:),Z(i,:)); -%Bz(i,:) = B_Ellip_Cryogenic_170('bz','asg_modified_201206',I,R(i,:),Z(i,:)); -%Br(i,:) = B_Ellip_Cryogenic_170('br','asg_modified_201206',I,R(i,:),Z(i,:)); end Aphi=res(:,:,1); Bz=res(:,:,2); Br=res(:,:,3); -% zphi=linspace(0.005,0.0195,100); -% rphi=linspace(0.06375,0.079,150); -% b=max(rphi); -% a=min(rphi); -% phib=0; -% phia=-30000; -% phi=((phib-phia)*log(rphi)+phia*log(b)-phib*log(a))/log(b/a); -% Er=((phib-phia)./rphi)/log(b/a); -% [Zphi,Rphi] = meshgrid(zphi,rphi); -% Phi=repmat(phi',1,length(zphi)); - - %% Define the individual boundaries -geom=importproto('../geometry_proto.data',[2,Inf]); +%read from input file +geom=importproto('geometry_proto_modif.data',[2,Inf]); geomcells={}; j=1; i=1; n=1; +% clean and separate the curves while i<=size(geom.Z,1) if isnan(geom.Z(i)) j=j+1; while isnan(geom.Z(i)) i=i+1; end n=1; end geomcells{j}.Z(n)=geom.Z(i); geomcells{j}.R(n)=geom.R(i); i=i+1; n=n+1; end + +geo_old = geomcells{2}; +%% Correcting the corners + +% Body ----------------------------------------- +N = 30; + +% First corner (lone corner) +start1 = 32; +end1 = 36; +R_corner = cat(2,linspace(geomcells{2}.R(start1),0.0869994,30),linspace(0.0869994,geomcells{2}.R(end1),30)); +Z_corner = cat(2,linspace(geomcells{2}.Z(start1),0.104844,30),linspace(0.104844,geomcells{2}.Z(end1),30)); +% Incomplete square +start2 =510; +R_square = cat(2,linspace(geomcells{2}.R(start2),0.0819994,N),linspace(0.0819994,0.0899991,N),linspace(0.08999911,0.0900011,N),linspace(0.0899911,0.085,N)); +Z_square = cat(2,linspace(geomcells{2}.Z(start2),0.0143061,N),linspace(0.0143061,0.0143061,N),linspace(0.0143061,0.00316498,N),linspace(0.00316498,0.00330665,N)); + +geomcells{2}.R = cat(2,geomcells{2}.R(1:start1),R_corner,geomcells{2}.R(end1:start2),R_square); +geomcells{2}.Z = cat(2,geomcells{2}.Z(1:start1),Z_corner,geomcells{2}.Z(end1:start2),Z_square); + +% Cathode -------------------------------------- + +%Incomplete square +end1 = 320; +R_square = cat(2,linspace(0.085,0.0900011,N),linspace(0.0900011,0.0899991,N),linspace(0.089991,0.08,N)); +Z_square = cat(2,linspace(-0.0829438,-0.0829483,N),linspace(-0.0829483,-0.0939433,N),linspace(-0.0939433,-0.094,N)); + +geomcells{1}.R = cat(2,R_square,geomcells{1}.R(end1:end)); +geomcells{1}.Z = cat(2,Z_square,geomcells{1}.Z(end1:end)); + +%% Define the ceramic boundaries + +geomcells{4}.Z = cat(2,geomcells{2}.Z,linspace(geomcells{2}.Z(end),-0.08294,50),geomcells{1}.Z(1:668),linspace(geomcells{1}.Z(668),geomcells{3}.Z(4),50),geomcells{3}.Z(4:end)); +geomcells{4}.R = cat(2,geomcells{2}.R,0.085*ones(1,50),geomcells{1}.R(1:668),0.085*ones(1,50),geomcells{3}.R(4:end)); + +% create the geometry structure geomcells{1}.name='Cathode'; geomcells{1}.Dval=-30000; geomcells{2}.name='Body'; geomcells{2}.Dval=0; geomcells{3}.name='Coaxial insert'; geomcells{3}.Dval=0; -geomcells{4}.name='insulator left'; +geomcells{4}.name='Ceramics'; geomcells{4}.Dval=0; -geomcells{5}.name='insulator right'; -geomcells{5}.Dval=0; -for k=1:length(geomcells) -geomcells{k}.order=3; +% geomcells{4}.name='insulator left'; +% geomcells{4}.Dval=0; +% geomcells{5}.name='insulator right'; +% geomcells{5}.Dval=0; + +for k=1:(length(geomcells)-1) +geomcells{k}.order=2; geomcells{k}.dim=2; geomcells{k}.epsce=1e-9; geomcells{k}.epsge=1e-9; geomcells{k}.type=0; geomcells{k}.periodic=0; end -geomcells{end-1}.type=2; -geomcells{end}.type=2; +% geomcells{end-1}.type=2; +% geomcells{end}.type=2; + +k = length(geomcells); +geomcells{k}.order=2; +geomcells{k}.dim=2; +geomcells{k}.epsce=1e-9; +geomcells{k}.epsge=1e-9; +geomcells{k}.type=2; +geomcells{k}.periodic=0; + %% Plots f=figure; for k=1:length(geomcells) plothandle=plot(geomcells{k}.Z, geomcells{k}.R,'k-','linewidth',1.5); hold on geomcells{k}.points=[geomcells{k}.Z; geomcells{k}.R]; order=geomcells{k}.order; knots=linspace(0,1,length(geomcells{k}.Z)-(order-2)); knots=augknt(knots, order); sizec=size(geomcells{k}.Z); order=length(knots)-sizec(end); coeffs=[geomcells{k}.Z; geomcells{k}.R]; pp=spmak(knots,coeffs); s=linspace(0,1,10000); fittedpos=fnval(pp,s); plot(fittedpos(1,:),fittedpos(2,:),'x-') end hold on %end %axis equal -[~,cont1]=contour(Z,R,R.*Aphi,15,'r:','linewidth',3); -%xlim([min(z) 0.04]) -%ylim([0.055 0.085]) -%[~,cont2]=contour(Zphi,Rphi,Phi,20,'b'); -%rectangle('Position',[-0.011, 0.06375, 0.032+0.011, 0.081-0.06375],'EdgeColor','magenta','Linestyle','--') +%Plot the magnetic field lines +%[~,cont1]=contour(Z,R,R.*Aphi,15,'r:','linewidth',3); -legend([plothandle,cont1],{'Gun geometry', 'Magnetic field lines','Wall approximation'},'location','southwest') +% legend([plothandle,cont1],{'Gun geometry', 'Magnetic field lines','Wall approximation'},'location','southwest') f.PaperUnits='centimeters'; f.PaperSize=[12,8]; xlabel('z [m]') ylabel('r [m]') print(f,'phiBprofile_protoasg','-dpdf','-fillpage') savefig(f,'phiBprofile_protoasg') hold off %% Save magnetic field and geometry to disk save=true; overwrite=true; if save -% idr=1:1:length(magnet.r); +% idr=1:1:length(magnet.r); % idz=1:1:length(magnet.z); % Aphi=zeros(length(idr),length(idz)); % Bz=zeros(length(idr),length(idz)); % Br=zeros(length(idr),length(idz)); % % for i=1:size(magnet.subcoils,1) % Aphi=Aphi+Icoil(i)*magnet.subcoils{i,1}(idr,idz); % Bz=Bz+Icoil(i)*magnet.subcoils{i,2}(idr,idz); % Br=Br+Icoil(i)*magnet.subcoils{i,3}(idr,idz); % end % savemagtoh5([name,'.h5'],magnet.r,magnet.z,Aphi,Br,Bz,overwrite); - savegeomtoh5('proto_geom.h5',geomcells,1e-2,overwrite); + savegeomtoh5('geom_crmcs.h5',geomcells,1e-2,overwrite); end diff --git a/matlab/Magnetic_Field_GLS_2020/Cryogenic/Draw_B_Ellip_170_geom_refurb.m b/matlab/Magnetic_Field_GLS_2020/Cryogenic/Draw_B_Ellip_170_geom_refurb.m index e15ba48..92d0c74 100644 --- a/matlab/Magnetic_Field_GLS_2020/Cryogenic/Draw_B_Ellip_170_geom_refurb.m +++ b/matlab/Magnetic_Field_GLS_2020/Cryogenic/Draw_B_Ellip_170_geom_refurb.m @@ -1,222 +1,244 @@ % creates the magnetic field h5 file necessary for fennecs % This uses the geometry of the refurbished 170GHz coaxial gyrotron gun % The magnetic field is the one created by the asg magnet % Both the magnet and geometry are definde in the Final report of the % Development of the european gyrotron CCCGDS6 %% Import and calculate the magnetic field -magnet=load('asg_magnet_red.mat'); +magnet=load('asg_magnet_ceramics.mat'); %% calculate magnetic field and magnetic vector idr=1:10:length(magnet.r); idz=1:10:length(magnet.z); % Define the coils currents %I = [ 9.7 7.7 96.2499 91.0124 87.4723]; % I = [0 0 4.5 1.967 88.2799]; % name='phiBprofile_refurbasg_4_5_red'; % I = [0 0 4.8 1.967 88.2799]; % name='phiBprofile_refurbasg_well_red'; % I = [0 0 5 1.967 88.2799]; % name='phiBprofile_refurbasg_5_red'; % I = [0 0 5.15 1.967 88.2799]; % name='phiBprofile_refurbasg_limdwn_red'; %I = [0 0 5.35 1.967 88.2799]; %name='phiBprofile_refurbasg_limup_red'; % I = [0 0 5.5 1.967 88.2799]; % name='phiBprofile_refurbasg_5_5_red'; %I = [0 0 5.6 1.967 88.2799]; %name='phiBprofile_refurbasg_5_6_red'; %I = [0 0 5.7 1.967 88.2799]; %name='phiBprofile_refurbasg_5_7_red'; % I = [0 0 5.75 1.967 88.2799]; % name='phiBprofile_refurbasg_5_75_red'; % I = [0 0 5.8 1.967 88.2799]; % name='phiBprofile_refurbasg_5_8_red'; % I = [0 0 5.85 1.967 88.2799]; % name='phiBprofile_refurbasg_5_85_red'; % I = [0 0 5.9 1.967 88.2799]; % name='phiBprofile_refurbasg_5_9_red'; % I = [0 0 6 1.967 88.2799]; % name='phiBprofile_refurbasg_6_red'; % I = [0 0 6.4 1.967 88.2799]; % name='phiBprofile_refurbasg_6_4_red'; % I = [0 0 6.79541 1.967 88.2799]; % name='phiBprofile_refurbasg_nominal_red'; % for Ib=linspace(4.6,6.8,12) %Ib=4.4; %I = [0 0 Ib 1.967 88.2799]; %name=sprintf('phiBprofile_refurbasg_%2d_red',floor(Ib*10)); I = 0.6*[0 0 6.79541 1.967 88.2799]; name='phiBprofile_refurbasg_nominal_0_6'; Icoil=[I(1) I(2) -I(3)-I(5) -I(3)-I(5) I(4)+I(5) I(4)+I(5) I(5) I(5)]; rmin = [0.13958 0.13958 0.13710 0.16952 0.13670 0.19101 0.13458 0.17197 ]; rmax = [0.142694 0.142694 0.16733 0.20005 0.18941 0.23666 0.16972 0.20412 ]; zmin = [0.07854 0.13854 0.18200 0.18200 0.31778 0.31748 0.50685 0.50685 ]; zmax = [0.09846 0.15846 0.28350 0.28350 0.48977 0.49017 0.71128 0.71128 ]; % Nturn = [ 218 218 2528.5 2626 7689 11517.5 6110 9656.5 ]; % na = 2. *[ 10 10 10 10 10 10 10 10 ]; % Number of subcoils (radial direction % nb = 2. *[ 10 10 10 10 10 10 10 10 ]; % Number of subcoils in longitudinal direction % Icoil=Icoil.*Nturn./(na.*nb); [Z,R] = meshgrid(magnet.z(idz),magnet.r(idr)); Aphi=zeros(length(idr),length(idz)); Bz=zeros(length(idr),length(idz)); Br=zeros(length(idr),length(idz)); for i=1:size(magnet.subcoils,1) Aphi=Aphi+Icoil(i)*magnet.subcoils{i,1}(idr,idz); Bz=Bz+Icoil(i)*magnet.subcoils{i,2}(idr,idz); Br=Br+Icoil(i)*magnet.subcoils{i,3}(idr,idz); end %% Define the individual boundaries -geom=importrefurb('../refurb_modif_improv.data'); +geom=importrefurb('Magnetic_Field_GLS_2020/refurb_modif_improv.data'); geomcells={}; j=1; i=1; n=1; while i<=size(geom.Z,1) if isnan(geom.Z(i)) j=j+1; while isnan(geom.Z(i)) i=i+1; end n=1; end geomcells{j}.Z(n)=geom.Z(i)/1e3; geomcells{j}.R(n)=geom.R(i)/1e3; i=i+1; n=n+1; end + +geomcells{4}.Z = cat(2,geomcells{2}.Z,geomcells{1}.Z,geomcells{3}.Z); +geomcells{4}.R = cat(2,geomcells{2}.R,geomcells{1}.R,geomcells{3}.R); + geomcells{1}.name='Cathode'; geomcells{1}.Dval=-30000; geomcells{2}.name='Body'; geomcells{2}.Dval=0; geomcells{3}.name='Coaxial insert'; geomcells{3}.Dval=0; -for k=1:length(geomcells) +geomcells{4}.name='Ceramics'; +geomcells{4}.Dval=0; + +for k=1:(length(geomcells)-1) geomcells{k}.order=3; geomcells{k}.dim=2; geomcells{k}.epsce=1e-9; geomcells{k}.epsge=1e-9; geomcells{k}.type=0; geomcells{k}.periodic=0; end +k = length(geomcells); +geomcells{k}.order=3; +geomcells{k}.dim=2; +geomcells{k}.epsce=1e-9; +geomcells{k}.epsge=1e-9; +geomcells{k}.type=2; +geomcells{k}.periodic=0; + %% Load the original boundary -geomorig=importrefurb('../refurb.data'); +geomorig=importrefurb('Magnetic_Field_GLS_2020/refurb.data'); geomcellsorig={}; j=1; i=1; n=1; while i<=size(geomorig.Z,1) if isnan(geomorig.Z(i)) j=j+1; while isnan(geomorig.Z(i)) i=i+1; end n=1; end geomcellsorig{j}.Z(n)=geomorig.Z(i)/1e3; geomcellsorig{j}.R(n)=geomorig.R(i)/1e3; i=i+1; n=n+1; end + +geomcellsorig{4}.Z = cat(2,geomcells{2}.Z,geomcells{1}.Z,geomcells{3}.Z); +geomcellsorig{4}.R = cat(2,geomcells{2}.R,geomcells{1}.R,geomcells{3}.R); + geomcellsorig{1}.name='Cathode'; geomcellsorig{1}.Dval=-30000; geomcellsorig{2}.name='Body'; geomcellsorig{2}.Dval=0; geomcellsorig{3}.name='Coaxial insert'; geomcellsorig{3}.Dval=0; +geomcellsorig{4}.name='Ceramics'; +geomcellsorig{4}.Dval=0; + for k=1:length(geomcellsorig) geomcellsorig{k}.order=3; geomcellsorig{k}.dim=2; geomcellsorig{k}.epsce=1e-9; geomcellsorig{k}.epsge=1e-9; end %% Plots f=figure; for k=1:length(geomcellsorig) plothandleorig=plot(geomcellsorig{k}.Z*1e3, geomcellsorig{k}.R*1e3,'r-','linewidth',2); hold on end for k=1:length(geomcells) %plothandle=plot(geomcells{k}.Z, geomcells{k}.R,'k-','linewidth',1.5); hold on geomcells{k}.points=[geomcells{k}.Z; geomcells{k}.R]; order=geomcells{k}.order; knots=linspace(0,1,length(geomcells{k}.Z)-(order-2)); knots=augknt(knots, order); sizec=size(geomcells{k}.Z); order=length(knots)-sizec(end); coeffs=[geomcells{k}.Z; geomcells{k}.R]; pp=spmak(knots,coeffs); s=linspace(0,1,500); fittedpos=fnval(pp,s); plot(fittedpos(1,:)*1e3,fittedpos(2,:)*1e3,'--','linewidth',2.2) end - +%% %axis equal raphi=R.*Aphi; lvls=logspace(-4,log10(max(raphi(:))),50); %[~,cont1]=contour(Z*1e3,R*1e3,raphi,lvls,'b:','linewidth',1.5); %xlim([min(magnet.z) max(magnet.z)]*1e3) %ylim([min(magnet.r) max(magnet.r)]*1e3) %[~,cont2]=contour(Zphi,Rphi,Phi,20,'b'); rectangle('Position',[-0.1, 0.045, 0.292, 0.092-0.045]*1e3,'EdgeColor','black','Linestyle','--','linewidth',2) %legend([plothandle,cont1],{'Gun geometry', 'Magnetic field lines'},'location','southwest') %f.PaperUnits='centimeters'; %f.PaperSize=[12,8]; xlabel('z [mm]') ylabel('r [mm]') axis equal grid on ylim([0,inf]) xlim([-300 750]) for i=1:length(rmin) rectangle('Position',[zmin(i) rmin(i) zmax(i)-zmin(i) rmax(i)-rmin(i)]*1e3,'edgecolor','r') text(zmin(i)*1e3-10,(rmin(i)+rmax(i))/2*1e3,sprintf('S%i',i),'fontsize',10,'color','r') end print(f,name,'-dpdf','-fillpage') savefig(f,name) set(f, 'Color', 'w'); export_fig(f,name,'-eps') hold off %% Save magnetic field and geometry to disk save=true; overwrite=true; if save idr=1:1:length(magnet.r); idz=1:1:length(magnet.z); Aphi=zeros(length(idr),length(idz)); Bz=zeros(length(idr),length(idz)); Br=zeros(length(idr),length(idz)); for i=1:size(magnet.subcoils,1) Aphi=Aphi+Icoil(i)*magnet.subcoils{i,1}(idr,idz); Bz=Bz+Icoil(i)*magnet.subcoils{i,2}(idr,idz); Br=Br+Icoil(i)*magnet.subcoils{i,3}(idr,idz); end savemagtoh5([name,'.h5'],magnet.r,magnet.z,Aphi,Br,Bz,overwrite); - savegeomtoh5('refurb_geom_improv.h5',geomcells,1e-2,overwrite); + savegeomtoh5('geom_ceramics.h5',geomcells,1e-2,overwrite); end % end diff --git a/matlab/Magnetic_Field_GLS_2020/Cryogenic/asg_magnet.m b/matlab/Magnetic_Field_GLS_2020/Cryogenic/asg_magnet.m index 6aa3496..9819f31 100644 --- a/matlab/Magnetic_Field_GLS_2020/Cryogenic/asg_magnet.m +++ b/matlab/Magnetic_Field_GLS_2020/Cryogenic/asg_magnet.m @@ -1,148 +1,148 @@ % Computes the magnetic field of each submagnet in the asg % superconducting magnet assembly % % The magnetic field is calculated for a unit current at positions % defined by r and z. % The fields are then saved in asg_magnet.mat % Evaluation grid space -z=linspace(-0.3,0.3,6000); -r=linspace(0.0,0.1,1000); +z=linspace(-0.2,0.35,5500); +r=linspace(0.008,0.095,1000); % z=linspace(-0.15,0.3,100); % r=linspace(0.0,0.1,100); [Z,R] = meshgrid(z,r); mode={'aphi','bz','br'}; % Transform array of evaluation points in a one-dimensional array % r1d = reshape(R,1,numel(R)); z1d = reshape(Z,1,numel(Z)); % % Geometry of the asg magnet as built with 2 layers removed on coils#3 and 3 on % coil #5 after the problems with the conductor % % Reference current % [0 0 6.79541 1.967 88.2799]=[I(1) I(2) I(3) I(4) I(5)] Ncoil = 8; Nturn = [ 218 218 2528.5 2626 7689 11517.5 6110 9656.5 ]; rmin = [0.13958 0.13958 0.13710 0.16952 0.13670 0.19101 0.13458 0.17197 ]; rmax = [0.142694 0.142694 0.16733 0.20005 0.18941 0.23666 0.16972 0.20412 ]; zmin = [0.07854 0.13854 0.18200 0.18200 0.31778 0.31748 0.50685 0.50685 ]; zmax = [0.09846 0.15846 0.28350 0.28350 0.48977 0.49017 0.71128 0.71128 ]; rcen = (rmin + rmax)/2; zcen = (zmin + zmax)/2; dr = rmax - rmin; dz = zmax - zmin; coilsurf= dr .* dz; a = rcen-dr/2; b = rcen+dr/2; l = dz; na = 2. *[ 10 10 10 10 10 10 10 10 ]; % Number of subcoils (radial direction nb = 2. *[ 10 10 10 10 10 10 10 10 ]; % Number of subcoils in longitudinal direction % Current = [ I(1) I(2) -I(3)-I(5) -I(3)-I(5) I(4)+I(5) I(4)+I(5) I(5) I(5) ]; % Current flowing in each power supply % JTot = Nturn; % Total current flowing in each coil % J = JTot./(dr.*dz); % Current density in each coil Ic1 = JTot ; % Total current flowing in each coil [A] %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Construct array of current loops %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % % Number of coils ncoil = length(rcen); if(~iscell(mode)) nbmodes=length({mode}); else nbmodes=length(mode); end subcoils=cell(ncoil,nbmodes); for icoil=1:ncoil for imode=1:nbmodes subcoils{icoil,imode}=zeros(length(r),length(z)); end end for icoil=1:ncoil %subcoils{icoil}=zeros(length(r),length(z),nbmodes); B=zeros(length(r),length(z),nbmodes); parfor ir=1:length(r) [Z,R] = meshgrid(z,r(ir)); r1d = reshape(R,1,numel(R)); z1d = reshape(Z,1,numel(Z)); % Number of subcoils ntot = sum(na(icoil).*nb(icoil)); % Preallocate space r2 = zeros(1,ntot); z2 = zeros(1,ntot); Ic = zeros(ntot,1); % isub is the index of subcoil % Positions of the coils subdivided isub = 0; for ib=1:nb(icoil) for ia=1:na(icoil) isub = isub +1; r2(isub) = (rcen(icoil) - dr(icoil)/2) + (dr(icoil)/nb(icoil))*(ib-0.5); z2(isub) = (zcen(icoil) - dz(icoil)/2) + (dz(icoil)/na(icoil))*(ia-0.5); Ic(isub) = Ic1(icoil) / (na(icoil)*nb(icoil)); end end % % Suppress locations where Green function diverges (i.e. on the sources) % for i=1:length(r2) [I] = find((r1d.*r1d -r2(i)*r2(i) == 0) & (z1d.*z1d - z2(i)^2)==0); if ~isempty(I) r1d(I)=NaN; z1d(I)=NaN; end end % Call greenem b = greenem_jph(mode,r1d,z1d,r2,z2); temp2=zeros(size(R,2),nbmodes); % b is a rectangular matrix (length(r1d) x length(r2)) which needs to be multiplied by the % column vector Ic to get the proper contribution of each source. for i=1:size(b,3) temp = b(:,:,i)*Ic; temp((temp==Inf)|(temp==-Inf))= 0; temp(isnan(temp)) = 0; temp2(:,i) = reshape(temp,size(R,1),size(R,2)); end B(ir,:,:)=temp2(:,:); end for imode=1:nbmodes subcoils{icoil,imode}=B(:,:,imode); end end -save('asg_magnet_red.mat','r','z','subcoils','mode','Ncoil','na','nb','Nturn','rmin','rmax','zmin','zmax','-v7.3') +save('asg_magnet_ceramics.mat','r','z','subcoils','mode','Ncoil','na','nb','Nturn','rmin','rmax','zmin','zmax','-v7.3')