classdef vertices2meq_test < meq_test
  % Tests of conversion of vertices-based passive structure descriptions to MEQ
  %
  % [+MEQ MatlabEQuilibrium Toolbox+] Swiss Plasma Center EPFL Lausanne 2022. All rights reserved.
  
  properties
    doplot = false;
  end
  
  properties(TestParameter)
    A = {1,2,2.7,3,10};         % aspect ratio
    orient = {0,pi/2,pi/4,pi/6} % orientation angle
    Resistance_method = {'Resistance','Resistivity'}; % how to specify resistance
    reshuffle = {true,false}; % reshuffle vertex order
    shape = {'Rectangle','Parallellogram','Trapezoid'};
  end
  
  methods (TestClassSetup)
    function matlab_version_check(testCase)
      testCase.assumeFalse(verLessThan('matlab','9.3'),...
        'can not run vertices2meq for this version, skipping test');
    end
  end
  
  methods ( Test,TestTags = {'Unit'} )
    function test_vertices2meq(testCase,A,orient,shape,Resistance_method,reshuffle)
      
      RZ = [...
        [0 0 1 1];
        [0 1 1 0]/A];
      RZ = RZ-mean(RZ,2);
      
      switch shape
        case 'Parallellogram'
          RZ(2,1:2) = RZ(2,1:2) + 1/A*sin(-orient*2); % shift left side
        case 'Rectangle' % rotate
          Rot = [cos(orient),sin(orient);-sin(orient),cos(orient)]; % rotation matrix
          RZ = Rot*RZ;
        case 'Trapezoid'
          RZ(2,1) = RZ(2,1) + 1/A*sin(-orient); % shift left side
        otherwise
          error('undefined %s',shape)
      end
      
      % Reshuffle vertex order
      if reshuffle
        rng(1); % set to fixed
        RZ = RZ(:,randperm(4));
        rng('shuffle'); % reset to random
      end
      
      V.Rvertex = RZ(1,:) - min(RZ(1,:)) + 1; % position offset
      V.Zvertex = RZ(2,:);
      
      switch Resistance_method
        case 'Resistance'
          V.Res = 1 ; V.Resy = []; 
        case 'Resistivity'
          V.Res = []; V.Resy = 1e-6; % resistivity in Ohm*m
        otherwise, error('untreated case %s',Resistance_method);
      end
      
      if testCase.doplot; clf; end
      [rv,zv,wv,hv,Tivs,Rs,Rv,dims,dimv] = vertices2meq(V,testCase.doplot); %#ok<ASGLU>
      if testCase.doplot; drawnow; end
      
      poly = polyshape(V.Rvertex,V.Zvertex);
      
      % check rv,zv are inside the polygon
      testCase.verifyTrue(all(poly.isinterior(rv,zv)),'Some filaments are not inside the polygon')
      % check area is maintained
      testCase.verifyEqual(poly.area,sum(wv.*hv),'RelTol',10*eps,'Total area not maintained')
      % check size of Tivs
      testCase.verifySize(Tivs,[numel(dims),numel(dimv)],'unexpected Tivs size');
      % Check that these give reasonable mutuals
      mut = greenem('mut',rv,zv,hv,zv);
      testCase.assertGreaterThan(mut,0)
      
      if ~isempty(V.Res)
        testCase.verifyEqual(Rs,V.Res,'RelTol',10*eps,...
          'Resistance not equal after splitting into elements');
      end
    end
  end
  
end
