Resolved inconsistencies in _density declaration, corrected minor error in test cases

Catch2 @ 255aa5f2

-Subproject commit 0e9bae1cdb5d7bcd967bcab10539805cc3bccb32
+Subproject commit 255aa5f2afe1a622c97422f65ace6ca915be0d8d

blaze @ 965fdd13

-Subproject commit 69e2fe554e72d90efc39ceec72d8cbddf88802a8
+Subproject commit 965fdd13a3e87c7c398054d39add0c14fe3779de

parmetis @ fba1fa6b

-Subproject commit e68c8e3947f9383e28cb8623d116ca996b3651a9
+Subproject commit fba1fa6b493d5030d2adfd4f04c98ff75bd5b8cf

spdlog @ cca004ef

-Subproject commit 83b9149930f392d7797b54fe97a66ab3f2120671
+Subproject commit cca004efe4e66136a5f9f37e007d28a23bb729e0

code/include/solvers/csdsnsolver.h

@@ -14,7 +14,7 @@ class CSDSNSolver : public SNSolver
     // Physics acess
     Vector _energies; /*! @brief: energy levels for CSD, lenght = _nEnergies */
     Vector _angle;    /*! @brief: angles for SN */
-    Vector _density;  /*! @brief: patient density for each grid cell */
+    std::vector<double> _density;  /*! @brief: patient density for each grid cell */
 
     std::vector<Matrix> _sigmaSE; /*!  @brief scattering cross section for all energies*/
     Vector _sigmaTE;              /*!  @brief total cross section for all energies*/

code/include/solvers/csdsnsolverfp.h

@@ -14,7 +14,7 @@ class CSDSNSolverFP : public SNSolver
     // Physics acess
     Vector _energies; /*! @brief: energy levels for CSD, lenght = _nEnergies */
     Vector _angle;    /*! @brief: angles for SN */
-    Vector _density;  /*! @brief: patient density for each grid cell */
+    std::vector<double> _density;  /*! @brief: patient density for each grid cell */
 
     std::vector<Matrix> _sigmaSE; /*!  @brief scattering cross section for all energies*/
     Vector _sigmaTE;              /*!  @brief total cross section for all energies*/

code/include/solvers/csdsnsolvernotrafo.h

@@ -14,7 +14,7 @@ class CSDSNSolverNoTrafo : public SNSolver
     // Physics acess
     Vector _energies; /*! @brief: energy levels for CSD, lenght = _nEnergies */
     Vector _angle;    /*! @brief: angles for SN */
-    Vector _density;  /*! @brief: patient density for each grid cell */
+    std::vector<double> _density;  /*! @brief: patient density for each grid cell */
 
     std::vector<Matrix> _sigmaSE; /*!  @brief scattering cross section for all energies*/
     Vector _sigmaTE;              /*!  @brief total cross section for all energies*/

code/include/solvers/csdsolvertrafofp.h

@@ -13,7 +13,8 @@ class CSDSolverTrafoFP : public SNSolver
 
     // Physics acess
     Vector _energies; /*! @brief: energy levels for CSD, lenght = _nEnergies */
-    Vector _density;  /*! @brief: patient density for each grid cell */
+    Vector _angle;    /*! @brief: angles for SN */
+    std::vector<double> _density;  /*! @brief: patient density for each grid cell */
 
     std::vector<Matrix> _sigmaSE; /*!  @brief scattering cross section for all energies*/
     Vector _sigmaTE;              /*!  @brief total cross section for all energies*/

code/src/problems/aircavity1d.cpp

@@ -25,6 +25,6 @@ std::vector<double> AirCavity1D::GetDensity( const VectorVector& cellMidPoints )
     std::vector<double> densities ( 4*cellMidPoints.size()/9, 1.0 );
     std::vector<double> air( 2*cellMidPoints.size()/9, 0.00125  );
     std::vector<double> water ( 3*cellMidPoints.size()/9, 1.0 );
-    densities.insert(densities.end(),air.begin(),air.end()-1);
-    densities.insert(densities.end(),water.begin(),water.end()-1);
+    densities.insert(densities.end(),air.begin(),air.end());
+    densities.insert(densities.end(),water.begin(),water.end());
     return densities; }

code/src/problems/musclebonelung.cpp

@@ -21,10 +21,12 @@ VectorVector MuscleBoneLung::SetupIC() {
     return psi;
 }
 
