nubspline_create.h File Reference

#include "nubspline_structs.h"

Go to the source code of this file.

Functions
NUgrid *	create_center_grid (double start, double end, double ratio, int num_points)
NUgrid *	create_general_grid (double *points, int num_points)
NUBspline_1d_c *	create_NUBspline_1d_c (NUgrid *x_grid, BCtype_c xBC, complex_float *data)
NUBspline_1d_d *	create_NUBspline_1d_d (NUgrid *x_grid, BCtype_d xBC, double *data)
NUBspline_1d_s *	create_NUBspline_1d_s (NUgrid *x_grid, BCtype_s xBC, float *data)
NUBspline_1d_z *	create_NUBspline_1d_z (NUgrid *x_grid, BCtype_z xBC, complex_double *data)
NUBspline_2d_c *	create_NUBspline_2d_c (NUgrid *x_grid, NUgrid *y_grid, BCtype_c xBC, BCtype_c yBC, complex_float *data)
NUBspline_2d_d *	create_NUBspline_2d_d (NUgrid *x_grid, NUgrid *y_grid, BCtype_d xBC, BCtype_d yBC, double *data)
NUBspline_2d_s *	create_NUBspline_2d_s (NUgrid *x_grid, NUgrid *y_grid, BCtype_s xBC, BCtype_s yBC, float *data)
NUBspline_2d_z *	create_NUBspline_2d_z (NUgrid *x_grid, NUgrid *restrict y_grid, BCtype_z xBC, BCtype_z yBC, complex_double *data)
NUBspline_3d_c *	create_NUBspline_3d_c (NUgrid *x_grid, NUgrid *y_grid, NUgrid *z_grid, BCtype_c xBC, BCtype_c yBC, BCtype_c zBC, complex_float *data)
NUBspline_3d_d *	create_NUBspline_3d_d (NUgrid *x_grid, NUgrid *y_grid, NUgrid *z_grid, BCtype_d xBC, BCtype_d yBC, BCtype_d zBC, double *data)
NUBspline_3d_s *	create_NUBspline_3d_s (NUgrid *x_grid, NUgrid *y_grid, NUgrid *z_grid, BCtype_s xBC, BCtype_s yBC, BCtype_s zBC, float *data)
NUBspline_3d_z *	create_NUBspline_3d_z (NUgrid *x_grid, NUgrid *y_grid, NUgrid *z_grid, BCtype_z xBC, BCtype_z yBC, BCtype_z zBC, complex_double *data)

Function Documentation

create_center_grid()

NUgrid * create_center_grid	(	double	start,
		double	end,
		double	ratio,
		int	num_points )

Definition at line 83 of file nugrid.cpp.

{
  center_grid *grid = new center_grid;
  if (grid != NULL) {
    assert (ratio > 1.0);
    grid->start       = start;
    grid->end         = end;
    grid->center      = 0.5*(start + end);
    grid->num_points  = num_points;
    grid->half_points = num_points/2;
    grid->odd = ((num_points % 2) == 1);
    grid->b = log(ratio) / (double)(grid->half_points-1);
    grid->bInv = 1.0/grid->b;
    grid->points = new double[num_points];
    if (grid->odd) {
      grid->even_half = 0.0;  
      grid->odd_one   = 1;
      grid->a = 0.5*(end-start)/expm1(grid->b*grid->half_points);
      grid->aInv = 1.0/grid->a;
      for (int i=-grid->half_points; i<=grid->half_points; i++) {
    double sign;
    if (i<0) 
      sign = -1.0;
    else
      sign =  1.0;
    grid->points[i+grid->half_points] = 
      sign * grid->a*expm1(grid->b*abs(i))+grid->center;
      }
    }
    else {
      grid->even_half = 0.5;  
      grid->odd_one   = 0;
      grid->a = 
    0.5*(end-start)/expm1(grid->b*(-0.5+grid->half_points));
      grid->aInv = 1.0/grid->a;
      for (int i=-grid->half_points; i<grid->half_points; i++) {
    double sign;
    if (i<0) sign = -1.0; 
    else     sign =  1.0;
    grid->points[i+grid->half_points] = 
      sign * grid->a*expm1(grid->b*fabs(0.5+i)) + grid->center;
      }
    }
    grid->reverse_map = center_grid_reverse_map;
    grid->code = CENTER;
  }
  return (NUgrid*) grid;
}

