Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2026-07-21 07:51:35

0001 // This file is part of the ACTS project.
0002 //
0003 // Copyright (C) 2016 CERN for the benefit of the ACTS project
0004 //
0005 // This Source Code Form is subject to the terms of the Mozilla Public
0006 // License, v. 2.0. If a copy of the MPL was not distributed with this
0007 // file, You can obtain one at https://mozilla.org/MPL/2.0/.
0008 
0009 #include <boost/test/unit_test.hpp>
0010 
0011 #include "Acts/Material/GridSurfaceMaterial.hpp"
0012 #include "Acts/Material/GridSurfaceMaterialFactory.hpp"
0013 #include "Acts/Material/Material.hpp"
0014 #include "Acts/Material/MaterialSlab.hpp"
0015 #include "Acts/Utilities/GridAxisGenerators.hpp"
0016 #include "Acts/Utilities/ProtoAxis.hpp"
0017 #include "Acts/Utilities/VectorHelpers.hpp"
0018 
0019 #include <numbers>
0020 #include <vector>
0021 
0022 using namespace Acts;
0023 
0024 // this is a global access to the x coordinate
0025 class GlobalAccessX final : public GridAccess::IGlobalToGridLocal {
0026  public:
0027   std::array<double, 1u> g2X(const Vector3& global) const {
0028     return {global.x()};
0029   }
0030 };
0031 
0032 class LocalAccessX final : public GridAccess::IBoundToGridLocal {
0033  public:
0034   std::array<double, 1u> l2X(const Vector2& local) const { return {local.x()}; }
0035 };
0036 
0037 class GlobalAccessPhi final : public GridAccess::IGlobalToGridLocal {
0038  public:
0039   std::array<double, 1u> g2Phi(const Vector3& global) const {
0040     return {std::atan2(global.y(), global.x())};
0041   }
0042 };
0043 
0044 class LocalAccessPhi final : public GridAccess::IBoundToGridLocal {
0045  public:
0046   std::array<double, 1u> l2Phi(const Vector2& local) const {
0047     return {std::atan2(local.y(), local.x())};
0048   }
0049 };
0050 
0051 class GlobalAccessXY final : public GridAccess::IGlobalToGridLocal {
0052  public:
0053   std::array<double, 2u> g2XY(const Vector3& global) const {
0054     return {global.x(), global.y()};
0055   }
0056 };
0057 
0058 class LocalAccessXY final : public GridAccess::IBoundToGridLocal {
0059  public:
0060   std::array<double, 2u> l2XY(const Vector2& local) const {
0061     return {local.x(), local.y()};
0062   }
0063 };
0064 
0065 class GlobalToZPhi final : public GridAccess::IGlobalToGridLocal {
0066  public:
0067   double zShift = 0.;
0068 
0069   explicit GlobalToZPhi(double shift) : zShift(shift) {}
0070 
0071   std::array<double, 2u> g2ZPhi(const Vector3& global) const {
0072     return {global.z() + zShift, VectorHelpers::phi(global)};
0073   }
0074 };
0075 
0076 // Local on cylinder surface is rPhi, z
0077 class LocalToZPhi final : public GridAccess::IBoundToGridLocal {
0078  public:
0079   double radius = 1.;
0080 
0081   explicit LocalToZPhi(double r) : radius(r) {}
0082 
0083   std::array<double, 2u> l2ZPhi(const Vector2& local) const {
0084     return {local[1u], local[0u] / radius};
0085   }
0086 };
0087 
0088 namespace ActsTests {
0089 
0090 BOOST_AUTO_TEST_SUITE(MaterialSuite)
0091 
0092 // This test covers some wrongly configured cases
0093 BOOST_AUTO_TEST_CASE(GridIndexedMaterial_invalid_bound2Grid_Unconnected) {
0094   std::vector<MaterialSlab> material;
0095 
0096   using EqBound = GridAxisGenerators::EqBound;
0097   using EqGrid = EqBound::grid_type<std::size_t>;
0098 
0099   EqBound eqBound{{0., 5.}, 5};
0100   EqGrid eqGrid{eqBound()};
0101 
0102   IndexedSurfaceMaterial<EqGrid>::BoundToGridLocalDelegate bToX;
0103 
0104   auto globalX = std::make_unique<const GlobalAccessX>();
0105   IndexedSurfaceMaterial<EqGrid>::GlobalToGridLocalDelegate gToX;
0106   gToX.connect<&GlobalAccessX::g2X>(std::move(globalX));
0107 
0108   BOOST_CHECK_THROW(
0109       auto ism = IndexedSurfaceMaterial<EqGrid>(
0110           std::move(eqGrid), IndexedMaterialAccessor{std::move(material)},
0111           std::move(bToX), std::move(gToX)),
0112       std::invalid_argument);
0113 }
0114 
0115 // This test covers some wrongly configured cases
0116 BOOST_AUTO_TEST_CASE(GridIndexedMaterial_invalid_global2Grid_Unconnected) {
0117   std::vector<MaterialSlab> material;
0118 
0119   using EqBound = GridAxisGenerators::EqBound;
0120   using EqGrid = EqBound::grid_type<std::size_t>;
0121 
0122   EqBound eqBound{{0., 5.}, 5};
0123   EqGrid eqGrid{eqBound()};
0124 
0125   auto localX = std::make_unique<const LocalAccessX>();
0126   IndexedSurfaceMaterial<EqGrid>::BoundToGridLocalDelegate bToX;
0127   bToX.connect<&LocalAccessX::l2X>(std::move(localX));
0128 
0129   IndexedSurfaceMaterial<EqGrid>::GlobalToGridLocalDelegate gToX;
0130 
0131   BOOST_CHECK_THROW(
0132       auto ism = IndexedSurfaceMaterial<EqGrid>(
0133           std::move(eqGrid), IndexedMaterialAccessor{std::move(material)},
0134           std::move(bToX), std::move(gToX)),
0135       std::invalid_argument);
0136 }
0137 
0138 // This test covers the locally indexed grid material in 1D
0139 BOOST_AUTO_TEST_CASE(GridMaterial1D) {
0140   std::vector<MaterialSlab> material;
0141   material.emplace_back(Material::Vacuum(), 0.0);  // vacuum
0142   material.emplace_back(Material::fromMolarDensity(1.0, 2.0, 3.0, 4.0, 5.0),
0143                         1.0);
0144   material.emplace_back(
0145       Material::fromMolarDensity(11.0, 12.0, 13.0, 14.0, 15.0), 2.0);
0146   material.emplace_back(
0147       Material::fromMolarDensity(21.0, 22.0, 23.0, 24.0, 25.0), 3.0);
0148   material.emplace_back(
0149       Material::fromMolarDensity(31.0, 32.0, 33.0, 34.0, 35.0), 4.0);
0150 
0151   // Bound, equidistant axis
0152   ProtoAxis pAxisX(AxisBoundaryType::Bound, 0.0, 5.0, 5);
0153 
0154   auto localX = std::make_unique<const LocalAccessX>();
0155   GridAccess::BoundToGridLocal1DimDelegate bToX;
0156   bToX.connect<&LocalAccessX::l2X>(std::move(localX));
0157 
0158   auto globalX = std::make_unique<const GlobalAccessX>();
0159   GridAccess::GlobalToGridLocal1DimDelegate gToX;
0160   gToX.connect<&GlobalAccessX::g2X>(std::move(globalX));
0161 
0162   auto ismX = GridSurfaceMaterialFactory::create(pAxisX, GridMaterialAccessor{},
0163                                                  std::move(bToX),
0164                                                  std::move(gToX), material);
0165 
0166   BOOST_CHECK(ismX != nullptr);
0167 
0168   // Local access test
0169   Vector2 l0(0.5, 0.);
0170   Vector2 l1(1.5, 0.);
0171   Vector2 l2(2.5, 0.);
0172   Vector2 l3(3.5, 0.);
0173   Vector2 l4(4.5, 0.);
0174 
0175   const MaterialSlab& ml0 = ismX->materialSlab(l0);
0176   const MaterialSlab& ml1 = ismX->materialSlab(l1);
0177   const MaterialSlab& ml2 = ismX->materialSlab(l2);
0178   const MaterialSlab& ml3 = ismX->materialSlab(l3);
0179   const MaterialSlab& ml4 = ismX->materialSlab(l4);
0180 
0181   BOOST_CHECK(ml0.material().isVacuum());
0182   BOOST_CHECK_EQUAL(ml1.material().X0(), 1.);
0183   BOOST_CHECK_EQUAL(ml2.material().X0(), 11.);
0184   BOOST_CHECK_EQUAL(ml3.material().X0(), 21.);
0185   BOOST_CHECK_EQUAL(ml4.material().X0(), 31.);
0186 
0187   // Try the same with Closed access
0188   // Bound, equidistant axis
0189   ProtoAxis pAxisPhi(AxisBoundaryType::Closed, -std::numbers::pi,
0190                      std::numbers::pi, 8);
0191 
0192   auto localPhi = std::make_unique<const LocalAccessPhi>();
0193   GridAccess::BoundToGridLocal1DimDelegate bToPhi;
0194   bToPhi.connect<&LocalAccessPhi::l2Phi>(std::move(localPhi));
0195 
0196   auto globalPhi = std::make_unique<const GlobalAccessPhi>();
0197   GridAccess::GlobalToGridLocal1DimDelegate gToPhi;
0198   gToPhi.connect<&GlobalAccessPhi::g2Phi>(std::move(globalPhi));
0199 
0200   std::vector<MaterialSlab> materialPhi;
0201   materialPhi.emplace_back(Material::Vacuum(), 0.0);  // vacuum
0202   materialPhi.emplace_back(Material::fromMolarDensity(1.0, 2.0, 3.0, 4.0, 5.0),
0203                            1.0);
0204   materialPhi.emplace_back(Material::Vacuum(), 0.0);  // vacuum
0205   materialPhi.emplace_back(
0206       Material::fromMolarDensity(11.0, 12.0, 13.0, 14.0, 15.0), 2.0);
0207   materialPhi.emplace_back(Material::Vacuum(), 0.0);  // vacuum
0208   materialPhi.emplace_back(
0209       Material::fromMolarDensity(21.0, 22.0, 23.0, 24.0, 25.0), 3.0);
0210   materialPhi.emplace_back(Material::Vacuum(), 0.0);  // vacuum
0211   materialPhi.emplace_back(
0212       Material::fromMolarDensity(31.0, 32.0, 33.0, 34.0, 35.0), 4.0);
0213 
0214   auto ismPhi = GridSurfaceMaterialFactory::create(
0215       pAxisPhi, GridMaterialAccessor{}, std::move(bToPhi), std::move(gToPhi),
0216       materialPhi);
0217 
0218   BOOST_CHECK(ismPhi != nullptr);
0219 
0220   for (std::size_t i = 0; i < 8; ++i) {
0221     double alpha = -std::numbers::pi + (i + 0.5) * std::numbers::pi / 4.;
0222     Vector2 query{std::cos(alpha), std::sin(alpha)};
0223     const MaterialSlab& m = ismPhi->materialSlab(query);
0224     if (i % 2 == 0) {
0225       BOOST_CHECK(m.material().isVacuum());
0226     } else {
0227       BOOST_CHECK_EQUAL(m.material().X0(), materialPhi[i].material().X0());
0228     }
0229   }
0230 }
0231 
0232 // This test covers the locally indexed grid material in 1D
0233 BOOST_AUTO_TEST_CASE(GridMaterial2D) {
0234   std::vector<std::vector<MaterialSlab>> material2x3;
0235   // This is a material matrix 2 bins in x and 3 bins in y
0236   std::vector<MaterialSlab> materialRow0;
0237   materialRow0.emplace_back(Material::fromMolarDensity(1.0, 2.0, 3.0, 4.0, 5.0),
0238                             1.0);
0239   materialRow0.emplace_back(
0240       Material::fromMolarDensity(11.0, 12.0, 13.0, 14.0, 15.0), 2.0);
0241   materialRow0.emplace_back(
0242       Material::fromMolarDensity(21.0, 22.0, 23.0, 24.0, 25.0), 3.0);
0243   std::vector<MaterialSlab> materialRow1;
0244   materialRow1.emplace_back(Material::fromMolarDensity(2.0, 2.0, 3.0, 4.0, 5.0),
0245                             1.0);
0246   materialRow1.emplace_back(
0247       Material::fromMolarDensity(12.0, 12.0, 13.0, 14.0, 15.0), 2.0);
0248   materialRow1.emplace_back(
0249       Material::fromMolarDensity(22.0, 22.0, 23.0, 24.0, 25.0), 3.0);
0250   // This gives a row major matrix
0251   material2x3.push_back(std::move(materialRow0));
0252   material2x3.push_back(std::move(materialRow1));
0253 
0254   BOOST_CHECK(material2x3[0][0].material().X0() == 1.);
0255   BOOST_CHECK(material2x3[0][1].material().X0() == 11.);
0256   BOOST_CHECK(material2x3[0][2].material().X0() == 21.);
0257   BOOST_CHECK(material2x3[1][0].material().X0() == 2.);
0258   BOOST_CHECK(material2x3[1][1].material().X0() == 12.);
0259   BOOST_CHECK(material2x3[1][2].material().X0() == 22.);
0260 
0261   ProtoAxis pAxisX(AxisBoundaryType::Bound, -1.0, 1.0, 2);
0262   ProtoAxis pAxisY(AxisBoundaryType::Bound, -1.5, 1.5, 3);
0263 
0264   std::vector<std::vector<MaterialSlab>> materialXY = material2x3;
0265 
0266   auto localXY = std::make_unique<const LocalAccessXY>();
0267   GridAccess::BoundToGridLocal2DimDelegate bToXY;
0268   bToXY.connect<&LocalAccessXY::l2XY>(std::move(localXY));
0269 
0270   auto globalXY = std::make_unique<const GlobalAccessXY>();
0271   GridAccess::GlobalToGridLocal2DimDelegate gToXY;
0272   gToXY.connect<&GlobalAccessXY::g2XY>(std::move(globalXY));
0273 
0274   auto ismXY = GridSurfaceMaterialFactory::create(
0275       pAxisX, pAxisY, GridMaterialAccessor{}, std::move(bToXY),
0276       std::move(gToXY), materialXY);
0277 
0278   BOOST_CHECK(ismXY != nullptr);
0279 
0280   // Local access test
0281   Vector2 l00(-0.5, -1.5);
0282   Vector2 l01(-0.5, 0.);
0283   Vector2 l02(-0.5, 1.5);
0284   Vector2 l10(0.5, -1.5);
0285   Vector2 l11(0.5, 0.);
0286   Vector2 l12(0.5, 1.5);
0287 
0288   const MaterialSlab& ml00 = ismXY->materialSlab(l00);
0289   const MaterialSlab& ml01 = ismXY->materialSlab(l01);
0290   const MaterialSlab& ml02 = ismXY->materialSlab(l02);
0291   const MaterialSlab& ml10 = ismXY->materialSlab(l10);
0292   const MaterialSlab& ml11 = ismXY->materialSlab(l11);
0293   const MaterialSlab& ml12 = ismXY->materialSlab(l12);
0294 
0295   BOOST_CHECK_EQUAL(ml00.material().X0(), 1.);
0296   BOOST_CHECK_EQUAL(ml01.material().X0(), 11.);
0297   BOOST_CHECK_EQUAL(ml02.material().X0(), 21.);
0298   BOOST_CHECK_EQUAL(ml10.material().X0(), 2.);
0299   BOOST_CHECK_EQUAL(ml11.material().X0(), 12.);
0300   BOOST_CHECK_EQUAL(ml12.material().X0(), 22.);
0301 
0302   // Let's try a ZPhi model as well
0303   auto materialZPhi = material2x3;
0304 
0305   ProtoAxis pAxisZ(AxisBoundaryType::Bound, -1.0, 1.0, 2);
0306   ProtoAxis pAxisPhi(AxisBoundaryType::Closed, -std::numbers::pi,
0307                      std::numbers::pi, 3);
0308 
0309   auto localZPhi = std::make_unique<const LocalToZPhi>(1.);
0310   GridAccess::BoundToGridLocal2DimDelegate bToZPhi;
0311   bToZPhi.connect<&LocalToZPhi::l2ZPhi>(std::move(localZPhi));
0312 
0313   auto globalZPhi = std::make_unique<const GlobalToZPhi>(0.);
0314   GridAccess::GlobalToGridLocal2DimDelegate gToZPhi;
0315   gToZPhi.connect<&GlobalToZPhi::g2ZPhi>(std::move(globalZPhi));
0316 
0317   auto ismZPhi = GridSurfaceMaterialFactory::create(
0318       pAxisZ, pAxisPhi, GridMaterialAccessor{}, std::move(bToZPhi),
0319       std::move(gToZPhi), materialZPhi);
0320 
0321   BOOST_CHECK(ismZPhi != nullptr);
0322 
0323   // Local access test - trick here is also to
0324   // see BoundtoGridLocal switches from r*phi, z -> phi,z
0325   Vector2 cl00(-0.5 * std::numbers::pi, -0.5);
0326   Vector2 cl01(0., -0.5);
0327   Vector2 cl02(0.5 * std::numbers::pi, -0.5);
0328   Vector2 cl10(-0.5 * std::numbers::pi, 0.5);
0329   Vector2 cl11(0., 0.5);
0330   Vector2 cl12(0.5 * std::numbers::pi, 0.5);
0331 
0332   const MaterialSlab& cml00 = ismZPhi->materialSlab(cl00);
0333   const MaterialSlab& cml01 = ismZPhi->materialSlab(cl01);
0334   const MaterialSlab& cml02 = ismZPhi->materialSlab(cl02);
0335   const MaterialSlab& cml10 = ismZPhi->materialSlab(cl10);
0336   const MaterialSlab& cml11 = ismZPhi->materialSlab(cl11);
0337   const MaterialSlab& cml12 = ismZPhi->materialSlab(cl12);
0338 
0339   BOOST_CHECK_EQUAL(cml00.material().X0(), 1.);
0340   BOOST_CHECK_EQUAL(cml01.material().X0(), 11.);
0341   BOOST_CHECK_EQUAL(cml02.material().X0(), 21.);
0342   BOOST_CHECK_EQUAL(cml10.material().X0(), 2.);
0343   BOOST_CHECK_EQUAL(cml11.material().X0(), 12.);
0344   BOOST_CHECK_EQUAL(cml12.material().X0(), 22.);
0345 
0346   // Test the closed character of the phi axis
0347   Vector2 cl03(1.05 * std::numbers::pi, -0.5);
0348   const MaterialSlab& cml03 = ismZPhi->materialSlab(cl03);
0349   BOOST_CHECK(cml03.material().X0() == cml00.material().X0());
0350 }
0351 
0352 // This test covers the locally indexed grid material in 1D
0353 BOOST_AUTO_TEST_CASE(GridIndexedMaterial1D) {
0354   std::vector<MaterialSlab> material;
0355   material.emplace_back(Material::Vacuum(), 0.0);  // vacuum
0356   material.emplace_back(Material::fromMolarDensity(1.0, 2.0, 3.0, 4.0, 5.0),
0357                         1.0);
0358   material.emplace_back(
0359       Material::fromMolarDensity(11.0, 12.0, 13.0, 14.0, 15.0), 2.0);
0360   material.emplace_back(
0361       Material::fromMolarDensity(21.0, 22.0, 23.0, 24.0, 25.0), 3.0);
0362 
0363   using EqBound = GridAxisGenerators::EqBound;
0364   using EqGrid = EqBound::grid_type<std::size_t>;
0365   using Point = EqGrid::point_t;
0366 
0367   EqBound eqBound{{0., 5.}, 5};
0368   EqGrid eqGrid{eqBound()};
0369 
0370   eqGrid.atPosition(Point{0.5}) = 1u;  // material 1
0371   eqGrid.atPosition(Point{1.5}) = 0u;  // vacuum
0372   eqGrid.atPosition(Point{2.5}) = 2u;  // material 2
0373   eqGrid.atPosition(Point{3.5}) = 2u;  // material 2
0374   eqGrid.atPosition(Point{4.5}) = 3u;  // material 3
0375 
0376   auto localX = std::make_unique<const LocalAccessX>();
0377   IndexedSurfaceMaterial<EqGrid>::BoundToGridLocalDelegate bToX;
0378   bToX.connect<&LocalAccessX::l2X>(std::move(localX));
0379 
0380   auto globalX = std::make_unique<const GlobalAccessX>();
0381   IndexedSurfaceMaterial<EqGrid>::GlobalToGridLocalDelegate gToX;
0382   gToX.connect<&GlobalAccessX::g2X>(std::move(globalX));
0383 
0384   IndexedSurfaceMaterial<EqGrid> ism(
0385       std::move(eqGrid), IndexedMaterialAccessor{std::move(material)},
0386       std::move(bToX), std::move(gToX));
0387 
0388   // Local access test
0389   Vector2 l0(0.5, 0.);
0390   Vector2 l1(1.5, 0.);
0391   Vector2 l2(2.5, 0.);
0392   Vector2 l3(3.5, 0.);
0393   Vector2 l4(4.5, 0.);
0394 
0395   const MaterialSlab& ml0 = ism.materialSlab(l0);
0396   const MaterialSlab& ml1 = ism.materialSlab(l1);
0397   const MaterialSlab& ml2 = ism.materialSlab(l2);
0398   const MaterialSlab& ml3 = ism.materialSlab(l3);
0399   const MaterialSlab& ml4 = ism.materialSlab(l4);
0400 
0401   BOOST_CHECK_EQUAL(ml0.material().X0(), 1.);
0402   BOOST_CHECK(ml1.material().isVacuum());
0403   BOOST_CHECK_EQUAL(ml2.material().X0(), 11.);
0404   BOOST_CHECK_EQUAL(ml3.material().X0(), 11.);
0405   BOOST_CHECK_EQUAL(ml4.material().X0(), 21.);
0406 
0407   // Now scale it - and access again
0408   ism.scale(2.);
0409   const MaterialSlab& sml0 = ism.materialSlab(l0);
0410   const MaterialSlab& sml1 = ism.materialSlab(l1);
0411   const MaterialSlab& sml2 = ism.materialSlab(l2);
0412   const MaterialSlab& sml3 = ism.materialSlab(l3);
0413   const MaterialSlab& sml4 = ism.materialSlab(l4);
0414 
0415   BOOST_CHECK_EQUAL(sml0.thickness(), 2.);
0416   BOOST_CHECK(sml1.material().isVacuum());
0417   BOOST_CHECK_EQUAL(sml2.thickness(), 4.);
0418   BOOST_CHECK_EQUAL(sml3.thickness(), 4.);
0419   BOOST_CHECK_EQUAL(sml4.thickness(), 6.);
0420 
0421   // Now test with the protoAxis creation method
0422 
0423   std::vector<MaterialSlab> materialStorage;
0424   materialStorage.emplace_back(Material::Vacuum(), 0.0);  // vacuum
0425   materialStorage.emplace_back(
0426       Material::fromMolarDensity(1.0, 2.0, 3.0, 4.0, 5.0), 1.0);
0427   materialStorage.emplace_back(
0428       Material::fromMolarDensity(11.0, 12.0, 13.0, 14.0, 15.0), 2.0);
0429   materialStorage.emplace_back(
0430       Material::fromMolarDensity(21.0, 22.0, 23.0, 24.0, 25.0), 3.0);
0431 
0432   std::vector<std::size_t> indexPayload = {0u, 1u, 2u, 3u, 0u, 3u, 2u, 1u, 0u};
0433 
0434   auto indexedAccessor = IndexedMaterialAccessor{std::move(materialStorage)};
0435 
0436   // An X proto axis
0437   ProtoAxis pAxisX(AxisBoundaryType::Bound, 0.0, 9.0, 9);
0438 
0439   auto localXidx = std::make_unique<const LocalAccessX>();
0440   IndexedSurfaceMaterial<EqGrid>::BoundToGridLocalDelegate bToXidx;
0441   bToXidx.connect<&LocalAccessX::l2X>(std::move(localXidx));
0442 
0443   auto globalXidx = std::make_unique<const GlobalAccessX>();
0444   IndexedSurfaceMaterial<EqGrid>::GlobalToGridLocalDelegate gToXidx;
0445   gToXidx.connect<&GlobalAccessX::g2X>(std::move(globalXidx));
0446 
0447   auto ismXidx = GridSurfaceMaterialFactory::create(
0448       pAxisX, std::move(indexedAccessor), std::move(bToXidx),
0449       std::move(gToXidx), indexPayload);
0450 
0451   // Check construction
0452   BOOST_CHECK(ismXidx != nullptr);
0453   // The vacuum (==0) indexed entries
0454   BOOST_CHECK(ismXidx->materialSlab(Vector2{0.5, 0.}).isVacuum());
0455   BOOST_CHECK(ismXidx->materialSlab(Vector2{4.5, 0.}).isVacuum());
0456   BOOST_CHECK(ismXidx->materialSlab(Vector2{8.5, 0.}).isVacuum());
0457   // The material 1 (==1) indexed entries
0458   BOOST_CHECK_EQUAL(ismXidx->materialSlab(Vector2{1.5, 0.}).material().X0(),
0459                     1.);
0460   BOOST_CHECK_EQUAL(ismXidx->materialSlab(Vector2{7.5, 0.}).material().X0(),
0461                     1.);
0462   // The material 2 (==2) indexed entries
0463   BOOST_CHECK_EQUAL(ismXidx->materialSlab(Vector2{2.5, 0.}).material().X0(),
0464                     11.);
0465   BOOST_CHECK_EQUAL(ismXidx->materialSlab(Vector2{6.5, 0.}).material().X0(),
0466                     11.);
0467   // The material 3 (==3) indexed entries
0468   BOOST_CHECK_EQUAL(ismXidx->materialSlab(Vector2{3.5, 0.}).material().X0(),
0469                     21.);
0470   BOOST_CHECK_EQUAL(ismXidx->materialSlab(Vector2{5.5, 0.}).material().X0(),
0471                     21.);
0472 }
0473 
0474 // This test covers the locally indexed grid material in 2D
0475 BOOST_AUTO_TEST_CASE(GridIndexedMaterial2D) {
0476   std::vector<MaterialSlab> material;
0477   material.emplace_back(Material::Vacuum(), 1.0);  // vacuum
0478   material.emplace_back(Material::fromMolarDensity(1.0, 2.0, 3.0, 4.0, 5.0),
0479                         1.0);
0480   material.emplace_back(
0481       Material::fromMolarDensity(11.0, 12.0, 13.0, 14.0, 15.0), 1.0);
0482   material.emplace_back(
0483       Material::fromMolarDensity(21.0, 22.0, 23.0, 24.0, 25.0), 1.0);
0484 
0485   //  Test (1) with explicit grid creation
0486   std::vector<MaterialSlab> materialT1 = material;
0487   using EqBoundEqClosed = GridAxisGenerators::EqBoundEqClosed;
0488   using EqEqGrid = EqBoundEqClosed::grid_type<std::size_t>;
0489   using Point = EqEqGrid::point_t;
0490 
0491   // 2 bins in z, 4 bins in phi
0492   EqBoundEqClosed eqeqBound{
0493       {-1., 1.}, 2, {-std::numbers::pi, std::numbers::pi}, 4};
0494   EqEqGrid eqeqGrid{eqeqBound()};
0495 
0496   eqeqGrid.atPosition(Point{-0.5, -std::numbers::pi * 0.75}) =
0497       1u;                                                          // material 1
0498   eqeqGrid.atPosition(Point{-0.5, -std::numbers::pi / 4.}) = 1u;   // material 1
0499   eqeqGrid.atPosition(Point{-0.5, std::numbers::pi / 4.}) = 0u;    // vacuum
0500   eqeqGrid.atPosition(Point{-0.5, std::numbers::pi * 0.75}) = 2u;  // material 2
0501 
0502   eqeqGrid.atPosition(Point{0.5, -std::numbers::pi * 0.75}) = 0u;  // vacuum
0503   eqeqGrid.atPosition(Point{0.5, -std::numbers::pi / 4.}) = 3u;    // material 3
0504   eqeqGrid.atPosition(Point{0.5, std::numbers::pi / 4.}) = 3u;     // material 3
0505   eqeqGrid.atPosition(Point{0.5, std::numbers::pi * 0.75}) = 0u;   // vacuum
0506 
0507   // With radius 20
0508   auto boundToGridT1 = std::make_unique<const LocalToZPhi>(20.);
0509   IndexedSurfaceMaterial<EqEqGrid>::BoundToGridLocalDelegate bToZPhiT1;
0510   bToZPhiT1.connect<&LocalToZPhi::l2ZPhi>(std::move(boundToGridT1));
0511 
0512   // With z shift 10
0513   auto globalToGridT1 = std::make_unique<const GlobalToZPhi>(10.);
0514   IndexedSurfaceMaterial<EqEqGrid>::GlobalToGridLocalDelegate gToZphiT1;
0515   gToZphiT1.connect<&GlobalToZPhi::g2ZPhi>(std::move(globalToGridT1));
0516 
0517   // Create the indexed material grid
0518   IndexedSurfaceMaterial<EqEqGrid> ism(
0519       std::move(eqeqGrid), IndexedMaterialAccessor{std::move(materialT1)},
0520       std::move(bToZPhiT1), std::move(gToZphiT1));
0521 
0522   // Test with proto grid greation
0523   auto materialT2 = material;
0524 
0525   auto boundToGridT2 = std::make_unique<const LocalToZPhi>(20.);
0526   GridAccess::BoundToGridLocal2DimDelegate bToZPhiT2;
0527   bToZPhiT2.connect<&LocalToZPhi::l2ZPhi>(std::move(boundToGridT2));
0528 
0529   auto globalToGridT2 = std::make_unique<const GlobalToZPhi>(10.);
0530   GridAccess::GlobalToGridLocal2DimDelegate gToZphiT2;
0531   gToZphiT2.connect<&GlobalToZPhi::g2ZPhi>(std::move(globalToGridT2));
0532 
0533   ProtoAxis pAxisZ(AxisBoundaryType::Bound, -1.0, 1.0, 2);
0534   ProtoAxis pAxisPhi(AxisBoundaryType::Closed, -std::numbers::pi,
0535                      std::numbers::pi, 4);
0536 
0537   std::vector<std::vector<std::size_t>> indexPayload = {
0538       std::vector<std::size_t>{1u, 1u, 0u, 2u},
0539       std::vector<std::size_t>{0u, 3u, 3u, 0u}};
0540 
0541   auto ismZPhi = GridSurfaceMaterialFactory::create(
0542       pAxisZ, pAxisPhi, IndexedMaterialAccessor{std::move(materialT2)},
0543       std::move(bToZPhiT2), std::move(gToZphiT2), indexPayload);
0544 
0545   // Local access test, both should give material 1
0546   Vector2 l0(-20 * std::numbers::pi * 0.75, -10.5);
0547   const MaterialSlab& ml0T1 = ism.materialSlab(l0);
0548   const MaterialSlab& ml0T2 = ismZPhi->materialSlab(l0);
0549   BOOST_CHECK_EQUAL(ml0T1.material().X0(), 1.);
0550   BOOST_CHECK_EQUAL(ml0T2.material().X0(), 1.);
0551 
0552   Vector2 l1(-20 * std::numbers::pi * 0.25,
0553              -11.5);  // checking out of bound access
0554 
0555   const MaterialSlab& mg1T1 = ism.materialSlab(l1);
0556   const MaterialSlab& mg1T2 = ismZPhi->materialSlab(l1);
0557   BOOST_CHECK_EQUAL(mg1T1.material().X0(), 1.);
0558   BOOST_CHECK_EQUAL(mg1T2.material().X0(), 1.);
0559 
0560   Vector2 l2(20 * std::numbers::pi * 0.25, -10.5);
0561   const MaterialSlab& mg2T1 = ism.materialSlab(l2);
0562   const MaterialSlab& mg2T2 = ismZPhi->materialSlab(l2);
0563   BOOST_CHECK(mg2T1.material().isVacuum());  // vacuum
0564   BOOST_CHECK(mg2T2.material().isVacuum());  // vacuum
0565 }
0566 
0567 // This test covers the globally indexed grid material with non-shared material
0568 BOOST_AUTO_TEST_CASE(GridGloballyIndexedMaterialNonShared) {
0569   auto material = std::make_shared<std::vector<MaterialSlab>>();
0570 
0571   material->emplace_back(Material::Vacuum(), 0.0);  // vacuum
0572   material->emplace_back(Material::fromMolarDensity(1.0, 2.0, 3.0, 4.0, 5.0),
0573                          1.0);
0574   material->emplace_back(
0575       Material::fromMolarDensity(11.0, 12.0, 13.0, 14.0, 15.0), 2.0);
0576   material->emplace_back(
0577       Material::fromMolarDensity(21.0, 22.0, 23.0, 24.0, 25.0), 3.0);
0578   material->emplace_back(
0579       Material::fromMolarDensity(31.0, 22.0, 23.0, 24.0, 25.0), 4.0);
0580 
0581   using EqBound = GridAxisGenerators::EqBound;
0582   using EqGrid = EqBound::grid_type<std::size_t>;
0583   using Point = EqGrid::point_t;
0584 
0585   EqBound eqBound{{0., 5.}, 5};
0586   EqGrid eqGrid{eqBound()};
0587 
0588   eqGrid.atPosition(Point{0.5}) = 1u;  // material 1
0589   eqGrid.atPosition(Point{1.5}) = 0u;  // vacuum
0590   eqGrid.atPosition(Point{2.5}) = 2u;  // material 2
0591   eqGrid.atPosition(Point{3.5}) = 2u;  // material 2
0592   eqGrid.atPosition(Point{4.5}) = 3u;  // material 3
0593 
0594   auto localX = std::make_unique<const LocalAccessX>();
0595   IndexedSurfaceMaterial<EqGrid>::BoundToGridLocalDelegate bToX;
0596   bToX.connect<&LocalAccessX::l2X>(std::move(localX));
0597 
0598   auto globalX = std::make_unique<const GlobalAccessX>();
0599   IndexedSurfaceMaterial<EqGrid>::GlobalToGridLocalDelegate gToX;
0600   gToX.connect<&GlobalAccessX::g2X>(std::move(globalX));
0601 
0602   GloballyIndexedSurfaceMaterial<EqGrid> ism(
0603       std::move(eqGrid), GloballyIndexedMaterialAccessor{material, false},
0604       std::move(bToX), std::move(gToX));
0605 
0606   // Local access test
0607   Vector2 l0(0.5, 0.);
0608   Vector2 l1(1.5, 0.);
0609   Vector2 l2(2.5, 0.);
0610   Vector2 l3(3.5, 0.);
0611   Vector2 l4(4.5, 0.);
0612 
0613   const MaterialSlab& ml0 = ism.materialSlab(l0);
0614   const MaterialSlab& ml1 = ism.materialSlab(l1);
0615   const MaterialSlab& ml2 = ism.materialSlab(l2);
0616   const MaterialSlab& ml3 = ism.materialSlab(l3);
0617   const MaterialSlab& ml4 = ism.materialSlab(l4);
0618 
0619   BOOST_CHECK_EQUAL(ml0.material().X0(), 1.);
0620   BOOST_CHECK(ml1.material().isVacuum());
0621   BOOST_CHECK_EQUAL(ml2.material().X0(), 11.);
0622   BOOST_CHECK_EQUAL(ml3.material().X0(), 11.);
0623   BOOST_CHECK_EQUAL(ml4.material().X0(), 21.);
0624 
0625   EqBound eqBound1{{0., 5.}, 1};
0626   EqGrid eqGrid1{eqBound1()};
0627 
0628   auto localX1 = std::make_unique<const LocalAccessX>();
0629   IndexedSurfaceMaterial<EqGrid>::BoundToGridLocalDelegate bToX1;
0630   bToX1.connect<&LocalAccessX::l2X>(std::move(localX1));
0631 
0632   auto globalX1 = std::make_unique<const GlobalAccessX>();
0633   IndexedSurfaceMaterial<EqGrid>::GlobalToGridLocalDelegate gToX1;
0634   gToX1.connect<&GlobalAccessX::g2X>(std::move(globalX1));
0635 
0636   eqGrid1.atPosition(Point{2.5}) = 4u;  // material 4
0637 
0638   GloballyIndexedSurfaceMaterial<EqGrid> ism1(
0639       std::move(eqGrid1), GloballyIndexedMaterialAccessor{material, false},
0640       std::move(bToX1), std::move(gToX1));
0641 
0642   Vector2 l0g1(2.5, 0.);
0643   const MaterialSlab& ml0g1 = ism1.materialSlab(l0g1);
0644   BOOST_CHECK_EQUAL(ml0g1.material().X0(), 31.);
0645 
0646   // Scale
0647   ism1.scale(2.);
0648   const MaterialSlab& sml0g1 = ism1.materialSlab(l0g1);
0649   BOOST_CHECK_EQUAL(sml0g1.thickness(), 8.);
0650 
0651   // First one stays unscaled
0652   const MaterialSlab& sml0 = ism.materialSlab(l0);
0653   BOOST_CHECK_EQUAL(sml0.thickness(), 1.);
0654 }
0655 
0656 // This test covers the globally indexed grid material with shared
0657 BOOST_AUTO_TEST_CASE(GridGloballyIndexedMaterialShared) {
0658   auto material = std::make_shared<std::vector<MaterialSlab>>();
0659 
0660   material->emplace_back(Material::Vacuum(), 0.0);  // vacuum
0661   material->emplace_back(Material::fromMolarDensity(1.0, 2.0, 3.0, 4.0, 5.0),
0662                          1.0);
0663 
0664   using EqBound = GridAxisGenerators::EqBound;
0665   using EqGrid = EqBound::grid_type<std::size_t>;
0666   using Point = EqGrid::point_t;
0667 
0668   EqBound eqBound0{{0., 5.}, 1};
0669   EqGrid eqGrid0{eqBound0()};
0670 
0671   eqGrid0.atPosition(Point{2.5}) = 1u;  // material 1
0672   auto localX0 = std::make_unique<const LocalAccessX>();
0673   IndexedSurfaceMaterial<EqGrid>::BoundToGridLocalDelegate bToX0;
0674   bToX0.connect<&LocalAccessX::l2X>(std::move(localX0));
0675 
0676   auto globalX0 = std::make_unique<const GlobalAccessX>();
0677   IndexedSurfaceMaterial<EqGrid>::GlobalToGridLocalDelegate gToX0;
0678   gToX0.connect<&GlobalAccessX::g2X>(std::move(globalX0));
0679 
0680   GloballyIndexedSurfaceMaterial<EqGrid> ism0(
0681       std::move(eqGrid0), GloballyIndexedMaterialAccessor{material, true},
0682       std::move(bToX0), std::move(gToX0));
0683 
0684   EqBound eqBound1{{0., 5.}, 1};
0685   EqGrid eqGrid1{eqBound1()};
0686 
0687   eqGrid1.atPosition(Point{2.5}) = 1u;  // material 1
0688   auto localX1 = std::make_unique<const LocalAccessX>();
0689   IndexedSurfaceMaterial<EqGrid>::BoundToGridLocalDelegate bToX1;
0690   bToX1.connect<&LocalAccessX::l2X>(std::move(localX1));
0691 
0692   auto globalX1 = std::make_unique<const GlobalAccessX>();
0693   IndexedSurfaceMaterial<EqGrid>::GlobalToGridLocalDelegate gToX1;
0694   gToX1.connect<&GlobalAccessX::g2X>(std::move(globalX1));
0695 
0696   GloballyIndexedSurfaceMaterial<EqGrid> ism1(
0697       std::move(eqGrid1), GloballyIndexedMaterialAccessor{material, true},
0698       std::move(bToX1), std::move(gToX1));
0699 
0700   Vector2 l0(2.5, 0.);
0701 
0702   // check grid material 0
0703   const MaterialSlab& ml0 = ism0.materialSlab(l0);
0704   BOOST_CHECK_EQUAL(ml0.material().X0(), 1.);
0705 
0706   const MaterialSlab& ml0g1 = ism1.materialSlab(l0);
0707   BOOST_CHECK_EQUAL(ml0g1.material().X0(), 1.);
0708 
0709   // scaling shared material should throw a std::invalid_argument
0710   BOOST_CHECK_THROW(ism1.scale(2.), std::invalid_argument);
0711 }
0712 
0713 // This test covers the grid material (non-indexed accessor)
0714 //
0715 // In this setup, the material is not indexed, but filled directly
0716 // into the grid structure.
0717 BOOST_AUTO_TEST_CASE(GridSurfaceMaterialTests) {
0718   using EqBound = GridAxisGenerators::EqBound;
0719   using EqGrid = EqBound::grid_type<MaterialSlab>;
0720   using Point = EqGrid::point_t;
0721 
0722   EqBound eqBound{{0., 5.}, 5};
0723   EqGrid eqGrid{eqBound()};
0724 
0725   eqGrid.atPosition(Point{0.5}) = MaterialSlab::Vacuum(0.0);
0726   eqGrid.atPosition(Point{1.5}) = MaterialSlab::Vacuum(1.0);
0727   eqGrid.atPosition(Point{2.5}) = MaterialSlab::Vacuum(2.0);
0728   eqGrid.atPosition(Point{3.5}) = MaterialSlab::Vacuum(3.0);
0729   eqGrid.atPosition(Point{4.5}) = MaterialSlab::Vacuum(4.0);
0730 
0731   auto localX = std::make_unique<const LocalAccessX>();
0732   GridSurfaceMaterial<EqGrid>::BoundToGridLocalDelegate bToX;
0733   bToX.connect<&LocalAccessX::l2X>(std::move(localX));
0734 
0735   auto globalX = std::make_unique<const GlobalAccessX>();
0736   GridSurfaceMaterial<EqGrid>::GlobalToGridLocalDelegate gToX;
0737   gToX.connect<&GlobalAccessX::g2X>(std::move(globalX));
0738 
0739   GridSurfaceMaterial<EqGrid> gsm(std::move(eqGrid), GridMaterialAccessor{},
0740                                   std::move(bToX), std::move(gToX));
0741 
0742   // Local access test
0743   Vector2 l0(0.5, 0.);
0744   Vector2 l1(1.5, 0.);
0745   Vector2 l2(2.5, 0.);
0746   Vector2 l3(3.5, 0.);
0747   Vector2 l4(4.5, 0.);
0748 
0749   const MaterialSlab& ml0 = gsm.materialSlab(l0);
0750   const MaterialSlab& ml1 = gsm.materialSlab(l1);
0751   const MaterialSlab& ml2 = gsm.materialSlab(l2);
0752   const MaterialSlab& ml3 = gsm.materialSlab(l3);
0753   const MaterialSlab& ml4 = gsm.materialSlab(l4);
0754 
0755   BOOST_CHECK_EQUAL(ml0.thickness(), 0.);
0756   BOOST_CHECK_EQUAL(ml1.thickness(), 1.);
0757   BOOST_CHECK_EQUAL(ml2.thickness(), 2.);
0758   BOOST_CHECK_EQUAL(ml3.thickness(), 3.);
0759   BOOST_CHECK_EQUAL(ml4.thickness(), 4.);
0760 
0761   // Now scale it - and access again
0762   gsm.scale(2.);
0763 
0764   const MaterialSlab& sml0 = gsm.materialSlab(l0);
0765   const MaterialSlab& sml1 = gsm.materialSlab(l1);
0766   const MaterialSlab& sml2 = gsm.materialSlab(l2);
0767   const MaterialSlab& sml3 = gsm.materialSlab(l3);
0768   const MaterialSlab& sml4 = gsm.materialSlab(l4);
0769 
0770   BOOST_CHECK_EQUAL(sml0.thickness(), 0.);
0771   BOOST_CHECK_EQUAL(sml1.thickness(), 2.);
0772   BOOST_CHECK_EQUAL(sml2.thickness(), 4.);
0773   BOOST_CHECK_EQUAL(sml3.thickness(), 6.);
0774   BOOST_CHECK_EQUAL(sml4.thickness(), 8.);
0775 }
0776 
0777 BOOST_AUTO_TEST_SUITE_END()
0778 
0779 }  // namespace ActsTests