move constants to nmie::

src/mesomie.hpp

@@ -0,0 +1,60 @@
+#ifndef SRC_NMIE_BASIC_HPP_
+#define SRC_NMIE_BASIC_HPP_
+//**********************************************************************************//
+//    Copyright (C) 2009-2018  Ovidio Pena <ovidio@bytesfall.com>                   //
+//    Copyright (C) 2013-2018  Konstantin Ladutenko <kostyfisik@gmail.com>          //
+//                                                                                  //
+//    This file is part of scattnlay                                                //
+//                                                                                  //
+//    This program is free software: you can redistribute it and/or modify          //
+//    it under the terms of the GNU General Public License as published by          //
+//    the Free Software Foundation, either version 3 of the License, or             //
+//    (at your option) any later version.                                           //
+//                                                                                  //
+//    This program is distributed in the hope that it will be useful,               //
+//    but WITHOUT ANY WARRANTY; without even the implied warranty of                //
+//    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the                 //
+//    GNU General Public License for more details.                                  //
+//                                                                                  //
+//    The only additional remark is that we expect that all publications            //
+//    describing work using this software, or all commercial products               //
+//    using it, cite at least one of the following references:                      //
+//    [1] O. Pena and U. Pal, "Scattering of electromagnetic radiation by           //
+//        a multilayered sphere," Computer Physics Communications,                  //
+//        vol. 180, Nov. 2009, pp. 2348-2354.                                       //
+//    [2] K. Ladutenko, U. Pal, A. Rivera, and O. Pena-Rodriguez, "Mie              //
+//        calculation of electromagnetic near-field for a multilayered              //
+//        sphere," Computer Physics Communications, vol. 214, May 2017,             //
+//        pp. 225-230.                                                              //
+//                                                                                  //
+//    You should have received a copy of the GNU General Public License             //
+//    along with this program.  If not, see <http://www.gnu.org/licenses/>.         //
+//**********************************************************************************//
+
+//**********************************************************************************//
+// This class implements the algorithm for a multilayered sphere described by:      //
+//    [1] W. Yang, "Improved recursive algorithm for light scattering by a          //
+//        multilayered sphere,” Applied Optics, vol. 42, Mar. 2003, pp. 1710-1720.  //
+//                                                                                  //
+// You can find the description of all the used equations in:                       //
+//    [2] O. Pena and U. Pal, "Scattering of electromagnetic radiation by           //
+//        a multilayered sphere," Computer Physics Communications,                  //
+//        vol. 180, Nov. 2009, pp. 2348-2354.                                       //
+//    [3] K. Ladutenko, U. Pal, A. Rivera, and O. Pena-Rodriguez, "Mie              //
+//        calculation of electromagnetic near-field for a multilayered              //
+//        sphere," Computer Physics Communications, vol. 214, May 2017,             //
+//        pp. 225-230.                                                              //
+//                                                                                  //
+// Hereinafter all equations numbers refer to [2]                                   //
+//**********************************************************************************//
+#include <iostream>
+#include <iomanip>
+#include <stdexcept>
+#include <vector>
+
+#include "special-functions-impl.hpp"
+#include "nmie.hpp"
+
+namespace nmie {
+}
+

src/nmie-applied-impl.hpp

@@ -178,11 +178,11 @@ namespace nmie {
     FloatType radius = 0.0;
     for (auto width : target_width_) {
       radius += width;
-      this->size_param_.push_back(2*this->PI_*radius/wavelength_);
+      this->size_param_.push_back(2*nmie::PI_*radius/wavelength_);
     }
     for (auto width : coating_width_) {
       radius += width;
-      this->size_param_.push_back(2*this->PI_*radius/wavelength_);
+      this->size_param_.push_back(2*nmie::PI_*radius/wavelength_);
     }
     this->total_radius_ = radius;
   }  // end of void MultiLayerMieApplied<FloatType>::GenerateSizeParameter();
