From 3bbe258eabab0d5f6022dd065c7841407d9a0fed Mon Sep 17 00:00:00 2001
From: Rowan Cockett <rowan@visiblegeology.com>
Date: Fri, 25 Oct 2013 11:14:21 -0700
Subject: [PATCH] beta schedualing

---
 SimPEG/inverse/BetaSchedual.py | 12 ++++++++++++
 SimPEG/inverse/Inversion.py    | 12 ++++++++++--
 SimPEG/inverse/Optimize.py     |  1 +
 SimPEG/inverse/__init__.py     |  1 +
 4 files changed, 24 insertions(+), 2 deletions(-)
 create mode 100644 SimPEG/inverse/BetaSchedual.py

diff --git a/SimPEG/inverse/BetaSchedual.py b/SimPEG/inverse/BetaSchedual.py
new file mode 100644
index 00000000..62fe8416
--- /dev/null
+++ b/SimPEG/inverse/BetaSchedual.py
@@ -0,0 +1,12 @@
+
+
+class Cooling(object):
+    """Simple Beta Schedual"""
+
+    beta0 = 1.e6
+    beta_coolingFactor = 5.
+
+    def getBeta(self):
+        if self._beta is None:
+            return beta0
+        return self._beta * beta_coolingFactor
diff --git a/SimPEG/inverse/Inversion.py b/SimPEG/inverse/Inversion.py
index cf96676f..b30b666e 100644
--- a/SimPEG/inverse/Inversion.py
+++ b/SimPEG/inverse/Inversion.py
@@ -1,8 +1,9 @@
 import numpy as np
 import scipy.sparse as sp
 from SimPEG.utils import sdiag, mkvc
+from SimPEG.inverse import BS
 
-class Inversion(object):
+class Inversion(object, BetaSchedual.Cooling):
     """docstring for Inversion"""
 
     maxIter = 10
@@ -42,6 +43,7 @@ class Inversion(object):
     def run(self, m0):
         m = m0
         self._iter = 0
+        self._beta = None
         while True:
             self._beta = self.getBeta()
             m = self.opt.minimize(self.evalFunction,m)
@@ -49,8 +51,13 @@ class Inversion(object):
             self._iter += 1
         return m
 
+    beta0 = 1.e6
+    beta_coolingFactor = 5.
+
     def getBeta(self):
-        return 1e-2
+        if self._beta is None:
+            return beta0
+        return self._beta * beta_coolingFactor
 
     def stoppingCriteria(self):
         self._STOP = np.zeros(2,dtype=bool)
@@ -191,3 +198,4 @@ class Inversion(object):
         dmisfit = self.prob.Jt(m, self.Wd * self.Wd * self.prob.J(m, v, u=u), u=u)
 
         return dmisfit
+
diff --git a/SimPEG/inverse/Optimize.py b/SimPEG/inverse/Optimize.py
index 4ac96d6a..ac0045c4 100644
--- a/SimPEG/inverse/Optimize.py
+++ b/SimPEG/inverse/Optimize.py
@@ -158,6 +158,7 @@ class GaussNewton(Minimize):
 class InexactGaussNewton(Minimize):
     name = 'InexactGaussNewton'
     def findSearchDirection(self):
+        # TODO: use BFGS as a preconditioner or gauss sidel of the WtW or solve WtW directly
         p, info = sp.linalg.cg(self.H, -self.g, tol=1e-05, maxiter=5)
         return p
 
diff --git a/SimPEG/inverse/__init__.py b/SimPEG/inverse/__init__.py
index c4ca335b..2ca95fb5 100644
--- a/SimPEG/inverse/__init__.py
+++ b/SimPEG/inverse/__init__.py
@@ -1,2 +1,3 @@
 from Optimize import *
 from Inversion import *
+import BetaSchedual