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
 Recommended replacement: mat_idx_to_electrode

 See also MAT_IDX_TO_ELECTRODE

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 % Recommended replacement: mat_idx_to_electrode
0015 %
0016 % See also MAT_IDX_TO_ELECTRODE
0017 
0018 % (C) 2017 Andy Adler and Bartek Grychtol. License: GPL v2 or v3. $Id: mk_hollow_electrode.m 6926 2024-05-31 15:34:13Z bgrychtol $
0019 
0020 if ischar(fmdl) && strcmp(fmdl,'UNIT_TEST'); do_unit_test; return; end
0021 
0022 if nargin==1; elec_idx = 1:length(fmdl.electrode); end
0023 elec_idx = elec_idx(:)'; % Row vector to index for loop
0024 
0025 elidx=[];
0026 for i = elec_idx
0027    elidx = [elidx, fmdl.electrode(i).nodes];
0028 end
0029 
0030 ELNODES = zeros(size(fmdl.nodes,1),1);
0031 ELNODES(elidx) = 1;
0032 % find all electrode elems
0033 ELELEM = sum(ELNODES(fmdl.elems),2)==3;
0034 % remove internal electrode elements
0035 fmdl.elems = fmdl.elems(~ELELEM,:);
0036 % used nodes
0037 unodes = unique(fmdl.elems);
0038 % remove internal electrode nodes
0039 for i = elec_idx
0040    idx = ismember(fmdl.electrode(i).nodes,unodes);
0041    fmdl.electrode(i).nodes = fmdl.electrode(i).nodes(idx);
0042 end
0043 
0044 nodemap = zeros(size(fmdl.nodes,1),1);
0045 nodemap(unodes) = 1:numel(unodes);
0046 % remove unused nodes
0047 fmdl.nodes = fmdl.nodes(unodes,:);
0048 %remap elems
0049 fmdl.elems = nodemap(fmdl.elems);
0050 %remap elecs
0051 for i = 1:length(fmdl.electrode)
0052    idx = ismember(fmdl.electrode(i).nodes,unodes);
0053    nidx =  nodemap( fmdl.electrode(i).nodes(idx) );
0054    if any(i==elec_idx);
0055       npts = fmdl.nodes(nidx,:);
0056       [npts_o, nnidx] = order_loop(npts);
0057       fmdl.electrode(i).nodes = nidx(nnidx);
0058    else % if not in specified index, don't reorder loop
0059       fmdl.electrode(i).nodes = nidx;
0060    end
0061 end
0062 
0063 
0064 fmdl.boundary = find_boundary(fmdl);
0065 
0066 [~,idx] = min(sum(fmdl.nodes.^2,2));
0067 fmdl.gnd_node = idx(1);
0068 
0069 % standard field order
0070 fmdl = eidors_obj('set', fmdl);
0071 
0072 function do_unit_test
0073     fmdl = unit_test_model1;
0074     subplot(221); show_fem(fmdl); axis([-.1,.3,0.4,0.8])
0075     title 'original model - filled electrodes';
0076     fmdl1= mk_hollow_electrode(fmdl,[1,2,5]);
0077     
0078     subplot(222); show_fem(fmdl1); axis([-.1,.3,0.4,0.8])
0079     title 'original model - hollow electrode #1';
0080 
0081     fmdl2= mk_hollow_electrode(fmdl);
0082     subplot(223); show_fem(fmdl2); axis([-.1,.3,0.4,0.8])
0083     title 'original model - all hollow electrode';
0084 
0085     fmdl = unit_test_model2;
0086     fmdl4= mk_hollow_electrode(fmdl, length(fmdl.electrode));
0087     subplot(224); show_fem(fmdl4); axis([-1,1,-1,1])
0088     title 'original model - one electrode';
0089     
0090 
0091 function fmdl = unit_test_model1
0092    R=0.65; Nel = 16;
0093    shape_str = [ ...
0094    'solid incyl  = cylinder (0,0,0; 0,0,1; 0.125) -maxh=0.015;\n', ...
0095    'solid farcyl = cylinder (0,0,0; 0,0,1;0.75) -maxh=0.55;\n' ...
0096    'solid pl1    =  plane(0,0,-0.09;0,0,-1);\n' ...
0097    'solid top    =  plane(0,0,0; 0,0,1) -maxh=0.13;\n' ...
0098    'solid mainobj= pl1 and top and farcyl and not incyl;\n'];
0099    elec_pos = zeros(Nel,6);
0100    Theta = 0:360/Nel:360-360/Nel;
0101    for i = 1:Nel
0102           elec_pos(i,1) = R*sind(Theta(i));
0103           elec_pos(i,2) = R*cosd(Theta(i));
0104           elec_pos(i,6) = 1;
0105           elec_obj(i) = {'top'};
0106    end
0107    elec_shape=[0.025];
0108    fmdl = ng_mk_gen_models(shape_str, elec_pos, elec_shape, elec_obj);
0109    fmdl = mdl2d_from3d(fmdl);
0110 
0111 function fmdl = unit_test_model2
0112    extra={'ball', ...
0113          ['solid cyls= cylinder(0.2,0.2,0;0.2,0.2,1;0.2) -maxh=0.05;' ...
0114           'solid ball= cyls and orthobrick(-1,-1,0;1,1,0.5);']};
0115    fmdl= ng_mk_cyl_models(0,[6],[0.1,0,0.05],extra); 
0116    eln   = find(elem_select(fmdl, '(x-0.2).^2+(y-0.2).^2<0.2^2'));
0117    eln   = unique(fmdl.elems(eln,:));
0118    fmdl.electrode(end+1) = struct('nodes', eln(:)', 'z_contact', .01);