LCOV - code coverage report
Current view: top level - test - testMagneticField.cpp (source / functions) Hit Total Coverage
Test: coverage.info.cleaned Lines: 138 138 100.0 %
Date: 2024-04-29 14:43:01 Functions: 16 16 100.0 %

          Line data    Source code
       1             : #include <stdexcept>
       2             : 
       3             : #include "crpropa/magneticField/MagneticFieldGrid.h"
       4             : #include "crpropa/magneticField/CMZField.h"
       5             : #include "crpropa/magneticField/PolarizedSingleModeMagneticField.h"
       6             : #include "crpropa/magneticField/GalacticMagneticField.h"
       7             : #include "crpropa/Grid.h"
       8             : #include "crpropa/Units.h"
       9             : #include "crpropa/Common.h"
      10             : 
      11             : #include "gtest/gtest.h"
      12             : 
      13             : using namespace crpropa;
      14             : 
      15           1 : TEST(testUniformMagneticField, SimpleTest) {
      16             :         UniformMagneticField B(Vector3d(-1, 5, 3));
      17             :         Vector3d b = B.getField(Vector3d(1, 0, 0));
      18           1 :         EXPECT_DOUBLE_EQ(b.getX(), -1);
      19           1 :         EXPECT_DOUBLE_EQ(b.getY(), 5);
      20           1 :         EXPECT_DOUBLE_EQ(b.getZ(), 3);
      21           1 : }
      22             : 
      23           2 : TEST(testMagneticDipoleField, SimpleTest) {
      24             :         // Test magnetic dipole
      25             :         // mu0 / (4*M_PI) * m / r^3 (2*cos(theta)*e_r + sin(theta)*e_theta)
      26             :         MagneticDipoleField B(Vector3d(0,0,0), Vector3d(0,0,1), 1);
      27           1 :         Vector3d b1 = B.getField(Vector3d(0, 0, 1)); // theta = 0
      28           1 :         Vector3d b2 = B.getField(Vector3d(1, 0, 0)); // theta = 0
      29           1 :         EXPECT_NEAR(b1.getX(), 0, 1E-8);
      30           1 :         EXPECT_NEAR(b1.getY(), 0, 1E-8);
      31           1 :         EXPECT_NEAR(b1.getZ(), mu0 / (4*M_PI) * 2, 1E-8);
      32           1 :         EXPECT_NEAR(b2.getX(), 0, 1E-8);
      33           1 :         EXPECT_NEAR(b2.getY(), 0, 1E-8);
      34           1 :         EXPECT_NEAR(b2.getZ(), -1 * mu0 / (4*M_PI), 1E-8);
      35           1 : }
      36             : 
      37             : #ifdef CRPROPA_HAVE_MUPARSER
      38           1 : TEST(testRenormalizeMagneticField, simpleTest) {
      39           1 :         ref_ptr<UniformMagneticField> field = new UniformMagneticField(Vector3d(2*nG, 0, 0));
      40           2 :         RenormalizeMagneticField modField(field, "B^2-1*nG");
      41           1 :         Vector3d b = modField.getField(Vector3d(5, 5, 5));
      42           1 :         EXPECT_NEAR(b.getR(), 3*nG, 0.001);
      43           2 : }
      44             : #endif
      45             : 
      46           2 : TEST(testMagneticFieldList, SimpleTest) {
      47             :         // Test a list of three magnetic fields
      48             :         MagneticFieldList B;
      49           1 :         B.addField(new UniformMagneticField(Vector3d(1, 0, 0)));
      50           1 :         B.addField(new UniformMagneticField(Vector3d(0, 2, 0)));
      51           2 :         B.addField(new UniformMagneticField(Vector3d(0, 0, 3)));
      52           1 :         Vector3d b = B.getField(Vector3d(0.));
      53           1 :         EXPECT_DOUBLE_EQ(b.x, 1);
      54           1 :         EXPECT_DOUBLE_EQ(b.y, 2);
      55           1 :         EXPECT_DOUBLE_EQ(b.z, 3);
      56           1 : }
      57             : 
      58           1 : TEST(testMagneticFieldEvolution, SimpleTest) {
      59             :         // Test if this decorator scales the underlying field as (1+z)^m
      60           1 :         ref_ptr<UniformMagneticField> B = new UniformMagneticField(Vector3d(1,0,0));
      61             :         double z = 1.2;
      62             :         double m = 3;
      63           2 :         MagneticFieldEvolution Bz(B, m);
      64             : 
      65             :         // scaled field
      66           1 :         Vector3d b = Bz.getField(Vector3d(0,0,0), z);
      67           1 :         EXPECT_DOUBLE_EQ(b.x, pow(1+z, m));
      68             : 
      69             :         // unscaled field
      70           1 :         b = Bz.getField(Vector3d(0,0,0));
      71           1 :         EXPECT_DOUBLE_EQ(b.x, 1);
      72           1 : }
      73             : 
      74           1 : class EchoMagneticField: public MagneticField {
      75             : public:
      76           3 :         Vector3d getField(const Vector3d &position) const {
      77           3 :                 return position;
      78             :         }
      79             : };
      80             : 
      81           1 : TEST(testPeriodicMagneticField, Exceptions) {
      82           1 :         ref_ptr<EchoMagneticField> f = new EchoMagneticField();
      83             :         ref_ptr<PeriodicMagneticField> p = new PeriodicMagneticField(f,
      84           1 :                         Vector3d(10000, 10000, 10000), Vector3d(1000, 1000, 1000), true);
      85             : 
      86             :         // box 0, 0, 0
      87           1 :         Vector3d v = p->getField(Vector3d(1000, 2000, 3000));
      88           1 :         EXPECT_DOUBLE_EQ(0, v.x);
      89           1 :         EXPECT_DOUBLE_EQ(1000, v.y);
      90           1 :         EXPECT_DOUBLE_EQ(2000, v.z);
      91             : 
      92             :         // box 1, 2, 3
      93           1 :         v = p->getField(Vector3d(12000, 23000, 34000));
      94           1 :         EXPECT_DOUBLE_EQ(9000, v.x);
      95           1 :         EXPECT_DOUBLE_EQ(2000, v.y);
      96           1 :         EXPECT_DOUBLE_EQ(7000, v.z);
      97             : 
      98             :         // box -1, -2, -3
      99           1 :         v = p->getField(Vector3d(0, -10000, -20000));
     100           1 :         EXPECT_DOUBLE_EQ(1000, v.x);
     101           1 :         EXPECT_DOUBLE_EQ(9000, v.y);
     102           1 :         EXPECT_DOUBLE_EQ(1000, v.z);
     103             : 
     104           1 : }
     105             : 
     106           1 : TEST(testCMZMagneticField, SimpleTest) {
     107           1 :         ref_ptr<CMZField> field = new CMZField();
     108             :         
     109             :         // check use-Values
     110           1 :         EXPECT_FALSE(field->getUseMCField());
     111           1 :         EXPECT_TRUE(field->getUseICField());
     112           1 :         EXPECT_FALSE(field->getUseNTFField());
     113           1 :         EXPECT_FALSE(field->getUseRadioArc());
     114             : 
     115             :         // check set function
     116           1 :         field->setUseMCField(true);
     117           1 :         EXPECT_TRUE(field->getUseMCField());
     118           1 :         field->setUseICField(false);
     119           1 :         EXPECT_FALSE(field->getUseICField());
     120           1 :         field->setUseNTFField(true);
     121           1 :         EXPECT_TRUE(field->getUseNTFField());
     122           1 :         field->setUseRadioArc(true);
     123           1 :         EXPECT_TRUE(field->getUseRadioArc());
     124           1 : }
     125             : 
     126           1 : TEST(testCMZMagneticField, TestICComponent) {
     127           1 :         ref_ptr<CMZField> field = new CMZField();
     128             :         Vector3d pos(10*pc,15*pc,-5*pc);        
     129             : 
     130             :         // check IC field at given position
     131           1 :         Vector3d bVec = field->getField(pos);
     132           1 :         EXPECT_NEAR(bVec.getR(), 10.501*muG, 1E-3*muG);
     133           1 :         EXPECT_NEAR(bVec.x, 0.225*muG, 1E-3*muG);
     134           1 :         EXPECT_NEAR(bVec.y, 0.524*muG, 1E-3*muG);
     135           1 :         EXPECT_NEAR(bVec.z, 10.486*muG, 1E-3*muG);
     136           1 : }
     137           1 : TEST(testCMZMagneticField, TestNTFField){
     138           1 :         ref_ptr<CMZField> field = new CMZField();
     139             :         Vector3d pos(10*pc,15*pc,-5*pc);        
     140             :         
     141             :         // check NFTField at given position
     142           1 :         Vector3d bVec = field->getNTFField(pos);
     143           1 :         EXPECT_NEAR(bVec.getR(),1.692*muG, 1e-3*muG);
     144           1 :         EXPECT_NEAR(bVec.x, -0.584*muG, 1e-3*muG);
     145           1 :         EXPECT_NEAR(bVec.y, -1.185*muG, 1e-3*muG);
     146           1 :         EXPECT_NEAR(bVec.z, 1.057*muG, 1e-3*muG);
     147           1 : }
     148           1 : TEST(testCMZMagneticField, TestRaidoArcField){
     149           1 :         ref_ptr<CMZField> field = new CMZField();
     150             :         Vector3d pos(10*pc,15*pc,-5*pc);        
     151             : 
     152             :         // check RadioArcField at given position
     153           1 :         Vector3d bVec = field->getRadioArcField(pos);
     154           1 :         EXPECT_NEAR(bVec.getR(), 31.616*muG, 1e-3*muG);
     155           1 :         EXPECT_NEAR(bVec.x, -4.671*muG, 1e-3*muG);
     156           1 :         EXPECT_NEAR(bVec.y, 5.465*muG, 1e-3*muG);
     157           1 :         EXPECT_NEAR(bVec.z, 30.788*muG, 1e-3*muG);
     158           1 : }
     159             : 
     160           1 : TEST(testCMZMagneticField, TestAzimutalComponent){
     161           1 :         ref_ptr<CMZField> field = new CMZField();
     162             :         Vector3d mid(12*pc, 9*pc, 20*pc);
     163             :         Vector3d pos(9*pc, 10*pc, 25*pc);       
     164             : 
     165             :         // simple Test for inner part
     166           1 :         Vector3d bVec = field->BAz(pos, mid, 100, 0.2, 60*pc);
     167           1 :         EXPECT_NEAR(bVec.x, 3939.782, 1e-3);
     168           1 :         EXPECT_NEAR(bVec.y, 14347.304, 1e-3);
     169           1 :         EXPECT_DOUBLE_EQ(bVec.z, 0);
     170             : 
     171             :         // simple Test for outer part
     172           1 :         bVec = field->BAz(pos, mid, 100, 0.2, 10*pc);
     173           1 :         EXPECT_NEAR(bVec.x, -164.659, 1e-3);
     174           1 :         EXPECT_NEAR(bVec.y, -1317.270, 1e-3);
     175           1 :         EXPECT_DOUBLE_EQ(bVec.z, 0);
     176             : 
     177             :         // test for molecular Clouds
     178           1 :         bVec = field->getMCField(pos);
     179           1 :         EXPECT_NEAR(bVec.x, -8.339*muG, 1e-3*muG);
     180           1 :         EXPECT_NEAR(bVec.y, -0.850*muG, 1e-3*muG);
     181           1 :         EXPECT_DOUBLE_EQ(bVec.z, 0);
     182           1 : }
     183             : 
     184           1 : TEST(testToroidalHaloField, SimpleTest) {
     185           1 :         ref_ptr<ToroidalHaloField> field = new ToroidalHaloField();
     186           1 :         Vector3d b = field->getField(Vector3d(0.));
     187           1 :         EXPECT_DOUBLE_EQ(b.x, 0);
     188           1 :         EXPECT_DOUBLE_EQ(b.y, 0);
     189           1 :         EXPECT_DOUBLE_EQ(b.z, 0);
     190             : 
     191           1 :         b = field->getField(Vector3d(1,0,0));
     192           1 :         EXPECT_DOUBLE_EQ(b.x, 0.5);
     193           1 :         EXPECT_DOUBLE_EQ(b.y, 0);
     194           1 :         EXPECT_DOUBLE_EQ(b.z, 0);
     195           1 : }
     196             : 
     197           1 : TEST(testLogarithmicSpiralField, SimpleTest) {
     198           1 :         ref_ptr<LogarithmicSpiralField> field = new LogarithmicSpiralField();
     199           1 :         Vector3d b = field->getField(Vector3d(8.5, 0, 0)*kpc);
     200           1 :         EXPECT_NEAR(b.x, -1., 1E-2);
     201           1 :         EXPECT_NEAR(b.y, 0, 1E-10);
     202           1 :         EXPECT_NEAR(b.z, 0, 1E-10);
     203           1 : }
     204             : 
     205           1 : TEST(testPolarizedSingleModeMagneticField, SimpleTest) {
     206           1 :         PolarizedSingleModeMagneticField B(2, 4, 0.5, Vector3d(1,1,1), Vector3d(0,1,0), Vector3d(1,0,0), "amplitude", "polarization", "elliptical");
     207           1 :         Vector3d b = B.getField(Vector3d(1,1,2));
     208           1 :         EXPECT_DOUBLE_EQ(b.x, 1);
     209           1 :         EXPECT_NEAR(b.y, 0, 1E-10);
     210           1 :         EXPECT_NEAR(b.z, 0, 1E-10);
     211           1 : }
     212             : 
     213           1 : int main(int argc, char **argv) {
     214           1 :         ::testing::InitGoogleTest(&argc, argv);
     215             :         return RUN_ALL_TESTS();
     216             : }

Generated by: LCOV version 1.14