Changeset 20601

Timestamp:

05/10/16 15:45:24 (9 years ago)

Author:

Eric.Larour

Message:

CHG: introducing threshold in linear solver.

File:

• issm/trunk-jpl/src/c/solutionsequences/solutionsequence_linear.cpp (modified) (2 diffs)

issm/trunk-jpl/src/c/solutionsequences/solutionsequence_linear.cpp

@@ -19 +19 @@
         Vector<IssmDouble>*  ys  = NULL;
         int  configuration_type;
+        
+        /*solver convergence test: */
+        IssmDouble nKUF;
+        IssmDouble nF;
+        IssmDouble solver_residue;
+        Vector<IssmDouble>* KU=NULL;
+        Vector<IssmDouble>* KUF=NULL;
+
+
 
         /*Recover parameters: */
@@ -27 +36 @@
         Reduceloadx(pf, Kfs, ys); delete Kfs;
         Solverx(&uf, Kff, pf, NULL, df, femmodel->parameters); 
+        
+        /*Check that the solver converged nicely: */
+                
+        //compute KUF = KU - F = K*U - F
+        KU=uf->Duplicate(); Kff->MatMult(uf,KU);
+        KUF=KU->Duplicate(); KU->Copy(KUF); KUF->AYPX(pf,-1.0);
+
+        //compute norm(KUF), norm(F) and residue
+        nKUF=KUF->Norm(NORM_TWO);
+        nF=pf->Norm(NORM_TWO);
+        solver_residue=nKUF/nF;
+        if(solver_residue>pow(10,-6))_error_("   solver residue too high!: norm(KU-F)/norm(F)=" << solver_residue << "\n");
+
+        //clean up
+        delete KU; delete KUF;
         delete Kff; delete pf; delete df;
-//#ifdef  _HAVE_ADOLC_
-//        for (int i =0; i<uf->svector->M; ++i) {
-//          ADOLC_DUMP_MACRO(uf->svector->vector[i]);
-//        }
-//#endif
+
+        //#ifdef  _HAVE_ADOLC_
+        //        for (int i =0; i<uf->svector->M; ++i) {
+        //          ADOLC_DUMP_MACRO(uf->svector->vector[i]);
+        //        }
+        //#endif
+        
         Mergesolutionfromftogx(&ug, uf,ys,femmodel->nodes,femmodel->parameters);delete uf; delete ys;
         InputUpdateFromSolutionx(femmodel,ug);