Refactor QAdder from SINGULAR integration routine to DEFAULT (#16)

* Fixed openMP issue qVS

* Refactor QAdder

* qstls refactoring

* Updated checkAdrFixed method to print which parameter is inconsistent

* Fixed formatting

* Added Vector3D method fill for specific index i

* Fixed formatting

* Fixed const

* Fixed tests

* Removed QSTLS bug fix

* QAdder float fix

* Removed unnecessary std

src/qvs.cpp

@@ -213,7 +213,7 @@ void QStlsCSR::computeAdr() {
 
 double QStlsCSR::getQAdder() const {
   Integrator1D itg1(ItgType::DEFAULT, in.getIntError());
-  Integrator2D itg2(ItgType::SINGULAR, ItgType::DEFAULT, in.getIntError());
+  Integrator2D itg2(ItgType::DEFAULT, ItgType::DEFAULT, in.getIntError());
   const bool segregatedItg = in.getInt2DScheme() == "segregated";
   const vector<double> itgGrid = (segregatedItg) ? wvg : vector<double>();
   const Interpolator1D ssfItp(wvg, ssf);
@@ -274,16 +274,16 @@ double QAdder::integrandNumerator1(const double q) const {
   const double w = itg2.getX();
   if (q == 0.0) { return 0.0; };
   double w2 = w * w;
-  double w3 = w2 * w;
+  double q2 = q * q;
   double logarg = (w + 2 * q) / (w - 2 * q);
   logarg = (logarg < 0.0) ? -logarg : logarg;
-  return q / (exp(q * q / Theta - mu) + 1.0) * q / w3 *
-         (q / w * log(logarg) - 1.0);
+  if (w == 0.0) { return 1.0 / (12.0 * (exp(q2 / Theta - mu) + 1.0)); };
+  return q2 / (exp(q2 / Theta - mu) + 1.0) * (q / w * log(logarg) - 1.0) / w2;
 }
 
 // Numerator integrand2
 double QAdder::integrandNumerator2(const double w) const {
-  return w * (ssf(w) - 1.0);
+  return (ssf(w) - 1.0);
 }
 
 // Denominator integral