vb_spherical_grad_F

 calculate gradient of spherical harmonics
 [usage]
   [grad, ltbl, Dim] = vb_spherical_grad_F(X, Q, N, Rmax)
 [input]
      X : coordinates of internal current (NP x 3)
      Q : normal vector of internal current (NP x 3)
      N : maximum order of spherical function
   Rmax : radius
 [output]
   grad : gradient F
   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_F(X, Q, N, Rmax)
0002 % calculate gradient of spherical harmonics
0003 % [usage]
0004 %   [grad, ltbl, Dim] = vb_spherical_grad_F(X, Q, N, Rmax)
0005 % [input]
0006 %      X : coordinates of internal current (NP x 3)
0007 %      Q : normal vector of internal current (NP x 3)
0008 %      N : maximum order of spherical function
0009 %   Rmax : radius
0010 % [output]
0011 %   grad : gradient F
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 [F, Itbl, dim, Rmax] = vb_spherical_harmo_F(X,N,Rmax);
0024 
0025 NP = size(F,1);
0026 D  = (N+1)^2;
0027 
0028 W1 = zeros(NP,1);
0029 W2 = zeros(NP,1);
0030 
0031 grad = zeros(NP,D);
0032 
0033 for n = 0:N
0034   for m = 0:n
0035       
0036     n_idx = n + 1;
0037     m_idx = m + 1;
0038 
0039     if n_idx < 1 | n_idx > size(Itbl,1) | m_idx < 1 | m_idx > size(Itbl,2)
0040       continue;
0041     end
0042     
0043     % cosine zone
0044     column_num = Itbl(n_idx,m_idx,1);
0045     if column_num < 1, continue; end;
0046       
0047     % X
0048     [W1, valid1] = vb_spherical_harmo_W(n,m,F,Itbl,'+',0);
0049     [W2, valid2] = vb_spherical_harmo_W(n,m,F,Itbl,'-',0);
0050     x = 0.5 * (W1 + W2);
0051         
0052     % Y
0053     [W1, valid1] = vb_spherical_harmo_W(n,m,F,Itbl,'+',1);
0054     [W2, valid2] = vb_spherical_harmo_W(n,m,F,Itbl,'-',1);
0055     y = 0.5 * (W1 - W2);
0056 
0057     % Z
0058     [W1, valid] = vb_spherical_harmo_W(n,m,F,Itbl,'0',0);
0059     z = W1;
0060 
0061     grad(:,column_num) = x .* Q(:,1) + y .* Q(:,2) + z .* Q(:,3);
0062 
0063     % sine zone
0064     column_num = Itbl(n_idx,m_idx,2);
0065     if column_num < 1, continue; end;
0066 
0067     % X
0068     [W1, valid1] = vb_spherical_harmo_W(n,m,F,Itbl,'+',1);
0069     [W2, valid2] = vb_spherical_harmo_W(n,m,F,Itbl,'-',1);
0070     x = 0.5 * (W1 + W2);
0071 
0072     % Y
0073     [W1, valid1] = vb_spherical_harmo_W(n,m,F,Itbl,'+',0);
0074     [W2, valid2] = vb_spherical_harmo_W(n,m,F,Itbl,'-',0);
0075     y = 0.5 * (W2 - W1);
0076 
0077     % Z
0078     [W1, valid] = vb_spherical_harmo_W(n,m,F,Itbl,'0',1);
0079     z = W1;
0080 
0081     grad(:,column_num) = x .* Q(:,1) + y .* Q(:,2) + z .* Q(:,3);
0082     
0083   end
0084 end
0085 
0086 grad = grad(:,2:D);
0087 Itbl = Itbl - 1;
0088 Dim  = dim(:,2:D);
0089 return;