Changeset 8033a6d in tol

Timestamp:

May 3, 2010, 6:43:55 AM (16 years ago)

Author:

Víctor de Buen <vdebuen@…>

Branches:

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

Children:

d5cd750

Parents:

cbeb848

Message:

refs #899

git-svn-id: ​https://www.tol-project.org/svn/tolp/trunk@2127 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 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);
+    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())"));
       BText method = BText("Golub_Reinsch_Mod");
-      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())))
+      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();
+      err = PIW_.T()*(PIW_*CI+A_);
+      test = err.MaxAbsNorm(); 
+    }
+    CI_isKnown = FrobeniusNormU(CI.Data()).IsKnown();
+    if(!CI_isKnown ||(test>Sqrt(BDat::Tolerance())))
     {
       Warning(I2("Cannot calculate values for interruptions due "
@@ -1463 +1452 @@
                "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
+            "\n  Max(|PIW'(PIW*CI+A_)|)="+test);
+    }
+    if(CI_isKnown)
     {
       for(j=0; j<I; j++) { interruptionValue_[interruptionIndex_[j]] = CI(j,0); }