References center_grid::a, center_grid::aInv, center_grid::b, center_grid::bInv, CENTER, center_grid::center, center_grid_reverse_map(), center_grid::code, center_grid::end, center_grid::even_half, center_grid::half_points, center_grid::num_points, center_grid::odd, center_grid::odd_one, center_grid::points, center_grid::reverse_map, sign(), and center_grid::start.

create_general_grid()

NUgrid * create_general_grid	(	double *	points,
		int	num_points )

Definition at line 155 of file nugrid.cpp.

{
  NUgrid* grid = new NUgrid;
  if (grid != NULL) {
    grid->reverse_map = general_grid_reverse_map;
    grid->code = GENERAL;
    grid->points = new double[num_points];
    grid->start = points[0];
    grid->end   = points[num_points-1];
    grid->num_points = num_points;
    for (int i=0; i<num_points; i++) 
      grid->points[i] = points[i];
    grid->code = GENERAL;
  }
  return grid;
}

References NUgrid::code, NUgrid::end, GENERAL, general_grid_reverse_map(), NUgrid::num_points, NUgrid::points, NUgrid::reverse_map, and NUgrid::start.

create_NUBspline_1d_c()

NUBspline_1d_c * create_NUBspline_1d_c	(	NUgrid *	x_grid,
		BCtype_c	xBC,
		complex_float *	data )

Definition at line 763 of file nubspline_create.cpp.

{
  // First, create the spline structure
  NUBspline_1d_c* spline = static_cast<NUBspline_1d_c*>(malloc(sizeof(NUBspline_1d_c)));
  if (spline == NULL)
    return spline;
  spline->sp_code = NU1D;
  spline->t_code  = SINGLE_COMPLEX;
 
  // Next, create the basis
  spline->x_basis = create_NUBasis (x_grid, xBC.lCode==PERIODIC);
  // M is the number of data points
//  int M;
//  if (xBC.lCode == PERIODIC) M = x_grid->num_points - 1;
//  else                       M = x_grid->num_points;
  int N = x_grid->num_points + 2;
 
  // Allocate coefficients and solve  
  spline->coefs = (complex_float*)malloc (sizeof(complex_float)*N);
  find_NUBcoefs_1d_c (spline->x_basis, xBC, data, 1, spline->coefs, 1);
    
  return spline;
}

References NUBspline_1d_c::coefs, create_NUBasis(), find_NUBcoefs_1d_c(), BCtype_c::lCode, NU1D, NUgrid::num_points, PERIODIC, SINGLE_COMPLEX, NUBspline_1d_c::sp_code, NUBspline_1d_c::t_code, and NUBspline_1d_c::x_basis.

create_NUBspline_1d_d()

NUBspline_1d_d * create_NUBspline_1d_d	(	NUgrid *	x_grid,
		BCtype_d	xBC,
		double *	data )

Definition at line 591 of file nubspline_create.cpp.

{
  // First, create the spline structure
  NUBspline_1d_d* spline = static_cast<NUBspline_1d_d*>(malloc(sizeof(NUBspline_1d_d)));
  if (spline == NULL)
    return spline;
  spline->sp_code = NU1D;
  spline->t_code  = DOUBLE_REAL;
 
  // Next, create the basis
  spline->x_basis = create_NUBasis (x_grid, xBC.lCode==PERIODIC);
  // M is the number of data points (set but unused)
//  int M;
//  if (xBC.lCode == PERIODIC) M = x_grid->num_points - 1;
//  else                       M = x_grid->num_points;
  int N = x_grid->num_points + 2;
 
  // Allocate coefficients and solve
  spline->coefs = (double*)malloc (sizeof(double)*N);
  find_NUBcoefs_1d_d (spline->x_basis, xBC, data, 1, spline->coefs, 1);
    
  return spline;
}

References NUBspline_1d_d::coefs, create_NUBasis(), DOUBLE_REAL, find_NUBcoefs_1d_d(), BCtype_d::lCode, NU1D, NUgrid::num_points, PERIODIC, NUBspline_1d_d::sp_code, NUBspline_1d_d::t_code, and NUBspline_1d_d::x_basis.

create_NUBspline_1d_s()

NUBspline_1d_s * create_NUBspline_1d_s	(	NUgrid *	x_grid,
		BCtype_s	xBC,
		float *	data )

Definition at line 247 of file nubspline_create.cpp.

