function Nprime = femlib_NGrad(xi, eta, zeta)
% function Nprime = femlib_NGrad(xi, eta, zeta)
%
% Calculates the derivative matrix of shape functions with respect to each
% natrual co-ordinate.
% 
% The trial and test functions are given by
%
% N1 = 0.5*xi*(1+zeta);
% N2 = 0.5*eta*(1+zeta);
% N3 = 0.5*(1-xi-eta)*(1+zeta);
% N4 = 0.5*xi*(1-zeta);
% N5 = 0.5*eta*(1-zeta);
% N6 = 0.5*(1-xi-eta)*(1-zeta);
%          _                                   _
% Nprime = |dN1/d_xi   dN2/d_xi ....  dN6/d_xi  |
%          |dN1/d_eta  dN2/d_eta      dN6/d_eta | 
%          |dN1/d_zeta dN2/d_zeta     dN6/d_zeta| 
%           -                                  -
%
% Dr. A. I. Hanna (2006).
Nprime = [(1+zeta)    0      -(1+zeta)  (1-zeta)    0      -(1-zeta);
             0     (1+zeta)  -(1+zeta)      0    (1-zeta)  -(1-zeta);
            xi       eta     (1-xi-eta)    -xi     -eta   -(1-xi-eta)];
Nprime = 0.5*Nprime;
return;