m_3d_fields.sci

用来实现三维阻抗及光学断层成像重建的matlab程序

function [v_f]=m_3d_fields(vtx,el_no,m_ind,E,tol,gnd_ind,v_f)
[nargout,nargin] = argn(0)
//function [v_f] = m_3d_fields(vtx,el_no,m_ind,E,tol,gnd_ind,v_f);
// 
//This function calculates the measurement fields using preconditioned conjugate gradients.
// 
// 
// 
//vtx     = The vertices
//el_no   = The total number of electrodes in the system
//m_ind   = The measurements matrix (indices of electrode pairs)
//E       = The full rank system matrix
//tol     = The tolerance in the forward solution 
//gnd_ind = The ground index
//v_f     = The measurements fields
 
 
[vr,vc] = size(vtx);
 
Is_supl = zeros(vr,size(m_ind,1));
//no of electrodes x no of measurements (now currents)!
 
MC = [];
 
for i = 1:size(m_ind,1)
   
  m_n = zeros(el_no,1);
   
  m_n(m_ind(i,1)) = 1
  m_n(m_ind(i,2)) = -1
   
  MC = [MC,m_n];
   
end
 
I = [Is_supl;MC];
I(gnd_ind,:) = 0;
 
if nargin<7 then
  v_f = zeros(size(E,1),size(I,2));
end
 
 
if isreal(E,0)==%t then
   
  //Preconditioner
   
  //!! Unknown function cholinc ,the original calling sequence is used
  M = cholinc(E,tol*10);
   
  for y = 1:size(MC,2)
    //Set this line to suit your approximation needs. ***************
    //for more details use help pcg on Matlab's command window.
     
    //! mtlb(relres) can be replaced by relres() or relres whether relres is an m-file or not
     
    //! mtlb(iter) can be replaced by iter() or iter whether iter is an m-file or not
     
    //! mtlb(resvec) can be replaced by resvec() or resvec whether resvec is an m-file or not
     
    //!! Unknown function pcg ,the original calling sequence is used
    [v_f(:,y),ones(8,1).*.[1 0 0;1 1 1;0 0 1;0 0 0],mtlb(relres),mtlb(iter),mtlb(resvec)]==pcg(E,I(:,y),tol*norm(I(:,y)),1000,M',M,v_f(:,y));
  end
end
//is real
 
 
if isreal(E,0)==%f then
   
   
  //!! Unknown function luinc ,the original calling sequence is used
  [L,U] = luinc(E,tol*10);
   
  for y = 1:size(MC,2)
     
     
    //! mtlb(relres) can be replaced by relres() or relres whether relres is an m-file or not
     
    //! mtlb(iter) can be replaced by iter() or iter whether iter is an m-file or not
     
    //! mtlb(resvec) can be replaced by resvec() or resvec whether resvec is an m-file or not
     
    //!! Unknown function bicgstab ,the original calling sequence is used
    [v_f(:,y),ones(8,1).*.[1 0 0;1 1 1;0 0 1;0 0 0],mtlb(relres),mtlb(iter),mtlb(resvec)]==bicgstab(E,I(:,y),tol*norm(I(:,y)),1000,L,U);
     
  end
end
//is complex
 
 
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
// This is part of the EIDORS suite.
// Copyright (c) N. Polydorides 2001
// Copying permitted under terms of GNU GPL
// See enclosed file gpl.html for details.
// EIDORS 3D version 1.0
// MATLAB version 5.3 R11
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%