{
  // First, create the spline structure
  NUBspline_1d_s* spline = static_cast<NUBspline_1d_s*>(malloc(sizeof(NUBspline_1d_s)));
  if (spline == NULL)
    return spline;
  spline->sp_code = NU1D;
  spline->t_code  = SINGLE_REAL;
 
  // Next, create the basis
  spline->x_basis = create_NUBasis (x_grid, xBC.lCode==PERIODIC);
  // M is the number of data points (but is unused)
//  int M;
//  if (xBC.lCode == PERIODIC) M = x_grid->num_points - 1;
//  else                       M = x_grid->num_points;
  int N = x_grid->num_points + 2;
 
  // Allocate coefficients and solve  
  spline->coefs = (float*)malloc (sizeof(float)*N);
  find_NUBcoefs_1d_s (spline->x_basis, xBC, data, 1, spline->coefs, 1);
    
  return spline;
}

References NUBspline_1d_s::coefs, create_NUBasis(), find_NUBcoefs_1d_s(), BCtype_s::lCode, NU1D, NUgrid::num_points, PERIODIC, SINGLE_REAL, NUBspline_1d_s::sp_code, NUBspline_1d_s::t_code, and NUBspline_1d_s::x_basis.

create_NUBspline_1d_z()

NUBspline_1d_z * create_NUBspline_1d_z	(	NUgrid *	x_grid,
		BCtype_z	xBC,
		complex_double *	data )

Definition at line 932 of file nubspline_create.cpp.

{
  // First, create the spline structure
  NUBspline_1d_z* spline = static_cast<NUBspline_1d_z*>(malloc(sizeof(NUBspline_1d_z)));
  if (spline == NULL)
    return spline;
  spline->sp_code = NU1D;
  spline->t_code  = DOUBLE_COMPLEX;
 
  // Next, create the basis
  spline->x_basis = create_NUBasis (x_grid, xBC.lCode==PERIODIC);
  // M is the number of data points
//  int M;
//  if (xBC.lCode == PERIODIC) M = x_grid->num_points - 1;
//  else                       M = x_grid->num_points;
  int N = x_grid->num_points + 2;
 
  // Allocate coefficients and solve  
  spline->coefs = (complex_double*)malloc (sizeof(complex_double)*N);
  find_NUBcoefs_1d_z (spline->x_basis, xBC, data, 1, spline->coefs, 1);
    
  return spline;
}

References NUBspline_1d_z::coefs, create_NUBasis(), DOUBLE_COMPLEX, find_NUBcoefs_1d_z(), BCtype_z::lCode, NU1D, NUgrid::num_points, PERIODIC, NUBspline_1d_z::sp_code, NUBspline_1d_z::t_code, and NUBspline_1d_z::x_basis.

create_NUBspline_2d_c()

NUBspline_2d_c * create_NUBspline_2d_c	(	NUgrid *	x_grid,
		NUgrid *	y_grid,
		BCtype_c	xBC,
		BCtype_c	yBC,
		complex_float *	data )

Definition at line 788 of file nubspline_create.cpp.

{
  // First, create the spline structure
  NUBspline_2d_c* spline = static_cast<NUBspline_2d_c*>(malloc(sizeof(NUBspline_2d_c)));
  if (spline == NULL)
    return spline;
  spline->sp_code = NU2D;
  spline->t_code  = SINGLE_COMPLEX;
 
  // Next, create the bases
  spline->x_basis = create_NUBasis (x_grid, xBC.lCode==PERIODIC);
  spline->y_basis = create_NUBasis (y_grid, yBC.lCode==PERIODIC);
//  int Mx,
  int My, Nx, Ny;
//  if (xBC.lCode == PERIODIC) Mx = x_grid->num_points - 1;
//  else                       Mx = x_grid->num_points;
  if (yBC.lCode == PERIODIC) My = y_grid->num_points - 1;
  else                       My = y_grid->num_points;
 
  Nx = x_grid->num_points + 2;
  Ny = y_grid->num_points + 2;
    
  spline->x_stride = Ny;
#ifndef HAVE_SSE2
  spline->coefs = (complex_float*)malloc (sizeof(complex_float)*Nx*Ny);
#else
  posix_memalign ((void**)&spline->coefs, 16, sizeof(complex_float)*Nx*Ny);
#endif
 
  // First, solve in the X-direction 
  for (int iy=0; iy<My; iy++) {
    int doffset = iy;
    int coffset = iy;
    find_NUBcoefs_1d_c (spline->x_basis, xBC, data+doffset, My,
            spline->coefs+coffset, Ny);
  }
  
  // Now, solve in the Y-direction
  for (int ix=0; ix<Nx; ix++) {
    int doffset = ix*Ny;
    int coffset = ix*Ny;
    find_NUBcoefs_1d_c (spline->y_basis, yBC, spline->coefs+doffset, 1, 
            spline->coefs+coffset, 1);
  }
    
  return spline;
}

