Added option to initialise solution

SpMat.h

@@ -137,7 +137,20 @@ public:
                                  double                                tol = 1e-4, 
                                  unsigned int                          miter = 200,
 				 boost::shared_ptr<Preconditioner<T> > C = boost::shared_ptr<Preconditioner<T> >()) const;
-   
+
+NEWMAT::ReturnMatrix SolveForx(const NEWMAT::ColumnVector&            b,
+			                 MatrixType                             type,
+			                 double                                 tol,
+			                 unsigned int                           miter,
+					 const NEWMAT::ColumnVector&            x_init) const;
+
+NEWMAT::ReturnMatrix SolveForx(const NEWMAT::ColumnVector&            b,
+			                 MatrixType                             type,
+			                 double                                 tol,
+			                 unsigned int                           miter,
+					 boost::shared_ptr<Preconditioner<T> >  C,
+					 const NEWMAT::ColumnVector&           x_init) const;
+
 private:
   unsigned int                              _m;
   unsigned int                              _n;
@@ -368,18 +381,48 @@ NEWMAT::ReturnMatrix SpMat<T>::AsNEWMAT() const
 //
 /////////////////////////////////////////////////////////////////////
 
+template<class T>
+NEWMAT::ReturnMatrix SpMat<T>::SolveForx(const NEWMAT::ColumnVector&            b,
+			                 MatrixType                             type,
+			                 double                                 tol,
+			                 unsigned int                           miter,
+					 const NEWMAT::ColumnVector&            x_init) const
+{
+  return this->SolveForx(b,type,tol,miter,boost::shared_ptr<Preconditioner<T> >(),x_init);
+}
+
 template<class T>
 NEWMAT::ReturnMatrix SpMat<T>::SolveForx(const NEWMAT::ColumnVector&            b,
 			                 MatrixType                             type,
 			                 double                                 tol,
 			                 unsigned int                           miter,
 					 boost::shared_ptr<Preconditioner<T> >  C) const
+{
+  NEWMAT::ColumnVector x_init;
+  return this->SolveForx(b,type,tol,miter,C,x_init);
+}
+
+template<class T>
+NEWMAT::ReturnMatrix SpMat<T>::SolveForx(const NEWMAT::ColumnVector&            b,
+			                 MatrixType                             type,
+			                 double                                 tol,
+			                 unsigned int                           miter,
+					 boost::shared_ptr<Preconditioner<T> >  C,
+					 const NEWMAT::ColumnVector&           x_init) const
 {
   if (_m != _n) throw SpMatException("SolveForx: Matrix must be square");
   if (int(_m) != b.Nrows()) throw SpMatException("SolveForx: Mismatch between matrix and vector");
 
   NEWMAT::ColumnVector  x(_n);
-  x = 0.0;
+  if (x.Nrows() == x_init.Nrows()) {
+    x = x_init;
+  } else {
+    if (x_init.Nrows()>0) {
+      throw SpMatException("SolveForx: initialisation vector has incorrect size");
+    } else {
+      x = 0.0;
+    }
+  }
   int                   status = 0;
   int                   liter = int(miter);
   double                ltol = tol;