0001
0002
0003
0004
0005 imdl= mk_common_model('a3cr',16);
0006 imdl= mk_common_model('n3r2',16);
0007 imdl.fwd_model.nodes = imdl.fwd_model.nodes*.25;
0008 img= calc_jacobian_bkgnd(imdl);
0009
0010 img.fwd_model = mdl_normalize(img.fwd_model,0);
0011
0012 img.fwd_model.jacobian= @np_calc_jacobian;
0013 img.fwd_model.system_mat= @np_calc_system_mat;
0014 img.fwd_model.solve= @np_fwd_solve;
0015 J_np= calc_jacobian( img );
0016
0017 img.fwd_model.jacobian= @jacobian_perturb;
0018 img.fwd_model.system_mat= @np_calc_system_mat;
0019 img.fwd_model.solve= @np_fwd_solve;
0020 J_np_p= calc_jacobian( img );
0021
0022 img.fwd_model.jacobian= @jacobian_adjoint;
0023 img.fwd_model.system_mat= @system_mat_1st_order;
0024 img.fwd_model.solve= @fwd_solve_1st_order;
0025 J_aa= calc_jacobian( img );
0026
0027 img.fwd_model.jacobian= @jacobian_perturb;
0028 img.fwd_model.system_mat= @system_mat_1st_order;
0029 img.fwd_model.solve= @fwd_solve_1st_order;
0030 J_aa_p= calc_jacobian( img );
0031
0032 norm(J_aa - J_aa_p,'fro')/norm(J_aa,'fro')
0033 norm(J_np - J_np_p,'fro')/norm(J_np,'fro')
0034 norm(J_np - J_aa_p,'fro')/norm(J_np,'fro')
0035 norm(J_np - J_aa ,'fro')/norm(J_np,'fro')
0036
0037
0038 imdl= mk_common_model('n3r2',16);
0039 img= calc_jacobian_bkgnd(imdl);
0040 for i=1:16; imdl.fwd_model.electrode(i).z_contact= .01;end
0041
0042
0043 img.fwd_model.jacobian= @np_calc_jacobian;
0044 img.fwd_model.system_mat= @np_calc_system_mat;
0045 img.fwd_model.solve= @np_fwd_solve;
0046 J_np= calc_jacobian( img );
0047
0048
0049 vh= fwd_solve(img);
0050
0051 for el= 1:50:size(img.fwd_model.elems,1);
0052 fprintf('\nelem#%03d: ',el);
0053 for delta= [1e-4,1e-6,1e-8]
0054 img_i = img;
0055 img_i.elem_data(el)= img_i.elem_data(el) + delta;
0056 vi= fwd_solve(img_i);
0057 J_p = (vi.meas - vh.meas)/delta;
0058 fprintf(' %8.6f', norm(J_p - J_np(:,el))/norm(J_np(:,el)));
0059 end
0060 end
0061