--- src/mg/apbs/vgrid.h.original	2004-11-10 07:49:43.000000000 -0800
+++ src/mg/apbs/vgrid.h	2007-01-31 14:38:54.910353800 -0800
@@ -56,7 +56,7 @@
 
 /** @brief Number of decimal places for comparisons and formatting
  *  @ingroup Vgrid */
-#define VGRID_DIGITS 6
+#define VGRID_DIGITS 5
 
 /**
  *  @ingroup Vgrid
--- src/mg/vmgrid.c.original	2004-01-19 13:01:38.000000000 -0800
+++ src/mg/vmgrid.c	2007-01-31 14:38:16.332228800 -0800
@@ -103,6 +103,8 @@
 VPUBLIC int Vmgrid_value(Vmgrid *thee, double pt[3], double *value) {
 
     int i, rc;
+    static int maxErrorPrint = 1000;
+    static int hasPrinted    = 0;
     double tvalue;
   
     VASSERT(thee != VNULL);
@@ -115,8 +117,62 @@
         }
     }
 
-    Vnm_print(2, "Vmgrid_value:  Point (%g, %g, %g) not found in \
-hiearchy!\n", pt[0], pt[1], pt[2]);
+    if( hasPrinted < maxErrorPrint ){
+
+       int nx, ny, nz, ihi, jhi, khi, ilo, jlo, klo, ii;
+       double hx, hy, hzed, xmin, ymin, zmin, ifloat, jfloat, kfloat;
+       double xmax, ymax, zmax;
+       double u, dx, dy, dz;
+
+       if( ++hasPrinted > 1 ){
+
+          double Vcompare = 1.0e-04;
+
+          Vnm_print(2, "Vmgrid_value: Grid info: ngrids=%d.\n", thee->ngrids );
+
+          for (ii=0; ii<thee->ngrids; ii++) {
+             Vnm_print(2, "\nGrid %d Ndim=[%d %d %d] Spacing=[[%.5f %.5f %.5f] \n", ii+1,
+                           thee->grids[ii]->nx, thee->grids[ii]->ny, thee->grids[ii]->nz,
+                           thee->grids[ii]->hx, thee->grids[ii]->hy, thee->grids[ii]->hzed );
+
+             Vnm_print(2, "         Rng=[%.5f %.5f %.5f] [%.5f %.5f %.5f]\n",
+                           thee->grids[ii]->xmin, thee->grids[ii]->ymin, thee->grids[ii]->zmin,
+                           thee->grids[ii]->xmax, thee->grids[ii]->ymax, thee->grids[ii]->zmax );
+
+             ifloat = (pt[0] - thee->grids[ii]->xmin)/thee->grids[ii]->hx;
+             jfloat = (pt[1] - thee->grids[ii]->ymin)/thee->grids[ii]->hy;
+             kfloat = (pt[2] - thee->grids[ii]->zmin)/thee->grids[ii]->hzed;
+
+             ihi = (int)ceil(ifloat);
+             jhi = (int)ceil(jfloat);
+             khi = (int)ceil(kfloat);
+             ilo = (int)floor(ifloat);
+             jlo = (int)floor(jfloat);
+             klo = (int)floor(kfloat);
+
+             if (VABS(pt[0] - thee->grids[ii]->xmin) < Vcompare) ilo = 0;
+             if (VABS(pt[1] - thee->grids[ii]->ymin) < Vcompare) jlo = 0;
+             if (VABS(pt[2] - thee->grids[ii]->zmin) < Vcompare) klo = 0;
+             if (VABS(pt[0] - thee->grids[ii]->xmax) < Vcompare) ihi = nx-1;
+             if (VABS(pt[1] - thee->grids[ii]->ymax) < Vcompare) jhi = ny-1;
+             if (VABS(pt[2] - thee->grids[ii]->zmax) < Vcompare) khi = nz-1;
+             /* See if we're on the mesh */
+             if ((ihi<nx) && (jhi<ny) && (khi<nz) &&
+                (ilo>=0) && (jlo>=0) && (klo>=0)) {
+
+                Vnm_print(2, "          On grid:");
+             } else {
+                Vnm_print(2, "         Off grid:");
+             }
+             Vnm_print(2, " Lo[%d %d %d] Hi[%d %d %d] float=[%.5e %.5e %.5e] Vcom=%.4e\n", ilo, jlo, klo, ihi, jhi, khi,
+                          ifloat, jfloat, kfloat, Vcompare );
+          }
+       }
+       if( hasPrinted == maxErrorPrint ){
+          Vnm_print(2, "Vmgrid_value: last time printing following error message.\n" );
+       }
+       Vnm_print(2, "Vmgrid_value:  Point (%.6f, %.6f, %.6f) not found in hierarchy!\n", pt[0], pt[1], pt[2]);
+    }
 
     return 0;
 }
