#include "lu.h"

__kernel void erhs1(__global double *g_frct,
                    __global double *g_rsd,
                    __global double *g_ce,
                    int nx,
                    int ny,
                    int nz)
{
  int i, j, k, m;
  double xi, eta, zeta;

  k = get_global_id(2);
  j = get_global_id(1);
  i = get_global_id(0);
  if (k >= nz || j >= ny || i >= nx) return;

  __global double (*frct)[ISIZ2/2*2+1][ISIZ1/2*2+1][5] = 
    (__global double (*)[ISIZ2/2*2+1][ISIZ1/2*2+1][5])g_frct;
  __global double (*rsd)[ISIZ2/2*2+1][ISIZ1/2*2+1][5] = 
    (__global double (*)[ISIZ2/2*2+1][ISIZ1/2*2+1][5])g_rsd;
  __global double (*ce)[13] = (__global double (*)[13])g_ce;

  zeta = ( (double)k ) / ( nz - 1 );
  eta = ( (double)j ) / ( ny0 - 1 );
  xi = ( (double)i ) / ( nx0 - 1 );
  for (m = 0; m < 5; m++) {
    frct[k][j][i][m] = 0.0;

    rsd[k][j][i][m] =  ce[m][0]
      + (ce[m][1]
      + (ce[m][4]
      + (ce[m][7]
      +  ce[m][10] * xi) * xi) * xi) * xi
      + (ce[m][2]
      + (ce[m][5]
      + (ce[m][8]
      +  ce[m][11] * eta) * eta) * eta) * eta
      + (ce[m][3]
      + (ce[m][6]
      + (ce[m][9]
      +  ce[m][12] * zeta) * zeta) * zeta) * zeta;
  }
}