References NUBspline_2d_c::coefs, create_NUBasis(), find_NUBcoefs_1d_c(), BCtype_c::lCode, NU2D, NUgrid::num_points, PERIODIC, posix_memalign(), SINGLE_COMPLEX, NUBspline_2d_c::sp_code, NUBspline_2d_c::t_code, NUBspline_2d_c::x_basis, NUBspline_2d_c::x_stride, and NUBspline_2d_c::y_basis.

create_NUBspline_2d_d()

NUBspline_2d_d * create_NUBspline_2d_d	(	NUgrid *	x_grid,
		NUgrid *	y_grid,
		BCtype_d	xBC,
		BCtype_d	yBC,
		double *	data )

Definition at line 616 of file nubspline_create.cpp.

{
  // First, create the spline structure
  NUBspline_2d_d* spline = static_cast<NUBspline_2d_d*>(malloc(sizeof(NUBspline_2d_d)));
  if (spline == NULL)
    return spline;
  spline->sp_code = NU2D;
  spline->t_code  = DOUBLE_REAL;
 
  // Next, create the bases
  spline->x_basis = create_NUBasis (x_grid, xBC.lCode==PERIODIC);
  spline->y_basis = create_NUBasis (y_grid, yBC.lCode==PERIODIC);
 
  // int Mx, (set but unused)
  int My, Nx, Ny;
//  if (xBC.lCode == PERIODIC) Mx = x_grid->num_points - 1;
//  else                       Mx = x_grid->num_points;
  if (yBC.lCode == PERIODIC) My = y_grid->num_points - 1;
  else                       My = y_grid->num_points;
 
  Nx = x_grid->num_points + 2;
  Ny = y_grid->num_points + 2;
  
  spline->x_stride = Ny;
#ifndef HAVE_SSE2
  spline->coefs = (double*)malloc (sizeof(double)*Nx*Ny);
#else
  posix_memalign ((void**)&spline->coefs, 16, sizeof(double)*Nx*Ny);
#endif
 
  // First, solve in the X-direction 
  for (int iy=0; iy<My; iy++) {
    int doffset = iy;
    int coffset = iy;
    find_NUBcoefs_1d_d (spline->x_basis, xBC, data+doffset, My,
            spline->coefs+coffset, Ny);
  }
  
  // Now, solve in the Y-direction
  for (int ix=0; ix<Nx; ix++) {
    int doffset = ix*Ny;
    int coffset = ix*Ny;
    find_NUBcoefs_1d_d (spline->y_basis, yBC, spline->coefs+doffset, 1, 
            spline->coefs+coffset, 1);
  }
    
  return spline;
}

References NUBspline_2d_d::coefs, create_NUBasis(), DOUBLE_REAL, find_NUBcoefs_1d_d(), BCtype_d::lCode, NU2D, NUgrid::num_points, PERIODIC, posix_memalign(), NUBspline_2d_d::sp_code, NUBspline_2d_d::t_code, NUBspline_2d_d::x_basis, NUBspline_2d_d::x_stride, and NUBspline_2d_d::y_basis.

create_NUBspline_2d_s()

NUBspline_2d_s * create_NUBspline_2d_s	(	NUgrid *	x_grid,
		NUgrid *	y_grid,
		BCtype_s	xBC,
		BCtype_s	yBC,
		float *	data )

Definition at line 272 of file nubspline_create.cpp.

