European Space Agency — European Space Agency¶

Royal Belgian Institute for Aeronomy¶

../_images/whiteSpace.png

T.06 - What is the precision of the field line tracing algorithm ?¶

Answer¶

The evaluation of the UNILIB library precision, is performed by two shell tracings:

with a centered aligned magnetic dipole;

with the same magnetic field as the Hassitt code RHOBAR.

Aligned dipole¶

For the centered aligned magnetic dipole, only the element g10 of the spherical harmonic coefficients is different from zero. In this case the value of g10 is chosen such that the moment of the magnetic dipole is equal to the value used by McIlwain, i.e. 0.311653 G Re^-3. In the library, the drift shell tracing is controlled by the arguments prop, stepx, stplst, epsrel, epsfl and epslon of the common blocks uc190 and uc192.

As a test case, the drift shell ( B=0.3, L=1.5) has been computed by the subroutine ud310() (with knfl set to 360) for two cases:

when the default values of the controlling parameters (prop = 0.2, stepx = 0.075, stplst = 15.0, epsrel = 6.0×10^-6, epsfl = 0.001 and epslon = 0.05) are used;

when the argument prop, which controls the step size along the field line, is divided by two (prop = 0.1).

In each case, the 360 magnetic field line segments are identical due to the axial symmetry. The geocentric spherical coordinates of the computed points for the magnetic field line segment at 90° of longitude are displayed in the table below for both cases. The segments computed with prop = 0.2 include 23 points while the segments computed when prop = 0.1 include 35 points. 19 values have been theoretically calculated for the same field line and are also given in the table.

prop = 0.2

prop = 0.1

theory

rho (Re)

theta (rad)

B (Gauss)

rho (Re)

theta (rad)

B (Gauss)

rho (Re)

theta (rad)

B (Gauss)

1.1139

1.03948

0.300000

1.1155

1.04071

0.298428

1.1354

1.05604

0.279835

1.1899

1.09964

0.235339

1.2392

1.14168

0.201867

1.2836

1.18227

0.176261

1.1139

1.03946

0.300000

1.3231

1.22153

0.156399

1.1155

1.04069

0.298428

1.3580

1.25948

0.140867

1.1465

1.06468

0.270055

1.3886

1.29635

0.128593

1.1141

1.03884

0.300000

1.2120

1.11810

0.219634

1.4150

1.33228

0.118861

1.1889

1.09794

0.236219

1.2761

1.17516

0.180326

1.4375

1.36741

0.111155

1.2575

1.15705

0.190964

1.3370

1.23608

0.150045

1.4563

1.40188

0.105101

1.3191

1.21616

0.158439

1.3922

1.30095

0.127222

1.4715

1.43583

0.100432

1.3728

1.27526

0.134933

1.4388

1.36959

0.110720

1.4833

1.46942

0.096955

1.4177

1.33437

0.118032

1.4718

1.43652

0.100341

1.4917

1.50278

0.094537

1.4533

1.39348

0.106159

1.4918

1.50282

0.094527

1.4968

1.53606

0.093098

1.4791

1.45258

0.098293

1.4987

1.56938

0.092587

1.4986

1.56938

0.092595

1.4948

1.51169

0.093802

1.4987

1.57080

0.092586

1.4986

1.57080

0.092594

1.5000

1.57080

0.092342

1.4987

1.57221

0.092587

1.4986

1.57221

0.092595

1.4948

1.62990

0.093802

1.4918

1.63877

0.094527

1.4968

1.60553

0.093098

1.4791

1.68901

0.098293

1.4718

1.70507

0.100341

1.4917

1.63881

0.094537

1.4533

1.74812

0.106159

1.4388

1.77200

0.110720

1.4833

1.67217

0.096955

1.4177

1.80722

0.118032

1.3922

1.84065

0.127222

1.4715

1.70576

0.100432

1.3728

1.86633

0.134933

1.3370

1.90551

0.150045

1.4563

1.73971

0.105101

1.3191

1.92543

0.158439

1.2761

1.96643

0.180326

