diff --git a/SimPEG/Mesh/DiffOperators.py b/SimPEG/Mesh/DiffOperators.py
index 8b431a65..7a965e0f 100644
--- a/SimPEG/Mesh/DiffOperators.py
+++ b/SimPEG/Mesh/DiffOperators.py
@@ -478,6 +478,23 @@ class DiffOperators(object):
                                                    kron3(av(n[2]), speye(n[1]), speye(n[0]))), format="csr")
         return self._aveF2CC
 
+
+    @property
+    def aveF2CCV(self):
+        "Construct the averaging operator on cell faces to cell centers."
+        if getattr(self, '_aveF2CCV', None) is None:
+            n = self.nCv
+            if(self.dim == 1):
+                self._aveF2CCV = av(n[0])
+            elif(self.dim == 2):
+                self._aveF2CCV = sp.block_diag((sp.kron(speye(n[1]), av(n[0])),
+                                               sp.kron(av(n[1]), speye(n[0]))), format="csr")
+            elif(self.dim == 3):
+                self._aveF2CCV = sp.block_diag((kron3(speye(n[2]), speye(n[1]), av(n[0])),
+                                               kron3(speye(n[2]), av(n[1]), speye(n[0])),
+                                               kron3(av(n[2]), speye(n[1]), speye(n[0]))), format="csr")
+        return self._aveF2CCV
+
     @property
     def aveCC2F(self):
         "Construct the averaging operator on cell cell centers to faces."
diff --git a/SimPEG/Tests/test_operators.py b/SimPEG/Tests/test_operators.py
index d1c38443..67741aa7 100644
--- a/SimPEG/Tests/test_operators.py
+++ b/SimPEG/Tests/test_operators.py
@@ -338,6 +338,15 @@ class TestAveraging2D(OrderTest):
         self.getAve = lambda M: M.aveF2CC
         self.orderTest()
 
+    def test_orderF2CCV(self):
+        self.name = "Averaging 2D: F2CCV"
+        funX = lambda x, y: (np.cos(x)+np.sin(y))
+        funY = lambda x, y: (np.cos(y)*np.sin(x))
+        self.getHere = lambda M: np.r_[call2(funX, M.gridFx), call2(funY, M.gridFy)]
+        self.getThere = lambda M: np.r_[call2(funX, M.gridCC), call2(funY, M.gridCC)]
+        self.getAve = lambda M: M.aveF2CCV
+        self.orderTest()
+
     def test_orderCC2F(self):
         self.name = "Averaging 2D: CC2F"
         fun = lambda x, y: (np.cos(x)+np.sin(y))
@@ -348,7 +357,6 @@ class TestAveraging2D(OrderTest):
         self.orderTest()
         self.expectedOrders = 2
 
-
     def test_orderE2CC(self):
         self.name = "Averaging 2D: E2CC"
         fun = lambda x, y: (np.cos(x)+np.sin(y))
@@ -357,6 +365,15 @@ class TestAveraging2D(OrderTest):
         self.getAve = lambda M: M.aveE2CC
         self.orderTest()
 
+    def test_orderE2CCV(self):
+        self.name = "Averaging 2D: E2CCV"
+        funX = lambda x, y: (np.cos(x)+np.sin(y))
+        funY = lambda x, y: (np.cos(y)*np.sin(x))
+        self.getHere  = lambda M: np.r_[call2(funX, M.gridEx), call2(funY, M.gridEy)]
+        self.getThere = lambda M: np.r_[call2(funX, M.gridCC), call2(funY, M.gridCC)]
+        self.getAve = lambda M: M.aveE2CCV
+        self.orderTest()
+
 
 class TestAveraging3D(OrderTest):
     name = "Averaging 3D"
@@ -400,6 +417,15 @@ class TestAveraging3D(OrderTest):
         self.getAve = lambda M: M.aveF2CC
         self.orderTest()
 
+    def test_orderF2CCV(self):
+        self.name = "Averaging 3D: F2CCV"
+        funX = lambda x, y, z: (np.cos(x)+np.sin(y)+np.exp(z))
+        funY = lambda x, y, z: (np.cos(x)+np.sin(y)*np.exp(z))
+        funZ = lambda x, y, z: (np.cos(x)*np.sin(y)+np.exp(z))
+        self.getHere = lambda M: np.r_[call3(funX, M.gridFx), call3(funY, M.gridFy), call3(funZ, M.gridFz)]
+        self.getThere = lambda M: np.r_[call3(funX, M.gridCC), call3(funY, M.gridCC), call3(funZ, M.gridCC)]
+        self.getAve = lambda M: M.aveF2CCV
+        self.orderTest()
 
     def test_orderE2CC(self):
         self.name = "Averaging 3D: E2CC"
@@ -409,6 +435,16 @@ class TestAveraging3D(OrderTest):
         self.getAve = lambda M: M.aveE2CC
         self.orderTest()
 
+    def test_orderE2CCV(self):
+        self.name = "Averaging 3D: E2CCV"
+        funX = lambda x, y, z: (np.cos(x)+np.sin(y)+np.exp(z))
+        funY = lambda x, y, z: (np.cos(x)+np.sin(y)*np.exp(z))
+        funZ = lambda x, y, z: (np.cos(x)*np.sin(y)+np.exp(z))
+        self.getHere = lambda M: np.r_[call3(funX, M.gridEx), call3(funY, M.gridEy), call3(funZ, M.gridEz)]
+        self.getThere = lambda M: np.r_[call3(funX, M.gridCC), call3(funY, M.gridCC), call3(funZ, M.gridCC)]
+        self.getAve = lambda M: M.aveE2CCV
+        self.orderTest()
+
     def test_orderCC2F(self):
         self.name = "Averaging 3D: CC2F"
         fun = lambda x, y, z: (np.cos(x)+np.sin(y)+np.exp(z))