{
  // First, create the spline structure
  NUBspline_2d_s* spline = static_cast<NUBspline_2d_s*>(malloc(sizeof(NUBspline_2d_s)));
  if (spline == NULL)
    return spline;
  spline->sp_code = NU2D;
  spline->t_code  = SINGLE_REAL;
 
  // Next, create the bases
  spline->x_basis = create_NUBasis (x_grid, xBC.lCode==PERIODIC);
  spline->y_basis = create_NUBasis (y_grid, yBC.lCode==PERIODIC);
  // set but unused
  //int Mx,
  int My, Nx, Ny;
//  if (xBC.lCode == PERIODIC) Mx = x_grid->num_points - 1;
//  else                       Mx = x_grid->num_points;
  if (yBC.lCode == PERIODIC) My = y_grid->num_points - 1;
  else                       My = y_grid->num_points;
 
  Nx = x_grid->num_points + 2;
  Ny = y_grid->num_points + 2;
    
  spline->x_stride = Ny;
#ifndef HAVE_SSE2
  spline->coefs = (float*)malloc (sizeof(float)*Nx*Ny);
#else
  posix_memalign ((void**)&spline->coefs, 16, sizeof(float)*Nx*Ny);
#endif
 
  // First, solve in the X-direction 
  for (int iy=0; iy<My; iy++) {
    int doffset = iy;
    int coffset = iy;
    find_NUBcoefs_1d_s (spline->x_basis, xBC, data+doffset, My,
            spline->coefs+coffset, Ny);
  }
  
  // Now, solve in the Y-direction
  for (int ix=0; ix<Nx; ix++) {
    int doffset = ix*Ny;
    int coffset = ix*Ny;
    find_NUBcoefs_1d_s (spline->y_basis, yBC, spline->coefs+doffset, 1, 
            spline->coefs+coffset, 1);
  }
    
  return spline;
}

References NUBspline_2d_s::coefs, create_NUBasis(), find_NUBcoefs_1d_s(), BCtype_s::lCode, NU2D, NUgrid::num_points, PERIODIC, posix_memalign(), SINGLE_REAL, NUBspline_2d_s::sp_code, NUBspline_2d_s::t_code, NUBspline_2d_s::x_basis, NUBspline_2d_s::x_stride, and NUBspline_2d_s::y_basis.

create_NUBspline_2d_z()

NUBspline_2d_z * create_NUBspline_2d_z	(	NUgrid *	x_grid,
		NUgrid *restrict	y_grid,
		BCtype_z	xBC,
		BCtype_z	yBC,
		complex_double *	data )

create_NUBspline_3d_c()

NUBspline_3d_c * create_NUBspline_3d_c	(	NUgrid *	x_grid,
		NUgrid *	y_grid,
		NUgrid *	z_grid,
		BCtype_c	xBC,
		BCtype_c	yBC,
		BCtype_c	zBC,
		complex_float *	data )

Definition at line 839 of file nubspline_create.cpp.

{
  // First, create the spline structure
  NUBspline_3d_c* spline = static_cast<NUBspline_3d_c*>(malloc(sizeof(NUBspline_3d_c)));
  if (spline == NULL)
    return spline;
  spline->sp_code = NU3D;
  spline->t_code  = SINGLE_COMPLEX;
 
  // Next, create the bases
  spline->x_basis = create_NUBasis (x_grid, xBC.lCode==PERIODIC);
  spline->y_basis = create_NUBasis (y_grid, yBC.lCode==PERIODIC);
  spline->z_basis = create_NUBasis (z_grid, zBC.lCode==PERIODIC);
//  int Mx,
  int My, Mz, Nx, Ny, Nz;
//  if (xBC.lCode == PERIODIC) Mx = x_grid->num_points - 1;
//  else                       Mx = x_grid->num_points;
  if (yBC.lCode == PERIODIC) My = y_grid->num_points - 1;
  else                       My = y_grid->num_points;
  if (zBC.lCode == PERIODIC) Mz = z_grid->num_points - 1;
  else                       Mz = z_grid->num_points;
 
  Nx = x_grid->num_points + 2;
  Ny = y_grid->num_points + 2;
  Nz = z_grid->num_points + 2;
 
  // Allocate coefficients and solve  
  spline->x_stride = Ny*Nz;
  spline->y_stride = Nz;
#ifndef HAVE_SSE2
  spline->coefs = (complex_float*)malloc (sizeof(complex_float)*Nx*Ny*Nz);
#else
  posix_memalign ((void**)&spline->coefs, 16, sizeof(complex_float)*Nx*Ny*Nz);
#endif
 
  // First, solve in the X-direction 
  for (int iy=0; iy<My; iy++) 
    for (int iz=0; iz<Mz; iz++) {
      int doffset = iy*Mz+iz;
      int coffset = iy*Nz+iz;
      find_NUBcoefs_1d_c (spline->x_basis, xBC, data+doffset, My*Mz,
              spline->coefs+coffset, Ny*Nz);
    }
  
  // Now, solve in the Y-direction
  for (int ix=0; ix<Nx; ix++) 
    for (int iz=0; iz<Nz; iz++) {
      int doffset = ix*Ny*Nz + iz;
      int coffset = ix*Ny*Nz + iz;
      find_NUBcoefs_1d_c (spline->y_basis, yBC, spline->coefs+doffset, Nz, 
              spline->coefs+coffset, Nz);
    }
 
  // Now, solve in the Z-direction
  for (int ix=0; ix<Nx; ix++) 
    for (int iy=0; iy<Ny; iy++) {
      int doffset = (ix*Ny+iy)*Nz;
      int coffset = (ix*Ny+iy)*Nz;
      find_NUBcoefs_1d_c (spline->z_basis, zBC, spline->coefs+doffset, 1, 
              spline->coefs+coffset, 1);
    }
  return spline;
}

