guard operator()

Paul-St-Young · Jun 25, 2018 · 01f7b95 · 01f7b95
1 parent f9f278f
commit 01f7b95
Show file tree

Hide file tree

Showing 2 changed files with 11 additions and 0 deletions.
diff --git a/src/QMCFiniteSize/NaturalSpline3D.cpp b/src/QMCFiniteSize/NaturalSpline3D.cpp
@@ -25,8 +25,18 @@ BCtype_d NaturalSpline3D::natural_boundary()
   bc.rVal = 1.0;
   return bc;
 }
+bool NaturalSpline3D::in_convex_hull(double x, double y, double z)
+{
+  bool inx, iny, inz;
+  inx = iny = inz = true;
+  if (x<grid3d_.x.start | x>grid3d_.x.end) inx = false;
+  if (y<grid3d_.y.start | y>grid3d_.y.end) iny = false;
+  if (z<grid3d_.z.start | z>grid3d_.z.end) inz = false;
+  return (inx & iny & inz);
+}
 double NaturalSpline3D::operator()(double x, double y, double z)
 {
+  if (not in_convex_hull(x, y, z)) return NAN;
   double val;
   eval_UBspline_3d_d(spline3d_, x, y, z, &val);
   return val;

diff --git a/src/QMCFiniteSize/NaturalSpline3D.h b/src/QMCFiniteSize/NaturalSpline3D.h
@@ -17,6 +17,7 @@ class NaturalSpline3D
   NaturalSpline3D(const Ugrid3D grid3d, double* vals);
   ~NaturalSpline3D();
   BCtype_d natural_boundary();
+  bool in_convex_hull(double x, double y, double z);
 
   // goal in life: evaluate Bspline at a point in 3D space
   double operator()(double x, double y, double z);