1.4375

1.77418

0.111155

1.2575

1.98454

0.190964

1.2120

2.02350

0.219634

1.4150

1.80931

0.118861

1.1889

2.04365

0.236219

1.1465

2.07691

0.270055

1.3886

1.84524

0.128593

1.1141

2.10275

0.300000

1.1155

2.10090

0.298428

1.3580

1.88212

0.140867

1.1139

2.10213

0.300000

1.3231

1.92006

0.156399

1.2836

1.95932

0.176261

1.2392

1.99991

0.201867

1.1899

2.04195

0.235339

1.1354

2.08556

0.279835

1.1155

2.10088

0.298428

1.1139

2.10211

0.300000

Results obtained with UNILIB v1.11

The discrepancies at the equatorial point reflect mainly the dimness in the evaluation of L. The second table below shows the values of the McIlwain parameter obtained by the library for field line segments limited by the magnetic field intensity B = 0.3 G as a function of the crossing distance at the equatorial plane. The results are obtained at 90° of longitude with the default controlling parameters (prop = 0.2) by using subroutines uf420(), ul230() and ul240().

Equatorial crossing

I

L

B equatorial

1.45000

6632.9

1.44949

0.102227

1.50000

7436.3

1.49975

0.092342

1.50025

7440.3

1.50000

0.092295

1.55000

8237.0

1.54950

0.083691

2.00000

15661.3

1.99984

0.038957

2.50000

24126.1

2.49998

0.019946

For each equatorial crossing distance, a slight systematic error appears in the computation of L. This error is more important at lower L. Nevertheless one should note that at L = 1.5, the error is only about 1.5 km.

The value of the Roederer shell parameter associated with the drift shell ( B = 0.3, L = 1.5) obtained by the library (subroutine ud330()) with the default controlling parameters is equal to 1.4986. This value has to be compared with the crossing distance at the equator, which is equal to 1.4987. We assume therefore that the slight systematic error observed in the tables above for the value of the McIlwain parameter is mainly due to the subroutines ul240() and ul242().

Comparison with RHOBAR¶

RHOBAR is a Fortran code developped by A. Hassitt and C.E. McIlwain in the sixties. This programme can be used to evaluate an average atmospheric densities over a magnetic drift shell. The programme is restricted to the Jensen & Cain geomagnetic field model. The strategies of the programme RHOBAR and of the library UNILIB to trace a magnetic drift shell widely differs. In UNILIB, the field line segments computed by subroutine ud310() are equally spaced in longitude and the number of segments is controlled by the user, e.g. 270. In the programme RHOBAR, the first few segments are equally spaced in longitude. The other segments are only traced at longitudes where the average densities is the more varying, i.e. at longitudes where the magnetic field lines reach lower altitudes. The programme RHOBAR makes also use of interpolation to determine the altitude of the mirror points. In consequence, the RHOBAR is better optimized to search the mirror point with the lowest altitude on a magnetic drift shell than the UNILIB library. On the other hand, the programme RHOBAR is restricted to the evaluation of the average density for a given magnetic field and atmospheric model, and, a part of its algorithm seems to be especially tuned for the Jensen and Cain geomagnetic field model.

To evaluate the precision of the drift shell tracing, the mirror points with the lowest altitude have been searched for different magnetic drift shells using the subroutine ud310() with 120, 270 or 360 field line segments, and using the programme RHOBAR. The results are displayed in the table below. The UNILIB library has been initialized with the same geomagnetic field model as used in the programme RHOBAR.

B

L

UD310, 120

UD310, 270

UD310, 360

RHOBAR

0.40

04.601

15.9

113.4

113.3

122.6

0.35

03.601

88.9

181.4

171.0

184.2

0.27

52.101

94.9

203.4

218.5

223.9

0.37

54.602

50.6

254.5

251.2

267.8

0.35

04.102

94.7

315.8

311.8

321.7

0.35

04.604

01.3

420.5

419.3

428.0

0.22

51.604

42.8

460.6

474.0

481.4

0.25

02.104

89.5