References NUBspline_3d_c::coefs, create_NUBasis(), find_NUBcoefs_1d_c(), BCtype_c::lCode, NU3D, NUgrid::num_points, PERIODIC, posix_memalign(), SINGLE_COMPLEX, NUBspline_3d_c::sp_code, NUBspline_3d_c::t_code, NUBspline_3d_c::x_basis, NUBspline_3d_c::x_stride, NUBspline_3d_c::y_basis, NUBspline_3d_c::y_stride, and NUBspline_3d_c::z_basis.

create_NUBspline_3d_d()

NUBspline_3d_d * create_NUBspline_3d_d	(	NUgrid *	x_grid,
		NUgrid *	y_grid,
		NUgrid *	z_grid,
		BCtype_d	xBC,
		BCtype_d	yBC,
		BCtype_d	zBC,
		double *	data )

Definition at line 668 of file nubspline_create.cpp.

{
  // First, create the spline structure
  NUBspline_3d_d* spline = static_cast<NUBspline_3d_d*>(malloc(sizeof(NUBspline_3d_d)));
  if (spline == NULL)
    return spline;
  spline->sp_code = NU3D;
  spline->t_code  = DOUBLE_REAL;
 
  // Next, create the bases
  spline->x_basis = create_NUBasis (x_grid, xBC.lCode==PERIODIC);
  spline->y_basis = create_NUBasis (y_grid, yBC.lCode==PERIODIC);
  spline->z_basis = create_NUBasis (z_grid, zBC.lCode==PERIODIC);
 
  // set but unused
  //int Mx,
  int My, Mz, Nx, Ny, Nz;
//  if (xBC.lCode == PERIODIC) Mx = x_grid->num_points - 1;
//  else                       Mx = x_grid->num_points;
  if (yBC.lCode == PERIODIC) My = y_grid->num_points - 1;
  else                       My = y_grid->num_points;
  if (zBC.lCode == PERIODIC) Mz = z_grid->num_points - 1;
  else                       Mz = z_grid->num_points;
 
  Nx = x_grid->num_points + 2;
  Ny = y_grid->num_points + 2;
  Nz = z_grid->num_points + 2;
  
  spline->x_stride = Ny*Nz;
  spline->y_stride = Nz;
#ifndef HAVE_SSE2
  spline->coefs = (double*)malloc (sizeof(double)*Nx*Ny*Nz);
#else
  posix_memalign ((void**)&spline->coefs, 16, sizeof(double)*Nx*Ny*Nz);
#endif
 
  // First, solve in the X-direction 
  for (int iy=0; iy<My; iy++) 
    for (int iz=0; iz<Mz; iz++) {
      int doffset = iy*Mz+iz;
      int coffset = iy*Nz+iz;
      find_NUBcoefs_1d_d (spline->x_basis, xBC, data+doffset, My*Mz,
              spline->coefs+coffset, Ny*Nz);
    }
  
  // Now, solve in the Y-direction
  for (int ix=0; ix<Nx; ix++) 
    for (int iz=0; iz<Nz; iz++) {
      int doffset = ix*Ny*Nz + iz;
      int coffset = ix*Ny*Nz + iz;
      find_NUBcoefs_1d_d (spline->y_basis, yBC, spline->coefs+doffset, Nz, 
              spline->coefs+coffset, Nz);
    }
 
  // Now, solve in the Z-direction
  for (int ix=0; ix<Nx; ix++) 
    for (int iy=0; iy<Ny; iy++) {
      int doffset = (ix*Ny+iy)*Nz;
      int coffset = (ix*Ny+iy)*Nz;
      find_NUBcoefs_1d_d (spline->z_basis, zBC, spline->coefs+doffset, 1, 
              spline->coefs+coffset, 1);
    }
  return spline;
}

