From 0590752c0b042c687fa82da619a204f13b51a13f Mon Sep 17 00:00:00 2001 From: Francesco Carpanese <francesco.carpanese@epfl.ch> Date: Sat, 7 Dec 2019 19:48:20 +0100 Subject: [PATCH] Clean up versio of the function to generate the discretized filemants --- IDS/gen_filament.m | 182 +++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 182 insertions(+) create mode 100644 IDS/gen_filament.m diff --git a/IDS/gen_filament.m b/IDS/gen_filament.m new file mode 100644 index 00000000..bee54bfb --- /dev/null +++ b/IDS/gen_filament.m @@ -0,0 +1,182 @@ +function [discretizationr, discretizationz, Twc] = eval_trap(vvdata, num, type, varargin ) +% Evaluate the trapezoid from vv data structure, and computes the +% discretization +% type = 2 -> rectangle description +% type = 3 -> oblique description +% type = 4 -> contour of the coils are specified varargin = [r, z ] +% contours +% vv data contains R, Z, dR, dZ for type = 2. +% vv data contains R, Z, dR, dZ, AC, AC2 for type = 3 where AC, AC2 are the +% angles of the parallelogram. Refer to GA document for specifications. + +if mod(num,1)~=0, warning('The number of turns must be an integer. Ceil operation on it'); end +num = ceil(num); + +R = vvdata(1, :); +Z = vvdata(2, :); +dR = vvdata(3, :); +dZ = vvdata(4, :); + +switch type + case 1 % outline option + disp('Not yet implemented. It will contain outline option IMAS') + return + case 2 % rectangle + AC = 0*ones(size(R)); + AC2 = 90*ones(size(R)); + + case 3 % oblique + AC = vvdata(5, :); + AC2 = vvdata(6, :); + + case 4 % contour of the polygon + + otherwise + disp('case not defined') + return +end + +% Initialize the extrema of the trapezoid +pointr = zeros(numel(AC),4); +pointz = zeros(numel(AC),4); + +discretizationr = []; +discretizationz = []; + +Twc = []; % Grouping matrix from windings to active coil +for ii = 1:numel(AC) % loop all the coils + if AC2(ii) ~= 0 && AC(ii) ==0 + DH = dZ(ii)/tand(AC2(ii)); + pointr(ii,1) = R(ii) -DH/2 -dR(ii)/2; + pointr(ii,2) = R(ii) +DH/2 -dR(ii)/2; + pointr(ii,3) = R(ii) +DH/2 +dR(ii)/2; + pointr(ii,4) = R(ii) -DH/2 +dR(ii)/2; + + pointz(ii,1) = Z(ii) -dZ(ii)/2; + pointz(ii,2) = Z(ii) +dZ(ii)/2; + pointz(ii,3) = Z(ii) +dZ(ii)/2; + pointz(ii,4) = Z(ii) -dZ(ii)/2; + + elseif AC(ii) ~= 0 && AC2(ii) ==0 + DH = dR(ii)*tand(AC(ii)); + + pointr(ii,1) = R(ii) - dR(ii)/2; + pointr(ii,2) = R(ii) - dR(ii)/2; + pointr(ii,3) = R(ii) + dR(ii)/2; + pointr(ii,4) = R(ii) + dR(ii)/2; + + pointz(ii,1) = Z(ii) - DH/2- dZ(ii)/2; + pointz(ii,2) = Z(ii) - DH/2 +dZ(ii)/2; + pointz(ii,3) = Z(ii) + DH/2 +dZ(ii)/2; + pointz(ii,4) = Z(ii) + DH/2 -dZ(ii)/2; + + elseif AC2(ii) == 0 && AC(ii) ==0 + pointr(ii,1) = R(ii) -dR(ii)/2; + pointr(ii,2) = R(ii) -dR(ii)/2; + pointr(ii,3) = R(ii) +dR(ii)/2; + pointr(ii,4) = R(ii) +dR(ii)/2; + + pointz(ii,1) = Z(ii) -dZ(ii)/2; + pointz(ii,2) = Z(ii) +dZ(ii)/2; + pointz(ii,3) = Z(ii) +dZ(ii)/2; + pointz(ii,4) = Z(ii) -dZ(ii)/2; + end + + + AD = [pointr(ii,4)-pointr(ii,1), pointz(ii,4)-pointz(ii,1)]; + AB = [pointr(ii,2)-pointr(ii,1), pointz(ii,2)-pointz(ii,1)]; + ex = [1,0]; + ey = [0,1]; + + + S = ([ex*AD', ex*AB'; ey*AD', ey*AB']); + + %fixed_area = 0.01; %Decide the area represented by each single coils + %area_of_element = abs(AB(1)*AD(2)-AB(2)*AD(1)); %Evaluate the closer area to the element + + %num = round(area_of_element/fixed_area); %Find the number of integer + %filaments closer to a fixed desired area for each element. + + + if numel(num) == 1 %selfsimilar equally spaced + + pointx = zeros(num*num, 1); + pointy = zeros(num*num, 1); + + % Mode with selfsimilar shape + for kk=1:num + for jj=1:num + pointx((kk-1)*num+jj) = (2*kk-1)/(2*num); %uniform grid in the deformed space defined by the parallelogram + pointy((kk-1)*num+jj) = (2*jj-1)/(2*num); + + tmp = S*[pointx((kk-1)*(num)+jj), pointy((kk-1)*(num)+jj)]'; % evaluate component in the r,z system of reference + discretizationr(end+1) = pointr(ii,1)+ tmp(1) ; % translate to the original point + discretizationz(end+1) = pointz(ii,1) + tmp(2); + end + end + + % Creating the grouping matrix + tmp = zeros(num*num, numel(AC)); + tmp(:,ii) = 1; + + Twc = [Twc ; tmp]; + + else %Depending on the turns and same similar aspect ration + ratio = sqrt(dot(AB,AB)/dot(AD,AD)); + + nx = ceil(sqrt(num(ii)/ratio)); + ny = ceil(nx*ratio); + + if ny>=2 + pointy = linspace(1/(ny-1)/2, 1- 1/(ny-1)/2, ny); + else + pointy = 0.5; + end + + if nx>=2 + pointx = linspace(1/(nx-1)/2, 1- 1/(nx-1)/2, nx); + else + pointx = 0.5; + end + + + [ pointxx, pointyy ] = meshgrid(pointx, pointy); + + for kk=1:ny + for jj=1:nx + + tmp = S*[pointxx(kk,jj); pointyy(kk,jj)]; % evaluate component in the r,z system of reference + discretizationr(end+1) = pointr(ii,1)+ tmp(1) ; % translate to the original point + discretizationz(end+1) = pointz(ii,1) + tmp(2); + end + end + + % Creating the grouping matrix + tmp = zeros(nx*ny, numel(num)); + tmp(:,ii) = 1; + Twc = [Twc ; tmp]; + + end + + + +end + + +%% Plot the coils discretization from the original data +%{ +figure +hold on +axis equal +for ii = 1:numel(AC) % loop all the coils + patch(pointr(ii,:),pointz(ii,:),'red'); +end +scatter(discretizationr, discretizationz, '.'); +%} +discretizationr = reshape(discretizationr, numel(discretizationr),1); +discretizationz = reshape(discretizationz, numel(discretizationz),1); + + +end + + -- GitLab