From 90a7d5a4672d8b4883fce5722939a29439b059d2 Mon Sep 17 00:00:00 2001
From: raphael couturier <couturie@extinction>
Date: Wed, 13 Aug 2014 19:35:33 +0200
Subject: [PATCH] modif de LSQR qui fonctionne quand meme mieux (un poil moins
 rapide que CGLS quand meme)

---
 code/ex15.c | 10 +++++-----
 code/ex29.c | 11 ++---------
 code/ex45.c |  8 --------
 3 files changed, 7 insertions(+), 22 deletions(-)

diff --git a/code/ex15.c b/code/ex15.c
index 0e0a8c4..94c9365 100644
--- a/code/ex15.c
+++ b/code/ex15.c
@@ -1,5 +1,5 @@
 
-// /home/couturie/work/petsc-3.5.1/arch-linux2-c-debug/bin/mpirun -np 4    ./ex15 -m 400 -n 400
+// /home/couturie/work/petsc-3.5.1/arch-linux2-c-debug/bin/mpirun -np 4    ./ex15 -m 400 -n 400 -ksp_type fgmres
 
 
 
@@ -232,7 +232,7 @@ int KrylovMinimizeLSQR(Mat A, Vec b, Vec x) {
 
   PetscScalar  alpha, beta;
   PetscReal norm=20, Eprecision=1e-8, tol=1e-40;     
-  PetscInt giter=0, ColS=12, col=0, Emaxiter=50000000, iter=0, iterations=15, Iiter=0;
+  PetscInt giter=0, ColS=12, col=0, Emaxiter=50000000, iter=0, iterations=20, Iiter=0;
   PetscErrorCode ierr;
   PetscScalar T1, T2;
   KSP ksp;
@@ -419,14 +419,14 @@ int KrylovMinimizeLSQR(Mat A, Vec b, Vec x) {
         VecAYPX(d,-thet,v);       //d = (v - thet * d);
         VecAYPX(d,1/rho,zero_short);     //d=d/ rho;
 
-
+        /*
         if (normar/(norma*normr) <= tol) { // check for convergence in min{|b-A*x|}
           break;
         }
         if (normr <= tolb) {           // check for convergence in A*x=b
           break;
         }
-
+         */
 
         VecAXPY(x_lsqr,phi,d);     // x_lsqr=x_lsqr+phi*d
         normr = abs(s) * normr;
@@ -630,7 +630,7 @@ int main(int argc,char **args)
     VecDuplicate(b,&x2);
     VecDuplicate(b,&sol);
     
-    KrylovMinimize(A, b, x2);
+    KrylovMinimizeLSQR(A, b, x2);
 
 
 
diff --git a/code/ex29.c b/code/ex29.c
index 5e21c65..ede2355 100644
--- a/code/ex29.c
+++ b/code/ex29.c
@@ -1,6 +1,6 @@
 
 
-//  /home/couturie/work/petsc-3.5.1/arch-linux2-c-debug/bin/mpirun -np 3 ex29 -da_grid_x 900 -da_grid_y 900
+//  /home/couturie/work/petsc-3.5.1/arch-linux2-c-debug/bin/mpirun -np 3 ex29 -da_grid_x 900 -da_grid_y 900 -ksp_type fgmres
 
 
 
@@ -240,7 +240,7 @@ int KrylovMinimizeLSQR(Mat A, Vec b, Vec x) {
 
   PetscScalar  alpha, beta;
   PetscReal norm=20, Eprecision=1e-5, tol=1e-40;     
-  PetscInt giter=0, ColS=12, col=0, Emaxiter=50000000, iter=0, iterations=15, Iiter=0;
+  PetscInt giter=0, ColS=12, col=0, Emaxiter=50000000, iter=0, iterations=20, Iiter=0;
   PetscErrorCode ierr;
   PetscScalar T1, T2;
   KSP ksp;
@@ -428,13 +428,6 @@ int KrylovMinimizeLSQR(Mat A, Vec b, Vec x) {
         VecAYPX(d,1/rho,zero_short);     //d=d/ rho;
 
 
-        if (normar/(norma*normr) <= tol) { // check for convergence in min{|b-A*x|}
-          break;
-        }
-        if (normr <= tolb) {           // check for convergence in A*x=b
-          break;
-        }
-
 
         VecAXPY(x_lsqr,phi,d);     // x_lsqr=x_lsqr+phi*d
         normr = abs(s) * normr;
diff --git a/code/ex45.c b/code/ex45.c
index a44a2d0..11fd5c2 100644
--- a/code/ex45.c
+++ b/code/ex45.c
@@ -408,14 +408,6 @@ int KrylovMinimizeLSQR(Mat A, Vec b, Vec x) {
         VecAYPX(d,1/rho,zero_short);     //d=d/ rho;
 
 
-        if (normar/(norma*normr) <= tol) { // check for convergence in min{|b-A*x|}
-          break;
-        }
-        if (normr <= tolb) {           // check for convergence in A*x=b
-          break;
-        }
-
-
         VecAXPY(x_lsqr,phi,d);     // x_lsqr=x_lsqr+phi*d
         normr = abs(s) * normr;
         MatMultTranspose(AS, u, vt); //vt=A'*u;
-- 
2.39.5