Line data    Source code

       1             : #include "crpropa/magneticField/ArchimedeanSpiralField.h"
       2             : 
       3             : namespace crpropa {
       4             : 
       5           0 : ArchimedeanSpiralField::ArchimedeanSpiralField(double B_0, double R_0, double Omega, double V_w) {
       6           0 :         setB0(B_0);
       7           0 :         setR0(R_0);
       8           0 :         setOmega(Omega);
       9           0 :         setVw(V_w);
      10           0 : }
      11             : 
      12           0 : Vector3d ArchimedeanSpiralField::getField(const Vector3d &pos) const {
      13             :         
      14             :         double r = pos.getR();
      15           0 :         double theta = pos.getTheta();
      16           0 :         double phi =pos.getPhi();
      17             :         
      18           0 :         double cos_phi = cos(phi);
      19           0 :         double sin_phi = sin(phi);
      20           0 :         double cos_theta = cos(theta);
      21           0 :         double sin_theta = sin(theta);
      22             : 
      23             :         Vector3d B(0.);
      24             : 
      25             :         // radial direction
      26           0 :         double C1 = R_0*R_0/r/r;
      27           0 :         B.x += C1 * cos_phi * sin_theta;
      28           0 :         B.y += C1 * sin_phi * sin_theta;
      29           0 :         B.z += C1 * cos_theta;
      30             :         
      31             :         // azimuthal direction  
      32           0 :         double C2 = - (Omega*R_0*R_0*sin_theta) / (r*V_w);
      33           0 :         B.x += C2 * (-sin_phi);
      34           0 :         B.y += C2 * cos_phi;
      35             : 
      36             :         // magnetic field switch at z = 0
      37           0 :         if (pos.z<0.) {
      38             :                 B *= -1;
      39             :         }
      40             : 
      41             :         // overall scaling
      42             :         B *= B_0;
      43             :         
      44             : 
      45           0 :         return B;
      46             : }
      47             : 
      48           0 : void ArchimedeanSpiralField::setR0(double R) {
      49           0 :         R_0 = R;
      50           0 :         return;
      51             : }
      52             : 
      53           0 : void ArchimedeanSpiralField::setB0(double B) {
      54           0 :         B_0 = B;
      55           0 :         return;
      56             : }
      57             : 
      58           0 : void ArchimedeanSpiralField::setOmega(double Om) {
      59           0 :         Omega = Om;
      60           0 :         return;
      61             : }
      62             : 
      63           0 : void ArchimedeanSpiralField::setVw(double v) {
      64           0 :         V_w = v;
      65           0 :         return;
      66             : }
      67             : 
      68             : 
      69           0 : double ArchimedeanSpiralField::getR0() const {
      70           0 :         return R_0;
      71             : }
      72             : 
      73           0 : double ArchimedeanSpiralField::getB0() const{
      74           0 :         return B_0;
      75             : }
      76             : 
      77           0 : double ArchimedeanSpiralField::getOmega() const{
      78           0 :         return Omega;
      79             : }
      80             : 
      81           0 : double ArchimedeanSpiralField::getVw() const {
      82           0 :         return V_w;
      83             : }
      84             : 
      85             : 
      86             : } //end namespace crpropa