From f8259456d03f2db29279730b432ae83c25f2d16c Mon Sep 17 00:00:00 2001
From: lheagy <lheagy@eos.ubc.ca>
Date: Thu, 5 Feb 2015 12:44:46 -0800
Subject: [PATCH] start of cyl mesh anisotropy

---
 SimPEG/Maps.py         | 24 ++++++++++++++++++++++++
 SimPEG/Mesh/CylMesh.py | 25 +++++++++++++++++++++++++
 2 files changed, 49 insertions(+)

diff --git a/SimPEG/Maps.py b/SimPEG/Maps.py
index 4ab1d855..e37d435a 100644
--- a/SimPEG/Maps.py
+++ b/SimPEG/Maps.py
@@ -640,3 +640,27 @@ class CircleMap(IdentityMap):
         g4 = a*(-Y + y)*(-sig1 + sig2)/(np.pi*(a**2*(-r + np.sqrt((X - x)**2 + (Y - y)**2))**2 + 1)*np.sqrt((X - x)**2 + (Y - y)**2))
         g5 = -a*(-sig1 + sig2)/(np.pi*(a**2*(-r + np.sqrt((X - x)**2 + (Y - y)**2))**2 + 1))
         return np.c_[g1,g2,g3,g4,g5]
+
+class AnisotropyMap(IdentityMap):
+
+    def __init__(self, mesh, **kwargs):
+        Utils.setKwargs(self, **kwargs)
+        self.mesh = mesh
+
+    @property
+    def nP(self):
+        """
+            :rtype: int
+            :return: number of parameters in the model
+        """
+        return self.mesh.nC*2
+
+    @property
+    def shape(self):
+        """
+           The default shape is (mesh.nC, nP).
+
+           :rtype: (int,int)
+           :return: shape of the operator as a tuple
+        """
+        return (self.mesh.nC*2, self.nP)
diff --git a/SimPEG/Mesh/CylMesh.py b/SimPEG/Mesh/CylMesh.py
index 93cd5451..d3a40dec 100644
--- a/SimPEG/Mesh/CylMesh.py
+++ b/SimPEG/Mesh/CylMesh.py
@@ -278,6 +278,19 @@ class CylMesh(BaseTensorMesh, InnerProducts, CylView):
                 #                                   kron3(speye(n[2]), av(n[1]), av(n[0]))), format="csr")
         return self._aveE2CC
 
+    @property
+    def aveE2CCV(self):
+        "Construct the averaging operator on cell edges to cell centers."
+        if getattr(self, '_aveE2CCV', None) is None:
+            # The number of cell centers in each direction
+            n = self.vnC
+            if self.isSymmetric:
+                self._aveE2CCV = sp.block_diag((sp.kron(av(n[1]), speye(n[0])),
+                                                sp.kron(speye(n[1]), av(n[0]))), format="csr")
+            else:
+                raise NotImplementedError('wrapping in the averaging is not yet implemented')
+        return self._aveE2CCV
+
 
     @property
     def aveF2CC(self):
@@ -295,6 +308,18 @@ class CylMesh(BaseTensorMesh, InnerProducts, CylView):
                 #                                    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.vnC
+            if self.isSymmetric:
+                self._aveF2CCV = sp.block_diag((sp.kron(speye(n[1]), av(n[0])),
+                                               sp.kron(av(n[1]), speye(n[0]))), format="csr")
+            else:
+                raise NotImplementedError('wrapping in the averaging is not yet implemented')
+        return self._aveF2CCV
+
 
 if __name__ == '__main__':