[cig-commits] r8226 - mc/3D/CitcomS/trunk/lib

[tan2 at geodynamics.org]

Author: tan2
Date: 2007-11-07 11:51:56 -0800 (Wed, 07 Nov 2007)
New Revision: 8226

Modified:
   mc/3D/CitcomS/trunk/lib/Stokes_flow_Incomp.c
Log:
Print a warning when matrix eqn solver not converging

Modified: mc/3D/CitcomS/trunk/lib/Stokes_flow_Incomp.c
===================================================================
--- mc/3D/CitcomS/trunk/lib/Stokes_flow_Incomp.c	2007-11-07 19:48:27 UTC (rev 8225)
+++ mc/3D/CitcomS/trunk/lib/Stokes_flow_Incomp.c	2007-11-07 19:51:56 UTC (rev 8226)
@@ -241,10 +241,9 @@
     }
 
 
-    valid = 1;
     r0dotz0 = 0;
 
-    while( (valid) && (count < *steps_max) &&
+    while( (count < *steps_max) &&
            (E->monitor.incompressibility >= E->control.tole_comp) &&
            (dpressure >= imp) && (dvelocity >= imp) )  {
 
@@ -277,6 +276,10 @@
         /* solve K*u1 = grad(s2) for u1 */
         assemble_grad_p(E, s2, F, lev);
         valid = solve_del2_u(E, E->u1, F, imp*v_res, lev);
+        if(!valid && (E->parallel.me==0)) {
+            fputs("Warning: solver not converging! 1\n", stderr);
+            fputs("Warning: solver not converging! 1\n", E->fp);
+        }
         strip_bcs_from_residual(E, E->u1, lev);
 
 
@@ -458,7 +461,7 @@
     valid = 1;
     r0dotrt = alpha = omega = 0;
 
-    while( (valid) && (count < *steps_max) &&
+    while( (count < *steps_max) &&
            ((E->monitor.incompressibility >= E->control.tole_comp) &&
             (dpressure >= imp) && (dvelocity >= imp)) )  {
 
@@ -498,7 +501,10 @@
         /* solve K*u0 = grad(pt) for u1 */
         assemble_grad_p(E, pt, F, lev);
         valid = solve_del2_u(E, u0, F, imp*v_res, lev);
-        if(!valid) fprintf(stderr, "not valid 1\n");
+        if(!valid && (E->parallel.me==0)) {
+            fputs("Warning: solver not converging! 1\n", stderr);
+            fputs("Warning: solver not converging! 1\n", E->fp);
+        }
         strip_bcs_from_residual(E, u0, lev);
 
 
@@ -525,7 +531,10 @@
         /* solve K*u1 = grad(st) for u1 */
         assemble_grad_p(E, st, F, lev);
         valid = solve_del2_u(E, E->u1, F, imp*v_res, lev);
-        if(!valid) fprintf(stderr, "not valid 2\n");
+        if(!valid && (E->parallel.me==0)) {
+            fputs("Warning: solver not converging! 2\n", stderr);
+            fputs("Warning: solver not converging! 2\n", E->fp);
+        }
         strip_bcs_from_residual(E, E->u1, lev);
 
 
@@ -714,7 +723,7 @@
     int neq = E->lmesh.neq;
     int gneq = E->mesh.neq;
     int lev = E->mesh.levmax;
-    int i, m;
+    int i, m, valid;
     double v_res;
 
     v_res = sqrt(global_vdot(E, F, F, lev) / gneq);
@@ -742,7 +751,11 @@
 
 
     /* solve K*u1 = F for u1 */
-    solve_del2_u(E, E->u1, F, imp*v_res, lev);
+    valid = solve_del2_u(E, E->u1, F, imp*v_res, lev);
+    if(!valid && (E->parallel.me==0)) {
+        fputs("Warning: solver not converging! 0\n", stderr);
+        fputs("Warning: solver not converging! 0\n", E->fp);
+    }
     strip_bcs_from_residual(E, E->u1, lev);