diff --git a/src/test/java/org/apache/commons/math3/analysis/interpolation/BicubicSplineInterpolatorTest.java b/src/test/java/org/apache/commons/math3/analysis/interpolation/BicubicSplineInterpolatorTest.java
index a5a8dc65e..c45f9ac56 100644
--- a/src/test/java/org/apache/commons/math3/analysis/interpolation/BicubicSplineInterpolatorTest.java
+++ b/src/test/java/org/apache/commons/math3/analysis/interpolation/BicubicSplineInterpolatorTest.java
@@ -19,6 +19,9 @@ package org.apache.commons.math3.analysis.interpolation;
 import org.apache.commons.math3.exception.DimensionMismatchException;
 import org.apache.commons.math3.exception.MathIllegalArgumentException;
 import org.apache.commons.math3.analysis.BivariateFunction;
+import org.apache.commons.math3.distribution.UniformRealDistribution;
+import org.apache.commons.math3.random.RandomGenerator;
+import org.apache.commons.math3.random.Well19937c;
 import org.junit.Assert;
 import org.junit.Test;
 import org.junit.Ignore;
@@ -76,22 +79,28 @@ public final class BicubicSplineInterpolatorTest {
     }
 
     /**
-     * Test of interpolator for a plane.
+     * Interpolating a plane.
      * <p>
      * z = 2 x - 3 y + 5
      */
-    @Ignore@Test
-    public void testPlane() {
+    @Test
+    public void testInterpolation1() {
+        final int sz = 21;
+        double[] xval = new double[sz];
+        double[] yval = new double[sz];
+        // Coordinate values
+        final double delta = 1d / (sz - 1);
+        for (int i = 0; i < sz; i++) {
+            xval[i] = -1 + 15 * i * delta;
+            yval[i] = -20 + 30 * i * delta;
+        }
+
+        // Function values
         BivariateFunction f = new BivariateFunction() {
                 public double value(double x, double y) {
                     return 2 * x - 3 * y + 5;
                 }
             };
-
-        BivariateGridInterpolator interpolator = new BicubicSplineInterpolator();
-
-        double[] xval = new double[] {3, 4, 5, 6.5};
-        double[] yval = new double[] {-4, -3, -1, 2, 2.5};
         double[][] zval = new double[xval.length][yval.length];
         for (int i = 0; i < xval.length; i++) {
             for (int j = 0; j < yval.length; j++) {
@@ -99,46 +108,52 @@ public final class BicubicSplineInterpolatorTest {
             }
         }
 
+        BivariateGridInterpolator interpolator = new BicubicSplineInterpolator();
         BivariateFunction p = interpolator.interpolate(xval, yval, zval);
         double x, y;
-        double expected, result;
-        
-        x = 4;
-        y = -3;
-        expected = f.value(x, y);
-        result = p.value(x, y);
-        Assert.assertEquals("On sample point", expected, result, 1e-15);
 
-        x = 4.5;
-        y = -1.5;
-        expected = f.value(x, y);
-        result = p.value(x, y);
-        Assert.assertEquals("half-way between sample points (middle of the patch)", expected, result, 0.3);
+        final RandomGenerator rng = new Well19937c(1234567L); // "tol" depends on the seed.
+        final UniformRealDistribution distX
+            = new UniformRealDistribution(rng, xval[0], xval[xval.length - 1]);
+        final UniformRealDistribution distY
+            = new UniformRealDistribution(rng, yval[0], yval[yval.length - 1]);
 
-        x = 3.5;
-        y = -3.5;
-        expected = f.value(x, y);
-        result = p.value(x, y);
-        Assert.assertEquals("half-way between sample points (border of the patch)", expected, result, 0.3);
+        final int numSamples = 50;
+        final double tol = 6;
+        for (int i = 0; i < numSamples; i++) {
+            x = distX.sample();
+            for (int j = 0; j < numSamples; j++) {
+                y = distY.sample();
+//                 System.out.println(x + " " + y + " " + f.value(x, y) + " " + p.value(x, y));
+                Assert.assertEquals(f.value(x, y),  p.value(x, y), tol);
+            }
+//             System.out.println();
+        }
     }
 
     /**
-     * Test of interpolator for a paraboloid.
+     * Interpolating a paraboloid.
      * <p>
      * z = 2 x<sup>2</sup> - 3 y<sup>2</sup> + 4 x y - 5
      */
     @Test
-    public void testParaboloid() {
+    public void testInterpolation2() {
+        final int sz = 21;
+        double[] xval = new double[sz];
+        double[] yval = new double[sz];
+        // Coordinate values
+        final double delta = 1d / (sz - 1);
+        for (int i = 0; i < sz; i++) {
+            xval[i] = -1 + 15 * i * delta;
+            yval[i] = -20 + 30 * i * delta;
+        }
+
+        // Function values
         BivariateFunction f = new BivariateFunction() {
                 public double value(double x, double y) {
                     return 2 * x * x - 3 * y * y + 4 * x * y - 5;
                 }
             };
-
-        BivariateGridInterpolator interpolator = new BicubicSplineInterpolator();
-
-        double[] xval = new double[] {3, 4, 5, 6.5};
-        double[] yval = new double[] {-4, -3, -2, -1, 0.5, 2.5};
         double[][] zval = new double[xval.length][yval.length];
         for (int i = 0; i < xval.length; i++) {
             for (int j = 0; j < yval.length; j++) {
@@ -146,26 +161,26 @@ public final class BicubicSplineInterpolatorTest {
             }
         }
 
+        BivariateGridInterpolator interpolator = new BicubicSplineInterpolator();
         BivariateFunction p = interpolator.interpolate(xval, yval, zval);
         double x, y;
-        double expected, result;
-        
-        x = 5;
-        y = 0.5;
-        expected = f.value(x, y);
-        result = p.value(x, y);
-        Assert.assertEquals("On sample point", expected, result, 1e-13);
 
-        x = 4.5;
-        y = -1.5;
-        expected = f.value(x, y);
-        result = p.value(x, y);
-        Assert.assertEquals("half-way between sample points (middle of the patch)", expected, result, 0.2);
+        final RandomGenerator rng = new Well19937c(1234567L); // "tol" depends on the seed.
+        final UniformRealDistribution distX
+            = new UniformRealDistribution(rng, xval[0], xval[xval.length - 1]);
+        final UniformRealDistribution distY
+            = new UniformRealDistribution(rng, yval[0], yval[yval.length - 1]);
 
-        x = 3.5;
-        y = -3.5;
-        expected = f.value(x, y);
-        result = p.value(x, y);
-        Assert.assertEquals("half-way between sample points (border of the patch)", expected, result, 0.2);
+        final int numSamples = 50;
+        final double tol = 251;
+        for (int i = 0; i < numSamples; i++) {
+            x = distX.sample();
+            for (int j = 0; j < numSamples; j++) {
+                y = distY.sample();
+//                 System.out.println(x + " " + y + " " + f.value(x, y) + " " + p.value(x, y));
+                Assert.assertEquals(f.value(x, y),  p.value(x, y), tol);
+            }
+//             System.out.println();
+        }
     }
 }