0001 function J= jacobian_elem2nodes( fwd_model, img)
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014 if nargin == 1
0015 img= fwd_model;
0016 elseif strcmp(getfield(warning('query','EIDORS:DeprecatedInterface'),'state'),'on')
0017 warning('EIDORS:DeprecatedInterface', ...
0018 ['Calling JACOBIAN_elem2nodes with two arguments is deprecated and will cause' ...
0019 ' an error in a future version. First argument ignored.']);
0020 end
0021 fwd_model= img.fwd_model;
0022
0023 fem_fmdl= fwd_model.jacobian_elem2nodes.fwd_model;
0024 EtoN = mapper_nodes_elems( fem_fmdl);
0025
0026 if ~isfield(img,'elem_data')
0027 img.elem_data= calc_elem_from_node_image(EtoN, img.node_data);
0028 end
0029 img.fwd_model = fem_fmdl;
0030 J= calc_jacobian(img)*EtoN';
0031
0032
0033 function e_d = calc_elem_from_node_image(EtoN, node_data);
0034 n_elems= size(EtoN,2);
0035 mu=1e-3;
0036
0037 [e_d, jnkflag] = lsqr([EtoN;mu*speye(n_elems)], ...
0038 [node_data;zeros(n_elems,1)], ...
0039 1e-8,1000);
0040
0041 return;
0042
0043
0044 [Nn,Ne]= size(EtoN);
0045 n1= ones(Nn,1); mu=1e-3;
0046 t=cputime;
0047 ed1= (EtoN'/(EtoN*EtoN'+mu^2*speye(Nn)))*n1;
0048 disp([cputime-t, std(ed1)]);
0049 t=cputime;
0050 ed2= ((EtoN'*EtoN+mu^2*speye(Ne))\EtoN')*n1;
0051 disp([cputime-t, std(ed2)]);
0052 t=cputime;
0053
0054 ed3= lsqr([EtoN;mu*speye(Ne)],[n1;zeros(Ne,1)],1e-8,1000);
0055 disp([cputime-t, std(ed3)]);