References NUBspline_3d_d::coefs, create_NUBasis(), DOUBLE_REAL, find_NUBcoefs_1d_d(), BCtype_d::lCode, NU3D, NUgrid::num_points, PERIODIC, posix_memalign(), NUBspline_3d_d::sp_code, NUBspline_3d_d::t_code, NUBspline_3d_d::x_basis, NUBspline_3d_d::x_stride, NUBspline_3d_d::y_basis, NUBspline_3d_d::y_stride, and NUBspline_3d_d::z_basis.

create_NUBspline_3d_s()

NUBspline_3d_s * create_NUBspline_3d_s	(	NUgrid *	x_grid,
		NUgrid *	y_grid,
		NUgrid *	z_grid,
		BCtype_s	xBC,
		BCtype_s	yBC,
		BCtype_s	zBC,
		float *	data )

Definition at line 324 of file nubspline_create.cpp.

{
  // First, create the spline structure
  NUBspline_3d_s* spline = static_cast<NUBspline_3d_s*>(malloc(sizeof(NUBspline_3d_s)));
  if (spline == NULL)
    return spline;
  spline->sp_code = NU3D;
  spline->t_code  = SINGLE_REAL;
 
  // Next, create the bases
  spline->x_basis = create_NUBasis (x_grid, xBC.lCode==PERIODIC);
  spline->y_basis = create_NUBasis (y_grid, yBC.lCode==PERIODIC);
  spline->z_basis = create_NUBasis (z_grid, zBC.lCode==PERIODIC);
  // set but unused
  // int Mx,
  int My, Mz, Nx, Ny, Nz;
//  if (xBC.lCode == PERIODIC) Mx = x_grid->num_points - 1;
//  else                       Mx = x_grid->num_points;
  if (yBC.lCode == PERIODIC) My = y_grid->num_points - 1;
  else                       My = y_grid->num_points;
  if (zBC.lCode == PERIODIC) Mz = z_grid->num_points - 1;
  else                       Mz = z_grid->num_points;
 
  Nx = x_grid->num_points + 2;
  Ny = y_grid->num_points + 2;
  Nz = z_grid->num_points + 2;
 
  // Allocate coefficients and solve  
  spline->x_stride = Ny*Nz;
  spline->y_stride = Nz;
#ifndef HAVE_SSE2
  spline->coefs = (float*)malloc (sizeof(float)*Nx*Ny*Nz);
#else
  posix_memalign ((void**)&spline->coefs, 16, sizeof(float)*Nx*Ny*Nz);
#endif
 
  // First, solve in the X-direction 
  for (int iy=0; iy<My; iy++) 
    for (int iz=0; iz<Mz; iz++) {
      int doffset = iy*Mz+iz;
      int coffset = iy*Nz+iz;
      find_NUBcoefs_1d_s (spline->x_basis, xBC, data+doffset, My*Mz,
              spline->coefs+coffset, Ny*Nz);
    }
  
  // Now, solve in the Y-direction
  for (int ix=0; ix<Nx; ix++) 
    for (int iz=0; iz<Nz; iz++) {
      int doffset = ix*Ny*Nz + iz;
      int coffset = ix*Ny*Nz + iz;
      find_NUBcoefs_1d_s (spline->y_basis, yBC, spline->coefs+doffset, Nz, 
              spline->coefs+coffset, Nz);
    }
 
  // Now, solve in the Z-direction
  for (int ix=0; ix<Nx; ix++) 
    for (int iy=0; iy<Ny; iy++) {
      int doffset = (ix*Ny+iy)*Nz;
      int coffset = (ix*Ny+iy)*Nz;
      find_NUBcoefs_1d_s (spline->z_basis, zBC, spline->coefs+doffset, 1, 
              spline->coefs+coffset, 1);
    }
  return spline;
}

References NUBspline_3d_s::coefs, create_NUBasis(), find_NUBcoefs_1d_s(), BCtype_s::lCode, NU3D, NUgrid::num_points, PERIODIC, posix_memalign(), SINGLE_REAL, NUBspline_3d_s::sp_code, NUBspline_3d_s::t_code, NUBspline_3d_s::x_basis, NUBspline_3d_s::x_stride, NUBspline_3d_s::y_basis, NUBspline_3d_s::y_stride, and NUBspline_3d_s::z_basis.

create_NUBspline_3d_z()