489.8

489.6

487.9

0.35

05.104

86.0

494.3

504.9

512.4

0.32

54.606

07.5

603.9

595.8

605.4

0.22

52.107

42.4

772.6

764.7

777.8

0.22

52.609

77.8

1002.7

1002.9

1003.2

0.20

02.101

088.5

1090.9

1099.7

1111.5

0.22

53.101

197.3

1175.5

1180.8

1195.2

The UNILIB results are obtained with the version 1.12 of the library. The discrepancies that remain sometimes between the results are mainly due to an incorrect use of the Jensen & Cain magnetic field model (Heynderickx et al., 1994) by the RHOBAR programme. This error gets its main effect at the mid latitudes.

References¶

Heynderickx, D., Lemaire, J., and Daly, E.J., Historical review of the different procedures used to compute the L parameter, Aeronomica Acta A-Nb 380 (1994), J. Nuclear Tracks Radiat. Meas. (1995)

Illustration¶

The subroutine faqt06 is used to initialize the UNILIB library with a centered magnetic field dipole aligned with the Earth’s axis.

      SUBROUTINE faqt06 (kunit, ifail)
C
C     Initialization of the UNILIB library with a
C     aligned centered magnetic dipole
C
      INTEGER*4      kunit, ifail
C
      include 'unilib.h'
C
      COMMON /UC140/  mint, mext, msun
C
      TYPE(zimf) mint
      TYPE(zsun) msun
      TYPE(zemf) mext
C
      COMMON /UC160/ pi, deg, re, gmagmo, eclipt, geoid, uma
C
      REAL*8         pi, deg, re
      REAL*8         gmagmo
      REAL*8         eclipt, geoid(3), uma(30)
C
      REAL*8         amjd, param(10)
      INTEGER*4      kext, kinit, iver
      CHARACTER*32   lbext
C
      DATA kinit, kext, amjd, param/ 1, 0, 0.0d0, 10*0.0d0/
C
      ifail       = 0
C
C     Initialize UNILIB
C
      CALL UT990 ( kunit, kinit, ifail )
      if( ifail.lt.0 )return
      iver        = ifail
C
C     Initialize the geomagnetic field model
C
      mint.kinner = -1
C                    any value less than 4 and not equal to 1
      mint.label  = 'Aligned Dipolar magnetic field  '
      mint.tzero  = 1960.0d0
      mint.epoch  = 1960.0d0
      mint.saarot = 0.0d0
      mint.elong  = 0.0d0
      mint.colat  = 0.0d0
C
      mint.norder = 2
      mint.gmmo   = gmagmo
      mint.coef(2,1) = 1.0D+5 * gmagmo
      mint.coef(1,2) = 0.0d0
      mint.coef(2,2) = 0.0d0
C
C     Simulate the outputs of UM510
C
      if( kunit.gt.0 ) write(kunit,1000) mint.label,
   :            mint.tzero, mint.norder, mint.epoch,
   :            mint.colat, mint.elong, mint.gmmo,
   :            mint.saarot
C
C     Initialize the external magnetic field model
C
      CALL UM520 (kext, amjd, param, lbext, kunit, ifail)
      if( ifail.lt.0 )return
C
C     Format
C
1000 format(/' --- Geomagnetic field model ---',
   : //6x,'Model (sp): ',a32,3x,'Epoch',f6.0,/29x,
   : 'Order + 1 =',i5,/21x,'Calculation epoch =',f7.1,5x,' year',/9x,
   :           'Colatitude of the dipole pole =',f8.2,4x,' deg',/10x,
   :            'Longitude of the dipole pole =',f8.2,4x,' deg',/19x,
   :                         'Earth dipole moment =',f12.6,' G/Re^3',
   :        /10x,'Correction for the SAA drift =',f8.2,4x,' deg')
C
      end

See also¶

T.04 - Why are Roederer’s L* and McIlwain’s L parameters different ?

Ex #2. Evaluation of the (B-L) coordinates

Ex #3. Evaluation of the third invariant

UNILIB/tags/v3.02