実装

形状関数を自然座標で評価する関数は 以下のようになる。 これは、全節点の形状関数について一度に計算される。

/* 形状関数 */
void WH_Fem__Shape3D__Hexa20N__N_arrayIn
(double whV_xi[3],
 double OUT__N_arrayIn[20])
/*
  入力引数:
    whV_xi はベクトルで、自然座標
  出力引数:
    OUT__N_arrayIn はスカラーの節点配列で、形状関数
*/
{
  double xi = whV_xi[0];
  double eta = whV_xi[1];
  double zeta = whV_xi[2];

  OUT__N_arrayIn[0] 
    = -0.125 * (1.0 - xi) * (1.0 - eta) * (1.0 - zeta) 
      * (2.0 + xi + eta + zeta);
  OUT__N_arrayIn[1] 
    = -0.125 * (1.0 + xi) * (1.0 - eta) * (1.0 - zeta) 
      * (2.0 - xi + eta + zeta);
  OUT__N_arrayIn[2] 
    = -0.125 * (1.0 + xi) * (1.0 + eta) * (1.0 - zeta) 
      * (2.0 - xi - eta + zeta);
  OUT__N_arrayIn[3] 
    = -0.125 * (1.0 - xi) * (1.0 + eta) * (1.0 - zeta) 
      * (2.0 + xi - eta + zeta);
  OUT__N_arrayIn[4] 
    = -0.125 * (1.0 - xi) * (1.0 - eta) * (1.0 + zeta) 
      * (2.0 + xi + eta - zeta);
  OUT__N_arrayIn[5] 
    = -0.125 * (1.0 + xi) * (1.0 - eta) * (1.0 + zeta) 
      * (2.0 - xi + eta - zeta);
  OUT__N_arrayIn[6] 
    = -0.125 * (1.0 + xi) * (1.0 + eta) * (1.0 + zeta) 
      * (2.0 - xi - eta - zeta);
  OUT__N_arrayIn[7] 
    = -0.125 * (1.0 - xi) * (1.0 + eta) * (1.0 + zeta) 
      * (2.0 + xi - eta - zeta);
  OUT__N_arrayIn[8]  
    = 0.25 * (1.0 - xi * xi) * (1.0 - eta) * (1.0 - zeta);
  OUT__N_arrayIn[9] 
    = 0.25 * (1.0 - eta * eta) * (1.0 - zeta) * (1.0 + xi);
  OUT__N_arrayIn[10] 
    = 0.25 * (1.0 - xi * xi) * (1.0 + eta) * (1.0 - zeta);
  OUT__N_arrayIn[11] 
    = 0.25 * (1.0 - eta * eta) * (1.0 - zeta) * (1.0 - xi);
  OUT__N_arrayIn[12] 
    = 0.25 * (1.0 - zeta * zeta) * (1.0 - xi) * (1.0 - eta);
  OUT__N_arrayIn[13] 
    = 0.25 * (1.0 - zeta * zeta) * (1.0 + xi) * (1.0 - eta);
  OUT__N_arrayIn[14] 
    = 0.25 * (1.0 - zeta * zeta) * (1.0 + xi) * (1.0 + eta);
  OUT__N_arrayIn[15] 
    = 0.25 * (1.0 - zeta * zeta) * (1.0 - xi) * (1.0 + eta);
  OUT__N_arrayIn[16] 
    = 0.25 * (1.0 - xi * xi) * (1.0 - eta) * (1.0 + zeta);
  OUT__N_arrayIn[17] 
    = 0.25 * (1.0 - eta * eta) * (1.0 + zeta) * (1.0 + xi);
  OUT__N_arrayIn[18] 
    = 0.25 * (1.0 - xi * xi) * (1.0 + eta) * (1.0 + zeta); 
  OUT__N_arrayIn[19] 
    = 0.25 * (1.0 - eta * eta) *(1.0 + zeta) * (1.0 - xi);
}