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 : }
|