ex_fwd3d_high_order

Make common model, and make an image

0001 %Make common model, and make an image
0002 imdl = mk_common_model('n3r2',16);
0003 fmdl = imdl.fwd_model;
0004 img0 = mk_image(fmdl,1);
0005 
0006 %Forward solve using eidors solver (linear)
0007 v0 = fwd_solve(img0); v0e=v0.meas;
0008 
0009 %Solvers (including mc_fem_modify and element type)
0010 fmdl.solve = @mc_fwd_solve;
0011 fmdl.system_mat = @mc_calc_system_mat;
0012 fmdl.jacobian   = @mc_calc_jacobian;
0013 fmdl.fem_modify = @mc_fem_modify;
0014 
0015 %Forward solve using mc solvers (linear, quadratic)
0016 
0017 %Linear FEM
0018 fmdl.mc_type    = 'tet4';
0019 img1 = mk_image(fmdl,1);
0020 img1.fwd_solve.get_all_meas = 1; %Internal voltage
0021 v1 = fwd_solve(img1); v1e=v1.meas;
0022 
0023 %Quadratic FEM
0024 fmdl.mc_type    = 'tet10';
0025 img2 = mk_image(fmdl,1);
0026 img2.fwd_solve.get_all_meas = 1; %Internal voltage
0027 v2 = fwd_solve(img2); v2e=v2.meas;
0028 
0029 %Electrode voltages
0030 figure; plot([v0e,v1e,v2e,[v0e-v1e,v2e-v0e]*2]);
0031 legend('0','1','2','1-0','2-0')
0032 xlim([1,100]);
0033 
0034 
0035 %Internal voltage for linear
0036 v1all = v1.volt; 
0037 img1n = rmfield(img1,'elem_data');
0038 img1n.node_data = v1all(1:size(fmdl.nodes,1),1);
0039 figure; show_slices(img1n,[inf,inf,2.5]);
0040 
0041 %Internal voltage or quadratic
0042 v2all = v2.volt; 
0043 img2n = rmfield(img2,'elem_data');
0044 img2n.node_data = v2all(1:size(fmdl.nodes,1),1);
0045 figure; show_slices(img2n,[inf,inf,2.5]);
0046 
0047 %Internal voltage difference of linear-quadratic
0048 img12n=img1n; img12n.node_data=v1all(1:size(fmdl.nodes,1),1)-v2all(1:size(fmdl.nodes,1),1);
0049 figure; show_slices(img12n,[inf,inf,2.5]);