Changeset e8916d6 in tol

Timestamp:

Nov 7, 2009, 4:24:13 PM (17 years ago)

Author:

Víctor de Buen <vdebuen@…>

Branches:

1941-docker-image-to-test, devel, master

Children:

4006126

Parents:

f7b99ac

Message:

Enhanced messages and error control in CleanInterruptions

git-svn-id: ​https://www.tol-project.org/svn/tolp/trunk@1622 be4846ae-f051-483a-831f-b64e7d08527f

File:

• tol/btol/bmodel/modcalc.cpp (modified) (2 diffs)

tol/btol/bmodel/modcalc.cpp

@@ -1428 +1428 @@
     CI = MinimumResidualsSolve(PIW_,-A_, CI);
     BBool CI_isKnown = FrobeniusNormU(CI.Data()).IsKnown();
-    BMat err = PIW_.T()*(PIW_*CI+A_);
-    BDat test = err.MaxAbsNorm(); 
-    if(!CI_isKnown ||(test>Sqrt(BDat::Tolerance())))
-    {
-    //Std(BText("\nBModel::CleanInterruptions test>Sqrt(BDat::Tolerance())"));
+    BMat CI_res = PIW_.T()*(PIW_*CI+A_);
+    BDat CI_maxAbsErr = CI_res.MaxAbsNorm(); 
+    BMat CI_U(N,I), CI_V(I,I);
+    BDiagMatrix<BDat> CI_D(I,1.0), CI_Dp;
+    if(!CI_isKnown ||(CI_maxAbsErr>Sqrt(BDat::Tolerance())))
+    {
+      Warning(I2("Cannot calculate values for interruptions due "
+                 "to a numerical error solving related linear programming "
+                 "problem by means of minimum residuals method.\n"
+                 "Singular value decomposition method will be tried."
+                 "This is a more robust but much slower method.",
+                 "No se pudieron calcular los valores de las interrupciones "
+                 "debido a un error de tipo numérico resolviendo el "
+                 "problema de programación lineal correspondiente "
+                 "mediante el método de mínimos residuos.\n"
+                 "Se intentará resolver mediante descomposición de valor "
+                 "singular que es más robusto pero mucho más lento.")+
+               "\n  Max(|PIW'(PIW*CI+A_)|)="<<CI_maxAbsErr);
       BText method = BText("Golub_Reinsch_Mod");
-      BMat U(N,I), V(I,I);
-      BDiagMatrix<BDat> D(I,1.0), Dp;
-      gsl_SingularValueDecomposition(PIW_,U,D,V, method);
-      Dp=D.P(Sqrt(DEpsilon()));
-      CI=V*Dp*(A_.T()*U).T();
-    }
-    CI_isKnown = FrobeniusNormU(CI.Data()).IsKnown();
-    if(CI_isKnown)
+      gsl_SingularValueDecomposition(PIW_,CI_U,CI_D,CI_V, method);
+      CI_Dp=CI_D.P(Sqrt(DEpsilon()));
+      CI=CI_V*CI_Dp*(A_.T()*CI_U).T();
+      CI_isKnown = FrobeniusNormU(CI.Data()).IsKnown();
+      CI_res = PIW_.T()*(PIW_*CI+A_);
+      CI_maxAbsErr = CI_res.MaxAbsNorm(); 
+    }
+    if(!CI_isKnown ||(CI_maxAbsErr>Sqrt(BDat::Tolerance())))
+    {
+      Error(I2("Cannot calculate values for interruptions due "
+               "to a numerical error solving related linear programming "
+               "problem by means of singular value decomposition method.",
+               "No se pudieron calcular los valores de las interrupciones "
+               "debido a un error de tipo numérico resolviendo el "
+               "problema de programación lineal correspondiente "
+               "mediante el método de descomposición de valor singular.")+
+            "\n  Max(|PIW'(PIW*CI+A_)|)="<<CI_maxAbsErr);
+      A_ = A*BDat::PosInf();
+    }
+    else
     {
       for(j=0; j<I; j++) { interruptionValue_[interruptionIndex_[j]] = CI(j,0); }
@@ -1460 +1485 @@
         }
       }
-    }
-    else
-    {
-      Warning(I2("Cannot calculate values for interruptions, probably due "
-                 "to a numerical problem in SVD of partial PI-weights matrix.",
-                 "No se pudieron calcular los valores de las interrupciones "
-                 "debido probablemente a un error de tipo numérico en la "
-                 "SVD de la matriz parcial de los PI-weights."));
-      A_ = A*BDat::PosInf();
     }
   }