Nutation, IAU 2000A model (MHB2000 luni-solar and planetary nutation with free core nutation omitted).
Status: canonical model.
The TT date DATE1+DATE2 is a Julian Date, apportioned in any convenient way between the two arguments. For example, JD(TT)=2450123.7 could be expressed in any of these ways, among others:
DATE1 DATE2
2450123.7D0 0D0 (JD method)
2451545D0 -1421.3D0 (J2000 method)
2400000.5D0 50123.2D0 (MJD method)
2450123.5D0 0.2D0 (date & time method)
The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
The nutation components in longitude and obliquity are in radians and with respect to the equinox and ecliptic of date. The obliquity at J2000.0 is assumed to be the Lieske et al. (1977) value of 84381.448 arcsec.
Both the luni-solar and planetary nutations are included. The latter are due to direct planetary nutations and the perturbations of the lunar and terrestrial orbits.
The routine computes the MHB2000 nutation series with the associated corrections for planetary nutations. It is an implementation of the nutation part of the IAU 2000A precession- nutation model, formally adopted by the IAU General Assembly in 2000, namely MHB2000 (Mathews et al. 2002), but with the free core nutation (FCN - see Note 4) omitted.
The full MHB2000 model also contains contributions to the nutations in longitude and obliquity due to the free-excitation of the free-core-nutation during the period 1979-2000. These FCN terms, which are time-dependent and unpredictable, are NOT included in the present routine and, if required, must be independently computed. With the FCN corrections included, the present routine delivers a pole which is at current epochs accurate to a few hundred microarcseconds. The omission of FCN introduces further errors of about that size.
The present routine provides classical nutation. The MHB2000 algorithm, from which it is adapted, deals also with (i) the offsets between the GCRS and mean poles and (ii) the adjustments in longitude and obliquity due to the changed precession rates. These additional functions, namely frame bias and precession adjustments, are supported by the SOFA routines BI00 and PR00.
The MHB2000 algorithm also provides "total" nutations, comprising the arithmetic sum of the frame bias, precession adjustments, luni-solar nutation and planetary nutation. These total nutations can be used in combination with an existing IAU 1976 precession implementation, such as PMAT76, to deliver GCRS-to-true predictions of sub-mas accuracy at current epochs. However, there are three shortcomings in the MHB2000 model that must be taken into account if more accurate or definitive results are required (see Wallace 2002):
(i) The MHB2000 total nutations are simply arithmetic sums, yet in reality the various components are successive Euler rotations. This slight lack of rigor leads to cross terms that exceed 1 mas after a century. The rigorous procedure is to form the GCRS-to-true rotation matrix by applying the bias, precession and nutation in that order.
(ii) Although the precession adjustments are stated to be with respect to Lieske et al. (1977), the MHB2000 model does not specify which set of Euler angles are to be used and how the adjustments are to be applied. The most literal and straightforward procedure is to adopt the 4-rotation epsilon_0, psi_A, omega_A, xi_A option, and to add DPSIPR to psi_A and DEPSPR to both omega_A and eps_A.
(iii) The MHB2000 model predates the determination by Chapront et al. (2002) of a 14.6 mas displacement between the J2000.0 mean equinox and the origin of the ICRS frame. It should, however, be noted that neglecting this displacement when calculating star coordinates does not lead to a 14.6 mas change in right ascension, only a small second-order distortion in the pattern of the precession-nutation effect.
For these reasons, the SOFA routines do not generate the "total nutations" directly, though they can of course easily be generated by calling BI00, PR00 and the present routine and adding the results.
The MHB2000 model contains 41 instances where the same frequency appears multiple times, of which 38 are duplicates and three are triplicates. To keep the present code close to the original MHB algorithm, this small inefficiency has not been corrected.
Chapront, J., Chapront-Touze, M. & Francou, G. 2002, Astron.Astrophys. 387, 700
Lieske, J.H., Lederle, T., Fricke, W. & Morando, B. 1977, Astron.Astrophys. 58, 1-16
Mathews, P.M., Herring, T.A., Buffet, B.A. 2002, J.Geophys.Res. 107, B4. The MHB_2000 code itself was obtained on 9th September 2002 from ftp//maia.usno.navy.mil/conv2000/chapter5/IAU2000A.
Simon, J.-L., Bretagnon, P., Chapront, J., Chapront-Touze, M., Francou, G., Laskar, J. 1994, Astron.Astrophys. 282, 663-683
Souchay, J., Loysel, B., Kinoshita, H., Folgueira, M. 1999, Astron.Astrophys.Supp.Ser. 135, 111
Wallace, P.T., "Software for Implementing the IAU 2000 Resolutions", in IERS Workshop 5.1 (2002)
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
real(kind=wp), | intent(in) | :: | date1 | TT as a 2-part Julian Date (Note 1) |
||
real(kind=wp), | intent(in) | :: | date2 | TT as a 2-part Julian Date (Note 1) |
||
real(kind=wp), | intent(out) | :: | dpsi | nutation, luni-solar + planetary (Note 2) |
||
real(kind=wp), | intent(out) | :: | deps | nutation, luni-solar + planetary (Note 2) |
subroutine NUT00A ( date1, date2, dpsi, deps )
implicit none
real(wp),intent(in) :: date1 !! TT as a 2-part Julian Date (Note 1)
real(wp),intent(in) :: date2 !! TT as a 2-part Julian Date (Note 1)
real(wp),intent(out) :: dpsi !! nutation, luni-solar + planetary (Note 2)
real(wp),intent(out) :: deps !! nutation, luni-solar + planetary (Note 2)
! Arcseconds in a full circle
real(wp),parameter :: turnas = 1296000.0_wp
! Units of 0.1 microarcsecond to radians
real(wp),parameter :: u2r = das2r/1.0e7_wp
! Miscellaneous
integer :: i, j
real(wp) :: t, el, elp, f, d, om, arg, dp, de, sarg, carg, &
al, alsu, af, ad, aom, alme, alve, alea, alma, &
alju, alsa, alur, alne, apa, dpsils, depsls, &
dpsipl, depspl
! -------------------------
! Luni-Solar nutation model
! -------------------------
! Number of terms in the luni-solar nutation model
integer,parameter :: nls = 678
! ---------------
! Planetary terms
! ---------------
! Number of terms in the planetary nutation model
integer,parameter :: npl = 687
! ----------------------------------------
! Tables of argument and term coefficients
! ----------------------------------------
!
! Luni-Solar argument multipliers
! L L' F D Om
! Coefficients for fundamental arguments
integer,dimension(5,nls),parameter :: nals = reshape([&
0, 0, 0, 0, 1, &
0, 0, 2, -2, 2, &
0, 0, 2, 0, 2, &
0, 0, 0, 0, 2, &
0, 1, 0, 0, 0, &
0, 1, 2, -2, 2, &
1, 0, 0, 0, 0, &
0, 0, 2, 0, 1, &
1, 0, 2, 0, 2, &
0, -1, 2, -2, 2, &
0, 0, 2, -2, 1, &
-1, 0, 2, 0, 2, &
-1, 0, 0, 2, 0, &
1, 0, 0, 0, 1, &
-1, 0, 0, 0, 1, &
-1, 0, 2, 2, 2, &
1, 0, 2, 0, 1, &
-2, 0, 2, 0, 1, &
0, 0, 0, 2, 0, &
0, 0, 2, 2, 2, &
0, -2, 2, -2, 2, &
-2, 0, 0, 2, 0, &
2, 0, 2, 0, 2, &
1, 0, 2, -2, 2, &
-1, 0, 2, 0, 1, &
2, 0, 0, 0, 0, &
0, 0, 2, 0, 0, &
0, 1, 0, 0, 1, &
-1, 0, 0, 2, 1, &
0, 2, 2, -2, 2, &
0, 0, -2, 2, 0, &
1, 0, 0, -2, 1, &
0, -1, 0, 0, 1, &
-1, 0, 2, 2, 1, &
0, 2, 0, 0, 0, &
1, 0, 2, 2, 2, &
-2, 0, 2, 0, 0, &
0, 1, 2, 0, 2, &
0, 0, 2, 2, 1, &
0, -1, 2, 0, 2, &
0, 0, 0, 2, 1, &
1, 0, 2, -2, 1, &
2, 0, 2, -2, 2, &
-2, 0, 0, 2, 1, &
2, 0, 2, 0, 1, &
0, -1, 2, -2, 1, &
0, 0, 0, -2, 1, &
-1, -1, 0, 2, 0, &
2, 0, 0, -2, 1, &
1, 0, 0, 2, 0, &
0, 1, 2, -2, 1, &
1, -1, 0, 0, 0, &
-2, 0, 2, 0, 2, &
3, 0, 2, 0, 2, &
0, -1, 0, 2, 0, &
1, -1, 2, 0, 2, &
0, 0, 0, 1, 0, &
-1, -1, 2, 2, 2, &
-1, 0, 2, 0, 0, &
0, -1, 2, 2, 2, &
-2, 0, 0, 0, 1, &
1, 1, 2, 0, 2, &
2, 0, 0, 0, 1, &
-1, 1, 0, 1, 0, &
1, 1, 0, 0, 0, &
1, 0, 2, 0, 0, &
-1, 0, 2, -2, 1, &
1, 0, 0, 0, 2, &
-1, 0, 0, 1, 0, &
0, 0, 2, 1, 2, &
-1, 0, 2, 4, 2, &
-1, 1, 0, 1, 1, &
0, -2, 2, -2, 1, &
1, 0, 2, 2, 1, &
-2, 0, 2, 2, 2, &
-1, 0, 0, 0, 2, &
1, 1, 2, -2, 2, &
-2, 0, 2, 4, 2, &
-1, 0, 4, 0, 2, &
2, 0, 2, -2, 1, &
2, 0, 2, 2, 2, &
1, 0, 0, 2, 1, &
3, 0, 0, 0, 0, &
3, 0, 2, -2, 2, &
0, 0, 4, -2, 2, &
0, 1, 2, 0, 1, &
0, 0, -2, 2, 1, &
0, 0, 2, -2, 3, &
-1, 0, 0, 4, 0, &
2, 0, -2, 0, 1, &
-2, 0, 0, 4, 0, &
-1, -1, 0, 2, 1, &
-1, 0, 0, 1, 1, &
0, 1, 0, 0, 2, &
0, 0, -2, 0, 1, &
0, -1, 2, 0, 1, &
0, 0, 2, -1, 2, &
0, 0, 2, 4, 2, &
-2, -1, 0, 2, 0, &
1, 1, 0, -2, 1, &
-1, 1, 0, 2, 0, &
-1, 1, 0, 1, 2, &
1, -1, 0, 0, 1, &
1, -1, 2, 2, 2, &
-1, 1, 2, 2, 2, &
3, 0, 2, 0, 1, &
0, 1, -2, 2, 0, &
-1, 0, 0, -2, 1, &
0, 1, 2, 2, 2, &
-1, -1, 2, 2, 1, &
0, -1, 0, 0, 2, &
1, 0, 2, -4, 1, &
-1, 0, -2, 2, 0, &
0, -1, 2, 2, 1, &
2, -1, 2, 0, 2, &
0, 0, 0, 2, 2, &
1, -1, 2, 0, 1, &
-1, 1, 2, 0, 2, &
0, 1, 0, 2, 0, &
0, -1, -2, 2, 0, &
0, 3, 2, -2, 2, &
0, 0, 0, 1, 1, &
-1, 0, 2, 2, 0, &
2, 1, 2, 0, 2, &
1, 1, 0, 0, 1, &
1, 1, 2, 0, 1, &
2, 0, 0, 2, 0, &
1, 0, -2, 2, 0, &
-1, 0, 0, 2, 2, &
0, 1, 0, 1, 0, &
0, 1, 0, -2, 1, &
-1, 0, 2, -2, 2, &
0, 0, 0, -1, 1, &
-1, 1, 0, 0, 1, &
1, 0, 2, -1, 2, &
1, -1, 0, 2, 0, &
0, 0, 0, 4, 0, &
1, 0, 2, 1, 2, &
0, 0, 2, 1, 1, &
1, 0, 0, -2, 2, &
-1, 0, 2, 4, 1, &
1, 0, -2, 0, 1, &
1, 1, 2, -2, 1, &
0, 0, 2, 2, 0, &
-1, 0, 2, -1, 1, &
-2, 0, 2, 2, 1, &
4, 0, 2, 0, 2, &
2, -1, 0, 0, 0, &
2, 1, 2, -2, 2, &
0, 1, 2, 1, 2, &
1, 0, 4, -2, 2, &
-1, -1, 0, 0, 1, &
0, 1, 0, 2, 1, &
-2, 0, 2, 4, 1, &
2, 0, 2, 0, 0, &
1, 0, 0, 1, 0, &
-1, 0, 0, 4, 1, &
-1, 0, 4, 0, 1, &
2, 0, 2, 2, 1, &
0, 0, 2, -3, 2, &
-1, -2, 0, 2, 0, &
2, 1, 0, 0, 0, &
0, 0, 4, 0, 2, &
0, 0, 0, 0, 3, &
0, 3, 0, 0, 0, &
0, 0, 2, -4, 1, &
0, -1, 0, 2, 1, &
0, 0, 0, 4, 1, &
-1, -1, 2, 4, 2, &
1, 0, 2, 4, 2, &
-2, 2, 0, 2, 0, &
-2, -1, 2, 0, 1, &
-2, 0, 0, 2, 2, &
-1, -1, 2, 0, 2, &
0, 0, 4, -2, 1, &
3, 0, 2, -2, 1, &
-2, -1, 0, 2, 1, &
1, 0, 0, -1, 1, &
0, -2, 0, 2, 0, &
-2, 0, 0, 4, 1, &
-3, 0, 0, 0, 1, &
1, 1, 2, 2, 2, &
0, 0, 2, 4, 1, &
3, 0, 2, 2, 2, &
-1, 1, 2, -2, 1, &
2, 0, 0, -4, 1, &
0, 0, 0, -2, 2, &
2, 0, 2, -4, 1, &
-1, 1, 0, 2, 1, &
0, 0, 2, -1, 1, &
0, -2, 2, 2, 2, &
2, 0, 0, 2, 1, &
4, 0, 2, -2, 2, &
2, 0, 0, -2, 2, &
0, 2, 0, 0, 1, &
1, 0, 0, -4, 1, &
0, 2, 2, -2, 1, &
-3, 0, 0, 4, 0, &
-1, 1, 2, 0, 1, &
-1, -1, 0, 4, 0, &
-1, -2, 2, 2, 2, &
-2, -1, 2, 4, 2, &
1, -1, 2, 2, 1, &
-2, 1, 0, 2, 0, &
-2, 1, 2, 0, 1, &
2, 1, 0, -2, 1, &
-3, 0, 2, 0, 1, &
-2, 0, 2, -2, 1, &
-1, 1, 0, 2, 2, &
0, -1, 2, -1, 2, &
-1, 0, 4, -2, 2, &
0, -2, 2, 0, 2, &
-1, 0, 2, 1, 2, &
2, 0, 0, 0, 2, &
0, 0, 2, 0, 3, &
-2, 0, 4, 0, 2, &
-1, 0, -2, 0, 1, &
-1, 1, 2, 2, 1, &
3, 0, 0, 0, 1, &
-1, 0, 2, 3, 2, &
2, -1, 2, 0, 1, &
0, 1, 2, 2, 1, &
0, -1, 2, 4, 2, &
2, -1, 2, 2, 2, &
0, 2, -2, 2, 0, &
-1, -1, 2, -1, 1, &
0, -2, 0, 0, 1, &
1, 0, 2, -4, 2, &
1, -1, 0, -2, 1, &
-1, -1, 2, 0, 1, &
1, -1, 2, -2, 2, &
-2, -1, 0, 4, 0, &
-1, 0, 0, 3, 0, &
-2, -1, 2, 2, 2, &
0, 2, 2, 0, 2, &
1, 1, 0, 2, 0, &
2, 0, 2, -1, 2, &
1, 0, 2, 1, 1, &
4, 0, 0, 0, 0, &
2, 1, 2, 0, 1, &
3, -1, 2, 0, 2, &
-2, 2, 0, 2, 1, &
1, 0, 2, -3, 1, &
1, 1, 2, -4, 1, &
-1, -1, 2, -2, 1, &
0, -1, 0, -1, 1, &
0, -1, 0, -2, 1, &
-2, 0, 0, 0, 2, &
-2, 0, -2, 2, 0, &
-1, 0, -2, 4, 0, &
1, -2, 0, 0, 0, &
0, 1, 0, 1, 1, &
-1, 2, 0, 2, 0, &
1, -1, 2, -2, 1, &
1, 2, 2, -2, 2, &
2, -1, 2, -2, 2, &
1, 0, 2, -1, 1, &
2, 1, 2, -2, 1, &
-2, 0, 0, -2, 1, &
1, -2, 2, 0, 2, &
0, 1, 2, 1, 1, &
1, 0, 4, -2, 1, &
-2, 0, 4, 2, 2, &
1, 1, 2, 1, 2, &
1, 0, 0, 4, 0, &
1, 0, 2, 2, 0, &
2, 0, 2, 1, 2, &
3, 1, 2, 0, 2, &
4, 0, 2, 0, 1, &
-2, -1, 2, 0, 0, &
0, 1, -2, 2, 1, &
1, 0, -2, 1, 0, &
0, -1, -2, 2, 1, &
2, -1, 0, -2, 1, &
-1, 0, 2, -1, 2, &
1, 0, 2, -3, 2, &
0, 1, 2, -2, 3, &
0, 0, 2, -3, 1, &
-1, 0, -2, 2, 1, &
0, 0, 2, -4, 2, &
-2, 1, 0, 0, 1, &
-1, 0, 0, -1, 1, &
2, 0, 2, -4, 2, &
0, 0, 4, -4, 4, &
0, 0, 4, -4, 2, &
-1, -2, 0, 2, 1, &
-2, 0, 0, 3, 0, &
1, 0, -2, 2, 1, &
-3, 0, 2, 2, 2, &
-3, 0, 2, 2, 1, &
-2, 0, 2, 2, 0, &
2, -1, 0, 0, 1, &
-2, 1, 2, 2, 2, &
1, 1, 0, 1, 0, &
0, 1, 4, -2, 2, &
-1, 1, 0, -2, 1, &
0, 0, 0, -4, 1, &
1, -1, 0, 2, 1, &
1, 1, 0, 2, 1, &
-1, 2, 2, 2, 2, &
3, 1, 2, -2, 2, &
0, -1, 0, 4, 0, &
2, -1, 0, 2, 0, &
0, 0, 4, 0, 1, &
2, 0, 4, -2, 2, &
-1, -1, 2, 4, 1, &
1, 0, 0, 4, 1, &
1, -2, 2, 2, 2, &
0, 0, 2, 3, 2, &
-1, 1, 2, 4, 2, &
3, 0, 0, 2, 0, &
-1, 0, 4, 2, 2, &
1, 1, 2, 2, 1, &
-2, 0, 2, 6, 2, &
2, 1, 2, 2, 2, &
-1, 0, 2, 6, 2, &
1, 0, 2, 4, 1, &
2, 0, 2, 4, 2, &
1, 1, -2, 1, 0, &
-3, 1, 2, 1, 2, &
2, 0, -2, 0, 2, &
-1, 0, 0, 1, 2, &
-4, 0, 2, 2, 1, &
-1, -1, 0, 1, 0, &
0, 0, -2, 2, 2, &
1, 0, 0, -1, 2, &
0, -1, 2, -2, 3, &
-2, 1, 2, 0, 0, &
0, 0, 2, -2, 4, &
-2, -2, 0, 2, 0, &
-2, 0, -2, 4, 0, &
0, -2, -2, 2, 0, &
1, 2, 0, -2, 1, &
3, 0, 0, -4, 1, &
-1, 1, 2, -2, 2, &
1, -1, 2, -4, 1, &
1, 1, 0, -2, 2, &
-3, 0, 2, 0, 0, &
-3, 0, 2, 0, 2, &
-2, 0, 0, 1, 0, &
0, 0, -2, 1, 0, &
-3, 0, 0, 2, 1, &
-1, -1, -2, 2, 0, &
0, 1, 2, -4, 1, &
2, 1, 0, -4, 1, &
0, 2, 0, -2, 1, &
1, 0, 0, -3, 1, &
-2, 0, 2, -2, 2, &
-2, -1, 0, 0, 1, &
-4, 0, 0, 2, 0, &
1, 1, 0, -4, 1, &
-1, 0, 2, -4, 1, &
0, 0, 4, -4, 1, &
0, 3, 2, -2, 2, &
-3, -1, 0, 4, 0, &
-3, 0, 0, 4, 1, &
1, -1, -2, 2, 0, &
-1, -1, 0, 2, 2, &
1, -2, 0, 0, 1, &
1, -1, 0, 0, 2, &
0, 0, 0, 1, 2, &
-1, -1, 2, 0, 0, &
1, -2, 2, -2, 2, &
0, -1, 2, -1, 1, &
-1, 0, 2, 0, 3, &
1, 1, 0, 0, 2, &
-1, 1, 2, 0, 0, &
1, 2, 0, 0, 0, &
-1, 2, 2, 0, 2, &
-1, 0, 4, -2, 1, &
3, 0, 2, -4, 2, &
1, 2, 2, -2, 1, &
1, 0, 4, -4, 2, &
-2, -1, 0, 4, 1, &
0, -1, 0, 2, 2, &
-2, 1, 0, 4, 0, &
-2, -1, 2, 2, 1, &
2, 0, -2, 2, 0, &
1, 0, 0, 1, 1, &
0, 1, 0, 2, 2, &
1, -1, 2, -1, 2, &
-2, 0, 4, 0, 1, &
2, 1, 0, 0, 1, &
0, 1, 2, 0, 0, &
0, -1, 4, -2, 2, &
0, 0, 4, -2, 4, &
0, 2, 2, 0, 1, &
-3, 0, 0, 6, 0, &
-1, -1, 0, 4, 1, &
1, -2, 0, 2, 0, &
-1, 0, 0, 4, 2, &
-1, -2, 2, 2, 1, &
-1, 0, 0, -2, 2, &
1, 0, -2, -2, 1, &
0, 0, -2, -2, 1, &
-2, 0, -2, 0, 1, &
0, 0, 0, 3, 1, &
0, 0, 0, 3, 0, &
-1, 1, 0, 4, 0, &
-1, -1, 2, 2, 0, &
-2, 0, 2, 3, 2, &
1, 0, 0, 2, 2, &
0, -1, 2, 1, 2, &
3, -1, 0, 0, 0, &
2, 0, 0, 1, 0, &
1, -1, 2, 0, 0, &
0, 0, 2, 1, 0, &
1, 0, 2, 0, 3, &
3, 1, 0, 0, 0, &
3, -1, 2, -2, 2, &
2, 0, 2, -1, 1, &
1, 1, 2, 0, 0, &
0, 0, 4, -1, 2, &
1, 2, 2, 0, 2, &
-2, 0, 0, 6, 0, &
0, -1, 0, 4, 1, &
-2, -1, 2, 4, 1, &
0, -2, 2, 2, 1, &
0, -1, 2, 2, 0, &
-1, 0, 2, 3, 1, &
-2, 1, 2, 4, 2, &
2, 0, 0, 2, 2, &
2, -2, 2, 0, 2, &
-1, 1, 2, 3, 2, &
3, 0, 2, -1, 2, &
4, 0, 2, -2, 1, &
-1, 0, 0, 6, 0, &
-1, -2, 2, 4, 2, &
-3, 0, 2, 6, 2, &
-1, 0, 2, 4, 0, &
3, 0, 0, 2, 1, &
3, -1, 2, 0, 1, &
3, 0, 2, 0, 0, &
1, 0, 4, 0, 2, &
5, 0, 2, -2, 2, &
0, -1, 2, 4, 1, &
2, -1, 2, 2, 1, &
0, 1, 2, 4, 2, &
1, -1, 2, 4, 2, &
3, -1, 2, 2, 2, &
3, 0, 2, 2, 1, &
5, 0, 2, 0, 2, &
0, 0, 2, 6, 2, &
4, 0, 2, 2, 2, &
0, -1, 1, -1, 1, &
-1, 0, 1, 0, 3, &
0, -2, 2, -2, 3, &
1, 0, -1, 0, 1, &
2, -2, 0, -2, 1, &
-1, 0, 1, 0, 2, &
-1, 0, 1, 0, 1, &
-1, -1, 2, -1, 2, &
-2, 2, 0, 2, 2, &
-1, 0, 1, 0, 0, &
-4, 1, 2, 2, 2, &
-3, 0, 2, 1, 1, &
-2, -1, 2, 0, 2, &
1, 0, -2, 1, 1, &
2, -1, -2, 0, 1, &
-4, 0, 2, 2, 0, &
-3, 1, 0, 3, 0, &
-1, 0, -1, 2, 0, &
0, -2, 0, 0, 2, &
0, -2, 0, 0, 2, &
-3, 0, 0, 3, 0, &
-2, -1, 0, 2, 2, &
-1, 0, -2, 3, 0, &
-4, 0, 0, 4, 0, &
2, 1, -2, 0, 1, &
2, -1, 0, -2, 2, &
0, 0, 1, -1, 0, &
-1, 2, 0, 1, 0, &
-2, 1, 2, 0, 2, &
1, 1, 0, -1, 1, &
1, 0, 1, -2, 1, &
0, 2, 0, 0, 2, &
1, -1, 2, -3, 1, &
-1, 1, 2, -1, 1, &
-2, 0, 4, -2, 2, &
-2, 0, 4, -2, 1, &
-2, -2, 0, 2, 1, &
-2, 0, -2, 4, 0, &
1, 2, 2, -4, 1, &
1, 1, 2, -4, 2, &
-1, 2, 2, -2, 1, &
2, 0, 0, -3, 1, &
-1, 2, 0, 0, 1, &
0, 0, 0, -2, 0, &
-1, -1, 2, -2, 2, &
-1, 1, 0, 0, 2, &
0, 0, 0, -1, 2, &
-2, 1, 0, 1, 0, &
1, -2, 0, -2, 1, &
1, 0, -2, 0, 2, &
-3, 1, 0, 2, 0, &
-1, 1, -2, 2, 0, &
-1, -1, 0, 0, 2, &
-3, 0, 0, 2, 0, &
-3, -1, 0, 2, 0, &
2, 0, 2, -6, 1, &
0, 1, 2, -4, 2, &
2, 0, 0, -4, 2, &
-2, 1, 2, -2, 1, &
0, -1, 2, -4, 1, &
0, 1, 0, -2, 2, &
-1, 0, 0, -2, 0, &
2, 0, -2, -2, 1, &
-4, 0, 2, 0, 1, &
-1, -1, 0, -1, 1, &
0, 0, -2, 0, 2, &
-3, 0, 0, 1, 0, &
-1, 0, -2, 1, 0, &
-2, 0, -2, 2, 1, &
0, 0, -4, 2, 0, &
-2, -1, -2, 2, 0, &
1, 0, 2, -6, 1, &
-1, 0, 2, -4, 2, &
1, 0, 0, -4, 2, &
2, 1, 2, -4, 2, &
2, 1, 2, -4, 1, &
0, 1, 4, -4, 4, &
0, 1, 4, -4, 2, &
-1, -1, -2, 4, 0, &
-1, -3, 0, 2, 0, &
-1, 0, -2, 4, 1, &
-2, -1, 0, 3, 0, &
0, 0, -2, 3, 0, &
-2, 0, 0, 3, 1, &
0, -1, 0, 1, 0, &
-3, 0, 2, 2, 0, &
1, 1, -2, 2, 0, &
-1, 1, 0, 2, 2, &
1, -2, 2, -2, 1, &
0, 0, 1, 0, 2, &
0, 0, 1, 0, 1, &
0, 0, 1, 0, 0, &
-1, 2, 0, 2, 1, &
0, 0, 2, 0, 2, &
-2, 0, 2, 0, 2, &
2, 0, 0, -1, 1, &
3, 0, 0, -2, 1, &
1, 0, 2, -2, 3, &
1, 2, 0, 0, 1, &
2, 0, 2, -3, 2, &
-1, 1, 4, -2, 2, &
-2, -2, 0, 4, 0, &
0, -3, 0, 2, 0, &
0, 0, -2, 4, 0, &
-1, -1, 0, 3, 0, &
-2, 0, 0, 4, 2, &
-1, 0, 0, 3, 1, &
2, -2, 0, 0, 0, &
1, -1, 0, 1, 0, &
-1, 0, 0, 2, 0, &
0, -2, 2, 0, 1, &
-1, 0, 1, 2, 1, &
-1, 1, 0, 3, 0, &
-1, -1, 2, 1, 2, &
0, -1, 2, 0, 0, &
-2, 1, 2, 2, 1, &
2, -2, 2, -2, 2, &
1, 1, 0, 1, 1, &
1, 0, 1, 0, 1, &
1, 0, 1, 0, 0, &
0, 2, 0, 2, 0, &
2, -1, 2, -2, 1, &
0, -1, 4, -2, 1, &
0, 0, 4, -2, 3, &
0, 1, 4, -2, 1, &
4, 0, 2, -4, 2, &
2, 2, 2, -2, 2, &
2, 0, 4, -4, 2, &
-1, -2, 0, 4, 0, &
-1, -3, 2, 2, 2, &
-3, 0, 2, 4, 2, &
-3, 0, 2, -2, 1, &
-1, -1, 0, -2, 1, &
-3, 0, 0, 0, 2, &
-3, 0, -2, 2, 0, &
0, 1, 0, -4, 1, &
-2, 1, 0, -2, 1, &
-4, 0, 0, 0, 1, &
-1, 0, 0, -4, 1, &
-3, 0, 0, -2, 1, &
0, 0, 0, 3, 2, &
-1, 1, 0, 4, 1, &
1, -2, 2, 0, 1, &
0, 1, 0, 3, 0, &
-1, 0, 2, 2, 3, &
0, 0, 2, 2, 2, &
-2, 0, 2, 2, 2, &
-1, 1, 2, 2, 0, &
3, 0, 0, 0, 2, &
2, 1, 0, 1, 0, &
2, -1, 2, -1, 2, &
0, 0, 2, 0, 1, &
0, 0, 3, 0, 3, &
0, 0, 3, 0, 2, &
-1, 2, 2, 2, 1, &
-1, 0, 4, 0, 0, &
1, 2, 2, 0, 1, &
3, 1, 2, -2, 1, &
1, 1, 4, -2, 2, &
-2, -1, 0, 6, 0, &
0, -2, 0, 4, 0, &
-2, 0, 0, 6, 1, &
-2, -2, 2, 4, 2, &
0, -3, 2, 2, 2, &
0, 0, 0, 4, 2, &
-1, -1, 2, 3, 2, &
-2, 0, 2, 4, 0, &
2, -1, 0, 2, 1, &
1, 0, 0, 3, 0, &
0, 1, 0, 4, 1, &
0, 1, 0, 4, 0, &
1, -1, 2, 1, 2, &
0, 0, 2, 2, 3, &
1, 0, 2, 2, 2, &
-1, 0, 2, 2, 2, &
-2, 0, 4, 2, 1, &
2, 1, 0, 2, 1, &
2, 1, 0, 2, 0, &
2, -1, 2, 0, 0, &
1, 0, 2, 1, 0, &
0, 1, 2, 2, 0, &
2, 0, 2, 0, 3, &
3, 0, 2, 0, 2, &
1, 0, 2, 0, 2, &
1, 0, 3, 0, 3, &
1, 1, 2, 1, 1, &
0, 2, 2, 2, 2, &
2, 1, 2, 0, 0, &
2, 0, 4, -2, 1, &
4, 1, 2, -2, 2, &
-1, -1, 0, 6, 0, &
-3, -1, 2, 6, 2, &
-1, 0, 0, 6, 1, &
-3, 0, 2, 6, 1, &
1, -1, 0, 4, 1, &
1, -1, 0, 4, 0, &
-2, 0, 2, 5, 2, &
1, -2, 2, 2, 1, &
3, -1, 0, 2, 0, &
1, -1, 2, 2, 0, &
0, 0, 2, 3, 1, &
-1, 1, 2, 4, 1, &
0, 1, 2, 3, 2, &
-1, 0, 4, 2, 1, &
2, 0, 2, 1, 1, &
5, 0, 0, 0, 0, &
2, 1, 2, 1, 2, &
1, 0, 4, 0, 1, &
3, 1, 2, 0, 1, &
3, 0, 4, -2, 2, &
-2, -1, 2, 6, 2, &
0, 0, 0, 6, 0, &
0, -2, 2, 4, 2, &
-2, 0, 2, 6, 1, &
2, 0, 0, 4, 1, &
2, 0, 0, 4, 0, &
2, -2, 2, 2, 2, &
0, 0, 2, 4, 0, &
1, 0, 2, 3, 2, &
4, 0, 0, 2, 0, &
2, 0, 2, 2, 0, &
0, 0, 4, 2, 2, &
4, -1, 2, 0, 2, &
3, 0, 2, 1, 2, &
2, 1, 2, 2, 1, &
4, 1, 2, 0, 2, &
-1, -1, 2, 6, 2, &
-1, 0, 2, 6, 1, &
1, -1, 2, 4, 1, &
1, 1, 2, 4, 2, &
3, 1, 2, 2, 2, &
5, 0, 2, 0, 1, &
2, -1, 2, 4, 2, &
2, 0, 2, 4, 1 ], [5,nls])
!
! Luni-Solar nutation coefficients, unit 1e-7 arcsec
! longitude (sin, t*sin, cos), obliquity (cos, t*cos, sin)
!
real(wp),dimension(6,nls),parameter :: cls = reshape([&
-172064161.0_wp, -174666.0_wp, 33386.0_wp, 92052331.0_wp, 9086.0_wp, 15377.0_wp, &
-13170906.0_wp, -1675.0_wp, -13696.0_wp, 5730336.0_wp, -3015.0_wp, -4587.0_wp, &
-2276413.0_wp, -234.0_wp, 2796.0_wp, 978459.0_wp, -485.0_wp, 1374.0_wp, &
2074554.0_wp, 207.0_wp, -698.0_wp, -897492.0_wp, 470.0_wp, -291.0_wp, &
1475877.0_wp, -3633.0_wp, 11817.0_wp, 73871.0_wp, -184.0_wp, -1924.0_wp, &
-516821.0_wp, 1226.0_wp, -524.0_wp, 224386.0_wp, -677.0_wp, -174.0_wp, &
711159.0_wp, 73.0_wp, -872.0_wp, -6750.0_wp, 0.0_wp, 358.0_wp, &
-387298.0_wp, -367.0_wp, 380.0_wp, 200728.0_wp, 18.0_wp, 318.0_wp, &
-301461.0_wp, -36.0_wp, 816.0_wp, 129025.0_wp, -63.0_wp, 367.0_wp, &
215829.0_wp, -494.0_wp, 111.0_wp, -95929.0_wp, 299.0_wp, 132.0_wp, &
128227.0_wp, 137.0_wp, 181.0_wp, -68982.0_wp, -9.0_wp, 39.0_wp, &
123457.0_wp, 11.0_wp, 19.0_wp, -53311.0_wp, 32.0_wp, -4.0_wp, &
156994.0_wp, 10.0_wp, -168.0_wp, -1235.0_wp, 0.0_wp, 82.0_wp, &
63110.0_wp, 63.0_wp, 27.0_wp, -33228.0_wp, 0.0_wp, -9.0_wp, &
-57976.0_wp, -63.0_wp, -189.0_wp, 31429.0_wp, 0.0_wp, -75.0_wp, &
-59641.0_wp, -11.0_wp, 149.0_wp, 25543.0_wp, -11.0_wp, 66.0_wp, &
-51613.0_wp, -42.0_wp, 129.0_wp, 26366.0_wp, 0.0_wp, 78.0_wp, &
45893.0_wp, 50.0_wp, 31.0_wp, -24236.0_wp, -10.0_wp, 20.0_wp, &
63384.0_wp, 11.0_wp, -150.0_wp, -1220.0_wp, 0.0_wp, 29.0_wp, &
-38571.0_wp, -1.0_wp, 158.0_wp, 16452.0_wp, -11.0_wp, 68.0_wp, &
32481.0_wp, 0.0_wp, 0.0_wp, -13870.0_wp, 0.0_wp, 0.0_wp, &
-47722.0_wp, 0.0_wp, -18.0_wp, 477.0_wp, 0.0_wp, -25.0_wp, &
-31046.0_wp, -1.0_wp, 131.0_wp, 13238.0_wp, -11.0_wp, 59.0_wp, &
28593.0_wp, 0.0_wp, -1.0_wp, -12338.0_wp, 10.0_wp, -3.0_wp, &
20441.0_wp, 21.0_wp, 10.0_wp, -10758.0_wp, 0.0_wp, -3.0_wp, &
29243.0_wp, 0.0_wp, -74.0_wp, -609.0_wp, 0.0_wp, 13.0_wp, &
25887.0_wp, 0.0_wp, -66.0_wp, -550.0_wp, 0.0_wp, 11.0_wp, &
-14053.0_wp, -25.0_wp, 79.0_wp, 8551.0_wp, -2.0_wp, -45.0_wp, &
15164.0_wp, 10.0_wp, 11.0_wp, -8001.0_wp, 0.0_wp, -1.0_wp, &
-15794.0_wp, 72.0_wp, -16.0_wp, 6850.0_wp, -42.0_wp, -5.0_wp, &
21783.0_wp, 0.0_wp, 13.0_wp, -167.0_wp, 0.0_wp, 13.0_wp, &
-12873.0_wp, -10.0_wp, -37.0_wp, 6953.0_wp, 0.0_wp, -14.0_wp, &
-12654.0_wp, 11.0_wp, 63.0_wp, 6415.0_wp, 0.0_wp, 26.0_wp, &
-10204.0_wp, 0.0_wp, 25.0_wp, 5222.0_wp, 0.0_wp, 15.0_wp, &
16707.0_wp, -85.0_wp, -10.0_wp, 168.0_wp, -1.0_wp, 10.0_wp, &
-7691.0_wp, 0.0_wp, 44.0_wp, 3268.0_wp, 0.0_wp, 19.0_wp, &
-11024.0_wp, 0.0_wp, -14.0_wp, 104.0_wp, 0.0_wp, 2.0_wp, &
7566.0_wp, -21.0_wp, -11.0_wp, -3250.0_wp, 0.0_wp, -5.0_wp, &
-6637.0_wp, -11.0_wp, 25.0_wp, 3353.0_wp, 0.0_wp, 14.0_wp, &
-7141.0_wp, 21.0_wp, 8.0_wp, 3070.0_wp, 0.0_wp, 4.0_wp, &
-6302.0_wp, -11.0_wp, 2.0_wp, 3272.0_wp, 0.0_wp, 4.0_wp, &
5800.0_wp, 10.0_wp, 2.0_wp, -3045.0_wp, 0.0_wp, -1.0_wp, &
6443.0_wp, 0.0_wp, -7.0_wp, -2768.0_wp, 0.0_wp, -4.0_wp, &
-5774.0_wp, -11.0_wp, -15.0_wp, 3041.0_wp, 0.0_wp, -5.0_wp, &
-5350.0_wp, 0.0_wp, 21.0_wp, 2695.0_wp, 0.0_wp, 12.0_wp, &
-4752.0_wp, -11.0_wp, -3.0_wp, 2719.0_wp, 0.0_wp, -3.0_wp, &
-4940.0_wp, -11.0_wp, -21.0_wp, 2720.0_wp, 0.0_wp, -9.0_wp, &
7350.0_wp, 0.0_wp, -8.0_wp, -51.0_wp, 0.0_wp, 4.0_wp, &
4065.0_wp, 0.0_wp, 6.0_wp, -2206.0_wp, 0.0_wp, 1.0_wp, &
6579.0_wp, 0.0_wp, -24.0_wp, -199.0_wp, 0.0_wp, 2.0_wp, &
3579.0_wp, 0.0_wp, 5.0_wp, -1900.0_wp, 0.0_wp, 1.0_wp, &
4725.0_wp, 0.0_wp, -6.0_wp, -41.0_wp, 0.0_wp, 3.0_wp, &
-3075.0_wp, 0.0_wp, -2.0_wp, 1313.0_wp, 0.0_wp, -1.0_wp, &
-2904.0_wp, 0.0_wp, 15.0_wp, 1233.0_wp, 0.0_wp, 7.0_wp, &
4348.0_wp, 0.0_wp, -10.0_wp, -81.0_wp, 0.0_wp, 2.0_wp, &
-2878.0_wp, 0.0_wp, 8.0_wp, 1232.0_wp, 0.0_wp, 4.0_wp, &
-4230.0_wp, 0.0_wp, 5.0_wp, -20.0_wp, 0.0_wp, -2.0_wp, &
-2819.0_wp, 0.0_wp, 7.0_wp, 1207.0_wp, 0.0_wp, 3.0_wp, &
-4056.0_wp, 0.0_wp, 5.0_wp, 40.0_wp, 0.0_wp, -2.0_wp, &
-2647.0_wp, 0.0_wp, 11.0_wp, 1129.0_wp, 0.0_wp, 5.0_wp, &
-2294.0_wp, 0.0_wp, -10.0_wp, 1266.0_wp, 0.0_wp, -4.0_wp, &
2481.0_wp, 0.0_wp, -7.0_wp, -1062.0_wp, 0.0_wp, -3.0_wp, &
2179.0_wp, 0.0_wp, -2.0_wp, -1129.0_wp, 0.0_wp, -2.0_wp, &
3276.0_wp, 0.0_wp, 1.0_wp, -9.0_wp, 0.0_wp, 0.0_wp, &
-3389.0_wp, 0.0_wp, 5.0_wp, 35.0_wp, 0.0_wp, -2.0_wp, &
3339.0_wp, 0.0_wp, -13.0_wp, -107.0_wp, 0.0_wp, 1.0_wp, &
-1987.0_wp, 0.0_wp, -6.0_wp, 1073.0_wp, 0.0_wp, -2.0_wp, &
-1981.0_wp, 0.0_wp, 0.0_wp, 854.0_wp, 0.0_wp, 0.0_wp, &
4026.0_wp, 0.0_wp, -353.0_wp, -553.0_wp, 0.0_wp, -139.0_wp, &
1660.0_wp, 0.0_wp, -5.0_wp, -710.0_wp, 0.0_wp, -2.0_wp, &
-1521.0_wp, 0.0_wp, 9.0_wp, 647.0_wp, 0.0_wp, 4.0_wp, &
1314.0_wp, 0.0_wp, 0.0_wp, -700.0_wp, 0.0_wp, 0.0_wp, &
-1283.0_wp, 0.0_wp, 0.0_wp, 672.0_wp, 0.0_wp, 0.0_wp, &
-1331.0_wp, 0.0_wp, 8.0_wp, 663.0_wp, 0.0_wp, 4.0_wp, &
1383.0_wp, 0.0_wp, -2.0_wp, -594.0_wp, 0.0_wp, -2.0_wp, &
1405.0_wp, 0.0_wp, 4.0_wp, -610.0_wp, 0.0_wp, 2.0_wp, &
1290.0_wp, 0.0_wp, 0.0_wp, -556.0_wp, 0.0_wp, 0.0_wp, &
-1214.0_wp, 0.0_wp, 5.0_wp, 518.0_wp, 0.0_wp, 2.0_wp, &
1146.0_wp, 0.0_wp, -3.0_wp, -490.0_wp, 0.0_wp, -1.0_wp, &
1019.0_wp, 0.0_wp, -1.0_wp, -527.0_wp, 0.0_wp, -1.0_wp, &
-1100.0_wp, 0.0_wp, 9.0_wp, 465.0_wp, 0.0_wp, 4.0_wp, &
-970.0_wp, 0.0_wp, 2.0_wp, 496.0_wp, 0.0_wp, 1.0_wp, &
1575.0_wp, 0.0_wp, -6.0_wp, -50.0_wp, 0.0_wp, 0.0_wp, &
934.0_wp, 0.0_wp, -3.0_wp, -399.0_wp, 0.0_wp, -1.0_wp, &
922.0_wp, 0.0_wp, -1.0_wp, -395.0_wp, 0.0_wp, -1.0_wp, &
815.0_wp, 0.0_wp, -1.0_wp, -422.0_wp, 0.0_wp, -1.0_wp, &
834.0_wp, 0.0_wp, 2.0_wp, -440.0_wp, 0.0_wp, 1.0_wp, &
1248.0_wp, 0.0_wp, 0.0_wp, -170.0_wp, 0.0_wp, 1.0_wp, &
1338.0_wp, 0.0_wp, -5.0_wp, -39.0_wp, 0.0_wp, 0.0_wp, &
716.0_wp, 0.0_wp, -2.0_wp, -389.0_wp, 0.0_wp, -1.0_wp, &
1282.0_wp, 0.0_wp, -3.0_wp, -23.0_wp, 0.0_wp, 1.0_wp, &
742.0_wp, 0.0_wp, 1.0_wp, -391.0_wp, 0.0_wp, 0.0_wp, &
1020.0_wp, 0.0_wp, -25.0_wp, -495.0_wp, 0.0_wp, -10.0_wp, &
715.0_wp, 0.0_wp, -4.0_wp, -326.0_wp, 0.0_wp, 2.0_wp, &
-666.0_wp, 0.0_wp, -3.0_wp, 369.0_wp, 0.0_wp, -1.0_wp, &
-667.0_wp, 0.0_wp, 1.0_wp, 346.0_wp, 0.0_wp, 1.0_wp, &
-704.0_wp, 0.0_wp, 0.0_wp, 304.0_wp, 0.0_wp, 0.0_wp, &
-694.0_wp, 0.0_wp, 5.0_wp, 294.0_wp, 0.0_wp, 2.0_wp, &
-1014.0_wp, 0.0_wp, -1.0_wp, 4.0_wp, 0.0_wp, -1.0_wp, &
-585.0_wp, 0.0_wp, -2.0_wp, 316.0_wp, 0.0_wp, -1.0_wp, &
-949.0_wp, 0.0_wp, 1.0_wp, 8.0_wp, 0.0_wp, -1.0_wp, &
-595.0_wp, 0.0_wp, 0.0_wp, 258.0_wp, 0.0_wp, 0.0_wp, &
528.0_wp, 0.0_wp, 0.0_wp, -279.0_wp, 0.0_wp, 0.0_wp, &
-590.0_wp, 0.0_wp, 4.0_wp, 252.0_wp, 0.0_wp, 2.0_wp, &
570.0_wp, 0.0_wp, -2.0_wp, -244.0_wp, 0.0_wp, -1.0_wp, &
-502.0_wp, 0.0_wp, 3.0_wp, 250.0_wp, 0.0_wp, 2.0_wp, &
-875.0_wp, 0.0_wp, 1.0_wp, 29.0_wp, 0.0_wp, 0.0_wp, &
-492.0_wp, 0.0_wp, -3.0_wp, 275.0_wp, 0.0_wp, -1.0_wp, &
535.0_wp, 0.0_wp, -2.0_wp, -228.0_wp, 0.0_wp, -1.0_wp, &
-467.0_wp, 0.0_wp, 1.0_wp, 240.0_wp, 0.0_wp, 1.0_wp, &
591.0_wp, 0.0_wp, 0.0_wp, -253.0_wp, 0.0_wp, 0.0_wp, &
-453.0_wp, 0.0_wp, -1.0_wp, 244.0_wp, 0.0_wp, -1.0_wp, &
766.0_wp, 0.0_wp, 1.0_wp, 9.0_wp, 0.0_wp, 0.0_wp, &
-446.0_wp, 0.0_wp, 2.0_wp, 225.0_wp, 0.0_wp, 1.0_wp, &
-488.0_wp, 0.0_wp, 2.0_wp, 207.0_wp, 0.0_wp, 1.0_wp, &
-468.0_wp, 0.0_wp, 0.0_wp, 201.0_wp, 0.0_wp, 0.0_wp, &
-421.0_wp, 0.0_wp, 1.0_wp, 216.0_wp, 0.0_wp, 1.0_wp, &
463.0_wp, 0.0_wp, 0.0_wp, -200.0_wp, 0.0_wp, 0.0_wp, &
-673.0_wp, 0.0_wp, 2.0_wp, 14.0_wp, 0.0_wp, 0.0_wp, &
658.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
-438.0_wp, 0.0_wp, 0.0_wp, 188.0_wp, 0.0_wp, 0.0_wp, &
-390.0_wp, 0.0_wp, 0.0_wp, 205.0_wp, 0.0_wp, 0.0_wp, &
639.0_wp, -11.0_wp, -2.0_wp, -19.0_wp, 0.0_wp, 0.0_wp, &
412.0_wp, 0.0_wp, -2.0_wp, -176.0_wp, 0.0_wp, -1.0_wp, &
-361.0_wp, 0.0_wp, 0.0_wp, 189.0_wp, 0.0_wp, 0.0_wp, &
360.0_wp, 0.0_wp, -1.0_wp, -185.0_wp, 0.0_wp, -1.0_wp, &
588.0_wp, 0.0_wp, -3.0_wp, -24.0_wp, 0.0_wp, 0.0_wp, &
-578.0_wp, 0.0_wp, 1.0_wp, 5.0_wp, 0.0_wp, 0.0_wp, &
-396.0_wp, 0.0_wp, 0.0_wp, 171.0_wp, 0.0_wp, 0.0_wp, &
565.0_wp, 0.0_wp, -1.0_wp, -6.0_wp, 0.0_wp, 0.0_wp, &
-335.0_wp, 0.0_wp, -1.0_wp, 184.0_wp, 0.0_wp, -1.0_wp, &
357.0_wp, 0.0_wp, 1.0_wp, -154.0_wp, 0.0_wp, 0.0_wp, &
321.0_wp, 0.0_wp, 1.0_wp, -174.0_wp, 0.0_wp, 0.0_wp, &
-301.0_wp, 0.0_wp, -1.0_wp, 162.0_wp, 0.0_wp, 0.0_wp, &
-334.0_wp, 0.0_wp, 0.0_wp, 144.0_wp, 0.0_wp, 0.0_wp, &
493.0_wp, 0.0_wp, -2.0_wp, -15.0_wp, 0.0_wp, 0.0_wp, &
494.0_wp, 0.0_wp, -2.0_wp, -19.0_wp, 0.0_wp, 0.0_wp, &
337.0_wp, 0.0_wp, -1.0_wp, -143.0_wp, 0.0_wp, -1.0_wp, &
280.0_wp, 0.0_wp, -1.0_wp, -144.0_wp, 0.0_wp, 0.0_wp, &
309.0_wp, 0.0_wp, 1.0_wp, -134.0_wp, 0.0_wp, 0.0_wp, &
-263.0_wp, 0.0_wp, 2.0_wp, 131.0_wp, 0.0_wp, 1.0_wp, &
253.0_wp, 0.0_wp, 1.0_wp, -138.0_wp, 0.0_wp, 0.0_wp, &
245.0_wp, 0.0_wp, 0.0_wp, -128.0_wp, 0.0_wp, 0.0_wp, &
416.0_wp, 0.0_wp, -2.0_wp, -17.0_wp, 0.0_wp, 0.0_wp, &
-229.0_wp, 0.0_wp, 0.0_wp, 128.0_wp, 0.0_wp, 0.0_wp, &
231.0_wp, 0.0_wp, 0.0_wp, -120.0_wp, 0.0_wp, 0.0_wp, &
-259.0_wp, 0.0_wp, 2.0_wp, 109.0_wp, 0.0_wp, 1.0_wp, &
375.0_wp, 0.0_wp, -1.0_wp, -8.0_wp, 0.0_wp, 0.0_wp, &
252.0_wp, 0.0_wp, 0.0_wp, -108.0_wp, 0.0_wp, 0.0_wp, &
-245.0_wp, 0.0_wp, 1.0_wp, 104.0_wp, 0.0_wp, 0.0_wp, &
243.0_wp, 0.0_wp, -1.0_wp, -104.0_wp, 0.0_wp, 0.0_wp, &
208.0_wp, 0.0_wp, 1.0_wp, -112.0_wp, 0.0_wp, 0.0_wp, &
199.0_wp, 0.0_wp, 0.0_wp, -102.0_wp, 0.0_wp, 0.0_wp, &
-208.0_wp, 0.0_wp, 1.0_wp, 105.0_wp, 0.0_wp, 0.0_wp, &
335.0_wp, 0.0_wp, -2.0_wp, -14.0_wp, 0.0_wp, 0.0_wp, &
-325.0_wp, 0.0_wp, 1.0_wp, 7.0_wp, 0.0_wp, 0.0_wp, &
-187.0_wp, 0.0_wp, 0.0_wp, 96.0_wp, 0.0_wp, 0.0_wp, &
197.0_wp, 0.0_wp, -1.0_wp, -100.0_wp, 0.0_wp, 0.0_wp, &
-192.0_wp, 0.0_wp, 2.0_wp, 94.0_wp, 0.0_wp, 1.0_wp, &
-188.0_wp, 0.0_wp, 0.0_wp, 83.0_wp, 0.0_wp, 0.0_wp, &
276.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
-286.0_wp, 0.0_wp, 1.0_wp, 6.0_wp, 0.0_wp, 0.0_wp, &
186.0_wp, 0.0_wp, -1.0_wp, -79.0_wp, 0.0_wp, 0.0_wp, &
-219.0_wp, 0.0_wp, 0.0_wp, 43.0_wp, 0.0_wp, 0.0_wp, &
276.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
-153.0_wp, 0.0_wp, -1.0_wp, 84.0_wp, 0.0_wp, 0.0_wp, &
-156.0_wp, 0.0_wp, 0.0_wp, 81.0_wp, 0.0_wp, 0.0_wp, &
-154.0_wp, 0.0_wp, 1.0_wp, 78.0_wp, 0.0_wp, 0.0_wp, &
-174.0_wp, 0.0_wp, 1.0_wp, 75.0_wp, 0.0_wp, 0.0_wp, &
-163.0_wp, 0.0_wp, 2.0_wp, 69.0_wp, 0.0_wp, 1.0_wp, &
-228.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 0.0_wp, 0.0_wp, &
91.0_wp, 0.0_wp, -4.0_wp, -54.0_wp, 0.0_wp, -2.0_wp, &
175.0_wp, 0.0_wp, 0.0_wp, -75.0_wp, 0.0_wp, 0.0_wp, &
-159.0_wp, 0.0_wp, 0.0_wp, 69.0_wp, 0.0_wp, 0.0_wp, &
141.0_wp, 0.0_wp, 0.0_wp, -72.0_wp, 0.0_wp, 0.0_wp, &
147.0_wp, 0.0_wp, 0.0_wp, -75.0_wp, 0.0_wp, 0.0_wp, &
-132.0_wp, 0.0_wp, 0.0_wp, 69.0_wp, 0.0_wp, 0.0_wp, &
159.0_wp, 0.0_wp, -28.0_wp, -54.0_wp, 0.0_wp, 11.0_wp, &
213.0_wp, 0.0_wp, 0.0_wp, -4.0_wp, 0.0_wp, 0.0_wp, &
123.0_wp, 0.0_wp, 0.0_wp, -64.0_wp, 0.0_wp, 0.0_wp, &
-118.0_wp, 0.0_wp, -1.0_wp, 66.0_wp, 0.0_wp, 0.0_wp, &
144.0_wp, 0.0_wp, -1.0_wp, -61.0_wp, 0.0_wp, 0.0_wp, &
-121.0_wp, 0.0_wp, 1.0_wp, 60.0_wp, 0.0_wp, 0.0_wp, &
-134.0_wp, 0.0_wp, 1.0_wp, 56.0_wp, 0.0_wp, 1.0_wp, &
-105.0_wp, 0.0_wp, 0.0_wp, 57.0_wp, 0.0_wp, 0.0_wp, &
-102.0_wp, 0.0_wp, 0.0_wp, 56.0_wp, 0.0_wp, 0.0_wp, &
120.0_wp, 0.0_wp, 0.0_wp, -52.0_wp, 0.0_wp, 0.0_wp, &
101.0_wp, 0.0_wp, 0.0_wp, -54.0_wp, 0.0_wp, 0.0_wp, &
-113.0_wp, 0.0_wp, 0.0_wp, 59.0_wp, 0.0_wp, 0.0_wp, &
-106.0_wp, 0.0_wp, 0.0_wp, 61.0_wp, 0.0_wp, 0.0_wp, &
-129.0_wp, 0.0_wp, 1.0_wp, 55.0_wp, 0.0_wp, 0.0_wp, &
-114.0_wp, 0.0_wp, 0.0_wp, 57.0_wp, 0.0_wp, 0.0_wp, &
113.0_wp, 0.0_wp, -1.0_wp, -49.0_wp, 0.0_wp, 0.0_wp, &
-102.0_wp, 0.0_wp, 0.0_wp, 44.0_wp, 0.0_wp, 0.0_wp, &
-94.0_wp, 0.0_wp, 0.0_wp, 51.0_wp, 0.0_wp, 0.0_wp, &
-100.0_wp, 0.0_wp, -1.0_wp, 56.0_wp, 0.0_wp, 0.0_wp, &
87.0_wp, 0.0_wp, 0.0_wp, -47.0_wp, 0.0_wp, 0.0_wp, &
161.0_wp, 0.0_wp, 0.0_wp, -1.0_wp, 0.0_wp, 0.0_wp, &
96.0_wp, 0.0_wp, 0.0_wp, -50.0_wp, 0.0_wp, 0.0_wp, &
151.0_wp, 0.0_wp, -1.0_wp, -5.0_wp, 0.0_wp, 0.0_wp, &
-104.0_wp, 0.0_wp, 0.0_wp, 44.0_wp, 0.0_wp, 0.0_wp, &
-110.0_wp, 0.0_wp, 0.0_wp, 48.0_wp, 0.0_wp, 0.0_wp, &
-100.0_wp, 0.0_wp, 1.0_wp, 50.0_wp, 0.0_wp, 0.0_wp, &
92.0_wp, 0.0_wp, -5.0_wp, 12.0_wp, 0.0_wp, -2.0_wp, &
82.0_wp, 0.0_wp, 0.0_wp, -45.0_wp, 0.0_wp, 0.0_wp, &
82.0_wp, 0.0_wp, 0.0_wp, -45.0_wp, 0.0_wp, 0.0_wp, &
-78.0_wp, 0.0_wp, 0.0_wp, 41.0_wp, 0.0_wp, 0.0_wp, &
-77.0_wp, 0.0_wp, 0.0_wp, 43.0_wp, 0.0_wp, 0.0_wp, &
2.0_wp, 0.0_wp, 0.0_wp, 54.0_wp, 0.0_wp, 0.0_wp, &
94.0_wp, 0.0_wp, 0.0_wp, -40.0_wp, 0.0_wp, 0.0_wp, &
-93.0_wp, 0.0_wp, 0.0_wp, 40.0_wp, 0.0_wp, 0.0_wp, &
-83.0_wp, 0.0_wp, 10.0_wp, 40.0_wp, 0.0_wp, -2.0_wp, &
83.0_wp, 0.0_wp, 0.0_wp, -36.0_wp, 0.0_wp, 0.0_wp, &
-91.0_wp, 0.0_wp, 0.0_wp, 39.0_wp, 0.0_wp, 0.0_wp, &
128.0_wp, 0.0_wp, 0.0_wp, -1.0_wp, 0.0_wp, 0.0_wp, &
-79.0_wp, 0.0_wp, 0.0_wp, 34.0_wp, 0.0_wp, 0.0_wp, &
-83.0_wp, 0.0_wp, 0.0_wp, 47.0_wp, 0.0_wp, 0.0_wp, &
84.0_wp, 0.0_wp, 0.0_wp, -44.0_wp, 0.0_wp, 0.0_wp, &
83.0_wp, 0.0_wp, 0.0_wp, -43.0_wp, 0.0_wp, 0.0_wp, &
91.0_wp, 0.0_wp, 0.0_wp, -39.0_wp, 0.0_wp, 0.0_wp, &
-77.0_wp, 0.0_wp, 0.0_wp, 39.0_wp, 0.0_wp, 0.0_wp, &
84.0_wp, 0.0_wp, 0.0_wp, -43.0_wp, 0.0_wp, 0.0_wp, &
-92.0_wp, 0.0_wp, 1.0_wp, 39.0_wp, 0.0_wp, 0.0_wp, &
-92.0_wp, 0.0_wp, 1.0_wp, 39.0_wp, 0.0_wp, 0.0_wp, &
-94.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
68.0_wp, 0.0_wp, 0.0_wp, -36.0_wp, 0.0_wp, 0.0_wp, &
-61.0_wp, 0.0_wp, 0.0_wp, 32.0_wp, 0.0_wp, 0.0_wp, &
71.0_wp, 0.0_wp, 0.0_wp, -31.0_wp, 0.0_wp, 0.0_wp, &
62.0_wp, 0.0_wp, 0.0_wp, -34.0_wp, 0.0_wp, 0.0_wp, &
-63.0_wp, 0.0_wp, 0.0_wp, 33.0_wp, 0.0_wp, 0.0_wp, &
-73.0_wp, 0.0_wp, 0.0_wp, 32.0_wp, 0.0_wp, 0.0_wp, &
115.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
-103.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
63.0_wp, 0.0_wp, 0.0_wp, -28.0_wp, 0.0_wp, 0.0_wp, &
74.0_wp, 0.0_wp, 0.0_wp, -32.0_wp, 0.0_wp, 0.0_wp, &
-103.0_wp, 0.0_wp, -3.0_wp, 3.0_wp, 0.0_wp, -1.0_wp, &
-69.0_wp, 0.0_wp, 0.0_wp, 30.0_wp, 0.0_wp, 0.0_wp, &
57.0_wp, 0.0_wp, 0.0_wp, -29.0_wp, 0.0_wp, 0.0_wp, &
94.0_wp, 0.0_wp, 0.0_wp, -4.0_wp, 0.0_wp, 0.0_wp, &
64.0_wp, 0.0_wp, 0.0_wp, -33.0_wp, 0.0_wp, 0.0_wp, &
-63.0_wp, 0.0_wp, 0.0_wp, 26.0_wp, 0.0_wp, 0.0_wp, &
-38.0_wp, 0.0_wp, 0.0_wp, 20.0_wp, 0.0_wp, 0.0_wp, &
-43.0_wp, 0.0_wp, 0.0_wp, 24.0_wp, 0.0_wp, 0.0_wp, &
-45.0_wp, 0.0_wp, 0.0_wp, 23.0_wp, 0.0_wp, 0.0_wp, &
47.0_wp, 0.0_wp, 0.0_wp, -24.0_wp, 0.0_wp, 0.0_wp, &
-48.0_wp, 0.0_wp, 0.0_wp, 25.0_wp, 0.0_wp, 0.0_wp, &
45.0_wp, 0.0_wp, 0.0_wp, -26.0_wp, 0.0_wp, 0.0_wp, &
56.0_wp, 0.0_wp, 0.0_wp, -25.0_wp, 0.0_wp, 0.0_wp, &
88.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
-75.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
85.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
49.0_wp, 0.0_wp, 0.0_wp, -26.0_wp, 0.0_wp, 0.0_wp, &
-74.0_wp, 0.0_wp, -3.0_wp, -1.0_wp, 0.0_wp, -1.0_wp, &
-39.0_wp, 0.0_wp, 0.0_wp, 21.0_wp, 0.0_wp, 0.0_wp, &
45.0_wp, 0.0_wp, 0.0_wp, -20.0_wp, 0.0_wp, 0.0_wp, &
51.0_wp, 0.0_wp, 0.0_wp, -22.0_wp, 0.0_wp, 0.0_wp, &
-40.0_wp, 0.0_wp, 0.0_wp, 21.0_wp, 0.0_wp, 0.0_wp, &
41.0_wp, 0.0_wp, 0.0_wp, -21.0_wp, 0.0_wp, 0.0_wp, &
-42.0_wp, 0.0_wp, 0.0_wp, 24.0_wp, 0.0_wp, 0.0_wp, &
-51.0_wp, 0.0_wp, 0.0_wp, 22.0_wp, 0.0_wp, 0.0_wp, &
-42.0_wp, 0.0_wp, 0.0_wp, 22.0_wp, 0.0_wp, 0.0_wp, &
39.0_wp, 0.0_wp, 0.0_wp, -21.0_wp, 0.0_wp, 0.0_wp, &
46.0_wp, 0.0_wp, 0.0_wp, -18.0_wp, 0.0_wp, 0.0_wp, &
-53.0_wp, 0.0_wp, 0.0_wp, 22.0_wp, 0.0_wp, 0.0_wp, &
82.0_wp, 0.0_wp, 0.0_wp, -4.0_wp, 0.0_wp, 0.0_wp, &
81.0_wp, 0.0_wp, -1.0_wp, -4.0_wp, 0.0_wp, 0.0_wp, &
47.0_wp, 0.0_wp, 0.0_wp, -19.0_wp, 0.0_wp, 0.0_wp, &
53.0_wp, 0.0_wp, 0.0_wp, -23.0_wp, 0.0_wp, 0.0_wp, &
-45.0_wp, 0.0_wp, 0.0_wp, 22.0_wp, 0.0_wp, 0.0_wp, &
-44.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
-33.0_wp, 0.0_wp, 0.0_wp, 16.0_wp, 0.0_wp, 0.0_wp, &
-61.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 0.0_wp, 0.0_wp, &
28.0_wp, 0.0_wp, 0.0_wp, -15.0_wp, 0.0_wp, 0.0_wp, &
-38.0_wp, 0.0_wp, 0.0_wp, 19.0_wp, 0.0_wp, 0.0_wp, &
-33.0_wp, 0.0_wp, 0.0_wp, 21.0_wp, 0.0_wp, 0.0_wp, &
-60.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
48.0_wp, 0.0_wp, 0.0_wp, -10.0_wp, 0.0_wp, 0.0_wp, &
27.0_wp, 0.0_wp, 0.0_wp, -14.0_wp, 0.0_wp, 0.0_wp, &
38.0_wp, 0.0_wp, 0.0_wp, -20.0_wp, 0.0_wp, 0.0_wp, &
31.0_wp, 0.0_wp, 0.0_wp, -13.0_wp, 0.0_wp, 0.0_wp, &
-29.0_wp, 0.0_wp, 0.0_wp, 15.0_wp, 0.0_wp, 0.0_wp, &
28.0_wp, 0.0_wp, 0.0_wp, -15.0_wp, 0.0_wp, 0.0_wp, &
-32.0_wp, 0.0_wp, 0.0_wp, 15.0_wp, 0.0_wp, 0.0_wp, &
45.0_wp, 0.0_wp, 0.0_wp, -8.0_wp, 0.0_wp, 0.0_wp, &
-44.0_wp, 0.0_wp, 0.0_wp, 19.0_wp, 0.0_wp, 0.0_wp, &
28.0_wp, 0.0_wp, 0.0_wp, -15.0_wp, 0.0_wp, 0.0_wp, &
-51.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-36.0_wp, 0.0_wp, 0.0_wp, 20.0_wp, 0.0_wp, 0.0_wp, &
44.0_wp, 0.0_wp, 0.0_wp, -19.0_wp, 0.0_wp, 0.0_wp, &
26.0_wp, 0.0_wp, 0.0_wp, -14.0_wp, 0.0_wp, 0.0_wp, &
-60.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
35.0_wp, 0.0_wp, 0.0_wp, -18.0_wp, 0.0_wp, 0.0_wp, &
-27.0_wp, 0.0_wp, 0.0_wp, 11.0_wp, 0.0_wp, 0.0_wp, &
47.0_wp, 0.0_wp, 0.0_wp, -1.0_wp, 0.0_wp, 0.0_wp, &
36.0_wp, 0.0_wp, 0.0_wp, -15.0_wp, 0.0_wp, 0.0_wp, &
-36.0_wp, 0.0_wp, 0.0_wp, 20.0_wp, 0.0_wp, 0.0_wp, &
-35.0_wp, 0.0_wp, 0.0_wp, 19.0_wp, 0.0_wp, 0.0_wp, &
-37.0_wp, 0.0_wp, 0.0_wp, 19.0_wp, 0.0_wp, 0.0_wp, &
32.0_wp, 0.0_wp, 0.0_wp, -16.0_wp, 0.0_wp, 0.0_wp, &
35.0_wp, 0.0_wp, 0.0_wp, -14.0_wp, 0.0_wp, 0.0_wp, &
32.0_wp, 0.0_wp, 0.0_wp, -13.0_wp, 0.0_wp, 0.0_wp, &
65.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
47.0_wp, 0.0_wp, 0.0_wp, -1.0_wp, 0.0_wp, 0.0_wp, &
32.0_wp, 0.0_wp, 0.0_wp, -16.0_wp, 0.0_wp, 0.0_wp, &
37.0_wp, 0.0_wp, 0.0_wp, -16.0_wp, 0.0_wp, 0.0_wp, &
-30.0_wp, 0.0_wp, 0.0_wp, 15.0_wp, 0.0_wp, 0.0_wp, &
-32.0_wp, 0.0_wp, 0.0_wp, 16.0_wp, 0.0_wp, 0.0_wp, &
-31.0_wp, 0.0_wp, 0.0_wp, 13.0_wp, 0.0_wp, 0.0_wp, &
37.0_wp, 0.0_wp, 0.0_wp, -16.0_wp, 0.0_wp, 0.0_wp, &
31.0_wp, 0.0_wp, 0.0_wp, -13.0_wp, 0.0_wp, 0.0_wp, &
49.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
32.0_wp, 0.0_wp, 0.0_wp, -13.0_wp, 0.0_wp, 0.0_wp, &
23.0_wp, 0.0_wp, 0.0_wp, -12.0_wp, 0.0_wp, 0.0_wp, &
-43.0_wp, 0.0_wp, 0.0_wp, 18.0_wp, 0.0_wp, 0.0_wp, &
26.0_wp, 0.0_wp, 0.0_wp, -11.0_wp, 0.0_wp, 0.0_wp, &
-32.0_wp, 0.0_wp, 0.0_wp, 14.0_wp, 0.0_wp, 0.0_wp, &
-29.0_wp, 0.0_wp, 0.0_wp, 14.0_wp, 0.0_wp, 0.0_wp, &
-27.0_wp, 0.0_wp, 0.0_wp, 12.0_wp, 0.0_wp, 0.0_wp, &
30.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-11.0_wp, 0.0_wp, 0.0_wp, 5.0_wp, 0.0_wp, 0.0_wp, &
-21.0_wp, 0.0_wp, 0.0_wp, 10.0_wp, 0.0_wp, 0.0_wp, &
-34.0_wp, 0.0_wp, 0.0_wp, 15.0_wp, 0.0_wp, 0.0_wp, &
-10.0_wp, 0.0_wp, 0.0_wp, 6.0_wp, 0.0_wp, 0.0_wp, &
-36.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-9.0_wp, 0.0_wp, 0.0_wp, 4.0_wp, 0.0_wp, 0.0_wp, &
-12.0_wp, 0.0_wp, 0.0_wp, 5.0_wp, 0.0_wp, 0.0_wp, &
-21.0_wp, 0.0_wp, 0.0_wp, 5.0_wp, 0.0_wp, 0.0_wp, &
-29.0_wp, 0.0_wp, 0.0_wp, -1.0_wp, 0.0_wp, 0.0_wp, &
-15.0_wp, 0.0_wp, 0.0_wp, 3.0_wp, 0.0_wp, 0.0_wp, &
-20.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
28.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, &
17.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-22.0_wp, 0.0_wp, 0.0_wp, 12.0_wp, 0.0_wp, 0.0_wp, &
-14.0_wp, 0.0_wp, 0.0_wp, 7.0_wp, 0.0_wp, 0.0_wp, &
24.0_wp, 0.0_wp, 0.0_wp, -11.0_wp, 0.0_wp, 0.0_wp, &
11.0_wp, 0.0_wp, 0.0_wp, -6.0_wp, 0.0_wp, 0.0_wp, &
14.0_wp, 0.0_wp, 0.0_wp, -6.0_wp, 0.0_wp, 0.0_wp, &
24.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
18.0_wp, 0.0_wp, 0.0_wp, -8.0_wp, 0.0_wp, 0.0_wp, &
-38.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-31.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-16.0_wp, 0.0_wp, 0.0_wp, 8.0_wp, 0.0_wp, 0.0_wp, &
29.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-18.0_wp, 0.0_wp, 0.0_wp, 10.0_wp, 0.0_wp, 0.0_wp, &
-10.0_wp, 0.0_wp, 0.0_wp, 5.0_wp, 0.0_wp, 0.0_wp, &
-17.0_wp, 0.0_wp, 0.0_wp, 10.0_wp, 0.0_wp, 0.0_wp, &
9.0_wp, 0.0_wp, 0.0_wp, -4.0_wp, 0.0_wp, 0.0_wp, &
16.0_wp, 0.0_wp, 0.0_wp, -6.0_wp, 0.0_wp, 0.0_wp, &
22.0_wp, 0.0_wp, 0.0_wp, -12.0_wp, 0.0_wp, 0.0_wp, &
20.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-13.0_wp, 0.0_wp, 0.0_wp, 6.0_wp, 0.0_wp, 0.0_wp, &
-17.0_wp, 0.0_wp, 0.0_wp, 9.0_wp, 0.0_wp, 0.0_wp, &
-14.0_wp, 0.0_wp, 0.0_wp, 8.0_wp, 0.0_wp, 0.0_wp, &
0.0_wp, 0.0_wp, 0.0_wp, -7.0_wp, 0.0_wp, 0.0_wp, &
14.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
19.0_wp, 0.0_wp, 0.0_wp, -10.0_wp, 0.0_wp, 0.0_wp, &
-34.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-20.0_wp, 0.0_wp, 0.0_wp, 8.0_wp, 0.0_wp, 0.0_wp, &
9.0_wp, 0.0_wp, 0.0_wp, -5.0_wp, 0.0_wp, 0.0_wp, &
-18.0_wp, 0.0_wp, 0.0_wp, 7.0_wp, 0.0_wp, 0.0_wp, &
13.0_wp, 0.0_wp, 0.0_wp, -6.0_wp, 0.0_wp, 0.0_wp, &
17.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-12.0_wp, 0.0_wp, 0.0_wp, 5.0_wp, 0.0_wp, 0.0_wp, &
15.0_wp, 0.0_wp, 0.0_wp, -8.0_wp, 0.0_wp, 0.0_wp, &
-11.0_wp, 0.0_wp, 0.0_wp, 3.0_wp, 0.0_wp, 0.0_wp, &
13.0_wp, 0.0_wp, 0.0_wp, -5.0_wp, 0.0_wp, 0.0_wp, &
-18.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-35.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
9.0_wp, 0.0_wp, 0.0_wp, -4.0_wp, 0.0_wp, 0.0_wp, &
-19.0_wp, 0.0_wp, 0.0_wp, 10.0_wp, 0.0_wp, 0.0_wp, &
-26.0_wp, 0.0_wp, 0.0_wp, 11.0_wp, 0.0_wp, 0.0_wp, &
8.0_wp, 0.0_wp, 0.0_wp, -4.0_wp, 0.0_wp, 0.0_wp, &
-10.0_wp, 0.0_wp, 0.0_wp, 4.0_wp, 0.0_wp, 0.0_wp, &
10.0_wp, 0.0_wp, 0.0_wp, -6.0_wp, 0.0_wp, 0.0_wp, &
-21.0_wp, 0.0_wp, 0.0_wp, 9.0_wp, 0.0_wp, 0.0_wp, &
-15.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
9.0_wp, 0.0_wp, 0.0_wp, -5.0_wp, 0.0_wp, 0.0_wp, &
-29.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-19.0_wp, 0.0_wp, 0.0_wp, 10.0_wp, 0.0_wp, 0.0_wp, &
12.0_wp, 0.0_wp, 0.0_wp, -5.0_wp, 0.0_wp, 0.0_wp, &
22.0_wp, 0.0_wp, 0.0_wp, -9.0_wp, 0.0_wp, 0.0_wp, &
-10.0_wp, 0.0_wp, 0.0_wp, 5.0_wp, 0.0_wp, 0.0_wp, &
-20.0_wp, 0.0_wp, 0.0_wp, 11.0_wp, 0.0_wp, 0.0_wp, &
-20.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-17.0_wp, 0.0_wp, 0.0_wp, 7.0_wp, 0.0_wp, 0.0_wp, &
15.0_wp, 0.0_wp, 0.0_wp, -3.0_wp, 0.0_wp, 0.0_wp, &
8.0_wp, 0.0_wp, 0.0_wp, -4.0_wp, 0.0_wp, 0.0_wp, &
14.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-12.0_wp, 0.0_wp, 0.0_wp, 6.0_wp, 0.0_wp, 0.0_wp, &
25.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-13.0_wp, 0.0_wp, 0.0_wp, 6.0_wp, 0.0_wp, 0.0_wp, &
-14.0_wp, 0.0_wp, 0.0_wp, 8.0_wp, 0.0_wp, 0.0_wp, &
13.0_wp, 0.0_wp, 0.0_wp, -5.0_wp, 0.0_wp, 0.0_wp, &
-17.0_wp, 0.0_wp, 0.0_wp, 9.0_wp, 0.0_wp, 0.0_wp, &
-12.0_wp, 0.0_wp, 0.0_wp, 6.0_wp, 0.0_wp, 0.0_wp, &
-10.0_wp, 0.0_wp, 0.0_wp, 5.0_wp, 0.0_wp, 0.0_wp, &
10.0_wp, 0.0_wp, 0.0_wp, -6.0_wp, 0.0_wp, 0.0_wp, &
-15.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-22.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
28.0_wp, 0.0_wp, 0.0_wp, -1.0_wp, 0.0_wp, 0.0_wp, &
15.0_wp, 0.0_wp, 0.0_wp, -7.0_wp, 0.0_wp, 0.0_wp, &
23.0_wp, 0.0_wp, 0.0_wp, -10.0_wp, 0.0_wp, 0.0_wp, &
12.0_wp, 0.0_wp, 0.0_wp, -5.0_wp, 0.0_wp, 0.0_wp, &
29.0_wp, 0.0_wp, 0.0_wp, -1.0_wp, 0.0_wp, 0.0_wp, &
-25.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 0.0_wp, 0.0_wp, &
22.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-18.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
15.0_wp, 0.0_wp, 0.0_wp, 3.0_wp, 0.0_wp, 0.0_wp, &
-23.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
12.0_wp, 0.0_wp, 0.0_wp, -5.0_wp, 0.0_wp, 0.0_wp, &
-8.0_wp, 0.0_wp, 0.0_wp, 4.0_wp, 0.0_wp, 0.0_wp, &
-19.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-10.0_wp, 0.0_wp, 0.0_wp, 4.0_wp, 0.0_wp, 0.0_wp, &
21.0_wp, 0.0_wp, 0.0_wp, -9.0_wp, 0.0_wp, 0.0_wp, &
23.0_wp, 0.0_wp, 0.0_wp, -1.0_wp, 0.0_wp, 0.0_wp, &
-16.0_wp, 0.0_wp, 0.0_wp, 8.0_wp, 0.0_wp, 0.0_wp, &
-19.0_wp, 0.0_wp, 0.0_wp, 9.0_wp, 0.0_wp, 0.0_wp, &
-22.0_wp, 0.0_wp, 0.0_wp, 10.0_wp, 0.0_wp, 0.0_wp, &
27.0_wp, 0.0_wp, 0.0_wp, -1.0_wp, 0.0_wp, 0.0_wp, &
16.0_wp, 0.0_wp, 0.0_wp, -8.0_wp, 0.0_wp, 0.0_wp, &
19.0_wp, 0.0_wp, 0.0_wp, -8.0_wp, 0.0_wp, 0.0_wp, &
9.0_wp, 0.0_wp, 0.0_wp, -4.0_wp, 0.0_wp, 0.0_wp, &
-9.0_wp, 0.0_wp, 0.0_wp, 4.0_wp, 0.0_wp, 0.0_wp, &
-9.0_wp, 0.0_wp, 0.0_wp, 4.0_wp, 0.0_wp, 0.0_wp, &
-8.0_wp, 0.0_wp, 0.0_wp, 4.0_wp, 0.0_wp, 0.0_wp, &
18.0_wp, 0.0_wp, 0.0_wp, -9.0_wp, 0.0_wp, 0.0_wp, &
16.0_wp, 0.0_wp, 0.0_wp, -1.0_wp, 0.0_wp, 0.0_wp, &
-10.0_wp, 0.0_wp, 0.0_wp, 4.0_wp, 0.0_wp, 0.0_wp, &
-23.0_wp, 0.0_wp, 0.0_wp, 9.0_wp, 0.0_wp, 0.0_wp, &
16.0_wp, 0.0_wp, 0.0_wp, -1.0_wp, 0.0_wp, 0.0_wp, &
-12.0_wp, 0.0_wp, 0.0_wp, 6.0_wp, 0.0_wp, 0.0_wp, &
-8.0_wp, 0.0_wp, 0.0_wp, 4.0_wp, 0.0_wp, 0.0_wp, &
30.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
24.0_wp, 0.0_wp, 0.0_wp, -10.0_wp, 0.0_wp, 0.0_wp, &
10.0_wp, 0.0_wp, 0.0_wp, -4.0_wp, 0.0_wp, 0.0_wp, &
-16.0_wp, 0.0_wp, 0.0_wp, 7.0_wp, 0.0_wp, 0.0_wp, &
-16.0_wp, 0.0_wp, 0.0_wp, 7.0_wp, 0.0_wp, 0.0_wp, &
17.0_wp, 0.0_wp, 0.0_wp, -7.0_wp, 0.0_wp, 0.0_wp, &
-24.0_wp, 0.0_wp, 0.0_wp, 10.0_wp, 0.0_wp, 0.0_wp, &
-12.0_wp, 0.0_wp, 0.0_wp, 5.0_wp, 0.0_wp, 0.0_wp, &
-24.0_wp, 0.0_wp, 0.0_wp, 11.0_wp, 0.0_wp, 0.0_wp, &
-23.0_wp, 0.0_wp, 0.0_wp, 9.0_wp, 0.0_wp, 0.0_wp, &
-13.0_wp, 0.0_wp, 0.0_wp, 5.0_wp, 0.0_wp, 0.0_wp, &
-15.0_wp, 0.0_wp, 0.0_wp, 7.0_wp, 0.0_wp, 0.0_wp, &
0.0_wp, 0.0_wp, -1988.0_wp, 0.0_wp, 0.0_wp, -1679.0_wp, &
0.0_wp, 0.0_wp, -63.0_wp, 0.0_wp, 0.0_wp, -27.0_wp, &
-4.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
0.0_wp, 0.0_wp, 5.0_wp, 0.0_wp, 0.0_wp, 4.0_wp, &
5.0_wp, 0.0_wp, 0.0_wp, -3.0_wp, 0.0_wp, 0.0_wp, &
0.0_wp, 0.0_wp, 364.0_wp, 0.0_wp, 0.0_wp, 176.0_wp, &
0.0_wp, 0.0_wp, -1044.0_wp, 0.0_wp, 0.0_wp, -891.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 0.0_wp, 0.0_wp, &
4.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
0.0_wp, 0.0_wp, 330.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
5.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 0.0_wp, 0.0_wp, &
-5.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, -1.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
0.0_wp, 0.0_wp, 5.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
0.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 0.0_wp, 0.0_wp, &
4.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
6.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
5.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
-7.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-12.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
5.0_wp, 0.0_wp, 0.0_wp, -3.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, -1.0_wp, 0.0_wp, 0.0_wp, &
-5.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-7.0_wp, 0.0_wp, 0.0_wp, 3.0_wp, 0.0_wp, 0.0_wp, &
7.0_wp, 0.0_wp, 0.0_wp, -4.0_wp, 0.0_wp, 0.0_wp, &
0.0_wp, 0.0_wp, -12.0_wp, 0.0_wp, 0.0_wp, -10.0_wp, &
4.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
-7.0_wp, 0.0_wp, 0.0_wp, 3.0_wp, 0.0_wp, 0.0_wp, &
-4.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 0.0_wp, 0.0_wp, &
0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 0.0_wp, 0.0_wp, &
7.0_wp, 0.0_wp, 0.0_wp, -3.0_wp, 0.0_wp, 0.0_wp, &
-4.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
4.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
-5.0_wp, 0.0_wp, 0.0_wp, 3.0_wp, 0.0_wp, 0.0_wp, &
5.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-5.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
5.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
-8.0_wp, 0.0_wp, 0.0_wp, 3.0_wp, 0.0_wp, 0.0_wp, &
9.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
6.0_wp, 0.0_wp, 0.0_wp, -3.0_wp, 0.0_wp, 0.0_wp, &
-5.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-7.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 0.0_wp, 0.0_wp, &
5.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
4.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, -1.0_wp, 0.0_wp, 0.0_wp, &
-5.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
4.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
9.0_wp, 0.0_wp, 0.0_wp, -3.0_wp, 0.0_wp, 0.0_wp, &
4.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
4.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
-4.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
9.0_wp, 0.0_wp, 0.0_wp, -3.0_wp, 0.0_wp, 0.0_wp, &
-4.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-4.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
8.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, -1.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, -1.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 0.0_wp, 0.0_wp, &
6.0_wp, 0.0_wp, 0.0_wp, -3.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 0.0_wp, 0.0_wp, &
-7.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
9.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-4.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-5.0_wp, 0.0_wp, 0.0_wp, 3.0_wp, 0.0_wp, 0.0_wp, &
-13.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-7.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
10.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, -1.0_wp, 0.0_wp, 0.0_wp, &
10.0_wp, 0.0_wp, 13.0_wp, 6.0_wp, 0.0_wp, -5.0_wp, &
0.0_wp, 0.0_wp, 30.0_wp, 0.0_wp, 0.0_wp, 14.0_wp, &
0.0_wp, 0.0_wp, -162.0_wp, 0.0_wp, 0.0_wp, -138.0_wp, &
0.0_wp, 0.0_wp, 75.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-7.0_wp, 0.0_wp, 0.0_wp, 4.0_wp, 0.0_wp, 0.0_wp, &
-4.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
4.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
5.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
5.0_wp, 0.0_wp, 0.0_wp, -3.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
-4.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
-5.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
6.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
9.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
5.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-7.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 0.0_wp, 0.0_wp, &
-4.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
7.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-4.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
4.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-6.0_wp, 0.0_wp, -3.0_wp, 3.0_wp, 0.0_wp, 1.0_wp, &
0.0_wp, 0.0_wp, -3.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, &
11.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, -1.0_wp, 0.0_wp, 0.0_wp, &
11.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
-1.0_wp, 0.0_wp, 3.0_wp, 3.0_wp, 0.0_wp, -1.0_wp, &
4.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
0.0_wp, 0.0_wp, -13.0_wp, 0.0_wp, 0.0_wp, -11.0_wp, &
3.0_wp, 0.0_wp, 6.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-7.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
5.0_wp, 0.0_wp, 0.0_wp, -3.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
5.0_wp, 0.0_wp, 0.0_wp, -3.0_wp, 0.0_wp, 0.0_wp, &
-7.0_wp, 0.0_wp, 0.0_wp, 3.0_wp, 0.0_wp, 0.0_wp, &
8.0_wp, 0.0_wp, 0.0_wp, -3.0_wp, 0.0_wp, 0.0_wp, &
-4.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
11.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, -1.0_wp, 0.0_wp, 0.0_wp, &
-4.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
8.0_wp, 0.0_wp, 0.0_wp, -4.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, -1.0_wp, 0.0_wp, 0.0_wp, &
11.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-6.0_wp, 0.0_wp, 0.0_wp, 3.0_wp, 0.0_wp, 0.0_wp, &
-4.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
-8.0_wp, 0.0_wp, 0.0_wp, 4.0_wp, 0.0_wp, 0.0_wp, &
-7.0_wp, 0.0_wp, 0.0_wp, 3.0_wp, 0.0_wp, 0.0_wp, &
-4.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, -1.0_wp, 0.0_wp, 0.0_wp, &
6.0_wp, 0.0_wp, 0.0_wp, -3.0_wp, 0.0_wp, 0.0_wp, &
-6.0_wp, 0.0_wp, 0.0_wp, 3.0_wp, 0.0_wp, 0.0_wp, &
6.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
6.0_wp, 0.0_wp, 0.0_wp, -1.0_wp, 0.0_wp, 0.0_wp, &
5.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
-5.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
-4.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-4.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
4.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
6.0_wp, 0.0_wp, 0.0_wp, -3.0_wp, 0.0_wp, 0.0_wp, &
-4.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
0.0_wp, 0.0_wp, -26.0_wp, 0.0_wp, 0.0_wp, -11.0_wp, &
0.0_wp, 0.0_wp, -10.0_wp, 0.0_wp, 0.0_wp, -5.0_wp, &
5.0_wp, 0.0_wp, 0.0_wp, -3.0_wp, 0.0_wp, 0.0_wp, &
-13.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
4.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
7.0_wp, 0.0_wp, 0.0_wp, -3.0_wp, 0.0_wp, 0.0_wp, &
4.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
5.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
-6.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
-5.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
-7.0_wp, 0.0_wp, 0.0_wp, 3.0_wp, 0.0_wp, 0.0_wp, &
5.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
13.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-4.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
5.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
-11.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
5.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
4.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
4.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
-4.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
6.0_wp, 0.0_wp, 0.0_wp, -3.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
-12.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
4.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-4.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, -1.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 0.0_wp, 0.0_wp, &
0.0_wp, 0.0_wp, -5.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, &
-7.0_wp, 0.0_wp, 0.0_wp, 4.0_wp, 0.0_wp, 0.0_wp, &
6.0_wp, 0.0_wp, 0.0_wp, -3.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
5.0_wp, 0.0_wp, 0.0_wp, -3.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, -1.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 0.0_wp, 0.0_wp, &
-5.0_wp, 0.0_wp, 0.0_wp, 3.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
12.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, -1.0_wp, 0.0_wp, 0.0_wp, &
-4.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
4.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
6.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
5.0_wp, 0.0_wp, 0.0_wp, -3.0_wp, 0.0_wp, 0.0_wp, &
4.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
-6.0_wp, 0.0_wp, 0.0_wp, 3.0_wp, 0.0_wp, 0.0_wp, &
4.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
6.0_wp, 0.0_wp, 0.0_wp, -3.0_wp, 0.0_wp, 0.0_wp, &
6.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-6.0_wp, 0.0_wp, 0.0_wp, 3.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
7.0_wp, 0.0_wp, 0.0_wp, -4.0_wp, 0.0_wp, 0.0_wp, &
4.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
-5.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
5.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-6.0_wp, 0.0_wp, 0.0_wp, 3.0_wp, 0.0_wp, 0.0_wp, &
-6.0_wp, 0.0_wp, 0.0_wp, 3.0_wp, 0.0_wp, 0.0_wp, &
-4.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
10.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
-4.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
7.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
7.0_wp, 0.0_wp, 0.0_wp, -3.0_wp, 0.0_wp, 0.0_wp, &
4.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
11.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, &
5.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
-6.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
4.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
5.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
-4.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
-4.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp, &
4.0_wp, 0.0_wp, 0.0_wp, -2.0_wp, 0.0_wp, 0.0_wp, &
3.0_wp, 0.0_wp, 0.0_wp, -1.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 0.0_wp, 0.0_wp, &
-3.0_wp, 0.0_wp, 0.0_wp, 2.0_wp, 0.0_wp, 0.0_wp], [6,nls])
!
! Planetary argument multipliers
! : L L' F D Om Me Ve E Ma Ju Sa Ur Ne pre
! Coefficients for fundamental arguments
integer,dimension(14,npl),parameter :: napl = reshape( [ &
0, 0, 0, 0, 0, 0, 0, 8,-16, 4, 5, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -8, 16, -4, -5, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 8,-16, 4, 5, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -1, 2, 2, &
0, 0, 0, 0, 0, 0, 0, -4, 8, -1, -5, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 4, -8, 3, 0, 0, 0, 1, &
0, 0, 1, -1, 1, 0, 0, 3, -8, 3, 0, 0, 0, 0, &
-1, 0, 0, 0, 0, 0, 10, -3, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 0, 0, -2, 6, -3, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 4, -8, 3, 0, 0, 0, 0, &
0, 0, 1, -1, 1, 0, 0, -5, 8, -3, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -4, 8, -3, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 4, -8, 1, 5, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -5, 6, 4, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 2, -5, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 2, -5, 0, 0, 1, &
0, 0, 1, -1, 1, 0, 0, -1, 0, 2, -5, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 2, -5, 0, 0, 0, &
0, 0, 1, -1, 1, 0, 0, -1, 0, -2, 5, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 0, 0, -2, 5, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 0, 0, -2, 5, 0, 0, 2, &
2, 0, -1, -1, 0, 0, 0, 3, -7, 0, 0, 0, 0, 0, &
1, 0, 0, -2, 0, 0, 19,-21, 3, 0, 0, 0, 0, 0, &
0, 0, 1, -1, 1, 0, 2, -4, 0, -3, 0, 0, 0, 0, &
1, 0, 0, -1, 1, 0, 0, -1, 0, 2, 0, 0, 0, 0, &
0, 0, 1, -1, 1, 0, 0, -1, 0, -4, 10, 0, 0, 0, &
-2, 0, 0, 2, 1, 0, 0, 2, 0, 0, -5, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 3, -7, 4, 0, 0, 0, 0, 0, &
0, 0, -1, 1, 0, 0, 0, 1, 0, 1, -1, 0, 0, 0, &
-2, 0, 0, 2, 1, 0, 0, 2, 0, -2, 0, 0, 0, 0, &
-1, 0, 0, 0, 0, 0, 18,-16, 0, 0, 0, 0, 0, 0, &
-2, 0, 1, 1, 2, 0, 0, 1, 0, -2, 0, 0, 0, 0, &
-1, 0, 1, -1, 1, 0, 18,-17, 0, 0, 0, 0, 0, 0, &
-1, 0, 0, 1, 1, 0, 0, 2, -2, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, -8, 13, 0, 0, 0, 0, 0, 2, &
0, 0, 2, -2, 2, 0, -8, 11, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, -8, 13, 0, 0, 0, 0, 0, 1, &
0, 0, 1, -1, 1, 0, -8, 12, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 8,-13, 0, 0, 0, 0, 0, 0, &
0, 0, 1, -1, 1, 0, 8,-14, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 8,-13, 0, 0, 0, 0, 0, 1, &
-2, 0, 0, 2, 1, 0, 0, 2, 0, -4, 5, 0, 0, 0, &
-2, 0, 0, 2, 2, 0, 3, -3, 0, 0, 0, 0, 0, 0, &
-2, 0, 0, 2, 0, 0, 0, 2, 0, -3, 1, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 3, -5, 0, 2, 0, 0, 0, 0, &
-2, 0, 0, 2, 0, 0, 0, 2, 0, -4, 3, 0, 0, 0, &
0, 0, -1, 1, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 0, -1, 2, 0, 0, 0, 0, 0, &
0, 0, 1, -1, 2, 0, 0, -2, 2, 0, 0, 0, 0, 0, &
-1, 0, 1, 0, 1, 0, 3, -5, 0, 0, 0, 0, 0, 0, &
-1, 0, 0, 1, 0, 0, 3, -4, 0, 0, 0, 0, 0, 0, &
-2, 0, 0, 2, 0, 0, 0, 2, 0, -2, -2, 0, 0, 0, &
-2, 0, 2, 0, 2, 0, 0, -5, 9, 0, 0, 0, 0, 0, &
0, 0, 1, -1, 1, 0, 0, -1, 0, 0, 0, -1, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, &
0, 0, 1, -1, 1, 0, 0, -1, 0, 0, 0, 0, 2, 0, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 2, 1, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 2, 2, &
-1, 0, 0, 1, 0, 0, 0, 3, -4, 0, 0, 0, 0, 0, &
0, 0, -1, 1, 0, 0, 0, 1, 0, 0, 2, 0, 0, 0, &
0, 0, 1, -1, 2, 0, 0, -1, 0, 0, 2, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 0, -9, 17, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 2, 0, -3, 5, 0, 0, 0, 0, 0, 0, &
0, 0, 1, -1, 1, 0, 0, -1, 0, -1, 2, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 1, -2, 0, 0, 0, &
1, 0, 0, -2, 0, 0, 17,-16, 0, -2, 0, 0, 0, 0, &
0, 0, 1, -1, 1, 0, 0, -1, 0, 1, -3, 0, 0, 0, &
-2, 0, 0, 2, 1, 0, 0, 5, -6, 0, 0, 0, 0, 0, &
0, 0, -2, 2, 0, 0, 0, 9,-13, 0, 0, 0, 0, 0, &
0, 0, 1, -1, 2, 0, 0, -1, 0, 0, 1, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 0, &
0, 0, -1, 1, 0, 0, 0, 1, 0, 0, 1, 0, 0, 0, &
0, 0, -2, 2, 0, 0, 5, -6, 0, 0, 0, 0, 0, 0, &
0, 0, -1, 1, 1, 0, 5, -7, 0, 0, 0, 0, 0, 0, &
-2, 0, 0, 2, 0, 0, 6, -8, 0, 0, 0, 0, 0, 0, &
2, 0, 1, -3, 1, 0, -6, 7, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 2, 0, 0, 0, 0, 1, 0, 0, 0, 0, &
0, 0, -1, 1, 1, 0, 0, 1, 0, 1, 0, 0, 0, 0, &
0, 0, 1, -1, 1, 0, 0, -1, 0, 0, 0, 2, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 2, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 2, 0, 2, &
0, 0, 0, 0, 0, 0, 0, -8, 15, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, -8, 15, 0, 0, 0, 0, 1, &
0, 0, 1, -1, 1, 0, 0, -9, 15, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 8,-15, 0, 0, 0, 0, 0, &
1, 0, -1, -1, 0, 0, 0, 8,-15, 0, 0, 0, 0, 0, &
2, 0, 0, -2, 0, 0, 2, -5, 0, 0, 0, 0, 0, 0, &
-2, 0, 0, 2, 0, 0, 0, 2, 0, -5, 5, 0, 0, 0, &
2, 0, 0, -2, 1, 0, 0, -6, 8, 0, 0, 0, 0, 0, &
2, 0, 0, -2, 1, 0, 0, -2, 0, 3, 0, 0, 0, 0, &
-2, 0, 1, 1, 0, 0, 0, 1, 0, -3, 0, 0, 0, 0, &
-2, 0, 1, 1, 1, 0, 0, 1, 0, -3, 0, 0, 0, 0, &
-2, 0, 0, 2, 0, 0, 0, 2, 0, -3, 0, 0, 0, 0, &
-2, 0, 0, 2, 0, 0, 0, 6, -8, 0, 0, 0, 0, 0, &
-2, 0, 0, 2, 0, 0, 0, 2, 0, -1, -5, 0, 0, 0, &
-1, 0, 0, 1, 0, 0, 0, 1, 0, -1, 0, 0, 0, 0, &
-1, 0, 1, 1, 1, 0,-20, 20, 0, 0, 0, 0, 0, 0, &
1, 0, 0, -2, 0, 0, 20,-21, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 0, 8,-15, 0, 0, 0, 0, 0, &
0, 0, 2, -2, 1, 0, 0,-10, 15, 0, 0, 0, 0, 0, &
0, 0, -1, 1, 0, 0, 0, 1, 0, 1, 0, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, &
0, 0, 1, -1, 2, 0, 0, -1, 0, 1, 0, 0, 0, 0, &
0, 0, 1, -1, 1, 0, 0, -1, 0, -2, 4, 0, 0, 0, &
2, 0, 0, -2, 1, 0, -6, 8, 0, 0, 0, 0, 0, 0, &
0, 0, -2, 2, 1, 0, 5, -6, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -1, 0, 0, 1, &
0, 0, 1, -1, 1, 0, 0, -1, 0, 0, -1, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, &
0, 0, 1, -1, 1, 0, 0, -1, 0, 0, 1, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 2, &
0, 0, 2, -2, 1, 0, 0, -9, 13, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 0, 7,-13, 0, 0, 0, 0, 0, &
-2, 0, 0, 2, 0, 0, 0, 5, -6, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 9,-17, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -9, 17, 0, 0, 0, 0, 2, &
1, 0, 0, -1, 1, 0, 0, -3, 4, 0, 0, 0, 0, 0, &
1, 0, 0, -1, 1, 0, -3, 4, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 2, 0, 0, -1, 2, 0, 0, 0, 0, 0, &
0, 0, -1, 1, 1, 0, 0, 0, 2, 0, 0, 0, 0, 0, &
0, 0, -2, 2, 0, 1, 0, -2, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 3, -5, 0, 2, 0, 0, 0, 0, &
-2, 0, 0, 2, 1, 0, 0, 2, 0, -3, 1, 0, 0, 0, &
-2, 0, 0, 2, 1, 0, 3, -3, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 8,-13, 0, 0, 0, 0, 0, 0, &
0, 0, -1, 1, 0, 0, 8,-12, 0, 0, 0, 0, 0, 0, &
0, 0, 2, -2, 1, 0, -8, 11, 0, 0, 0, 0, 0, 0, &
-1, 0, 0, 1, 0, 0, 0, 2, -2, 0, 0, 0, 0, 0, &
-1, 0, 0, 0, 1, 0, 18,-16, 0, 0, 0, 0, 0, 0, &
0, 0, 1, -1, 1, 0, 0, -1, 0, -1, 1, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 3, -7, 4, 0, 0, 0, 0, 0, &
-2, 0, 1, 1, 1, 0, 0, -3, 7, 0, 0, 0, 0, 0, &
0, 0, 1, -1, 2, 0, 0, -1, 0, -2, 5, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 0, 0, 0, -2, 5, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 0, -4, 8, -3, 0, 0, 0, 0, &
1, 0, 0, 0, 1, 0,-10, 3, 0, 0, 0, 0, 0, 0, &
0, 0, 2, -2, 1, 0, 0, -2, 0, 0, 0, 0, 0, 0, &
-1, 0, 0, 0, 1, 0, 10, -3, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 0, 4, -8, 3, 0, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 0, 0, 0, 2, -5, 0, 0, 0, &
0, 0, -1, 1, 0, 0, 0, 1, 0, 2, -5, 0, 0, 0, &
2, 0, -1, -1, 1, 0, 0, 3, -7, 0, 0, 0, 0, 0, &
-2, 0, 0, 2, 0, 0, 0, 2, 0, 0, -5, 0, 0, 0, &
0, 0, 0, 0, 1, 0, -3, 7, -4, 0, 0, 0, 0, 0, &
-2, 0, 0, 2, 0, 0, 0, 2, 0, -2, 0, 0, 0, 0, &
1, 0, 0, 0, 1, 0,-18, 16, 0, 0, 0, 0, 0, 0, &
-2, 0, 1, 1, 1, 0, 0, 1, 0, -2, 0, 0, 0, 0, &
0, 0, 1, -1, 2, 0, -8, 12, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 1, 0, -8, 13, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 1, -2, 0, 0, 0, 0, 1, &
0, 0, 1, -1, 1, 0, 0, 0, -2, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 1, -2, 0, 0, 0, 0, 0, &
0, 0, 1, -1, 1, 0, 0, -2, 2, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -1, 2, 0, 0, 0, 0, 1, &
-1, 0, 0, 1, 1, 0, 3, -4, 0, 0, 0, 0, 0, 0, &
-1, 0, 0, 1, 1, 0, 0, 3, -4, 0, 0, 0, 0, 0, &
0, 0, 1, -1, 1, 0, 0, -1, 0, 0, -2, 0, 0, 0, &
0, 0, 1, -1, 1, 0, 0, -1, 0, 0, 2, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 2, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 2, 0, 0, 2, &
0, 0, 1, -1, 0, 0, 3, -6, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 1, 0, -3, 5, 0, 0, 0, 0, 0, 0, &
0, 0, 1, -1, 2, 0, -3, 4, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 0, -2, 4, 0, 0, 0, 0, 0, &
0, 0, 2, -2, 1, 0, -5, 6, 0, 0, 0, 0, 0, 0, &
0, 0, -1, 1, 0, 0, 5, -7, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 5, -8, 0, 0, 0, 0, 0, 0, &
-2, 0, 0, 2, 1, 0, 6, -8, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 0, -8, 15, 0, 0, 0, 0, 0, &
-2, 0, 0, 2, 1, 0, 0, 2, 0, -3, 0, 0, 0, 0, &
-2, 0, 0, 2, 1, 0, 0, 6, -8, 0, 0, 0, 0, 0, &
1, 0, 0, -1, 1, 0, 0, -1, 0, 1, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 3, -5, 0, 0, 0, &
0, 0, 1, -1, 1, 0, 0, -1, 0, -1, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 0, 0, -1, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 1, &
0, 0, 1, -1, 1, 0, 0, -1, 0, 1, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 2, &
0, 0, 1, -1, 2, 0, 0, -1, 0, 0, -1, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 0, 0, 0, 0, -1, 0, 0, 0, &
0, 0, -1, 1, 0, 0, 0, 1, 0, 0, -1, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -7, 13, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 7,-13, 0, 0, 0, 0, 0, &
2, 0, 0, -2, 1, 0, 0, -5, 6, 0, 0, 0, 0, 0, &
0, 0, 2, -2, 1, 0, 0, -8, 11, 0, 0, 0, 0, 0, &
0, 0, 2, -2, 1, -1, 0, 2, 0, 0, 0, 0, 0, 0, &
-2, 0, 0, 2, 0, 0, 0, 4, -4, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 2, -2, 0, 0, 0, &
0, 0, 1, -1, 1, 0, 0, -1, 0, 0, 3, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 3, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 3, 0, 0, 2, &
-2, 0, 0, 2, 0, 0, 3, -3, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 2, 0, 0, -4, 8, -3, 0, 0, 0, 0, &
0, 0, 0, 0, 2, 0, 0, 4, -8, 3, 0, 0, 0, 0, &
2, 0, 0, -2, 1, 0, 0, -2, 0, 2, 0, 0, 0, 0, &
0, 0, 1, -1, 2, 0, 0, -1, 0, 2, 0, 0, 0, 0, &
0, 0, 1, -1, 2, 0, 0, 0, -2, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 0, 1, -2, 0, 0, 0, 0, 0, &
0, 0, -1, 1, 0, 0, 0, 2, -2, 0, 0, 0, 0, 0, &
0, 0, -1, 1, 0, 0, 0, 1, 0, 0, -2, 0, 0, 0, &
0, 0, 2, -2, 1, 0, 0, -2, 0, 0, 2, 0, 0, 0, &
0, 0, 1, -1, 1, 0, 3, -6, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 3, -5, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 3, -5, 0, 0, 0, 0, 0, 0, &
0, 0, 1, -1, 1, 0, -3, 4, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, -3, 5, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, -3, 5, 0, 0, 0, 0, 0, 2, &
0, 0, 2, -2, 2, 0, -3, 3, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, -3, 5, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 2, -4, 0, 0, 0, 0, 1, &
0, 0, 1, -1, 1, 0, 0, 1, -4, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 2, -4, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -2, 4, 0, 0, 0, 0, 1, &
0, 0, 1, -1, 1, 0, 0, -3, 4, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -2, 4, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, -2, 4, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -5, 8, 0, 0, 0, 0, 0, 2, &
0, 0, 2, -2, 2, 0, -5, 6, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, -5, 8, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -5, 8, 0, 0, 0, 0, 0, 1, &
0, 0, 1, -1, 1, 0, -5, 7, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, -5, 8, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 5, -8, 0, 0, 0, 0, 0, 0, &
0, 0, 1, -1, 2, 0, 0, -1, 0, -1, 0, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 0, 0, 0, -1, 0, 0, 0, 0, &
0, 0, -1, 1, 0, 0, 0, 1, 0, -1, 0, 0, 0, 0, &
0, 0, 2, -2, 1, 0, 0, -2, 0, 1, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -6, 11, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 6,-11, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, -1, 0, 4, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 1, 0, -4, 0, 0, 0, 0, 0, 0, &
2, 0, 0, -2, 1, 0, -3, 3, 0, 0, 0, 0, 0, 0, &
-2, 0, 0, 2, 0, 0, 0, 2, 0, 0, -2, 0, 0, 0, &
0, 0, 2, -2, 1, 0, 0, -7, 9, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 4, -5, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 2, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 2, 0, 0, 0, 1, &
0, 0, 1, -1, 1, 0, 0, -1, 0, 2, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 2, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 2, 0, 0, 0, 2, &
0, 0, 2, -2, 2, 0, 0, -2, 0, 2, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 5, 0, 0, 2, &
0, 0, 0, 0, 1, 0, 3, -5, 0, 0, 0, 0, 0, 0, &
0, 0, -1, 1, 0, 0, 3, -4, 0, 0, 0, 0, 0, 0, &
0, 0, 2, -2, 1, 0, -3, 3, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 0, 2, -4, 0, 0, 0, 0, 0, &
0, 0, 2, -2, 1, 0, 0, -4, 4, 0, 0, 0, 0, 0, &
0, 0, 1, -1, 2, 0, -5, 7, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 3, -6, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -3, 6, 0, 0, 0, 0, 1, &
0, 0, 1, -1, 1, 0, 0, -4, 6, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -3, 6, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, -3, 6, 0, 0, 0, 0, 2, &
0, 0, -1, 1, 0, 0, 2, -2, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 2, -3, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -5, 9, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, -5, 9, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 5, -9, 0, 0, 0, 0, 0, &
0, 0, -1, 1, 0, 0, 0, 1, 0, -2, 0, 0, 0, 0, &
0, 0, 2, -2, 1, 0, 0, -2, 0, 2, 0, 0, 0, 0, &
-2, 0, 1, 1, 1, 0, 0, 1, 0, 0, 0, 0, 0, 0, &
0, 0, -2, 2, 0, 0, 3, -3, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, -6, 10, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, -6, 10, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -2, 3, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -2, 3, 0, 0, 0, 0, 0, 1, &
0, 0, 1, -1, 1, 0, -2, 2, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 2, -3, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 2, -3, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 3, 0, 0, 0, 1, &
0, 0, 1, -1, 1, 0, 0, -1, 0, 3, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 3, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 3, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 4, -8, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -4, 8, 0, 0, 0, 0, 2, &
0, 0, -2, 2, 0, 0, 0, 2, 0, -2, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -4, 7, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, -4, 7, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 4, -7, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 1, 0, -2, 3, 0, 0, 0, 0, 0, 0, &
0, 0, 2, -2, 1, 0, 0, -2, 0, 3, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -5, 10, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 1, 0, -1, 2, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 4, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, -3, 5, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, -3, 5, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 3, -5, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 1, -2, 0, 0, 0, 0, 0, 1, &
0, 0, 1, -1, 1, 0, 1, -3, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 1, -2, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, -1, 2, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, -1, 2, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -7, 11, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -7, 11, 0, 0, 0, 0, 0, 1, &
0, 0, -2, 2, 0, 0, 4, -4, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 2, -3, 0, 0, 0, 0, 0, &
0, 0, 2, -2, 1, 0, -4, 4, 0, 0, 0, 0, 0, 0, &
0, 0, -1, 1, 0, 0, 4, -5, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 1, -1, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, -4, 7, 0, 0, 0, 0, 0, 1, &
0, 0, 1, -1, 1, 0, -4, 6, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, -4, 7, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -4, 6, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -4, 6, 0, 0, 0, 0, 0, 1, &
0, 0, 1, -1, 1, 0, -4, 5, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, -4, 6, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 4, -6, 0, 0, 0, 0, 0, 0, &
-2, 0, 0, 2, 0, 0, 2, -2, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, &
0, 0, -1, 1, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 1, -1, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -1, 0, 5, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 1, -3, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -1, 3, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, -7, 12, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -1, 1, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -1, 1, 0, 0, 0, 0, 0, 1, &
0, 0, 1, -1, 1, 0, -1, 0, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 1, -1, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 1, -1, 0, 0, 0, 0, 0, 1, &
0, 0, 1, -1, 1, 0, 1, -2, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -2, 5, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, -1, 0, 4, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 1, 0, -4, 0, 0, 0, 0, &
0, 0, 0, 0, 1, 0, -1, 1, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -6, 10, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, -6, 10, 0, 0, 0, 0, 0, &
0, 0, 2, -2, 1, 0, 0, -3, 0, 3, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -3, 7, 0, 0, 0, 0, 2, &
-2, 0, 0, 2, 0, 0, 4, -4, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -5, 8, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 5, -8, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -1, 0, 3, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, -1, 0, 3, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 1, 0, -3, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 2, -4, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, -2, 4, 0, 0, 0, 0, 0, 1, &
0, 0, 1, -1, 1, 0, -2, 3, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, -2, 4, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -6, 9, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -6, 9, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 6, -9, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 0, 1, 0, -2, 0, 0, 0, 0, &
0, 0, 2, -2, 1, 0, -2, 2, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -4, 6, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 4, -6, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 3, -4, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -1, 0, 2, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 1, 0, -2, 0, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 0, 1, 0, -1, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, -5, 9, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 3, -4, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, -3, 4, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -3, 4, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 3, -4, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 3, -4, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 1, 0, 0, 2, -2, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 0, -1, 0, 2, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 1, 0, 0, -3, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 1, 0, 1, -5, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -1, 0, 1, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 1, 0, -1, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 1, 0, -1, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 1, 0, -3, 5, 0, 0, 0, &
0, 0, 0, 0, 1, 0, -3, 4, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 1, 0, 0, -2, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 2, -2, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 1, 0, 0, -1, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 0, -1, 0, 1, 0, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 0, -2, 2, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, -8, 14, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 1, 0, 2, -5, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 5, -8, 3, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 5, -8, 3, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, -1, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 3, -8, 3, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -3, 8, -3, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 1, 0, -2, 5, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -8, 12, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -8, 12, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 1, 0, 1, -2, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 1, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 2, 0, 0, 2, &
0, 0, 2, -2, 1, 0, -5, 5, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 1, 0, 1, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 1, 0, 1, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 1, 0, 1, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 3, -6, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, -3, 6, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, -3, 6, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, -1, 4, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -5, 7, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -5, 7, 0, 0, 0, 0, 0, 1, &
0, 0, 1, -1, 1, 0, -5, 6, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 5, -7, 0, 0, 0, 0, 0, 0, &
0, 0, 2, -2, 1, 0, 0, -1, 0, 1, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -1, 0, 1, 0, 0, 0, 0, &
0, 0, 0, 0, 0, -1, 0, 3, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 1, 0, 2, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, -2, 6, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 1, 0, 2, -2, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -6, 9, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 6, -9, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, -2, 2, 0, 0, 0, 0, 0, 1, &
0, 0, 1, -1, 1, 0, -2, 1, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 2, -2, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 2, -2, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 1, 0, 3, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, -5, 7, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 5, -7, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 1, 0, -2, 2, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 4, -5, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 1, -3, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, -1, 3, 0, 0, 0, 0, 0, 1, &
0, 0, 1, -1, 1, 0, -1, 2, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, -1, 3, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -7, 10, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -7, 10, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 3, -3, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, -4, 8, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -4, 5, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -4, 5, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 4, -5, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 1, 1, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, -2, 0, 5, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 0, 3, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -9, 13, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, -1, 5, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, -2, 0, 4, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 2, 0, -4, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -2, 7, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 2, 0, -3, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, -2, 5, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, -2, 5, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -6, 8, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -6, 8, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 6, -8, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 1, 0, 0, 2, 0, -2, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -3, 9, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 5, -6, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 5, -6, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 2, 0, -2, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 2, 0, -2, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 2, 0, -2, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -5, 10, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 4, -4, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 4, -4, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -3, 3, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 3, -3, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 3, -3, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 3, -3, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 2, 0, 0, -3, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, -5, 13, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 2, 0, -1, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 2, 0, -1, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 2, 0, 0, -2, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 2, 0, 0, -2, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 3, -2, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 3, -2, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 2, 0, 0, -1, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, -6, 15, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -8, 15, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -3, 9, -4, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 2, 0, 2, -5, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, -2, 8, -1, -5, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 6, -8, 3, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 1, &
0, 0, 1, -1, 1, 0, 0, 1, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, -6, 16, -4, -5, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, -2, 8, -3, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, -2, 8, -3, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 6, -8, 1, 5, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 2, 0, -2, 5, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 3, -5, 4, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -8, 11, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -8, 11, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, -8, 11, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 11, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 2, 0, 0, 1, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 3, -3, 0, 2, 0, 0, 0, 2, &
0, 0, 2, -2, 1, 0, 0, 4, -8, 3, 0, 0, 0, 0, &
0, 0, 1, -1, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, &
0, 0, 2, -2, 1, 0, 0, -4, 8, -3, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 1, 2, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 2, 0, 1, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -3, 7, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 0, 4, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -5, 6, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -5, 6, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 5, -6, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 5, -6, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 2, 0, 2, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, -1, 6, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 7, -9, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 2, -1, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 2, -1, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 6, -7, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 5, -5, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -1, 4, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, -1, 4, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -7, 9, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -7, 9, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 4, -3, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 3, -1, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -4, 4, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 4, -4, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 4, -4, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 4, -4, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 2, 1, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, -3, 0, 5, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 1, 1, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 1, 1, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -9, 12, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 3, 0, -4, 0, 0, 0, 0, &
0, 0, 2, -2, 1, 0, 1, -1, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 7, -8, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 3, 0, -3, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 3, 0, -3, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -2, 6, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -6, 7, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 6, -7, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 6, -6, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 3, 0, -2, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 3, 0, -2, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 5, -4, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 3, -2, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 3, -2, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 3, 0, -1, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 3, 0, -1, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 3, 0, 0, -2, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 4, -2, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 3, 0, 0, -1, 0, 0, 2, &
0, 0, 2, -2, 1, 0, 0, 1, 0, -1, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, -8, 16, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 3, 0, 2, -5, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 7, -8, 3, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, -5, 16, -4, -5, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 3, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, -1, 8, -3, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -8, 10, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -8, 10, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, -8, 10, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 2, 2, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 3, 0, 1, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -3, 8, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -5, 5, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 5, -5, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 5, -5, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 5, -5, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 7, -7, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 7, -7, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 6, -5, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 7, -8, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 5, -3, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 4, -3, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 1, 2, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -9, 11, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -9, 11, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 4, 0, -4, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 4, 0, -3, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -6, 6, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 6, -6, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 6, -6, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 4, 0, -2, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 6, -4, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 3, -1, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 3, -1, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 3, -1, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 4, 0, -1, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 4, 0, 0, -2, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 5, -2, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 4, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 8, -9, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 5, -4, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 2, 1, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 2, 1, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 2, 1, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, -7, 7, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 7, -7, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 4, -2, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 4, -2, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 4, -2, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 4, -2, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 5, 0, -4, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 5, 0, -3, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 5, 0, -2, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 3, 0, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -8, 8, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 8, -8, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 5, -3, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 5, -3, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, -9, 9, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, -9, 9, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, -9, 9, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 9, -9, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 6, -4, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 6, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 6, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 6, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 6, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 6, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 6, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 6, 0, 0, 0, 0, 0, 1, &
0, 0, 0, 0, 0, 0, 0, 6, 0, 0, 0, 0, 0, 2, &
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 2, &
1, 0, 0, -2, 0, 0, 0, 2, 0, -2, 0, 0, 0, 0, &
1, 0, 0, -2, 0, 0, 2, -2, 0, 0, 0, 0, 0, 0, &
1, 0, 0, -2, 0, 0, 0, 1, 0, -1, 0, 0, 0, 0, &
1, 0, 0, -2, 0, 0, 1, -1, 0, 0, 0, 0, 0, 0, &
-1, 0, 0, 0, 0, 0, 3, -3, 0, 0, 0, 0, 0, 0, &
-1, 0, 0, 0, 0, 0, 0, 2, 0, -2, 0, 0, 0, 0, &
-1, 0, 0, 2, 0, 0, 0, 4, -8, 3, 0, 0, 0, 0, &
1, 0, 0, -2, 0, 0, 0, 4, -8, 3, 0, 0, 0, 0, &
-2, 0, 0, 2, 0, 0, 0, 4, -8, 3, 0, 0, 0, 0, &
-1, 0, 0, 0, 0, 0, 0, 2, 0, -3, 0, 0, 0, 0, &
-1, 0, 0, 0, 0, 0, 0, 1, 0, -1, 0, 0, 0, 0, &
-1, 0, 0, 0, 0, 0, 1, -1, 0, 0, 0, 0, 0, 0, &
-1, 0, 0, 2, 0, 0, 2, -2, 0, 0, 0, 0, 0, 0, &
1, 0, -1, 1, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, &
-1, 0, 0, 2, 0, 0, 0, 2, 0, -3, 0, 0, 0, 0, &
-2, 0, 0, 0, 0, 0, 0, 2, 0, -3, 0, 0, 0, 0, &
1, 0, 0, 0, 0, 0, 0, 4, -8, 3, 0, 0, 0, 0, &
-1, 0, 1, -1, 1, 0, 0, -1, 0, 0, 0, 0, 0, 0, &
1, 0, 1, -1, 1, 0, 0, -1, 0, 0, 0, 0, 0, 0, &
-1, 0, 0, 0, 0, 0, 0, 4, -8, 3, 0, 0, 0, 0, &
-1, 0, 0, 2, 1, 0, 0, 2, 0, -2, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 2, 0, -2, 0, 0, 0, 0, &
-1, 0, 0, 2, 0, 0, 0, 2, 0, -2, 0, 0, 0, 0, &
-1, 0, 0, 2, 0, 0, 3, -3, 0, 0, 0, 0, 0, 0, &
1, 0, 0, -2, 1, 0, 0, -2, 0, 2, 0, 0, 0, 0, &
1, 0, 2, -2, 2, 0, -3, 3, 0, 0, 0, 0, 0, 0, &
1, 0, 2, -2, 2, 0, 0, -2, 0, 2, 0, 0, 0, 0, &
1, 0, 0, 0, 0, 0, 1, -1, 0, 0, 0, 0, 0, 0, &
1, 0, 0, 0, 0, 0, 0, 1, 0, -1, 0, 0, 0, 0, &
0, 0, 0, -2, 0, 0, 2, -2, 0, 0, 0, 0, 0, 0, &
0, 0, 0, -2, 0, 0, 0, 1, 0, -1, 0, 0, 0, 0, &
0, 0, 2, 0, 2, 0, -2, 2, 0, 0, 0, 0, 0, 0, &
0, 0, 2, 0, 2, 0, 0, -1, 0, 1, 0, 0, 0, 0, &
0, 0, 2, 0, 2, 0, -1, 1, 0, 0, 0, 0, 0, 0, &
0, 0, 2, 0, 2, 0, -2, 3, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 2, 0, 0, 0, 2, 0, -2, 0, 0, 0, 0, &
0, 0, 1, 1, 2, 0, 0, 1, 0, 0, 0, 0, 0, 0, &
1, 0, 2, 0, 2, 0, 0, 1, 0, 0, 0, 0, 0, 0, &
-1, 0, 2, 0, 2, 0, 10, -3, 0, 0, 0, 0, 0, 0, &
0, 0, 1, 1, 1, 0, 0, 1, 0, 0, 0, 0, 0, 0, &
1, 0, 2, 0, 2, 0, 0, 1, 0, 0, 0, 0, 0, 0, &
0, 0, 2, 0, 2, 0, 0, 4, -8, 3, 0, 0, 0, 0, &
0, 0, 2, 0, 2, 0, 0, -4, 8, -3, 0, 0, 0, 0, &
-1, 0, 2, 0, 2, 0, 0, -4, 8, -3, 0, 0, 0, 0, &
2, 0, 2, -2, 2, 0, 0, -2, 0, 3, 0, 0, 0, 0, &
1, 0, 2, 0, 1, 0, 0, -2, 0, 3, 0, 0, 0, 0, &
0, 0, 1, 1, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, &
-1, 0, 2, 0, 1, 0, 0, 1, 0, 0, 0, 0, 0, 0, &
-2, 0, 2, 2, 2, 0, 0, 2, 0, -2, 0, 0, 0, 0, &
0, 0, 2, 0, 2, 0, 2, -3, 0, 0, 0, 0, 0, 0, &
0, 0, 2, 0, 2, 0, 1, -1, 0, 0, 0, 0, 0, 0, &
0, 0, 2, 0, 2, 0, 0, 1, 0, -1, 0, 0, 0, 0, &
0, 0, 2, 0, 2, 0, 2, -2, 0, 0, 0, 0, 0, 0, &
-1, 0, 2, 2, 2, 0, 0, -1, 0, 1, 0, 0, 0, 0, &
1, 0, 2, 0, 2, 0, -1, 1, 0, 0, 0, 0, 0, 0, &
-1, 0, 2, 2, 2, 0, 0, 2, 0, -3, 0, 0, 0, 0, &
2, 0, 2, 0, 2, 0, 0, 2, 0, -3, 0, 0, 0, 0, &
1, 0, 2, 0, 2, 0, 0, -4, 8, -3, 0, 0, 0, 0, &
1, 0, 2, 0, 2, 0, 0, 4, -8, 3, 0, 0, 0, 0, &
1, 0, 1, 1, 1, 0, 0, 1, 0, 0, 0, 0, 0, 0, &
0, 0, 2, 0, 2, 0, 0, 1, 0, 0, 0, 0, 0, 0, &
2, 0, 2, 0, 1, 0, 0, 1, 0, 0, 0, 0, 0, 0, &
-1, 0, 2, 2, 2, 0, 0, 2, 0, -2, 0, 0, 0, 0, &
-1, 0, 2, 2, 2, 0, 3, -3, 0, 0, 0, 0, 0, 0, &
1, 0, 2, 0, 2, 0, 1, -1, 0, 0, 0, 0, 0, 0, &
0, 0, 2, 2, 2, 0, 0, 2, 0, -2, 0, 0, 0, 0 ] , [14,npl])
!
! Planetary nutation coefficients, unit 1e-7 arcsec
! longitude (sin, cos), obliquity (sin, cos)
!
! Longitude and obliquity coefficients
integer,dimension(4,npl),parameter :: icpl = reshape([ &
1440, 0, 0, 0, &
56, -117, -42, -40, &
125, -43, 0, -54, &
0, 5, 0, 0, &
3, -7, -3, 0, &
3, 0, 0, -2, &
-114, 0, 0, 61, &
-219, 89, 0, 0, &
-3, 0, 0, 0, &
-462, 1604, 0, 0, &
99, 0, 0, -53, &
-3, 0, 0, 2, &
0, 6, 2, 0, &
3, 0, 0, 0, &
-12, 0, 0, 0, &
14, -218, 117, 8, &
31, -481, -257, -17, &
-491, 128, 0, 0, &
-3084, 5123, 2735, 1647, &
-1444, 2409, -1286, -771, &
11, -24, -11, -9, &
26, -9, 0, 0, &
103, -60, 0, 0, &
0, -13, -7, 0, &
-26, -29, -16, 14, &
9, -27, -14, -5, &
12, 0, 0, -6, &
-7, 0, 0, 0, &
0, 24, 0, 0, &
284, 0, 0, -151, &
226, 101, 0, 0, &
0, -8, -2, 0, &
0, -6, -3, 0, &
5, 0, 0, -3, &
-41, 175, 76, 17, &
0, 15, 6, 0, &
425, 212, -133, 269, &
1200, 598, 319, -641, &
235, 334, 0, 0, &
11, -12, -7, -6, &
5, -6, 3, 3, &
-5, 0, 0, 3, &
6, 0, 0, -3, &
15, 0, 0, 0, &
13, 0, 0, -7, &
-6, -9, 0, 0, &
266, -78, 0, 0, &
-460, -435, -232, 246, &
0, 15, 7, 0, &
-3, 0, 0, 2, &
0, 131, 0, 0, &
4, 0, 0, 0, &
0, 3, 0, 0, &
0, 4, 2, 0, &
0, 3, 0, 0, &
-17, -19, -10, 9, &
-9, -11, 6, -5, &
-6, 0, 0, 3, &
-16, 8, 0, 0, &
0, 3, 0, 0, &
11, 24, 11, -5, &
-3, -4, -2, 1, &
3, 0, 0, -1, &
0, -8, -4, 0, &
0, 3, 0, 0, &
0, 5, 0, 0, &
0, 3, 2, 0, &
-6, 4, 2, 3, &
-3, -5, 0, 0, &
-5, 0, 0, 2, &
4, 24, 13, -2, &
-42, 20, 0, 0, &
-10, 233, 0, 0, &
-3, 0, 0, 1, &
78, -18, 0, 0, &
0, 3, 1, 0, &
0, -3, -1, 0, &
0, -4, -2, 1, &
0, -8, -4, -1, &
0, -5, 3, 0, &
-7, 0, 0, 3, &
-14, 8, 3, 6, &
0, 8, -4, 0, &
0, 19, 10, 0, &
45, -22, 0, 0, &
-3, 0, 0, 0, &
0, -3, 0, 0, &
0, 3, 0, 0, &
3, 5, 3, -2, &
89, -16, -9, -48, &
0, 3, 0, 0, &
-3, 7, 4, 2, &
-349, -62, 0, 0, &
-15, 22, 0, 0, &
-3, 0, 0, 0, &
-53, 0, 0, 0, &
5, 0, 0, -3, &
0, -8, 0, 0, &
15, -7, -4, -8, &
-3, 0, 0, 1, &
-21, -78, 0, 0, &
20, -70, -37, -11, &
0, 6, 3, 0, &
5, 3, 2, -2, &
-17, -4, -2, 9, &
0, 6, 3, 0, &
32, 15, -8, 17, &
174, 84, 45, -93, &
11, 56, 0, 0, &
-66, -12, -6, 35, &
47, 8, 4, -25, &
0, 8, 4, 0, &
10, -22, -12, -5, &
-3, 0, 0, 2, &
-24, 12, 0, 0, &
5, -6, 0, 0, &
3, 0, 0, -2, &
4, 3, 1, -2, &
0, 29, 15, 0, &
-5, -4, -2, 2, &
8, -3, -1, -5, &
0, -3, 0, 0, &
10, 0, 0, 0, &
3, 0, 0, -2, &
-5, 0, 0, 3, &
46, 66, 35, -25, &
-14, 7, 0, 0, &
0, 3, 2, 0, &
-5, 0, 0, 0, &
-68, -34, -18, 36, &
0, 14, 7, 0, &
10, -6, -3, -5, &
-5, -4, -2, 3, &
-3, 5, 2, 1, &
76, 17, 9, -41, &
84, 298, 159, -45, &
3, 0, 0, -1, &
-3, 0, 0, 2, &
-3, 0, 0, 1, &
-82, 292, 156, 44, &
-73, 17, 9, 39, &
-9, -16, 0, 0, &
3, 0, -1, -2, &
-3, 0, 0, 0, &
-9, -5, -3, 5, &
-439, 0, 0, 0, &
57, -28, -15, -30, &
0, -6, -3, 0, &
-4, 0, 0, 2, &
-40, 57, 30, 21, &
23, 7, 3, -13, &
273, 80, 43, -146, &
-449, 430, 0, 0, &
-8, -47, -25, 4, &
6, 47, 25, -3, &
0, 23, 13, 0, &
-3, 0, 0, 2, &
3, -4, -2, -2, &
-48, -110, -59, 26, &
51, 114, 61, -27, &
-133, 0, 0, 57, &
0, 4, 0, 0, &
-21, -6, -3, 11, &
0, -3, -1, 0, &
-11, -21, -11, 6, &
-18, -436, -233, 9, &
35, -7, 0, 0, &
0, 5, 3, 0, &
11, -3, -1, -6, &
-5, -3, -1, 3, &
-53, -9, -5, 28, &
0, 3, 2, 1, &
4, 0, 0, -2, &
0, -4, 0, 0, &
-50, 194, 103, 27, &
-13, 52, 28, 7, &
-91, 248, 0, 0, &
6, 49, 26, -3, &
-6, -47, -25, 3, &
0, 5, 3, 0, &
52, 23, 10, -23, &
-3, 0, 0, 1, &
0, 5, 3, 0, &
-4, 0, 0, 0, &
-4, 8, 3, 2, &
10, 0, 0, 0, &
3, 0, 0, -2, &
0, 8, 4, 0, &
0, 8, 4, 1, &
-4, 0, 0, 0, &
-4, 0, 0, 0, &
-8, 4, 2, 4, &
8, -4, -2, -4, &
0, 15, 7, 0, &
-138, 0, 0, 0, &
0, -7, -3, 0, &
0, -7, -3, 0, &
54, 0, 0, -29, &
0, 10, 4, 0, &
-7, 0, 0, 3, &
-37, 35, 19, 20, &
0, 4, 0, 0, &
-4, 9, 0, 0, &
8, 0, 0, -4, &
-9, -14, -8, 5, &
-3, -9, -5, 3, &
-145, 47, 0, 0, &
-10, 40, 21, 5, &
11, -49, -26, -7, &
-2150, 0, 0, 932, &
-12, 0, 0, 5, &
85, 0, 0, -37, &
4, 0, 0, -2, &
3, 0, 0, -2, &
-86, 153, 0, 0, &
-6, 9, 5, 3, &
9, -13, -7, -5, &
-8, 12, 6, 4, &
-51, 0, 0, 22, &
-11, -268, -116, 5, &
0, 12, 5, 0, &
0, 7, 3, 0, &
31, 6, 3, -17, &
140, 27, 14, -75, &
57, 11, 6, -30, &
-14, -39, 0, 0, &
0, -6, -2, 0, &
4, 15, 8, -2, &
0, 4, 0, 0, &
-3, 0, 0, 1, &
0, 11, 5, 0, &
9, 6, 0, 0, &
-4, 10, 4, 2, &
5, 3, 0, 0, &
16, 0, 0, -9, &
-3, 0, 0, 0, &
0, 3, 2, -1, &
7, 0, 0, -3, &
-25, 22, 0, 0, &
42, 223, 119, -22, &
-27, -143, -77, 14, &
9, 49, 26, -5, &
-1166, 0, 0, 505, &
-5, 0, 0, 2, &
-6, 0, 0, 3, &
-8, 0, 1, 4, &
0, -4, 0, 0, &
117, 0, 0, -63, &
-4, 8, 4, 2, &
3, 0, 0, -2, &
-5, 0, 0, 2, &
0, 31, 0, 0, &
-5, 0, 1, 3, &
4, 0, 0, -2, &
-4, 0, 0, 2, &
-24, -13, -6, 10, &
3, 0, 0, 0, &
0, -32, -17, 0, &
8, 12, 5, -3, &
3, 0, 0, -1, &
7, 13, 0, 0, &
-3, 16, 0, 0, &
50, 0, 0, -27, &
0, -5, -3, 0, &
13, 0, 0, 0, &
0, 5, 3, 1, &
24, 5, 2, -11, &
5, -11, -5, -2, &
30, -3, -2, -16, &
18, 0, 0, -9, &
8, 614, 0, 0, &
3, -3, -1, -2, &
6, 17, 9, -3, &
-3, -9, -5, 2, &
0, 6, 3, -1, &
-127, 21, 9, 55, &
3, 5, 0, 0, &
-6, -10, -4, 3, &
5, 0, 0, 0, &
16, 9, 4, -7, &
3, 0, 0, -2, &
0, 22, 0, 0, &
0, 19, 10, 0, &
7, 0, 0, -4, &
0, -5, -2, 0, &
0, 3, 1, 0, &
-9, 3, 1, 4, &
17, 0, 0, -7, &
0, -3, -2, -1, &
-20, 34, 0, 0, &
-10, 0, 1, 5, &
-4, 0, 0, 2, &
22, -87, 0, 0, &
-4, 0, 0, 2, &
-3, -6, -2, 1, &
-16, -3, -1, 7, &
0, -3, -2, 0, &
4, 0, 0, 0, &
-68, 39, 0, 0, &
27, 0, 0, -14, &
0, -4, 0, 0, &
-25, 0, 0, 0, &
-12, -3, -2, 6, &
3, 0, 0, -1, &
3, 66, 29, -1, &
490, 0, 0, -213, &
-22, 93, 49, 12, &
-7, 28, 15, 4, &
-3, 13, 7, 2, &
-46, 14, 0, 0, &
-5, 0, 0, 0, &
2, 1, 0, 0, &
0, -3, 0, 0, &
-28, 0, 0, 15, &
5, 0, 0, -2, &
0, 3, 0, 0, &
-11, 0, 0, 5, &
0, 3, 1, 0, &
-3, 0, 0, 1, &
25, 106, 57, -13, &
5, 21, 11, -3, &
1485, 0, 0, 0, &
-7, -32, -17, 4, &
0, 5, 3, 0, &
-6, -3, -2, 3, &
30, -6, -2, -13, &
-4, 4, 0, 0, &
-19, 0, 0, 10, &
0, 4, 2, -1, &
0, 3, 0, 0, &
4, 0, 0, -2, &
0, -3, -1, 0, &
-3, 0, 0, 0, &
5, 3, 1, -2, &
0, 11, 0, 0, &
118, 0, 0, -52, &
0, -5, -3, 0, &
-28, 36, 0, 0, &
5, -5, 0, 0, &
14, -59, -31, -8, &
0, 9, 5, 1, &
-458, 0, 0, 198, &
0, -45, -20, 0, &
9, 0, 0, -5, &
0, -3, 0, 0, &
0, -4, -2, -1, &
11, 0, 0, -6, &
6, 0, 0, -2, &
-16, 23, 0, 0, &
0, -4, -2, 0, &
-5, 0, 0, 2, &
-166, 269, 0, 0, &
15, 0, 0, -8, &
10, 0, 0, -4, &
-78, 45, 0, 0, &
0, -5, -2, 0, &
7, 0, 0, -4, &
-5, 328, 0, 0, &
3, 0, 0, -2, &
5, 0, 0, -2, &
0, 3, 1, 0, &
-3, 0, 0, 0, &
-3, 0, 0, 0, &
0, -4, -2, 0, &
-1223, -26, 0, 0, &
0, 7, 3, 0, &
3, 0, 0, 0, &
0, 3, 2, 0, &
-6, 20, 0, 0, &
-368, 0, 0, 0, &
-75, 0, 0, 0, &
11, 0, 0, -6, &
3, 0, 0, -2, &
-3, 0, 0, 1, &
-13, -30, 0, 0, &
21, 3, 0, 0, &
-3, 0, 0, 1, &
-4, 0, 0, 2, &
8, -27, 0, 0, &
-19, -11, 0, 0, &
-4, 0, 0, 2, &
0, 5, 2, 0, &
-6, 0, 0, 2, &
-8, 0, 0, 0, &
-1, 0, 0, 0, &
-14, 0, 0, 6, &
6, 0, 0, 0, &
-74, 0, 0, 32, &
0, -3, -1, 0, &
4, 0, 0, -2, &
8, 11, 0, 0, &
0, 3, 2, 0, &
-262, 0, 0, 114, &
0, -4, 0, 0, &
-7, 0, 0, 4, &
0, -27, -12, 0, &
-19, -8, -4, 8, &
202, 0, 0, -87, &
-8, 35, 19, 5, &
0, 4, 2, 0, &
16, -5, 0, 0, &
5, 0, 0, -3, &
0, -3, 0, 0, &
1, 0, 0, 0, &
-35, -48, -21, 15, &
-3, -5, -2, 1, &
6, 0, 0, -3, &
3, 0, 0, -1, &
0, -5, 0, 0, &
12, 55, 29, -6, &
0, 5, 3, 0, &
-598, 0, 0, 0, &
-3, -13, -7, 1, &
-5, -7, -3, 2, &
3, 0, 0, -1, &
5, -7, 0, 0, &
4, 0, 0, -2, &
16, -6, 0, 0, &
8, -3, 0, 0, &
8, -31, -16, -4, &
0, 3, 1, 0, &
113, 0, 0, -49, &
0, -24, -10, 0, &
4, 0, 0, -2, &
27, 0, 0, 0, &
-3, 0, 0, 1, &
0, -4, -2, 0, &
5, 0, 0, -2, &
0, -3, 0, 0, &
-13, 0, 0, 6, &
5, 0, 0, -2, &
-18, -10, -4, 8, &
-4, -28, 0, 0, &
-5, 6, 3, 2, &
-3, 0, 0, 1, &
-5, -9, -4, 2, &
17, 0, 0, -7, &
11, 4, 0, 0, &
0, -6, -2, 0, &
83, 15, 0, 0, &
-4, 0, 0, 2, &
0, -114, -49, 0, &
117, 0, 0, -51, &
-5, 19, 10, 2, &
-3, 0, 0, 0, &
-3, 0, 0, 2, &
0, -3, -1, 0, &
3, 0, 0, 0, &
0, -6, -2, 0, &
393, 3, 0, 0, &
-4, 21, 11, 2, &
-6, 0, -1, 3, &
-3, 8, 4, 1, &
8, 0, 0, 0, &
18, -29, -13, -8, &
8, 34, 18, -4, &
89, 0, 0, 0, &
3, 12, 6, -1, &
54, -15, -7, -24, &
0, 3, 0, 0, &
3, 0, 0, -1, &
0, 35, 0, 0, &
-154, -30, -13, 67, &
15, 0, 0, 0, &
0, 4, 2, 0, &
0, 9, 0, 0, &
80, -71, -31, -35, &
0, -20, -9, 0, &
11, 5, 2, -5, &
61, -96, -42, -27, &
14, 9, 4, -6, &
-11, -6, -3, 5, &
0, -3, -1, 0, &
123, -415, -180, -53, &
0, 0, 0, -35, &
-5, 0, 0, 0, &
7, -32, -17, -4, &
0, -9, -5, 0, &
0, -4, 2, 0, &
-89, 0, 0, 38, &
0, -86, -19, -6, &
0, 0, -19, 6, &
-123, -416, -180, 53, &
0, -3, -1, 0, &
12, -6, -3, -5, &
-13, 9, 4, 6, &
0, -15, -7, 0, &
3, 0, 0, -1, &
-62, -97, -42, 27, &
-11, 5, 2, 5, &
0, -19, -8, 0, &
-3, 0, 0, 1, &
0, 4, 2, 0, &
0, 3, 0, 0, &
0, 4, 2, 0, &
-85, -70, -31, 37, &
163, -12, -5, -72, &
-63, -16, -7, 28, &
-21, -32, -14, 9, &
0, -3, -1, 0, &
3, 0, 0, -2, &
0, 8, 0, 0, &
3, 10, 4, -1, &
3, 0, 0, -1, &
0, -7, -3, 0, &
0, -4, -2, 0, &
6, 19, 0, 0, &
5, -173, -75, -2, &
0, -7, -3, 0, &
7, -12, -5, -3, &
-3, 0, 0, 2, &
3, -4, -2, -1, &
74, 0, 0, -32, &
-3, 12, 6, 2, &
26, -14, -6, -11, &
19, 0, 0, -8, &
6, 24, 13, -3, &
83, 0, 0, 0, &
0, -10, -5, 0, &
11, -3, -1, -5, &
3, 0, 1, -1, &
3, 0, 0, -1, &
-4, 0, 0, 0, &
5, -23, -12, -3, &
-339, 0, 0, 147, &
0, -10, -5, 0, &
5, 0, 0, 0, &
3, 0, 0, -1, &
0, -4, -2, 0, &
18, -3, 0, 0, &
9, -11, -5, -4, &
-8, 0, 0, 4, &
3, 0, 0, -1, &
0, 9, 0, 0, &
6, -9, -4, -2, &
-4, -12, 0, 0, &
67, -91, -39, -29, &
30, -18, -8, -13, &
0, 0, 0, 0, &
0, -114, -50, 0, &
0, 0, 0, 23, &
517, 16, 7, -224, &
0, -7, -3, 0, &
143, -3, -1, -62, &
29, 0, 0, -13, &
-4, 0, 0, 2, &
-6, 0, 0, 3, &
5, 12, 5, -2, &
-25, 0, 0, 11, &
-3, 0, 0, 1, &
0, 4, 2, 0, &
-22, 12, 5, 10, &
50, 0, 0, -22, &
0, 7, 4, 0, &
0, 3, 1, 0, &
-4, 4, 2, 2, &
-5, -11, -5, 2, &
0, 4, 2, 0, &
4, 17, 9, -2, &
59, 0, 0, 0, &
0, -4, -2, 0, &
-8, 0, 0, 4, &
-3, 0, 0, 0, &
4, -15, -8, -2, &
370, -8, 0, -160, &
0, 0, -3, 0, &
0, 3, 1, 0, &
-6, 3, 1, 3, &
0, 6, 0, 0, &
-10, 0, 0, 4, &
0, 9, 4, 0, &
4, 17, 7, -2, &
34, 0, 0, -15, &
0, 5, 3, 0, &
-5, 0, 0, 2, &
-37, -7, -3, 16, &
3, 13, 7, -2, &
40, 0, 0, 0, &
0, -3, -2, 0, &
-184, -3, -1, 80, &
-3, 0, 0, 1, &
-3, 0, 0, 0, &
0, -10, -6, -1, &
31, -6, 0, -13, &
-3, -32, -14, 1, &
-7, 0, 0, 3, &
0, -8, -4, 0, &
3, -4, 0, 0, &
0, 4, 0, 0, &
0, 3, 1, 0, &
19, -23, -10, 2, &
0, 0, 0, -10, &
0, 3, 2, 0, &
0, 9, 5, -1, &
28, 0, 0, 0, &
0, -7, -4, 0, &
8, -4, 0, -4, &
0, 0, -2, 0, &
0, 3, 0, 0, &
-3, 0, 0, 1, &
-9, 0, 1, 4, &
3, 12, 5, -1, &
17, -3, -1, 0, &
0, 7, 4, 0, &
19, 0, 0, 0, &
0, -5, -3, 0, &
14, -3, 0, -1, &
0, 0, -1, 0, &
0, 0, 0, -5, &
0, 5, 3, 0, &
13, 0, 0, 0, &
0, -3, -2, 0, &
2, 9, 4, 3, &
0, 0, 0, -4, &
8, 0, 0, 0, &
0, 4, 2, 0, &
6, 0, 0, -3, &
6, 0, 0, 0, &
0, 3, 1, 0, &
5, 0, 0, -2, &
3, 0, 0, -1, &
-3, 0, 0, 0, &
6, 0, 0, 0, &
7, 0, 0, 0, &
-4, 0, 0, 0, &
4, 0, 0, 0, &
6, 0, 0, 0, &
0, -4, 0, 0, &
0, -4, 0, 0, &
5, 0, 0, 0, &
-3, 0, 0, 0, &
4, 0, 0, 0, &
-5, 0, 0, 0, &
4, 0, 0, 0, &
0, 3, 0, 0, &
13, 0, 0, 0, &
21, 11, 0, 0, &
0, -5, 0, 0, &
0, -5, -2, 0, &
0, 5, 3, 0, &
0, -5, 0, 0, &
-3, 0, 0, 2, &
20, 10, 0, 0, &
-34, 0, 0, 0, &
-19, 0, 0, 0, &
3, 0, 0, -2, &
-3, 0, 0, 1, &
-6, 0, 0, 3, &
-4, 0, 0, 0, &
3, 0, 0, 0, &
3, 0, 0, 0, &
4, 0, 0, 0, &
3, 0, 0, -1, &
6, 0, 0, -3, &
-8, 0, 0, 3, &
0, 3, 1, 0, &
-3, 0, 0, 0, &
0, -3, -2, 0, &
126, -63, -27, -55, &
-5, 0, 1, 2, &
-3, 28, 15, 2, &
5, 0, 1, -2, &
0, 9, 4, 1, &
0, 9, 4, -1, &
-126, -63, -27, 55, &
3, 0, 0, -1, &
21, -11, -6, -11, &
0, -4, 0, 0, &
-21, -11, -6, 11, &
-3, 0, 0, 1, &
0, 3, 1, 0, &
8, 0, 0, -4, &
-6, 0, 0, 3, &
-3, 0, 0, 1, &
3, 0, 0, -1, &
-3, 0, 0, 1, &
-5, 0, 0, 2, &
24, -12, -5, -11, &
0, 3, 1, 0, &
0, 3, 1, 0, &
0, 3, 2, 0, &
-24, -12, -5, 10, &
4, 0, -1, -2, &
13, 0, 0, -6, &
7, 0, 0, -3, &
3, 0, 0, -1, &
3, 0, 0, -1 ], [4,npl])
! Interval between fundamental date J2000.0 and given date (JC).
t = ( ( date1-dj00 ) + date2 ) / djc
! -------------------
! LUNI-SOLAR NUTATION
! -------------------
!
! Fundamental (Delaunay) arguments
!
! Mean anomaly of the Moon (IERS 2003).
el = FAL03 ( t )
! Mean anomaly of the Sun (MHB2000).
elp = mod ( 1287104.79305_wp + &
t*( 129596581.0481_wp + &
t*( - 0.5532_wp + &
t*( 0.000136_wp + &
t*( - 0.00001149_wp )))), turnas ) * das2r
! Mean longitude of the Moon minus that of the ascending node
! (IERS 2003.
f = FAF03 ( t )
! Mean elongation of the Moon from the Sun (MHB2000).
d = mod ( 1072260.70369_wp + &
t*( 1602961601.2090_wp + &
t*( - 6.3706_wp + &
t*( 0.006593_wp + &
t*( - 0.00003169_wp )))), turnas ) * das2r
! Mean longitude of the ascending node of the Moon (IERS 2003).
om = FAOM03 ( t )
! Initialize the nutation values.
dp = 0.0_wp
de = 0.0_wp
! Summation of luni-solar nutation series (in reverse order).
do i = nls, 1, -1
! Argument and functions.
arg = mod ( real ( nals(1,i), wp ) * el + &
real ( nals(2,i), wp ) * elp + &
real ( nals(3,i), wp ) * f + &
real ( nals(4,i), wp ) * d + &
real ( nals(5,i), wp ) * om, d2pi )
sarg = sin(arg)
carg = cos(arg)
! Term.
dp = dp + ( cls(1,i) + cls(2,i) * t ) * sarg &
+ cls(3,i) * carg
de = de + ( cls(4,i) + cls(5,i) * t ) * carg &
+ cls(6,i) * sarg
end do
! Convert from 0.1 microarcsec units to radians.
dpsils = dp * u2r
depsls = de * u2r
! ------------------
! PLANETARY NUTATION
! ------------------
! n.b. The MHB2000 code computes the luni-solar and planetary nutation
! in different routines, using slightly different Delaunay
! arguments in the two cases. This behaviour is faithfully
! reproduced here. Use of the IERS 2003 expressions for both
! cases leads to negligible changes, well below
! 0.1 microarcsecond.
! Mean anomaly of the Moon (MHB2000).
al = mod ( 2.35555598_wp + 8328.6914269554_wp * t, d2pi )
! Mean anomaly of the Sun (MHB2000).
alsu = mod ( 6.24006013_wp + 628.301955_wp * t, d2pi )
! Mean longitude of the Moon minus that of the ascending node
! (MHB2000).
af = mod ( 1.627905234_wp + 8433.466158131_wp * t, d2pi )
! Mean elongation of the Moon from the Sun (MHB2000).
ad = mod ( 5.198466741_wp + 7771.3771468121_wp * t, d2pi )
! Mean longitude of the ascending node of the Moon (MHB2000).
aom = mod ( 2.18243920_wp - 33.757045_wp * t, d2pi )
! General accumulated precession in longitude (IERS 2003).
apa = FAPA03 ( t )
! Planetary longitudes, Mercury through Uranus (IERS 2003).
alme = FAME03 ( t )
alve = FAVE03 ( t )
alea = FAE03 ( t )
alma = FAMA03 ( t )
alju = FAJU03 ( t )
alsa = FASA03 ( t )
alur = FAUR03 ( t )
! Neptune longitude (MHB2000).
alne = mod ( 5.321159000_wp + 3.8127774000_wp * t, d2pi )
! Initialize the nutation values.
dp = 0.0_wp
de = 0.0_wp
! Summation of planetary nutation series (in reverse order).
do i = npl, 1, -1
! Argument and functions.
arg = mod ( real ( napl( 1,i), wp ) * al + &
real ( napl( 2,i), wp ) * alsu + &
real ( napl( 3,i), wp ) * af + &
real ( napl( 4,i), wp ) * ad + &
real ( napl( 5,i), wp ) * aom + &
real ( napl( 6,i), wp ) * alme + &
real ( napl( 7,i), wp ) * alve + &
real ( napl( 8,i), wp ) * alea + &
real ( napl( 9,i), wp ) * alma + &
real ( napl(10,i), wp ) * alju + &
real ( napl(11,i), wp ) * alsa + &
real ( napl(12,i), wp ) * alur + &
real ( napl(13,i), wp ) * alne + &
real ( napl(14,i), wp ) * apa, d2pi )
sarg = sin(arg)
carg = cos(arg)
! Term.
dp = dp + real(icpl(1,i), wp) * sarg + real(icpl(2,i), wp) * carg
de = de + real(icpl(3,i), wp) * sarg + real(icpl(4,i), wp) * carg
end do
! Convert from 0.1 microarcsec units to radians.
dpsipl = dp * u2r
depspl = de * u2r
! -------
! RESULTS
! -------
! Add luni-solar and planetary components.
dpsi = dpsils + dpsipl
deps = depsls + depspl
end subroutine NUT00A