Repository: xgcm/xgcm (python)

Issue:
Globally set `fill_value` set on grid object is not propagated.
I just discovered what I think is a bug:

Consider this simple example:
```
from xgcm import Grid
import xarray as xr
import numpy as np

z = np.arange(3)
z_right = z+0.5
ds = xr.DataArray(np.ones_like(z)*1e5, dims=['z'], coords={'z':z}).to_dataset(name='data')
ds = ds.assign_coords({'z_right':z_right})

grid = Grid(ds, coords={'Z':{'center':'z', 'right':'z_right'}},
            boundary='fill', fill_value=np.nan, periodic=False)
grid.interp(ds.data, 'Z')
```
This results in:

![image](https://user-images.githubusercontent.com/14314623/92794864-f24fd080-f37d-11ea-85cf-485dd378c971.jpeg)

Which is not what I would expect (the edge value should be padded with nan) and is not equivalent to the 'old-school' manual way:
```
grid.interp(ds.data, 'Z',boundary='fill', fill_value=np.nan)
```
![image](https://user-images.githubusercontent.com/14314623/92795073-190e0700-f37e-11ea-90b1-0ae3ba25de54.jpeg)

I think the issue is that `fill_value=0.0` is hardcoded in all axis methods as default. See [here](https://github.com/xgcm/xgcm/blob/9c83d5368193237d80f7282f94a76561e0711c51/xgcm/grid.py#L614) for example.


We should add a test that checks that results are the same if you specify `boundary` and/or `fill_value` via grid object or in the axis method. 

Pinging @dcherian 

Pull request:
Bugfix 224
<!-- Feel free to remove check-list items aren't relevant to your change -->

 - [x] Closes #224 
 - [x] Tests added
 - [x] Passes `black . `
 - [ ] User visible changes (including notable bug fixes) are documented in `whats-new.rst`

Patch:
diff --git a/xgcm/grid.py b/xgcm/grid.py
index 42ee589..e8b1c68 100644
--- a/xgcm/grid.py
+++ b/xgcm/grid.py
@@ -581,7 +581,7 @@ class Axis:
         da,
         to=None,
         boundary=None,
-        fill_value=0.0,
+        fill_value=None,
         boundary_discontinuity=None,
         vector_partner=None,
         keep_coords=False,
@@ -619,7 +619,7 @@ class Axis:
         da,
         to=None,
         boundary=None,
-        fill_value=0.0,
+        fill_value=None,
         boundary_discontinuity=None,
         vector_partner=None,
         keep_coords=False,
@@ -702,7 +702,7 @@ class Axis:
         da,
         to=None,
         boundary=None,
-        fill_value=0.0,
+        fill_value=None,
         boundary_discontinuity=None,
         vector_partner=None,
         keep_coords=False,
@@ -737,7 +737,7 @@ class Axis:
         da,
         to=None,
         boundary=None,
-        fill_value=0.0,
+        fill_value=None,
         boundary_discontinuity=None,
         vector_partner=None,
         keep_coords=False,


Test patch:
diff --git a/xgcm/test/test_metrics_ops_single_axis.py b/xgcm/test/test_metrics_ops_single_axis.py
index 5835dca..6d1ffd9 100644
--- a/xgcm/test/test_metrics_ops_single_axis.py
+++ b/xgcm/test/test_metrics_ops_single_axis.py
@@ -71,6 +71,43 @@ class TestParametrized:
             assert new.equals(expected)
 
 
+@pytest.mark.parametrize(
+    "funcname",
+    ["interp", "diff", "min", "max", "cumsum", "derivative", "cumint"],
+)
+@pytest.mark.parametrize("boundary", ["fill", "extend"])
+@pytest.mark.parametrize("fill_value", [0, 10, None])
+def test_boundary_global_input(funcname, boundary, fill_value):
+    """Test that globally defined boundary values result in
+    the same output as when the parameters are defined on either
+    the grid or axis methods
+    """
+    ds, coords, metrics = datasets_grid_metric("C")
+    axis = "X"
+
+    # Test results by globally specifying fill value/boundary on grid object
+    grid_global = Grid(
+        ds,
+        coords=coords,
+        metrics=metrics,
+        periodic=False,
+        boundary=boundary,
+        fill_value=fill_value,
+    )
+    func_global = getattr(grid_global, funcname)
+    global_result = func_global(ds.tracer, axis)
+
+    # Test results by manually specifying fill value/boundary on grid method
+    grid_manual = Grid(
+        ds, coords=coords, metrics=metrics, periodic=False, boundary=boundary
+    )
+    func_manual = getattr(grid_manual, funcname)
+    manual_result = func_manual(
+        ds.tracer, axis, boundary=boundary, fill_value=fill_value
+    )
+    xr.testing.assert_allclose(global_result, manual_result)
+
+
 def test_derivative_uniform_grid():
     # this is a uniform grid
     # a non-uniform grid would provide a more rigorous test