NUBspline_3d_z * create_NUBspline_3d_z	(	NUgrid *	x_grid,
		NUgrid *	y_grid,
		NUgrid *	z_grid,
		BCtype_z	xBC,
		BCtype_z	yBC,
		BCtype_z	zBC,
		complex_double *	data )

Definition at line 1008 of file nubspline_create.cpp.

{
  // First, create the spline structure
  NUBspline_3d_z* spline = static_cast<NUBspline_3d_z*>(malloc(sizeof(NUBspline_3d_z)));
  if (spline == NULL)
    return spline;
  spline->sp_code = NU3D;
  spline->t_code  = DOUBLE_COMPLEX;
  spline->x_grid = x_grid;
  spline->y_grid = y_grid;
  spline->z_grid = z_grid;
 
  // Next, create the bases
  spline->x_basis = create_NUBasis (x_grid, xBC.lCode==PERIODIC);
  spline->y_basis = create_NUBasis (y_grid, yBC.lCode==PERIODIC);
  spline->z_basis = create_NUBasis (z_grid, zBC.lCode==PERIODIC);
//  int Mx,
  int My, Mz, Nx, Ny, Nz;
//  if (xBC.lCode == PERIODIC) Mx = x_grid->num_points - 1;
//  else                       Mx = x_grid->num_points;
  if (yBC.lCode == PERIODIC) My = y_grid->num_points - 1;
  else                       My = y_grid->num_points;
  if (zBC.lCode == PERIODIC) Mz = z_grid->num_points - 1;
  else                       Mz = z_grid->num_points;
 
  Nx = x_grid->num_points + 2;
  Ny = y_grid->num_points + 2;
  Nz = z_grid->num_points + 2;
 
  // Allocate coefficients and solve  
  spline->x_stride = Ny*Nz;
  spline->y_stride = Nz;
#ifndef HAVE_SSE2
  spline->coefs = (complex_double*)malloc (sizeof(complex_double)*Nx*Ny*Nz);
#else
  posix_memalign ((void**)&spline->coefs, 16, sizeof(complex_double)*Nx*Ny*Nz);
#endif
 
  // First, solve in the X-direction 
  for (int iy=0; iy<My; iy++) 
    for (int iz=0; iz<Mz; iz++) {
      int doffset = iy*Mz+iz;
      int coffset = iy*Nz+iz;
      find_NUBcoefs_1d_z (spline->x_basis, xBC, data+doffset, My*Mz,
              spline->coefs+coffset, Ny*Nz);
      /* for (int ix=0; ix<Nx; ix++) {
    complex_double z = spline->coefs[coffset+ix*spline->x_stride];
    if (isnan(creal(z)))
      fprintf (stderr, "NAN encountered in create_NUBspline_3d_z at real part of (%d,%d,%d)\n",
           ix,iy,iz);
    if (isnan(cimag(z)))
      fprintf (stderr, "NAN encountered in create_NUBspline_3d_z at imag part of (%d,%d,%d)\n",
           ix,iy,iz);
       } */
    }
  
  // Now, solve in the Y-direction
  for (int ix=0; ix<Nx; ix++) 
    for (int iz=0; iz<Nz; iz++) {
      int doffset = ix*Ny*Nz + iz;
      int coffset = ix*Ny*Nz + iz;
      find_NUBcoefs_1d_z (spline->y_basis, yBC, spline->coefs+doffset, Nz, 
              spline->coefs+coffset, Nz);
    }
 
  // Now, solve in the Z-direction
  for (int ix=0; ix<Nx; ix++) 
    for (int iy=0; iy<Ny; iy++) {
      int doffset = (ix*Ny+iy)*Nz;
      int coffset = (ix*Ny+iy)*Nz;
      find_NUBcoefs_1d_z (spline->z_basis, zBC, spline->coefs+doffset, 1, 
              spline->coefs+coffset, 1);
    }
  return spline;
}

References NUBspline_3d_z::coefs, create_NUBasis(), DOUBLE_COMPLEX, find_NUBcoefs_1d_z(), BCtype_z::lCode, NU3D, NUgrid::num_points, PERIODIC, posix_memalign(), NUBspline_3d_z::sp_code, NUBspline_3d_z::t_code, NUBspline_3d_z::x_basis, NUBspline_3d_z::x_grid, NUBspline_3d_z::x_stride, NUBspline_3d_z::y_basis, NUBspline_3d_z::y_grid, NUBspline_3d_z::y_stride, NUBspline_3d_z::z_basis, and NUBspline_3d_z::z_grid.