-std::vector<double> MuscleBoneLung::GetDensity( const VectorVector& cellMidPoints ) { 
-    std::vector<double> densities ( 1.5*cellMidPoints.size()/9, 1.05 ); //muscle layer
-    std::vector<double> bone( 1.5*cellMidPoints.size()/9, 1.92  );
-    std::vector<double> lung ( 6*cellMidPoints.size()/9, 0.26 );
-    densities.insert(densities.end(),bone.begin(),bone.end()-1);
-    densities.insert(densities.end(),lung.begin(),lung.end()-1);
+std::vector<double>  MuscleBoneLung::GetDensity( const VectorVector& cellMidPoints ) { 
+    std::cout<<"Length of mesh "<< cellMidPoints.size() << std::endl;
+    std::vector<double> densities ( 167, 1.05 ); //muscle layer
+    std::vector<double> bone( 167, 1.92  );
+    std::vector<double> lung ( 665, 0.26 );
+    densities.insert(densities.end(),bone.begin(),bone.end());
+    densities.insert(densities.end(),lung.begin(),lung.end());
+    std::cout<<"Length of densities "<< densities.size() << std::endl;
     return densities; }

code/src/problems/waterphantom.cpp

@@ -13,7 +13,6 @@ std::vector<VectorVector> WaterPhantom::GetExternalSource( const Vector& energie
 VectorVector WaterPhantom::SetupIC() {
     VectorVector psi( _mesh->GetNumCells(), Vector( _settings->GetNQuadPoints(), 1e-10 ) );
     auto cellMids = _mesh->GetCellMidPoints();
-    std::cout<<"Length of mesh "<< cellMids.size() << std::endl;
     double s      = 0.1;
     for( unsigned j = 0; j < cellMids.size(); ++j ) {
         double x                                = cellMids[j][0];

code/src/solvers/csdsnsolver.cpp

@@ -118,7 +118,7 @@ CSDSNSolver::CSDSNSolver( Config* settings ) : SNSolver( settings ) {
         // exit( EXIT_FAILURE );
     */
     // Get patient density
-    _density = Vector( _nCells, 1.0 );
+    _density = std::vector<double>( _nCells, 1.0 );
 }
 
 void CSDSNSolver::Solve() {

code/src/solvers/csdsnsolverfp.cpp

@@ -107,7 +107,8 @@ CSDSNSolverFP::CSDSNSolverFP( Config* settings ) : SNSolver( settings ) {
         _scatteringKernel( p, p ) = _weights[p];
     }
 
-    _density = Vector( _nCells, 1.0 );
+
+    _density = std::vector<double>( _nCells, 1.0 );
     //exit(EXIT_SUCCESS);
 }
 
@@ -118,6 +119,8 @@ void CSDSNSolverFP::Solve() {
     auto log      = spdlog::get( "event" );
     auto cellMids = _mesh->GetCellMidPoints();
 
+    _density = _problem->GetDensity(cellMids);
+
     // setup IC and incoming BC on left
     // auto cellMids = _settings->GetCellMidPoints();
     _sol = std::vector<Vector>( _nCells, Vector( _nq, 0.0 ) );

code/src/solvers/csdsnsolvernotrafo.cpp

@@ -86,7 +86,7 @@ CSDSNSolverNoTrafo::CSDSNSolverNoTrafo( Config* settings ) : SNSolver( settings
         _scatteringKernel( p, p ) = _weights[p];
     }
     // Get patient density
-    _density = Vector( _nCells, 1.0 );
+    _density = std::vector<double>( _nCells, 1.0 );
 }
 
 void CSDSNSolverNoTrafo::Solve() {

code/src/solvers/csdsolvertrafofp.cpp

@@ -98,6 +98,7 @@ CSDSolverTrafoFP::CSDSolverTrafoFP( Config* settings ) : SNSolver( settings ) {
                 }
                 std::cout<<"];"<<std::endl;*/
     }
+    
 
     // recompute scattering kernel. TODO: add this to kernel function
     for( unsigned p = 0; p < _nq; ++p ) {
@@ -107,8 +108,13 @@ CSDSolverTrafoFP::CSDSolverTrafoFP( Config* settings ) : SNSolver( settings ) {
         _scatteringKernel( p, p ) = _weights[p];
     }
 
+<<<<<<< HEAD
     _density = Vector( _nCells, 1.0 );
     // exit(EXIT_SUCCESS);
+=======
+    _density = std::vector<double> ( _nCells, 1.0 );
+    //exit(EXIT_SUCCESS);
+>>>>>>> Resolved inconsistencies in _density declaration, corrected minor error in test cases
 }
 
 void CSDSolverTrafoFP::Solve() {