From c4a5e26da2f1b15fe51420fa1678a74d251c7cde Mon Sep 17 00:00:00 2001
From: Amaury Johnan <amjohnen@gmail.com>
Date: Wed, 22 Oct 2014 14:05:53 +0000
Subject: [PATCH] fix for Coverty Scan warnings

---
 Numeric/MetricBasis.cpp | 47 +++++++++++++++++++++++++++++------------
 Numeric/MetricBasis.h   |  2 ++
 Numeric/bezierBasis.cpp | 23 ++++++++------------
 3 files changed, 44 insertions(+), 28 deletions(-)

diff --git a/Numeric/MetricBasis.cpp b/Numeric/MetricBasis.cpp
index 13b11cfbd9..5a66ba2638 100644
--- a/Numeric/MetricBasis.cpp
+++ b/Numeric/MetricBasis.cpp
@@ -8,7 +8,6 @@
 #include "MetricBasis.h"
 #include "BasisFactory.h"
 #include "pointsGenerators.h"
-#include "BasisFactory.h"
 #include "OS.h"
 #include <queue>
 #include <sstream>
@@ -55,8 +54,10 @@ namespace {
   }
 }
 
-MetricBasis::MetricBasis(int tag) : _type(ElementType::ParentTypeFromTag(tag)),
-    _dim(ElementType::DimensionFromTag(tag))
+MetricBasis::MetricBasis(int tag) :
+    _type(ElementType::ParentTypeFromTag(tag)),
+    _dim(ElementType::DimensionFromTag(tag)),
+    _jacobian(NULL)
 {
   const bool serendip = false;
   const int metOrder = metricOrder(tag);
@@ -76,25 +77,20 @@ MetricBasis::MetricBasis(int tag) : _type(ElementType::ParentTypeFromTag(tag)),
 
   // get jacobian
   FuncSpaceData *jacSpace = NULL;
-  if (_type == TYPE_TRI || _type == TYPE_QUA) {
-    jacSpace = NULL;
-  }
-  else if (_type == TYPE_TET || _type == TYPE_HEX || _type == TYPE_PRI) {
+  if (_type == TYPE_TET || _type == TYPE_HEX || _type == TYPE_PRI) {
     jacSpace = new FuncSpaceData(true, tag, jacOrder, &serendip);
   }
   else if (_type == TYPE_PYR) {
     jacSpace = new FuncSpaceData(true, tag, false, jacOrder+3,
                                  jacOrder, &serendip);
   }
-  else
+  else if (_type != TYPE_TRI || _type != TYPE_QUA)
     Msg::Fatal("metric not implemented for element tag %d", tag);
 
   if (jacSpace) {
     _jacobian = BasisFactory::getJacobianBasis(*jacSpace);
     delete jacSpace;
   }
-  else
-    _jacobian = NULL;
 
   _fillInequalities(metOrder);
 }
@@ -173,6 +169,10 @@ double MetricBasis::getMinSampledR(MElement *el, int deg) const
 
 double MetricBasis::minSampledRnew(MElement *el, int deg)
 {
+  static double tmEIG1 = 0;
+  static double tmEIG2 = 0;
+  static double tmOTH = 0;
+  double T = Cpu();
   int dim = el->getDim();
   if (dim < 2) return 1;
 
@@ -231,10 +231,19 @@ double MetricBasis::minSampledRnew(MElement *el, int deg)
     }
     fullVector<double> valReal(dim), valImag(dim);
     fullMatrix<double> vecLeft(dim, dim), vecRight(dim, dim);
+
+    tmOTH += Cpu()-T;
+    T = Cpu();
     metric.eig(valReal, valImag, vecLeft, vecRight, true);
+    tmEIG1 += Cpu()-T;
+    T = Cpu();
+    //metric.eigvalSym(valReal, true);
+    tmEIG2 += Cpu()-T;
+    T = Cpu();
     minR = std::min(minR, std::sqrt(valReal(0) / valReal(dim-1)));
   }
 
+  Msg::Warning("eig vs other %g %g %g", tmEIG1, tmEIG2, tmOTH);
   return minR;
 }
 
@@ -631,7 +640,7 @@ bool MetricBasis::validateBezierForMetricAndJacobian()
                                   dim > 2 ? nodes(i, 2) + symRand(range) : 0);
       }
       MElement *el = MElement::createElement(tag, vertices);
