mk_hollow_electrode

 MK_HOLLOW_ELECTRODE: remove nodes with indicated electrdoes
 fmdl_out = mk_hollow_electrode(fmdl, elec_idx)
   fmdl:     input model
   elec_idx: electrodes for which we remove internal nodes
             if not provided, remove for all electrodes
   fmdl_out: output model

 This function is useful for internal electrodes created by
   meshing software in which we don't want to calculate current
   flow within the electrode

 Limitations: currently works for 2D fems only

0001 function fmdl = mk_hollow_electrode(fmdl, elec_idx)
0002 % MK_HOLLOW_ELECTRODE: remove nodes with indicated electrdoes
0003 % fmdl_out = mk_hollow_electrode(fmdl, elec_idx)
0004 %   fmdl:     input model
0005 %   elec_idx: electrodes for which we remove internal nodes
0006 %             if not provided, remove for all electrodes
0007 %   fmdl_out: output model
0008 %
0009 % This function is useful for internal electrodes created by
0010 %   meshing software in which we don't want to calculate current
0011 %   flow within the electrode
0012 %
0013 % Limitations: currently works for 2D fems only
0014 
0015 % (C) 2017 Andy Adler and Bartek Grychtol. License: GPL v2 or v3. $Id: mk_hollow_electrode.m 5424 2017-04-25 17:45:19Z aadler $
0016 
0017 if ischar(fmdl) && strcmp(fmdl,'UNIT_TEST'); do_unit_test; return; end
0018 
0019 if nargin==1; elec_idx = 1:length(fmdl.electrode); end
0020 elec_idx = elec_idx(:)'; % Row vector to index for loop
0021 
0022 elidx=[];
0023 for i = elec_idx
0024    elidx = [elidx, fmdl.electrode(i).nodes];
0025 end
0026 
0027 ELNODES = zeros(size(fmdl.nodes,1),1);
0028 ELNODES(elidx) = 1;
0029 % find all electrode elems
0030 ELELEM = sum(ELNODES(fmdl.elems),2)==3;
0031 % remove internal electrode elements
0032 fmdl.elems = fmdl.elems(~ELELEM,:);
0033 % used nodes
0034 unodes = unique(fmdl.elems);
0035 % remove internal electrode nodes
0036 for i = elec_idx
0037    idx = ismember(fmdl.electrode(i).nodes,unodes);
0038    fmdl.electrode(i).nodes = fmdl.electrode(i).nodes(idx);
0039 end
0040 
0041 nodemap = zeros(size(fmdl.nodes,1),1);
0042 nodemap(unodes) = 1:numel(unodes);
0043 % remove unused nodes
0044 fmdl.nodes = fmdl.nodes(unodes,:);
0045 %remap elems
0046 fmdl.elems = nodemap(fmdl.elems);
0047 %remap elecs
0048 for i = 1:length(fmdl.electrode)
0049    idx = ismember(fmdl.electrode(i).nodes,unodes);
0050    nidx =  nodemap( fmdl.electrode(i).nodes(idx) );
0051    if any(i==elec_idx);
0052       npts = fmdl.nodes(nidx,:);
0053       [npts_o, nnidx] = order_loop(npts);
0054       fmdl.electrode(i).nodes = nidx(nnidx);
0055    else % if not in specified index, don't reorder loop
0056       fmdl.electrode(i).nodes = nidx;
0057    end
0058 end
0059 
0060 
0061 fmdl.boundary = find_boundary(fmdl);
0062 
0063 [~,idx] = min(sum(fmdl.nodes.^2,2));
0064 fmdl.gnd_node = idx(1);
0065 
0066 
0067 function do_unit_test
0068     fmdl = unit_test_model1;
0069     subplot(221); show_fem(fmdl); axis([-.1,.3,0.4,0.8])
0070     title 'original model - filled electrodes';
0071     fmdl1= mk_hollow_electrode(fmdl,[1,2,5]);
0072     
0073     subplot(222); show_fem(fmdl1); axis([-.1,.3,0.4,0.8])
0074     title 'original model - hollow electrode #1';
0075 
0076     fmdl2= mk_hollow_electrode(fmdl);
0077     subplot(223); show_fem(fmdl2); axis([-.1,.3,0.4,0.8])
0078     title 'original model - all hollow electrode';
0079 
0080     fmdl = unit_test_model2;
0081     fmdl4= mk_hollow_electrode(fmdl, length(fmdl.electrode));
0082     subplot(224); show_fem(fmdl4); axis([-1,1,-1,1])
0083     title 'original model - one electrode';
0084     
0085 
0086 function fmdl = unit_test_model1
0087    R=0.65; Nel = 16;
0088    shape_str = [ ...
0089    'solid incyl  = cylinder (0,0,0; 0,0,1; 0.125) -maxh=0.015;\n', ...
0090    'solid farcyl = cylinder (0,0,0; 0,0,1;0.75) -maxh=0.55;\n' ...
0091    'solid pl1    =  plane(0,0,-0.09;0,0,-1);\n' ...
0092    'solid top    =  plane(0,0,0; 0,0,1) -maxh=0.13;\n' ...
0093    'solid mainobj= pl1 and top and farcyl and not incyl;\n'];
0094    elec_pos = zeros(Nel,6);
0095    Theta = 0:360/Nel:360-360/Nel;
0096    for i = 1:Nel
0097           elec_pos(i,1) = R*sind(Theta(i));
0098           elec_pos(i,2) = R*cosd(Theta(i));
0099           elec_pos(i,6) = 1;
0100           elec_obj(i) = {'top'};
0101    end
0102    elec_shape=[0.025];
0103    fmdl = ng_mk_gen_models(shape_str, elec_pos, elec_shape, elec_obj);
0104    fmdl = mdl2d_from3d(fmdl);
0105 
0106 function fmdl = unit_test_model2
0107    extra={'ball', ...
0108          ['solid cyls= cylinder(0.2,0.2,0;0.2,0.2,1;0.2) -maxh=0.05;' ...
0109           'solid ball= cyls and orthobrick(-1,-1,0;1,1,0.5);']};
0110    fmdl= ng_mk_cyl_models(0,[6],[0.1,0,0.05],extra); 
0111    eln   = find(elem_select(fmdl, '(x-0.2).^2+(y-0.2).^2<0.2^2'));
0112    eln   = unique(fmdl.elems(eln,:));
0113    fmdl.electrode(end+1) = struct('nodes', eln(:)', 'z_contact', .01);