Contents
% This code is used to calculate the elegtromagnetic of cylindrical wires. % close all windows: clc; clear all; close all;
%%PART I - Calculate integral parameter of wire for different frequencies %In PART I the specific resistance and inductance parameters are calculated % for a couple of frequencies for a copper wire. As well the skindepth and % the ratios of different resistancemodels are also calculated. The values % are written in a cellarray, which is the result [tab] of this function. % Input data for a copper wire of diameter 3.3 mm % Gauge: 8 AWG. Average Wire Diameter: 0.1285 in (3.264) mm %copper wire info: https://en.wikipedia.org/wiki/American_wire_gauge %Define global variables global u0 u0=4*pi*1e-7; %permeability of free space constant in [Vs/Am] ur = 1; % relative Permeability of Material e.g. Copper %Define coil parameters rho = 1.72e-8; % Resistivity in [Ohm*m], e.g. Copper D_wire_m = 3.3e-3; % Wire-Diameter 0.002677 m f = [1e3, 1e4, 1e5, 1e6, 1e7, 1e8, 13.5e6, 1e9, 2e9, 3e9, 4e9, 10e9, 20e9]'; % frequency vectorCalculate additonial parameter
R_m = D_wire_m/2; % Wire-Radius in [m²] area = R_m^2*pi; % Circle_area in [m²] sigma = 1/rho; % Conductivity, [S/m] R_DC = 1/(sigma*area); % DC-Resistance load per unit length in [Ohm/m] omega = 2*pi*f; % angular frequency in [Hz] delta = 1./sqrt(omega.*sigma.*u0.*ur./2); % skindepth in [m] % calculate impedance of this wire; z_math = Z_wire (omega, R_m, sigma, u0, ur); % R_DC_Ratio by an tubemodel-calculation (DC-Resistance in depence of skindepth and wire-radius) R_dc_tube_ratio = tube_model_func(delta ,R_m ,sigma); % *Collect the data for the table* tab_header = {'frequency', 'skindepth', 'L lpul', 'R_DC lpul', 'R_DC-ratio ', 'R_DC-tubemodel_ratio'; ... 'Hz', 'm', 'H/m', 'Ohm/m' , '[1]', '[1]' }; tab2 = num2cell([ f, ... % 1. column frequency delta, ... % 2. column skindepth imag(z_math)./omega, ... % 3. column impedance (H/m) real(z_math), ... % 4. column resistance (Ohm/m) real(z_math)/R_DC,... % 5. column resistance-Ratio [1] R_dc_tube_ratio ]); % 6. column resistance tube_model ratio [1] tab = [tab_header; cellfun(@num2str,(tab2), 'UniformOutput',0)]; % lpul: load per unit length %%PART II - Plot currentdensity on wire for 13.5e6 Hz % In PART II the currentdensity is calculated as a function of radius of the % wire for a frequency of 13.56 MHz and a current of 1 Amps. The current % distribution is plotted as 3D-surface plot as well as over the radius. % f_circle_Hz = 13.5e6; % frequency in Hz current_A = 1; % current in wire in Amps dr = 200; % number of steps (dicretization % create a grid for that circle [X,Y] = meshgrid(-R_m: R_m/(dr-1) :R_m , -R_m: R_m/(dr-1) :R_m ); % Calculate radius for each Point: r = sqrt(X.^2+Y.^2); % calculation of currentdensity: A = sqrt(-1j*2*pi*f_circle_Hz*sigma*u0*ur); % Nomalized argument for Bessel function [J0] = besselj ( 0 , A.* r ); % first bessel-function of zeroth order [J1] = besselj ( 1 , A.*R_m ); % first bessel-function of first order J_vec = A.*current_A.*1./(2.*pi.*R_m).*J0./J1; % Recalc only the wire with radius J_vec((X.^2+Y.^2)>=R_m^2) = NaN; X(isnan(J_vec)) = NaN; Y(isnan(J_vec)) = NaN; % Plot resultsPlot the Current density
figure (1); surf(X,Y,real(J_vec),'LineStyle','none'); colorbar; axis square; view([0 0 90]); title (['Current density in wire (R=', num2str(R_m),' m, I=',num2str(current_A),' A, f=', num2str(f_circle_Hz), ' Hz)' ]); %%Plot the Current density over radius figure (2); plot(Y(X==0), J_vec(X==0)); grid on; xlabel('radius in m'); ylabel('J in A/m²'); title (['Current density distribution over radius in wire (R=', num2str(R_m)]) %%PART III - Calculate the electrical field strength and the voltage drop for the parameters of PART II % In PART III the electric field strength for the wire of a length of 1000m % and the voltage drop are calculated. The parameter for resistance and % inductance are also computed. length_m = 1000; % lenghth of wire [m]; E_Vpm = current_A .* A./(2 .* pi.* R_m.* sigma) .* besselj ( 0 , A.* R_m )./ besselj ( 1 , A.*R_m ) % electrical field_strength in Volt per m deltaU_V = E_Vpm.*length_m % complex voltage drop in Volt absU_V = abs(deltaU_V) % absolute voltage drop in Volt R_Ohm = real(deltaU_V/current_A) % Resistance in Ohm Li_H = imag(deltaU_V/(2*pi*f_circle_Hz*current_A)) % inductance in Henr I0=current_A; %Coil current in Amps a=.15; %Coil radius in m %Define coordinates of coil center point x_p=0; y_p=0; z_p=0; %%calculate the magnetic field at a single point in space %Input test point x=0; y=0; z=.1; [Bx,By,Bz] = magnetic_field_current_loop(x,y,z,x_p,y_p,z_p,a,I0) %Input vector of points x=0; y=0; z=linspace(0,.25,100); %These default coordinates calculates the magnetic field along the z axis through the center of the coil [Bx,By,Bz] = magnetic_field_current_loop(x,y,z,x_p,y_p,z_p,a,I0); figure(3); plot(z,Bz); xlabel('z [m]'); ylabel('Bz [T]'); title('1D magnetic field tests'); %%input mesh of points in 2D plane x=0; [y,z]=meshgrid(linspace(-.05,.05,25),linspace(0,.1,25)); %this is a 2d plane over the x=0 plane that extends away from the coils in the yz plane. [Bx,By,Bz] = magnetic_field_current_loop(x,y,z,x_p,y_p,z_p,a,I0); figure (5); surf(y,z,Bz); xlabel('y [m]'); ylabel('z [m]'); zlabel('Bz [T]'); title('2D magnetic field tests'); colorbar; %add colorbar shading flat; %Removes black lines from the meshWarning: Imaginary parts of complex X and/or Y arguments ignored E_Vpm = 0.0929 + 0.0924i deltaU_V = 92.8570 +92.3501i absU_V = 130.9617 R_Ohm = 92.8570 Li_H = 1.0887e-06 Bx = 0 By = 0 Bz = 2.4129e-06
No comments:
Post a Comment