-      if (!el) {
+      if (el != NULL) {
         Msg::Error("MElement was unable to create element for tag %d", tag);
         ++numError;
       }
@@ -648,13 +657,23 @@ bool MetricBasis::validateBezierForMetricAndJacobian()
 
       /*{
         static int deg = 2;
+        BasisFactory::getCondNumBasis(el->getTypeForMSH(), deg);
         MetricBasis::minSampledR(el, deg);
         MetricBasis::minSampledRnew(el, deg);
 
         double time = Cpu();
         double minR1;
-        for (int cn = 0; cn < 100; ++cn)
-          minR1 = MetricBasis::minSampledR(el, deg);
+        for (int cn = 0; cn < 100; ++cn) {
+          const CondNumBasis *cnb =
+              BasisFactory::getCondNumBasis(el->getTypeForMSH(), deg);
+          fullMatrix<double> nodesXYZ(el->getNumVertices(), 3);
+          el->getNodesCoord(nodesXYZ);
+          fullVector<double> invCond(cnb->getNumCondNumNodes());
+          cnb->getInvCondNum(nodesXYZ, invCond);
+          minR1 = 1;
+          for (int i = 0; i < invCond.size(); ++i)
+            minR1 = std::min(minR1, invCond(i));
+        }
         double tm1 = Cpu()-time;
 
         time = Cpu();
@@ -805,7 +824,7 @@ void MetricBasis::interpolateAfterNSubdivisions(
   const int nSub = _bezier->getNumDivision();
   const int numCoeff = md->_metcoeffs->size2();
   const int numMetPnts = md->_metcoeffs->size1();
-  const int numJacPnts = md->_jaccoeffs ? md->_jaccoeffs->size() : 0;
+  const int numJacPnts = md->haveJac() ? md->_jaccoeffs->size() : 0;
   const int numNodPnts = nodes.size1();
 
   const bezierBasis *bezierMapping;
diff --git a/Numeric/MetricBasis.h b/Numeric/MetricBasis.h
index 0817458ace..635b7fa711 100644
--- a/Numeric/MetricBasis.h
+++ b/Numeric/MetricBasis.h
@@ -51,6 +51,8 @@ private:
       delete _metcoeffs;
       delete _jaccoeffs;
     }
+
+    bool haveJac() {return _jaccoeffs != NULL;}
   };
 
 public:
diff --git a/Numeric/bezierBasis.cpp b/Numeric/bezierBasis.cpp
index 2a56445b35..78df24305d 100644
--- a/Numeric/bezierBasis.cpp
+++ b/Numeric/bezierBasis.cpp
@@ -584,61 +584,56 @@ void bezierBasis::_construct()
       subPoints.push_back(gmshGeneratePointsLine(0));
       break;
     case TYPE_LIN : {
-      _numLagCoeff = order == 0 ? 1 : 2;
+      _numLagCoeff = order ? 2 : 1;
       _dimSimplex = 0;
       _exponents = gmshGenerateMonomialsLine(order);
       subPoints = generateSubPointsLine(order);
       break;
     }
     case TYPE_TRI : {
-      _numLagCoeff = order == 0 ? 1 : 3;
+      _numLagCoeff = order ? 3 : 1;
       _dimSimplex = 2;
       _exponents = gmshGenerateMonomialsTriangle(order);
       subPoints = generateSubPointsTriangle(order);
       break;
     }
     case TYPE_QUA : {
-      _numLagCoeff = order == 0 ? 1 : 4;
+      _numLagCoeff = order ? 4 : 1;
       _dimSimplex = 0;
       _exponents = gmshGenerateMonomialsQuadrangle(order);
       subPoints = generateSubPointsQuad(order);
       break;
     }
     case TYPE_TET : {
-      _numLagCoeff = order == 0 ? 1 : 4;
+      _numLagCoeff = order ? 4 : 1;
       _dimSimplex = 3;
       _exponents = gmshGenerateMonomialsTetrahedron(order);
       subPoints = generateSubPointsTetrahedron(order);
       break;
     }
     case TYPE_PRI : {
-      _numLagCoeff = order == 0 ? 1 : 6;
+      _numLagCoeff = order ? 6 : 1;
       _dimSimplex = 2;
       _exponents = gmshGenerateMonomialsPrism(order);
       subPoints = generateSubPointsPrism(order);
       break;
     }
     case TYPE_HEX : {
-      _numLagCoeff = order == 0 ? 1 : 8;
+      _numLagCoeff = order ? 8 : 1;
       _dimSimplex = 0;
       _exponents = gmshGenerateMonomialsHexahedron(order);
       subPoints = generateSubPointsHex(order);
       break;
     }
     default : {
-      Msg::Error("Unknown function space of parentType %d : "
-          "reverting to TET_1", _data.elementType());
-      _numLagCoeff = 4;
-      _dimSimplex = 3;
-      _exponents = gmshGenerateMonomialsTetrahedron(0);
-      subPoints = generateSubPointsTetrahedron(0);
-      break;
+      Msg::Fatal("Unknown function space for parentType %d", _data.elementType());
+      return;
     }
   }
   _numDivisions = static_cast<int>(subPoints.size());
 
   fullMatrix<double> bezierPoints = _exponents;
-  bezierPoints.scale(1./order);
+  if (order) bezierPoints.scale(1./order);
 
   matrixBez2Lag = generateBez2LagMatrix(_exponents, bezierPoints, order, _dimSimplex);
   matrixBez2Lag.invert(matrixLag2Bez);
-- 
GitLab