@@ -306,7 +306,7 @@ namespace nmie {
                                                                  const int isIgnoreAvailableNmax) {
     ConvertToSP(); // Converts to size parameter units only the particle design,
     // so we need to convert input parameters too...
-    const FloatType a = 2*this->PI_/wavelength_;
+    const FloatType a = 2*nmie::PI_/wavelength_;
     this->MultiLayerMie<FloatType>::RunFieldCalculationPolar(outer_arc_points, radius_points, a*from_Rho, a*to_Rho,
                                                                from_Theta, to_Theta, from_Phi, to_Phi,
                                                                isIgnoreAvailableNmax == 0 ? false : true);

src/nmie-nearfield.hpp

@@ -394,7 +394,7 @@ namespace nmie {
     }
 
     // magnetic field
-    std::complex<evalType> hffact = ml/static_cast<evalType>(cc_*mu_);
+    std::complex<evalType> hffact = ml/static_cast<evalType>(nmie::cc_*nmie::mu_);
     for (int i = 0; i < 3; i++) {
       H[i] = hffact*H[i];
     }

src/special-functions-impl.hpp

@@ -94,6 +94,8 @@ int getNStar(int nmax, std::complex<FloatType> z, const int valid_digits) {
 
 // Custom implementation of complex cot function to avoid overflow
 // if Im(z) < 0, then it evaluates cot(z) as conj(cot(conj(z)))
+// see Eqs. 10-12 of [1] for details.
+// [1]H. Du, Mie-Scattering Calculation, Appl. Opt. 43, 1951 (2004).
 template <typename FloatType>
 std::complex<FloatType> complex_cot(const std::complex<FloatType> z) {
   auto Remx = z.real();

tests/mpmath_special_functions_test_generator.py

@@ -204,7 +204,7 @@ def main():
                                                  # output_dps=7, max_num_elements_of_nlist=51)
                                                  # output_dps=5, max_num_elements_of_nlist=3)
     # sf_evals.run_test(mrb.D1, 'D1')
-    # sf_evals.run_test(mrb.D2, 'D2')
+    #
     # sf_evals.run_test(mrb.D3, 'D3')
     # sf_evals.run_test(mrb.psi, 'psi', is_only_x=True)
     # sf_evals.run_test(mrb.xi, 'xi', is_only_x=True)

tests/test_bulk_sphere.cc

@@ -12,15 +12,15 @@ TEST(BulkSphere, ArgPi) {
   nmie::MultiLayerMie<nmie::FloatType> nmie;
   nmie.SetLayersIndex({std::complex<double>(4,0)});
   for (auto WL:WLs) {
-    nmie.SetLayersSize({2*nmie.PI_*host_index*core_radius/(WL+delta)});
+    nmie.SetLayersSize({2*nmie::PI_*host_index*core_radius/(WL+delta)});
     nmie.RunMieCalculation();
     double Qabs_p = std::abs(static_cast<double>(nmie.GetQabs()));
 
-    nmie.SetLayersSize({2*nmie.PI_*host_index*core_radius/(WL-delta)});
+    nmie.SetLayersSize({2*nmie::PI_*host_index*core_radius/(WL-delta)});
     nmie.RunMieCalculation();
     double Qabs_m = std::abs(static_cast<double>(nmie.GetQabs()));
 
-    nmie.SetLayersSize({2*nmie.PI_*host_index*core_radius/(WL)});
+    nmie.SetLayersSize({2*nmie::PI_*host_index*core_radius/(WL)});
     nmie.RunMieCalculation();
     double Qabs = std::abs(static_cast<double>(nmie.GetQabs()));
     EXPECT_GT(Qabs_p+Qabs_m, Qabs);

tests/test_near_field.cc

@@ -8,7 +8,7 @@ TEST(RunFieldCalculationCartesian, HandlesInput) {
   nmie::MultiLayerMie<nmie::FloatType> nmie;
 //  EXPECT_THROW(nmie.RunFieldCalculationPolar(0), std::invalid_argument);
 //  EXPECT_THROW(nmie.RunFieldCalculationPolar(1,1,10,5), std::invalid_argument);
-  nmie::FloatType total_r = 2*nmie.PI_*1000/532;
+  nmie::FloatType total_r = 2*nmie::PI_*1000/532;
 //  double r = 1500;
   nmie.SetLayersSize({total_r/2, total_r});
   nmie.SetLayersIndex({ {1.330,0}, {1.33,0}});
@@ -58,8 +58,8 @@ TEST(RunFieldCalculationCartesian, HandlesInput) {
 TEST(LargeBubbleSpectrum, DISABLED_HandlesInput) {
 //TEST(LargeBubbleSpectrum, HandlesInput) { // TODO fix fail...
   nmie::MultiLayerMie<nmie::FloatType> nmie;
-  nmie::FloatType core_r  = 2*nmie.PI_*100;
-  nmie::FloatType shell_r = 2*nmie.PI_*(100+0.1);
+  nmie::FloatType core_r  = 2*nmie::PI_*100;
+  nmie::FloatType shell_r = 2*nmie::PI_*(100+0.1);
   nmie.SetLayersIndex({ {1,0}, {1.33,0}});
   double central_WL = 0.7007;
   double relative_distance = 1e-10;
@@ -121,12 +121,12 @@ TEST(BulkSphere, HandlesInput) {
     nmie.SetMaxTerms(-1);
 //    nmie.RunMieCalculation();
 //    std::cout<<" test case: "<<std::get<2>(data)<<" Qsca="<<nmie.GetQsca()<<std::endl;
-    nmie.RunFieldCalculationPolar(4,3,x,x*3, 0, static_cast<double>(nmie.PI_), 0, static_cast<double>(nmie.PI_),true, -1);
+    nmie.RunFieldCalculationPolar(4,3,x,x*3, 0, static_cast<double>(nmie::PI_), 0, static_cast<double>(nmie::PI_),true, -1);
     auto Eabs = nmie.GetFieldEabs();
     for (auto &E:Eabs) E=nmie::pow2(E);
 //    print(Eabs)
     EXPECT_TRUE(nmie.GetFieldConvergence())<<"Outside test for x="<<x<<" m="<<m<<" test case: "<<std::get<2>(data)<<std::endl;
-    nmie.RunFieldCalculationPolar(4,10,x*0.01,x, 0, static_cast<double>(nmie.PI_), 0, static_cast<double>(nmie.PI_),true, -1);
+    nmie.RunFieldCalculationPolar(4,10,x*0.01,x, 0, static_cast<double>(nmie::PI_), 0, static_cast<double>(nmie::PI_),true, -1);
     EXPECT_TRUE(nmie.GetFieldConvergence())<<"Inside test for x="<<x<<" m="<<m<<" test case: "<<std::get<2>(data)<<std::endl;
   }
 }