mdl2d_from3d

 mdl2d_from3d: Create 2D mdl from z=0 plane of 3d model

0001 function mdl2 = mdl2d_from3d(mdl3);
0002 % mdl2d_from3d: Create 2D mdl from z=0 plane of 3d model
0003 
0004 % (C) Andy Adler, Alistair Boyle 2013. Licenced under GPL v2 or v3
0005 % $Id: mdl2d_from3d.m 5908 2019-02-12 07:23:19Z aadler $
0006    % set name
0007    mdl2 = eidors_obj('fwd_model',sprintf('%s 2D',mdl3.name));
0008 
0009    % set nodes
0010    [bdy,idx] = find_boundary(mdl3.elems);
0011    vtx = mdl3.nodes;
0012    z_vtx = reshape(vtx(bdy,3), size(bdy) );
0013    lay0  = find( all(z_vtx==0,2) );
0014    bdy0  = bdy( lay0, :);
0015    
0016    vtx0  = unique(bdy0(:));
0017    mdl2.nodes = vtx(vtx0,1:2);
0018    if isempty(mdl2.nodes)
0019       error('mdl2d_from3d: Something went wrong; mdl2 has no nodes');
0020    end
0021 
0022    % set elems
0023    nmap  = zeros(size(vtx,1),1); nmap(vtx0) = 1:length(vtx0);
0024    bdy0  = reshape(nmap(bdy0), size(bdy0) ); % renumber to new scheme
0025    mdl2.elems = bdy0;
0026 
0027    % set boundary
0028    mdl2.boundary = find_boundary( mdl2.elems);
0029 
0030    % set gnd_node
0031    mdl2.gnd_node = nmap(mdl3.gnd_node);
0032    if mdl2.gnd_node == 0 % we've just killed it
0033       ctr = mean(mdl2.nodes);
0034       d = bsxfun(@minus, mdl2.nodes, ctr).^2;
0035       [jnk, mdl2.gnd_node] = min(d);
0036       mdl2.gnd_node = mdl2.gnd_node(1);
0037    end
0038 
0039    % set material indices
0040    % TODO: vectorize code
0041    if isfield(mdl3,'mat_idx');
0042    mdl2.mat_idx = {};
0043    idx0  = idx( lay0, :);
0044    for i=1:size(mdl3.mat_idx,2)
0045      mdl2.mat_idx{i} = [];
0046      ii = 1;
0047      for j=1:size(mdl3.mat_idx{i},1)
0048          idx_tmp = find( idx0==mdl3.mat_idx{i}(j) );
0049          if not(isempty(idx_tmp))
0050            mdl2.mat_idx{i}(ii,1) = idx_tmp(1,1);
0051            ii = ii + 1;
0052          end
0053      end
0054    end
0055    end %isfield
0056 
0057    
0058    % set electrode
0059    if isfield(mdl3,'electrode')
0060      mdl2.electrode = mdl3.electrode;
0061      for i=1:length(mdl2.electrode);
0062         enodes = nmap( mdl2.electrode(i).nodes );
0063         enodes(enodes==0) = []; % Remove 3D layers
0064         mdl2.electrode(i).nodes = enodes(:)';
0065      end
0066    end
0067 
0068    % copy other fields
0069    if isfield(mdl3,'stimulation'); mdl2.stimulation= mdl3.stimulation; end
0070    try   
0071        mdl2.solve      = mdl3.solve;
0072    catch
0073        mdl2.solve      = 'eidors_default';end
0074    try   
0075        mdl2.jacobian   = mdl3.jacobian;
0076    catch
0077        mdl2.jacobian   = 'eidors_default';end
0078    try   
0079        mdl2.system_mat = mdl3.system_mat;  
0080    catch
0081        mdl2.system_mat = 'eidors_default'; end;
0082    try   
0083        mdl2.normalize_measurements = mdl3.normalize_measurements;  
0084    catch
0085        mdl2.normalize_measurements = 0; end;
0086 
0087    % update cache
0088    mdl2 = eidors_obj('fwd_model',mdl2);
0089