From f86ebe12dbdcd83aa217cc2b1c61939c57363e7b Mon Sep 17 00:00:00 2001
From: Rowan Cockett <rowan@visiblegeology.com>
Date: Mon, 22 Jul 2013 13:40:18 -0700
Subject: [PATCH] Added Poisson Equation Tests.

Forwards and Backwards.
---
 SimPEG/tests/OrderTest.py           |  5 ++--
 SimPEG/tests/test_operatorOrders.py | 37 ++++++++++++++++++++++++++---
 2 files changed, 37 insertions(+), 5 deletions(-)

diff --git a/SimPEG/tests/OrderTest.py b/SimPEG/tests/OrderTest.py
index 68837f7a..652b224a 100644
--- a/SimPEG/tests/OrderTest.py
+++ b/SimPEG/tests/OrderTest.py
@@ -34,15 +34,16 @@ class OrderTest(unittest.TestCase):
             if ii == 0:
                 print ''
                 print 'Testing order of:  ' + self.name
+                print '__________________________________________'
                 print '   h  |   inf norm  |  ratio   |  order'
-                print '------------------------------------------'
+                print '~~~~~~|~~~~~~~~~~~~~|~~~~~~~~~~|~~~~~~~~~~'
                 print '%4i  |  %8.2e   |' % (nc, err)
             else:
                 order.append(np.log(err/err_old)/np.log(float(nc_old)/float(nc)))
                 print '%4i  |  %8.2e   |  %6.4f  |  %6.4f' % (nc, err, err_old/err, order[-1])
             err_old = err
             nc_old = nc
-
+        print '------------------------------------------'
         self.assertTrue(len(np.where(np.array(order) > 0.9*self.expectedOrder)[0]) > np.floor(0.75*len(order)))
 
 
diff --git a/SimPEG/tests/test_operatorOrders.py b/SimPEG/tests/test_operatorOrders.py
index 9bb9fa03..9fdebc10 100644
--- a/SimPEG/tests/test_operatorOrders.py
+++ b/SimPEG/tests/test_operatorOrders.py
@@ -1,10 +1,10 @@
 import numpy as np
 from OrderTest import OrderTest
 import unittest
+from scipy.sparse.linalg import dsolve
 
 
 class TestCurl(OrderTest):
-
     name = "Curl"
 
     def getError(self):
@@ -33,7 +33,6 @@ class TestCurl(OrderTest):
 
 
 class TestFaceDiv(OrderTest):
-
     name = "Face Divergence"
 
     def getError(self):
@@ -59,7 +58,6 @@ class TestFaceDiv(OrderTest):
 
 
 class TestNodalGrad(OrderTest):
-
     name = "Nodal Gradient"
 
     def getError(self):
@@ -84,5 +82,38 @@ class TestNodalGrad(OrderTest):
         self.orderTest()
 
 
+class TestPoissonEqn(OrderTest):
+    name = "Poisson Equation"
+    meshSizes = [16, 20, 24]
+
+    def getError(self):
+        # Create some functions to integrate
+        fun = lambda x: np.sin(2*np.pi*x[:, 0])*np.sin(2*np.pi*x[:, 1])*np.sin(2*np.pi*x[:, 2])
+        sol = lambda x: -3.*((2*np.pi)**2)*fun(x)
+
+        self.M.setCellGradBC('dirichlet')
+
+        D = self.M.faceDiv
+        G = self.M.cellGrad
+        if self.forward:
+            sA = sol(self.M.gridCC)
+            sN = D*G*fun(self.M.gridCC)
+            err = np.linalg.norm((sA - sN), np.inf)
+        else:
+            fA = fun(self.M.gridCC)
+            fN = dsolve.spsolve(D*G, sol(self.M.gridCC))
+            err = np.linalg.norm((fA - fN), np.inf)
+        return err
+
+    def test_orderForward(self):
+        self.name = "Poisson Equation - Forward"
+        self.forward = True
+        self.orderTest()
+
+    def test_orderBackward(self):
+        self.name = "Poisson Equation - Backward"
+        self.forward = False
+        self.orderTest()
+
 if __name__ == '__main__':
     unittest.main()