vb_spherical_grad_G

 spherical harmonics gradient G external current
 [usage]
   [grad, ltbl, Dim] = vb_spherical_grad_G(X, Q, N, Rmax)
 [input]
      X : coordinates of external current (NP x 3)
      Q : normal vector of external current (NP x 3)
      N : maximum order of spherical function
   Rmax : radial
 [output]
   grad : gradient G
   ltbl : legendre index table
    Dim : dimension
 [note]
  
 [history]
   2006.06.30 (Sako) initial version


 Copyright (C) 2011, ATR All Rights Reserved.
 License : New BSD License(see VBMEG_LICENSE.txt)

0001 function [grad, Itbl, Dim] = vb_spherical_grad_G(X, Q, N, Rmax)
0002 % spherical harmonics gradient G external current
0003 % [usage]
0004 %   [grad, ltbl, Dim] = vb_spherical_grad_G(X, Q, N, Rmax)
0005 % [input]
0006 %      X : coordinates of external current (NP x 3)
0007 %      Q : normal vector of external current (NP x 3)
0008 %      N : maximum order of spherical function
0009 %   Rmax : radial
0010 % [output]
0011 %   grad : gradient G
0012 %   ltbl : legendre index table
0013 %    Dim : dimension
0014 % [note]
0015 %
0016 % [history]
0017 %   2006.06.30 (Sako) initial version
0018 %
0019 %
0020 % Copyright (C) 2011, ATR All Rights Reserved.
0021 % License : New BSD License(see VBMEG_LICENSE.txt)
0022 
0023 % (N+1) order is necessary
0024 [G, Itbl, dim, Rmax] = vb_spherical_harmo_G(X,N+1,Rmax);
0025 
0026 NP = size(G,1);
0027 D  = (N+1)^2;
0028 
0029 V1 = zeros(NP,1);
0030 V2 = zeros(NP,1);
0031 
0032 grad = zeros(NP,D);
0033 
0034 for n = 0:N
0035   for m = 0:n
0036       
0037     n_idx = n + 1;
0038     m_idx = m + 1;
0039       
0040     if n_idx < 1 | n_idx > size(Itbl,1) | m_idx < 1 | m_idx > size(Itbl,2)
0041       continue;
0042     end
0043 
0044     % cosine zone
0045     column_num = Itbl(n_idx,m_idx,1);
0046     if column_num < 1, continue; end;
0047 
0048     % X
0049     [V1, valid1] = vb_spherical_harmo_V(n,m,G,Itbl,'+',0);
0050     [V2, valid2] = vb_spherical_harmo_V(n,m,G,Itbl,'-',0);
0051     x = 0.5 * (V1 + V2);
0052         
0053     % Y
0054     [V1, valid1] = vb_spherical_harmo_V(n,m,G,Itbl,'+',1);
0055     [V2, valid2] = vb_spherical_harmo_V(n,m,G,Itbl,'-',1);
0056     y = 0.5 * (V1 - V2);
0057 
0058     % Z
0059     [V1, valid] = vb_spherical_harmo_V(n,m,G,Itbl,'0',0);
0060     z = V1;
0061 
0062     grad(:,column_num) = x .* Q(:,1) + y .* Q(:,2) + z .* Q(:,3);
0063 
0064     % sine zone
0065     column_num = Itbl(n_idx,m_idx,2);
0066     if column_num < 1, continue; end;
0067       
0068     % X
0069     [V1, valid1] = vb_spherical_harmo_V(n,m,G,Itbl,'+',1);
0070     [V2, valid2] = vb_spherical_harmo_V(n,m,G,Itbl,'-',1);
0071     x = 0.5 * (V1 + V2);
0072 
0073     % Y
0074     [V1, valid1] = vb_spherical_harmo_V(n,m,G,Itbl,'+',0);
0075     [V2, valid2] = vb_spherical_harmo_V(n,m,G,Itbl,'-',0);
0076     y = 0.5 * (V2 - V1);
0077 
0078     % Z
0079     [V1, valid] = vb_spherical_harmo_V(n,m,G,Itbl,'0',1);
0080     z = V1;
0081 
0082     grad(:,column_num) = x .* Q(:,1) + y .* Q(:,2) + z .* Q(:,3);
0083   
0084   end
0085 end
0086 
0087 %grad = grad(:,2:D);
0088 %Itbl = Itbl(1:N+1,1:N+1,:) - 1;
0089 %Dim  = dim(:,2:D);
0090 Itbl = Itbl(1:N+1,1:N+1,:);
0091 Dim  = dim(:,1:D);
0092 return;