Overview
-This module deals with Cosmology in general relativity, Earth gravity,
geodesic distances and different time scales ( UT , TDB , TCB ).
-Ellipsoidal models of the Earth are used, spheroidal and scalene ellipsoids.
-A Luni-Solar calendar may also be used... just for fun !
XROM | Function | Desciption |
16,00
16,01 16,02 16,03 16,04 16,05 16,06 16,07 16,08 16,09 16,10 16,11 16,12 16,13 16,14 16,15 16,16 16,17 16,18 16,19 16,20 16,21 16,22 16,23 16,24 16,25 16,26 16,27 16,28 16,29 16,30 16,31 16,32 16,33 16,34 16,35 16,36 16,37 16,38 16,39 16,40 16,41 16,42 16,43 16,44 16,45 16,46 16,47 16,48 16,49 16,50 16,51 16,52 16,53 16,54 16,55 16,56 16,57 16,58 16,59 16,60 16,61 16,62 16,63 |
-GRVITME 1B
ABC AF "AUM" "AUR" "BCT" "BDT" "COSMO" "COSMO+" "CSM" "DT" "DTTM" "ETP" "EUCL" "FWD" "GAB" "GAB2" "GEO" "GEO+" "GR-LS" "GRV" "GRV+" "GRV3" "GRV3+" "GRV4" "GRV4+" HUBBLE "J0" J2 "J2TM" "LS-GR" "SPOL" "SPOL+" "TGD" "TGD2" "TGD3" "TOL" "TT-UT" -AUX FNS 1/Z CDAY CBRT CMD CROSS ELLOX+ DOT FMD GW GX JDAY NORM "P2" "P3" "P4" R-S RFZ RJZ S-R ST<>A HMT LEG #GR X/E3 X^3 |
Section Header
3 semi-axes of Triaxial Ellipsoid model of the Earth Equatorial radius & 1 / flattening of the Earth Age of the Universe ( Radiation + Dust ) Age of a Universe ( Radiation without Dust ) Barycentric Coordinate Time Barycentric Dynamical Time Redshift -> Several Distances Idem + prompts for inputs Redshift -> Several Distances in an Empty Universe Number of Days since 2000/01/01 -> Date ( YMD ) Idem with TIME module Einstein's Twin Paradox Redshift -> Distance in an Euclidean Universe Geodesic Distance: Forward Problem An Example of Subroutine to get the Speed of Light 2nd Example for the Speed of Light Geocentric Coordinates Idem + Prompts Example of a Luni-Solar Calendar Gravity on a Scalene Ellipsoid Idem + Prompts Gravity on a Spheroidal Earth Idem + Prompts Gravity from the Geopotential ( only a few terms ) Stores the constants used by "GRV4" Hubble's Constant Number of Days since 2000/01/01 + Day of the week Idem for Gregorian dates ( without DOW ) Idem with TIME Module Example of a Luni-Solar Calendar Speed of Light in General Relativity / coord. time Idem + Prompt Terrestrial Geodesic Distance ( Andoyer's 2nd order ) Terrestrial Geodesic Distance ( Vincenty's formula ) Geodesic Distance - Triaxial Ellipsoidal Earth Redshift -> Distance in Tolman Universes Terrestrial Time minus Universal Time Section Header Inverse of a complex number Julian Day -> Date Cube Root Ceil ( Y Z / X ) ( M-code routine used by LuniSol cal. ) Cross Product ( X Y Z ) x ( M N O ) Loxodromics on an Ellipoidal Earth Dot Product ( X Y Z ) . ( M N O ) Floor ( Y Z / X ) ( M-code routine used by LuniSol cal. ) M-Code for Gauss-Legendre 6-point formula Idem Date -> Julian Day Norm of a 3D-vector ( X Y Z ) Quadratic Equation a.x^2+b.x+c = 0 Cubic Equation a.x^3+b.x^2+c.x+d = 0 Quartic Equationx x^4+a.x^3+b.x^2+c.x+d = 0 Rectangular-Spherical Conversion Carlson Elliptic Integral 1st kind RF(x,y+iz,y-iz) Carlson Elliptic Integral 3rd kind RJ(x,y+iz,y-iz,t>0) Spherical-Rectangular Conversion Exchanges X-Y-Z <> M-N-O Hermite Polynomial Legendre Polynomial An M-Code Subroutine called by "GRV3" X divided by 1000 Cube |
Three Semi-Axes of a Scalene Ellipsoidal Model of the
Earth
ABC places
a = 6378.171274 km in X-register
b = 6378.101946 km in Y-register
c = 6356.751868 km in Z-register
-The previous contents of registers X Y Z are overwritten.
-The last decimals are somewhat meaningless !
Equatorial Radius of the Earth & Inverse of the
Earth Flattening
AF returns
a = 6378.13661 km in X-register
1/f = 298.256421 in Y-register
-The previous contents of X & Y are saved in Z & T
-These values seem the most recent estimations.
Gravity on a Scalene Ellipsoid
-"GRV" employs the formulae given in reference [6] and the semi-axes
above.
STACK | INPUTS | OUTPUTS |
Z | L ( ° ' " ) | g(0) |
Y | B ( ° ' " ) | g(0) |
X | h ( m ) | g(h) |
Example1: L = 77°03'56"
W , B = 38°55'17"2 N , h = 67m
( U.S. Naval Observatory , Washington D.C. )
-77.0356 ENTER^
38.55172 ENTER^
67
XEQ "GRV" >>>> g(67) = 9.800516275 m/s2
RDN g(0) = 9.800723034 m/s2
Example2: L = 116°51'50"4 W , B = 33°21'22"4 N h = 1706 m ( Mount Palomar Observatory )
-116.51504 ENTER^
33.21224
ENTER^
1706
R/S >>>> g(1709) = 9.790660012
m/s2
RDN g(0) = 9.795923287
m/s2
Notes:
-The longitudes are positive East.
-Latitudes & longitudes are both supposed geodetic, but using a
geocentric longitude will not change the result significantly.
-The results are very accurate provided the altitude remains relatively
small
-Unfortunately, the Earth gravitational field is much more complex
than that of a triaxial ellipsoid,
so this program only gives slightly better results than the classical
Somigliana's formula.
-XEQ "GRV+" to get PROMPTs for the inputs.
-The gravity g above an ellipsoid of revolution may be computed by rigorous
formulae ( cf reference [7] ):
STACK | INPUTS | OUTPUTS |
Z | / | B ( deg ) |
Y | b ( ° ' " ) | R ( m ) |
X | h ( m ) | g ( m/s2 ) |
where R = distance from the center of the Earth
and B = geocentric
latitude
Example1: b = 38°55'17"2 h = 23456 m
38.55172 ENTER^
23456 XEQ "GRV3"
>>>> g = 9.728751603 m/s2
RDN R = 6393194.708 m
RDN B = 38°73415862
Example2: b = 38°55'17"2 h = 12345678 m
38.55172 ENTER^
12345678 XEQ "GRV3" >>>> g = 1.078713339
m/s2
RDN R = 18715394.22 m
RDN B = 38°85746754
Note:
XEQ "GRV3+" to get PROMPTs for the inputs
Gravity from the Geopotential ( a few terms )
-"GRV4" uses the Earth's gravitational potential U to compute the gravity
( spherical harmonics )
-Thousands - or at least hundreds - of terms are needed to get full
accuracy.
-This routine employs only a few terms of the geopotential which must
be stored in R19 thru R54 by "GRV4+" before executing "GRV4"
STACK | INPUTS | OUTPUTS |
Z | B (deg) | / |
Y | L (deg) | / |
X | R (m) | g ( m/s2) |
where R , L , B are the geocentric coordinates.
R = distance from the center of the Earth
L = longitude,
B = geocentric latitude ( not
geographic latitude )
-"GEO" may be used to calculate R & B from the altitude and the geographic latitude.
Example: R = 6369806 m L = -77.065556° B = 38.733471°
XEQ "GRV4+" to get the constants in the proper registers, then the HP-41 displays
38.733471 ENTER^
-77.065556 ENTER^
6369806 R/S
>>>> g = 9.800330901 m/s2
---Execution time = 79s---
gR = R01
= -9.800278379 m/s2
gL = R02
= 0.000091758 m/s2
gB = R03
= -0.032085221 m/s2
Note:
-The constants in R19 to R54 are unchanged by this program,
so you can XEQ "GRV4" instead of "GRV4+" to get another g-value
at another location
HUBBLE returns in X-register Hubble's constant: 1 / H0 = 1.377172143 E10 years
-Most of the decimals are of course doubtful
-They have been chosen to give H0 = 71 km / s / Mpc
-The relation betxeen H0 = x km/s/Mpc &
1 / H0 = y 109 years is x
y = 977.79222143
Age of the Universe ( Dust + Radiation )
-There are many subcases and the following program only works if the cosmological constant is positive
and if the quartic: (Omega)lambda.y4 + ( 1-(Omega)tot ).y2 + (Omega)mat.y + (Omega)rad
has 2 distinct non-positive roots and a pair of complex conjugate roots.
-These conditions are precisely satisfied by our Universe !
STACK | INPUTS | OUTPUTS |
Z | ¶0 | k |
Y | L0 | R0 |
X | (Omega)mat | t0 |
t0 = age of the Universe ( in years )
Where R0 = Radius of the
Universe ( in light-years )
k = +1 , 0 , -1 for Spherical , Euclidean , Hyperbolic spaces
Example1: With (Omega)mat = 0.044 , (Omega)lambda = 0.519 , (Omega)rad = 0.000049
0.000049 ENTER^
0.519 ENTER^
0.044 XEQ "AUM"
>>>> t0 = 1.5587728 E10 years
R0 = 2.0833961 E10 light-years
k = -1 ( hyperbolic space )
Example2: With (Omega)mat = 0.271 , (Omega)lambda = 0.732 , (Omega)rad = 0.000049
( These Omega-values are close to those suggested by WMAP = NASA's Wilkinson Microwave Anisotropy Probe )
0.000049 ENTER^
0.732
ENTER^
0.271
XEQ "AUM" >>>> t0 = 1.3667886 E10 years
R0 = 2.4940750 E11 light-years
k = +1 ( spherical space )
Example3: With (Omega)mat = 0.2679 , (Omega)lambda = 0.732 , (Omega)rad = 0.0001
.0001 ENTER^
.732 ENTER^
.2679 XEQ "AUM" >>>> t0
= 1.3693298 E10 years
R0 = 1.3771721 E10 light-years
R0 = c/H0 when the Universe is
Euclidean ( conventionally )
k = 0 ( Euclidean space )
Note:
-"AUM" calls "P4" ( quartic equations ) and the Carlson elliptic integrals
are computed with the M-Code routines RFZ & RJZ
Age of the Universe ( Radiation without Dust )
-These kinds of Universe are also called Tolman Universes.
-The calculations are much easier and do not require elliptic integrals
-Here, (Omega)mat = 0
STACK | INPUTS | OUTPUTS |
Z | / | k |
Y | ¶0 | R0 |
X | L0 | t0 |
t0 = age of the Universe ( in years )
Where R0 = Radius of the
Universe ( in light-years )
k = +1 , 0 , -1 for Spherical , Euclidean , Hyperbolic spaces
Example: ¶0 = (Omega)rad = 49 10 -6 , L0 = (Omega)lambda = 0.61
49 E-6 ENTER^
0.61
XEQ "AUR" >>>> t0 = 1.8236300 E10
years
RDN R0 = 2.2053788 E10 light-years
RDN k = -1 ( hyperbolic Universe )
Redshift >>> Distance of a Galaxy
-This problem involves the same kind of integrals as "AUM"
-So, quartic equation + Carlson elliptic integrals
-The same restrictions apply
Data Registers: R00 thru R32
R10 = z + 1
R13 = (Omega)rad
R11 = (Omega)mat
R14 = (Omega)tot - 1
R12 = (Omega)lambda
R15 thru R32: temp
-when the program stops:
R00 = k ( +1 = Spherical
, 0 = Euclidean , -1 = Hyperbolic )
R01 = D , R02 = D0 , R03 = DL , R04 = VR
, R05 = t0 , R06 = R0 ( current scale factor )
Flags: F00 set flag F00
if you don't want to re-calculate the age of the Universe t0
F01 & F02 are used by "P4"
STACK | INPUTS | OUTPUTS |
T | ¶0 | VR |
Z | L0 | DL |
Y | (Omega)mat | D0 |
X | z | D |
L | / | t0 |
D
= light-travel time distance ( in light-years ) = light-travel time
( in years )
D0 = comoving
distance ( in light-years )
DL = luminosity-distance
( in light-years )
VR = recessional
velocity ( speed of light = 1 )
t0 =
age of the Universe ( in years )
Example: With (Omega)mat = 0.044 , (Omega)lambda = 0.519 , (Omega)rad = 0.000049 & z = 7
0.000049 ENTER^
0.519
ENTER^
0.044
ENTER^
7 XEQ "COSMO"
>>>> D = 1.4119088 E10 l-y
= R01
RDN D0 = 3.4190016 E10 l-y
= R02
RDN DL = 4.1392272 E11 l-y
= R03
RDN VR = 2.482624772
= R04 ( unit c = 1 )
LASTX t0 = 1.5587728 E10 years
= R05
and R00 = k = -1 ( hyperbolic space )
R06 =
R0 = 2.0833961 E10 light-years
Notes:
-XEQ "COSMO+" to get PROMPTs for the 4 inputs
-Execution time is about 1 minute.
Redshift >>> Distance of a Galaxy in an Empty Universe
-The problem is of course much simpler in an empty universe.
-"CSM" only works with a cosmological constant such that 0 <
(Omega)lambda < 1
STACK | INPUTS | OUTPUTS |
T | / | VR |
Z | / | DL |
Y | L0 | D0 |
X | z | D |
L | / | t0 |
D
= light-travel time distance ( in light-years ) = light-travel time
( in years )
D0 = comoving
distance ( in light-years )
DL = luminosity-distance
( in light-years )
VR = recessional
velocity ( speed of light = 1 )
t0 =
age of the Universe ( in years )
Example: (Omega)lambda = 0.519 & z = 7
0.519
ENTER^
7 XEQ "CSM"
>>>> D = 1.4892172 E10 l-y
RDN D0 = 3.7410987 E10 l-y
RDN DL = 5.1055493 E11 l-y
RDN VR = 2.716507697
( unit c = 1 )
LASTX t0 = 1.7367386 E10 years
Note:
k = -1 for these models
Redshift >>> Distance of a Galaxy in a Euclidean Universe
-Given the cosmological parameter and the redshift z
"EUCL" returns the light-travel time distance D and the age
of the Universe.
STACK | INPUTS | OUTPUTS |
Y | L0 | t0 |
X | z | D |
D = light-travel
time distance ( in light-years ) = light-travel time ( in years )
t0 = age of the Universe
( in years )
Example: L0 = 0.732 & z = 7
0.732 ENTER^
7
XEQ "EUCL" >>>> D = 1.2915891 E10 ly
X<>Y t0 = 1.3698976 E10 years
Redshift >>> Distance of a Galaxy in Tolman Universes
-Like with"AUR", (Omega)mat = 0
-Such models have been studied by Tolman.
-There are many cases, and "TOL" only works for Big Bang Universes with
a positive Cosmological Constant.
-The comoving radial distance involves elliptic integrals and "TOL"
calculates the light-travel time Distance D,
the age of the Universe t0 , the current radius of
the Universe R0 and the constant k
( D & t0 may be expressed in terms of elementary
functions )
STACK | INPUTS | OUTPUTS |
T | / | k |
Z | ¶0 | R0 |
Y | L0 | t0 |
X | z | D |
Example: ¶0 = (Omega)rad = 49 10 -6 , L0 = (Omega)lambda = 0.61 , z = 7
49 E-6 ENTER^
0.61
ENTER^
7
XEQ "TOL" >>>> D = 1.5726496 1010
light-years
RDN t0 = 1.8236300 1010
years
RDN R0 = 2.2053788 1010 light-years
RDN k = -1 ( hyperbolic Universe )
-This program calculates TDB-TT. It employs the first terms of the series
given in reference [8],
those whose amplitude > 0.4 second over the time-span
[1000,3000]
-The errors should not exceed 1 or 2 microseconds over this interval
of time.
STACK | INPUTS | OUTPUTS |
Y | YYYY.MNDD | / |
X | HH.MNSS | TDB-TT (s) |
Example1: 2012/08/29
16h41m37s ( TT or TDB )
2012.0829 ENTER^
16.4137
XEQ "BDT" >>>> TDB-TT = - 0.0013227 s
Example2: 3000/04/04 1h23m45s ( TT or TDB )
3000.0404 ENTER^
1.2345
XEQ "BDT" ( or RTN R/S ) >>>>
TDB-TT = + 0.0015386 s
-Anothe time scale may also be used: TCB = Temps Coordonné
Barycentrique = Barycentric Coordinate Time.
-Whereas TDB - TT remains smaller than 0.002 second
at least several millenia aroud J2000, the difference TCB - TDB
is much larger:
-It was almost 0 in 1977 but in 2012 , TCB -TDB is about 17 seconds.
Formula:
TCB = TDB + 1.550519768 E-8 ( JDTCB - 2443144.5003725 ) 86400 + 6.55 E-5 seconds
where JDTCB is the Julian Date in the TCB scale.
-"BCT" uses an iteration to compute TCB - TDB for a given
date expressed in TDB
STACK | INPUTS | OUTPUTS |
Y | YYYY.MNDD | / |
X | HH.MNSS | TCB-TDB (s) |
Where the input time is expressed in TDB
Example: 3000/04/04 1h23m45s TDB
3000.0404 ENTER^
1.2345
XEQ "TCB" >>>> 500.6752393 seconds
Terrestrial Time minus Universal Time
-This routine uses the polynomials given in reference [9]
STACK | INPUTS | OUTPUTS |
X | year | TT-UT (sec) |
Examples:
-2000 XEQ "TT-UT" >>>> 46675.68
s
400 XEQ "TT-UT" >>>>
6699.22 s
1200 XEQ "TT-UT" >>>> 736.44
s
1680 XEQ "TT-UT" >>>>
15.31 s
1760 XEQ "TT-UT" >>>>
14.87 s
1841 XEQ "TT-UT" >>>>
5.52 s
1880 XEQ "TT-UT" >>>> -5.01
s
1906 XEQ "TT-UT" >>>>
5.10 s
1934 XEQ "TT-UT" >>>>
23.86 s
1951 XEQ "TT-UT" >>>>
29.47 s
1984 XEQ "TT-UT" >>>>
53.73 s
2000 XEQ "TT-UT" >>>>
63.86 s
2041 XEQ "TT-UT" >>>>
85.52 s
2100 XEQ "TT-UT" >>>> 202.74
s
3000 XEQ "TT-UT" >>>> 4435.68
s
Note:
-In most cases, the program does not stop at the END
-So, press XEQ "TT-UT" or RTN R/S for another date
Date <> Number of Days since 2000/01/01
-"J0" computes the number of days since 2000/01/01 0h ( i-e the
Julian Day number minus 2451544.5 ) in X-register and the day of the week
in T-register.
-Registers Y and Z are saved.
-The Gregorian calendar reform is taken into account: The day following
1582/10/04 is 1582/10/15
-"J0" ( and "DT" = the reverse operation ) are valid at least since
4713 B.C. up to ... the next calendar reform!
STACK | INPUTS | OUTPUTS |
T | / | dow |
Z | Z | Z |
Y | Y | Y |
X | YYYY.MNDDdd | N |
L | / | YYYY.MNDDdd |
( dow = day of week = 0 for Sunday, 1 for Monday, ... , 6 for Saturday )
-N= the number of days between a date YYYY.MMDDdd and 2000/01/01 at 0h = the Julian Day number minus 2451544.5
Examples:
April 4th 2134 at 6h AM 2134.040425 XEQ "J0" >>> N = 49036.25 RDN RDN dow = 0 = Sunday
1234.0428 XEQ "J0"
>>> N = -279651 RDN RDN dow = 5 = Friday
-4123.0707 XEQ "J0" >>>
N = -2236225 RDN RDN dow = 1 = Monday
-"DT" performs the reverse transformation
STACK | INPUTS | OUTPUTS |
Z | Z | Z |
Y | Y | Y |
X | N | YYYY.MNDDdd |
L | / | N |
Examples:
N = 49036.25 XEQ "DT"
>>> 2134.040425
N = -279651
XEQ "DT" >>> 1234.0428
N = -2236225 XEQ "DT"
>>> -4123.0707
Notes:
J2 is an M-Code routine that takes a Gregorian date ( even before
1582 ) YYYY.MNDD and returns the number of days since 2000/01/01
J2TM does the same thing but with the TIME module functions
Like "DT", DTTM takes the number of days since 2000/01/01
and returns the Gregorian date YYYY.MNDD
-"ETP" uses Einstein's special relativity to deal with uniformly accelerated
linear
motion of a spacecraft.
-In special relativity, the speed is also v = dx/dt , but acceleration
is a = ( dv/dt ) ( 1 - v2/c2 )-3/2
where t = time passed on Earth
-The proper Time T ( passed on board ship ) verifies dT = dt
( 1 - v2/c2 )1/2
-We use the following units:
-Times ( t , T ) are expressed in Julian years (
365.25 days )
-Distance D is expressed in light-years (
one light-year = 9.460,730,472,580,8 1015 m. )
-The unit of accelerations ( a ) is the "standard"
gravity on Earth g = 9.80665 m/s2
-The velocity of light = c = 299,792,458 m/s
( from the definition of 1 meter )
*********************************************************************************************
Return Journey
<<<<<Deceleration "a" <<<<< |
<<<<<Acceleration "-a" <<<<<<<
Earth=O-------------------------------D/2------------------------------------D=(Star)--------------->
x
>>>>>Acceleration "a" >>>>> | >>>>>Deceleration "-a" >>>>>>>
Outward Journey
-The graphic representation of x(t) looks like this:
x
|
| D
*
|
*
*
|
*
*
| D/2
*
*
|
*
*
|
*
*
|*---------------|----------------|-----------------|-----------------*-----
t
Takeoff
Spacecraft reaches the star
Landing
*************************************************************************************************
STACK | INPUTS | OUTPUTS |
T | / | 1-vmax/c |
Z | / | vmax/c |
Y | / | tE |
X | / | TS |
TS = time passed on
board ship
tE = time passed
on Earth.
vmax = maximum speed
( when x = D/2 )
Example1: On 2007/12/25 immediate
boarding for the first interstellar flight in the direction of Toliman
( alpha Centauri ) Distance D = 4.343 light-years.
Fox Mulder remains on Earth but Dana Scully's flying saucer takes off at
10h13m TT , acceleration = 0.7g
XEQ "ETP" >>>>
4.343 ENTER^
0.7 R/S
>>>> TS =
8.837390060 years
RDN tE
= 13.09998313 years
RDN vmax/c = 0.9211383858
vmax = 92.11383858% * 299792458 m/s = 276150341
m/s
RDN 1-vmax/c = 0.07886161418
-When she steps from the ship, Dana's watch shows 2016/10/26 , 6h46m41s
wheras according to Mulder, it's 4h40m08s o'clock and we are
on 2021/01/30
Example2: With a = 1g and D = 2,200,000 light-years ( approximately the distance of the Andromeda galaxy ) "ETP" returns:
TS
= 56.71150062 years
tE
= 4,400,003.874 years
vmax/c =
1
1-vmax/c = 3.877715920 E-13 whence
vmax/c was actually 0.999,999,999,999,612,228,408
Note:
-Caculating 1-vmax/c is only useful if the maximum
v/c is very close to 1
Speed of Light in General Relativity / coordinate time
-If the clocks are synchronized along the path of the light and if we
use the proper time, the speed of light is constant.
-In this case, no program is needed to calculate it: c = 299792458
m/s
-The following routine employs the coordinate-time t to measure the
velocity of light V = dL/dt
-We may have V < c or V > c or V = c
-The propagation of light is characterized by null geodesics:
ds2 = gab dxa dxb = 0 a , b = 0 , 1 , 2 , 3 where x0 = c.t & x1 , x2 , x3 are spatial coordinates
-This may be re-written
ds2 = g00 dx0 dx0 + 2 g0i
dx0 dxi + gij dxi dxj
i , j = 1 , 2 , 3
= [ sqrt(g00) dx0 + g0i dxi
/ sqrt(g00) ]2 - ( g0i g0j
/ g00 ) dxi dxj + gij
dxi dxj
= g00 ( dx0 )2 [ 1 + ( g0i
/ g00) ( dxi / dx0 ) ]2
- dL2
where dL2 = hij dxi dxj is the spatial metric with hij = - gij + ( g0i g0j / g00 )
[ The spatial tensor is usually denoted gij ( the geek letter gamma ), but I've prefered to use hij in case your browser ignores greek symbols ]
-The coordinates of the velocity of light are Vi = dxi / dt and its modulus is V = dL/dt = sqrt ( hij Vi Vj )
whence ds2 / dt2 = 0 =
g00 c2 [ 1 + ( g0i / g00) (
1/c ) ( dxi / dt ) ]2 - V2
= g00 c2 [ 1 + ( g0i / g00)
Vi / c ]2 - V2
-The direction of propagation is defined by a 3-vector ki = dxi / dL = ( dxi / dt ) ( dt / dL ) = Vi / V ( hij ki kj = 1 )
so, V = c Sqrt(g00) + ( g0i / sqrt(g00) ) ki V
-Finally, it yields: V/c
= sqrt (g00) / [ 1 - ki g0i / sqrt(g00)
]
Data Registers: R00 = function name ( Registers R00 thru R07 are to be initialized before executing "SPOL" )
R01 = x1
R05 = k1
R08 = g00 R12 = g01
R15 = g12 R17 = g23
R02 = x2
R06 = k2
R09 = g11 R13 = g02
R16 = g13
R03 = x3
R07 = k3
R10 = g22 R14 = g03
R04 = x0 = c.t
R11 = g33
R09 = h11 R13 = h02
R16 = h13
R10 = h22 R14 = h03
R11 = h33
Flags: /
Subroutine: A program that takes
x1 , x2 , x3 ,
x0 in registers R01-R02-R03-R04 and calculates
and stores the 10 components gab
into R08 to R17 as shown above.
STACK | INPUT | OUTPUT |
X | / | v/c |
where v is the speed of light measured with the time coordinate
Example1: ds2 = [ 1 - 1 / (4r) + r2 / 1000 ] c2 dt2 - ( dx2 + dy2 + dz2 ) / [ 1 - 1/(4r) + r2/1000 ] - ( 2 / r3 ) ( y.dx - x.dy ) c dt
with r = sqrt ( x2 + y2 + z2 )
Evaluate the speed of light at a point P(1,2,3,0)
a) in the
direction defined by the vector k(4,5,6)
b) in
the direction defined by the vector k(0,1,0)
-We have g00 = 1 - 1/(4r) + r2/1000 = -
1 / g11 = - 1 / g22 = - 1 / g33
, g01 = - y / r3 , g02
= x / r3 , the other components of the metric tensor =
0
>>> The subroutine "GAB" is already in the module and stores
the components of this metric tensor in the proper registers
"GAB" STO 00 ( use XRQ "SPOL+" and you will get PROMPTs for all the inputs )
1
STO 01
2
STO 02
3
STO 03
0
STO 04
a) In the 1st direction
4
STO 05
5
STO 06
6
STO 07
XEQ "SPOL" >>>> V/c = 0.966923596
b) In the 2nd direction
0
STO 05
1
STO 06
0
STO 07
XEQ "SPOL" or XEQ 10 ( faster )
>>>> V/c = 0.992171327
Evaluate the speed of light at a point P(7,8,9,0) in the direction defined by the vector k(2,3,4)
7 STO 01
2 STO 05
8 STO 02
3 STO 06
9 STO 03
4 STO 07
0 STO 04
XEQ "SPOL" >>>> V/c = 1.084831634
Notes:
-XEQ "SPOL+" to get PROMPTs for the inputs
-The components of k may be multiplied by a positive constant without changing the results.
-Line 48 LBL 10 is only useful to get faster results at the same point.
-The metric tensor of this example is not very realistic though it could
represent approximately the gravitational field of a rotating sphere
in a universe with a ( very ) large negative cosmological constant.
Example2: ds2 = [ 1 - r2 / 10000 ] c2 dt2 - ( dx2 + dy2 + dz2 ) - ( 0.02 ) ( x.dy - y.dx ) c dt
with r = sqrt ( x2 + y2 + z2 )
Evaluate the speed of light at the point P(1,0,0,0) in the directions k(0,1,0) and k(0,-1,0)
-We have g00 = 1 - r2/10000
, g11 = g22 = g33 = -1
, g01 = 0.01 y , g02 = -0.01
x , the other components of the metric tensor = 0
>>> The subroutine "GAB2" is already in the module and stores
the components of this metric tensor in the proper registers
"GAB2" STO 00 ( use XRQ "SPOL+" and you will get PROMPTs for all the inputs )
1
STO 01 0 STO 05
0
STO 02 1 STO 06
0
STO 03 0 STO 07
0
STO 04
XEQ "SPOL" >>>> V/c = 0.990049504
1 CHS STO 06 R/S
or XEQ 10 >>>> V'/c = 1.010050504
Notes:
-This metric corresponds to a rotating frame of reference.
-The difference between V and V' explains the Sagnac effect ( interference
)
-More generally, if the angular velocity = w , the metric may be written in cylindrical coordinates:
ds2 = ( 1 - r2w2/c2 ) c2 dt2 - dr2 - r2 df2 - dz2 - 2 w r2 df dt
-And the spatial metric is
dL2 = dr2 + dz2 + r2 df2 / ( 1 - r2w2/c2 )
-The ratio circumference / radius = 2.p / ( 1 - r2w2/c2 ) 1/2 > 2.p
-As for the speed of light ( if its path is the circumference ) V/c = 1 ± r.w/c ( if we neglect the terms of order 2 and higher )
-Note that we've just gotten a more accurate result in the example above
with r.w/c = 0.01
Remarks:
-"SPOL" is essentially useful if there is at least an index i
for which g0i # 0
-Otherwise, the formula is much simpler:
V/c = Sqrt(g00)
-In this case, V does not depend on the direction k and calculating
the other gab is unuseful.
-For instance, with the Schwarzschild metric:
ds2 = [ 1 - 2 GM / ( c2 r ) - L r2 / 3 ] c2 dt2 - ( dx2 + dy2 + dz2 ) / [ 1 - 2 GM / ( c2 r ) - L r2 / 3 ] where L is the cosmological constant,
G = gravitational constant = 6.673 E-11 m3/ kg / s2 and M = mass of the central body ( for example, the Sun M = 1.989 E30 kg )
-The speed of light is V/c = [ 1 - 2 GM / ( c2 r ) - L r2 / 3 ] 1/2
-Neglecting the cosmological constant, the speed of light near the Sun ( at a distance = Sun's radius = 6.96 E8 m ) is
V/c ~ 0.999997878 whence
V ~ 299791822 m / s
Example of a Luni-Solar Calendar
"GR-LS" & "LS-GR" perform the conversion between Gregorian
dates and a hypothetic Luni-Solar calendar
that employs the following approximations:
1 Lunar month = 334995 / 11344 days ~ 29.53058886
days
1 tropical year / 1 lunar month = 774439 / 62615 ~ 12.36826639
which leads to 1 tropical year ~ 365.2421896 days
STACK | INPUT1 | OUTPUT1 | INPUT2 | OUTPUT2 |
X | YYYY.MNDDGr | yyyy.mnddLS | yyyy.mnddLS | YYYY.MNDDGr |
-3101.0123 XEQ "GR-LS" >>>>
0.0101 XEQ "LS-GR" ( or simply
R/S ) >>>> -3101.0123
1979.0716 XEQ "GR-LS" >>>>
5080.0721 XEQ "LS-GR" ( or simply R/S )
>>>> 1979.0716
2012.1221 XEQ "GR-LS" >>>>
5113.1308 XEQ "LS-GR" ( or simply R/S )
>>>> 2012.1221
2015.1014 XEQ "GR-LS" >>>>
5116.1001 XEQ "LS-GR" ( or simply R/S )
>>>> 2015.1014
Notes:
-Here, the years yyyy may also be negative ( otherwise, delete lines
48-40-39-38 ).
-These programs take the onset of the Kali-Yuga as origin i-e
-3101.0123 ( Gregorian Calendar ) -3101 = 3102 B.C.
Loxodromics on an Ellipoidal Earth
Formulae:
L2 - L1 = ( Tan µ
).{ Ln [ ( Tan ( 45° + b2/2 ) ) / (Tan ( 45° +
b1/2 ) ) ] - (e/2).Ln [ ( 1 + e.sin b2 ).( 1 - e.sin
b1 )/( 1 - e.sin b2 )/( 1 + e.sin b1 )
] }
D = [ a.( 1 - e2 )/ cos µ ] §b1b2
( 1 - e2 sin2 t ) -3/2 dt
( § = Integral )
-One could use any integration formula but "ELLOX" evaluates this integral
by a series expansion up to the terms in e6
( the first neglected term is proportional to e8
)
D
= [ a.( 1 - e2 )/ cos µ ] [ ( 1 + 3e2/4 + 45e4/64
+ 175e6/256 ).( b2 - b1 ) - ( 3e2/8
+ 15e4/32 + 525e6/1024 ).( sin 2b2 - sin
2b1 )
+ ( 15e4/256 + 105e6/1024 ).( sin 4b2
- sin 4b1 ) - ( 35e6/3072 ).( sin 6b2
- sin 6b1 ) ] + .......................
-However, we cannot use this formula if cos µ = 0 , in this case ( i-e if b1 = b2 ) we have:
D = a.(
cos b ).( L2 - L1 ) / ( 1 - e2 sin2
b ) -1/2 the loxodromic line is the parallel
of latitude b = b1 = b2
STACK | INPUTS | OUTPUTS |
T | L1 ( ° ' " ) | / |
Z | b1 ( ° ' " ) | / |
Y | L2 ( ° ' " ) | µ ( ° ' " ) |
X | b2 ( ° ' " ) | D ( km ) |
Example1:
XEQ "ELLOX+" L1^B.1=?
-77.03560 ENTER^
38.55172 R/S
L2^B.2=?
2.20138 ENTER^
48.50112 R/S
>>>> D = 6453.389603 km
X<>Y µ = 80°10'15"31
-Compare this value with the length of the corresponding geodesic:
6181.621427 km
-Rhumb-sailing is easier but slower than orthodromy.
Example2:
XEQ "ELLOX+" L1^B.1=?
0
ENTER^
30
R/S
L2^B.2=?
120 ENTER^
30
R/S >>>> D = 11578.35295
km
X<>Y µ = 90°
Notes:
-Due to roundoff errors, there is a loss of significant digits if the
2 latitudes are very close but not equal.
-Otherwise, the accuracy is of the order of a few centimeters.
-In order to avoid a DATA ERROR & an OUT OF RANGE if the latitudes
are +/- 90°
key in 89°59'59"9999 instead of 90° , the difference
is negligible.
Terrestrial Geodesic Distance ( Andoyer's 2nd order
)
-Andoyer's formula gives a good accuracy.
-Moreover, the method in non-iterative and therefore, fast !
-It employs a spheroidal model of the Erath ( Ellipsoid of revolution
)
STACK | INPUTS | OUTPUTS |
T | L1 ( ° ' " ) | / |
Z | b1 ( ° ' " ) | / |
Y | L2 ( ° ' " ) | / |
X | b2 ( ° ' " ) | D ( km ) |
Example:
77.0356 ENTER^
38.55172 ENTER^
-2.20138 ENTER^
48.50112 XEQ "TGD" >>>>
D = 6181.621449 km
-The error smaller than 15 mm !
-The accuracy is very good ( less than 1 meter ), provided D is not
too large
-So, avoid nearly antipodal points...
Terrestrial Geodesic Distance ( Vincenty's formula
)
-Vincenty's formula gives even more accurate results but, on the other
hand, it's an iterative method.
-The forward and backward azimuths are also calculated.
-The longitudes are positive East
STACK | INPUTS | OUTPUTS |
T | L1 ( ° ' " ) | Az2 ( ° ' " ) |
Z | b1 ( ° ' " ) | Az2 ( ° ' " ) |
Y | L2 ( ° ' " ) | Az1 ( ° ' " ) |
X | b2 ( ° ' " ) | D ( km ) |
where Az1 = forward azimuth
( from P1 to P2 )
and Az2 =
back azimuth ( from P2 to P1 )
Example:
-77.0356 ENTER^
38.55172 ENTER^
2.20138 ENTER^
48.50112 XEQ "TGD2" >>>>
the successive Ln-values are displayed and finally
D = 6181.621427 km
RDN Az1 = 51°47'36"8132
RDN Az2 = -68°09'58"9654
-We also have R06 = µ = azimuth at the equator = 37.74571892°
Notes:
-The accuracy is excellent, provided the 2 points are not
nearly antipodal.
-Otherwise, the algorithm does not converge ! It happens when | Ln
| exceeds 180°
Geodesic Distance: Forward Problem
-Given the coordinates ( L1 , b1 ) of the first
point P1 , the forward azimuth Az1 and the
length D of the geodesic,
"FWD" computes the longitude and the latitude ( L2
, b2 ) of the second point P2
STACK | INPUTS | OUTPUTS |
T | L1 ( ° ' " ) | b2 ( ° ' " ) |
Z | b1 ( ° ' " ) | b2 ( ° ' " ) |
Y | Az1 ( ° ' " ) | b2 ( ° ' " ) |
X | D ( km ) | L2 ( ° ' " ) |
where Az1 is the forward azimuth in the direction from the first point P1 towards the second point P2
Example: L1 = 10°30' b1 = 49°41' Az1 = 12°24' D = 16000 km Calculate L2 & b2
10.30 ENTER^
49.41 ENTER^
12.24 ENTER^
16000 XEQ "FWD" >>> the successive sigma-approximations
are displayed and eventually, it yields:
L2 = -177°03'08"003 X<>Y b2 = -14°06'40"7483
-This method works for antipodal points too and - more generally - for
any length D !
Geodesic Distance - Triaxial Ellipsoidal Earth
"TGD3" works for a triaxial ellipsoid
-The geodesic distance involves integrals which are approximated by
Gauss-Legendre 6-point formula
-The method calculates 1 ( SF 01 ) or 2 terms ( CF 01 ) of Fourier
series.
-Flags F03 & F04 are also used as shown below:
Flags: F01-F03-F04
CF 01 2 terms in the Fourier series
SF 01 only 1 term is found
SF 04 geocentric longitudes and latitudes
CF 04 & CF 03 geocentric longitudes
& geodetic latitudes
CF 04 & SF 03 geodetic
longitudes and latitudes
STACK | INPUTS | OUTPUTS |
T | l1 ( ° ' " ) | / |
Z | b1 ( ° ' " ) | / |
Y | l2 ( ° ' " ) | D1 ( km ) |
X | b2 ( ° ' " ) | D2 or D1 ( km ) |
L | / | D0 ( km ) |
CF 01 ---Execution
time = 3mn41s---
SF 01
---Execution time = 2mn14s---
D0 is computed without Fourier series
D1 is computed with 1 term in the Fourier series
( SF 01 )
D2 is computed with 2 terms in the Fourier series
( CF 01 )
Example1: Geodesic distance between the U.S. Naval Observatory at Washington (D.C.) and the "Observatoire de Paris"
l1 = 77°03'56"0
W , b1 = 38°55'17"2 N
l2 = 2°20'13"8
E , b2 =
48°50'11"2 N
CF 01 CF 03 CF 04
-77.0356 ENTER^
38.55172 ENTER^
2.20138 ENTER^
48.50112 XEQ "TGD3" >>>>
D2 = 6181.625506 km = R01
X<>Y D1 = 6181.625507 km
LASTX D0 = 6181.628094 km = R21
>>> Compared to the value given by the rigorous Jacobi's method ( 6181.625476
), D2 is in error by only 3 centimeters !
Example2: Geodesic distance between the U.S. Naval Observatory at Washington (D.C.)
l1 = 77°03'56"0
W , b1 = 38°55'17"2 N
and Cape Town Observatory ( South Africa )
l2 = 18°28'35"7
E , b2 = 33°56'02"5
S
CF 01 CF 03 CF 04
-77.0356 ENTER^
38.55172 ENTER^
18.28357 ENTER^
-33.56025 XEQ "TGD3" >>>>
D2 = 12709.56468 km = R01
X<>Y D1 = 12709.56609 km
LASTX D0 = 12709.56636 km = R21
Example3: With ( -40° , -40° ) ( +120° , +50° ) it yields:
If CF 01 CF 03 CF 04 geocentric longitudes, geodetic latitudes
D2 = 18095.48412 km = R01
X<>Y D1 = 18095.48948
km
LASTX D0 = 18095.54448 km
= R21
If CF 01 SF 03 CF 04 geodetic longitudes and latitudes
D2 = 18095.49331 km = R01
X<>Y D1 = 18095.49868
km
LASTX D0 = 18095.55361 km
= R21
If CF 01 SF 04 geocentric longitudes and latitudes
D2 = 18099.69052 km = R01
X<>Y D1 = 18099.69587
km
LASTX D0 = 18099.75037 km
= R21
Notes:
-Longitudes are positive East.
-The accuracy of D0 is about 60 meters provided the points
are not nearly antipodal.
-If flag F01 is set, D2 is not computed ( faster but less
accurate in general ).
-The results should verify: D2 <= D1
<= D0 and D2 should
be the most accurate result.
-However, the formulae don't make the differences between
maximum and minimum: they just return an approximate extremum with
D1 & D2
-Unfortunately, the surface of the Earth is much more comlex than that
of a scalene ellipsoid.
-If the latitude b and the height h of an observer O are known, "GEO"
calculates
the geocentric latitude b' and the distance rho = OC where
C is the center of the Earth.
Formulae:
(rho) sin b' = a(1-f) sin u + h
sin b
a = 6378.13661 km
(rho) cos b' = a cos u + h cos
b
f = 1/298.256421
where tan u = (1-f) tan b
STACK | INPUTS | OUTPUTS |
Y | h (km) | rho (km) |
X | b ( ° ' " ) | b' ( ° ' " ) |
Example: Find the geocentric coordinates
of the Mount Palomar Observatory b = 33°21'22"4
h = 1706 m = 1.706 km
1.706 ENTER^
33.21224 XEQ "GEOC" >>>>
b' = 33°10'47"124
X<>Y rho = 6373.415625 km
Note:
XEQ "GEO+" to get PROMPTs for the inputs ... and the outputs !
JDAY and CDAY are mutually reciprocal date functions to convert a given
calendar date into the Julian day number and back to the calendar date.
Use flag 00 to select either Julian or Gregorian calendars in the conversions.
The date format is also MM.DDYYYY regardless of the time module settings
if theres one.
STACK | INPUT/OUTPUT | OUTPUT/INPUT |
X | MN.DDYYYY | Julian Day Nb |
Example: the date May 21st, 2014 corresponds to
2,456,799 (Gregorian calendar) or 2,456,812 (Julian calendar) day numbers.
You can also use JDAY to calculate the elapsed number of days between two dates simply converting both to their Julian day numbers and subtracting them. If you do that youll notice a small discrepancy (18 days) between this approach and the results from DAYS which leads me to believe that DAYS has some different convention perhaps it assumes a 30-day month for financial calculations - but unfortunately it appears to be a stealth function, as there is no documentation for it at all in the Securities Pac manual.
These two functions are based on the PPC routines JC and CJ ported
to an all-MCODE implementation to make it faster.
The formulas are given in the PPC ROM manual.
Rectangular<>Spherical Conversion
x = r cos b cos l
y = r cos b sin l
z = r sin b
where x , y , z = rectangular coordinates,
r = ( x2 + y2 + z2 )1/2
, l = longitude , b = latitude
STACK | INPUTS | OUTPUTS |
T | T | T |
Z | z | b |
Y | y | l |
X | x | r |
L | / | (x2+y2)1/2 |
Example: x = 3 ; y = 4 ; z =
-7 find the spherical coordinates ( in DEG mode )
-7 ENTER^
4 ENTER^
3 XEQ "R-S" r
= 8.602325267
RDN l = 53.13010235°
RDN b = -54.46232221°
-The reciprocal function work similarly
STACK | INPUTS | OUTPUTS |
T | T | T |
Z | b | z |
Y | l | y |
X | r | x |
L | / | (x2+y2)1/2 |
Example: r = 10 ; l =
124° ; b = 37° find x ; y
; z ( in DEG mode )
37 ENTER^
124 ENTER^
10 XEQ "S-R" x = -4.465913097
RDN y = 6.620988446
RDN z = 6.018150232
Note:
R-S and S-R work in all angular modes
Formulae: Hn(x)
= 2x.Hn-1(x) - 2(n-1).Hn-2(x) ; H0(x)
= 1 ; H1(x) = 2x
STACK | INPUTS | OUTPUTS |
Y | n | n |
X | x | Hn(x) |
L | / | x |
Example: Calculate H7(4.9)
7 ENTER^
3.14 XEQ "HMT" >>>> H7 (3.14)
= 73726.24332
Formulae: n.Pn(x)
= (2n-1).x.Pn-1(x) - (n-1).Pn-2(x) ;
P0(x) = 1 ; P1(x) = x
STACK | INPUTS | OUTPUTS |
Y | n | n |
X | x | Pn(x) |
L | / | x |
Example: Calculate P7(4.9)
7 ENTER^
4.9 XEQ "LEG" gives P7(4.9)
= 1698444.018
"P2" solves the equation: a.x2+b.x+c = 0
with a # 0
STACK | INPUTS | OUTPUTS |
Z | a # 0 | c/a |
Y | b | y |
X | c | x |
-If CF 00 the 2 solutions are x , y
-If SF 00 ------------------ x+i.y , x-i.y
Examples: Solve 2.x2 + 3.x - 4 = 0 and 2.x2 + 3.x + 4 = 0
2 ENTER^
3 ENTER^
-4
XEQ "P2" >>>> -2.350781059
X<>Y
0.850781060 Flag F00 is clear,
therefore the 2 solutions are -2.350781059 and
0.850781060
2
ENTER^
3 ENTER^
4 R/S
>>>> -0.75
X<>Y
1.198957881 Flag F00 is set,
therefore the 2 solutions are -0.75 + 1.198957881.i
and -0.75 - 1.198957881.i
"P3" finds the 3 roots of a.x3+b.x2+c.x+d
by the Cardan's ( or Tartaglia's ) formulae ( with a # 0 )
STACK | INPUTS | OUTPUTS |
T | a # 0 | / |
Z | b | z |
Y | c | y |
X | d | x |
-If CF 01 the 3 solutions are x ; y ; z
-If SF 01 ------------------ x ; y+i.z ; y-i.z
Example: Solve 2.x3-5.x2-7.x+1 = 0 and 2.x3-5.x2+7.x+1 = 0
2 ENTER^
-5 ENTER^
-7 ENTER^
1 XEQ "P3" >>>>
3.467727065 RDN 0.131206041
RDN -1.098933107 which are the 3 real solutions
since flag F01 is clear.
2 ENTER^
-5 ENTER^
7 ENTER^
1
R/S >>>> -0.130131632
RDN 1.315065816 RDN 1.453569820
-Flag F01 is set , therefore the 3 solutions are -0.130131633
; 1.315065816 - 1.453569820.i ; 1.315065816
+ 1.453569820.i
"P4" solves the equation x4+a.x3+b.x2+c.x+d
= 0 ( if the leading coefficient is not 1 , divide all the
equation by this coefficient )
STACK | INPUTS | OUTPUTS |
T | a | t |
Z | b | z |
Y | c | y |
X | d | x |
-If CF 01 & CF 02 the 4 solutions are
x ; y ; z ; t
-If SF 01 & CF 02 ------------------
x+i.y ; x-i.y ; z ; t
-If SF 01 & SF 02 ------------------
x+i.y ; x-i.y ; z+i.t ; z-i.t
Example1: Solve x4-2.x3-35.x2+36.x+180 = 0
-2 ENTER^
-35 ENTER^
36 ENTER^
180 XEQ "P4" >>>>
6 RDN 3 RDN -2 RDN
-5 Flags F01 & F02 are clear , the 4 solutions
are 6 ; 3 ; -2 ; -5
Example2: Solve x4-5.x3+11.x2-189.x+522 = 0
-5 ENTER^
11 ENTER^
-189 ENTER^
522
R/S >>>> -2 RDN
5 RDN 3 RDN 6 Flag
F01 is set & F02 is clear , the 4 solutions are -2+5.i ;
-2-5.i ; 3 ; 6
Example3: Solve x4-8.x3+26.x2-168.x+1305 = 0
-8 ENTER^
26 ENTER^
-168 ENTER^
1305 R/S
>>>> -2 RDN 5 RDN
6 RDN 3 Flags F01 & F02
are set , the 4 solutions are -2+5.i ; -2-5.i ; 6+3.i ; 6-3.i
Carlson Elliptic Integral 1st kind: RF(x,y+iz,y-iz)
STACK | INPUTS | OUTPUTS |
Z | z | / |
Y | y | / |
X | x | RF(x,y+iz,y-iz) |
Example:
4 ENTER^
2 ENTER^
1 XEQ "RFZ" >>>> RF
( 1 , 2+4i , 2-4i ) = 0.631488738
Carlson Elliptic Integral 3rd kind: RJ(x,y+iz,y-iz,p)
RJZ computes RJ ( x , y+i.z , y-i.z , p
) with p > 0
STACK | INPUTS | OUTPUTS |
T | p>0 | / |
Z | z | / |
Y | y | / |
X | x | RJ(x,y+iz,y-iz,p) |
Example:
4 ENTER^
3 ENTER^
2 ENTER^
1 XEQ "RJZ" >>>> RJ
( 1 , 2+3i , 2-3i , 4 ) = 0.205644141
References:
[1] Stamatia Mavridès - "L'Univers relativiste"
- Masson ISBN 2-225-36080-7 ( in French )
[2] Jean Heidmann - "Introduction à la cosmologie"
- PUF ( in French )
[3] Landau & Lifshitz - "Classical Theory of
Fields" - Pergamon Press
[4] Feige - "Elliptic Integrals for Cosmological
Constant Cosmologies" - Astron. Nachr 313 (1992) 3, 139-163
[5] Kantowski, Kao, Thomas - "Distance-Redshift Relations
in Inhomogeneous Friedmann-Lemaitre-Robertson-Walker Cosmology" -
The Astrophysical Journal,
545: 549-560
[6] Michèle Caputo - "Some Space Gravity Formulas
and the Dimensions and the Mass of the Earth"
[7] http://earth-info.nga.mil/GandG/publications/tr8350.2/wgs84fin.pdf
[8] Fairhead & Bretagnon - "An Analytica Formula
for the Time Transformation TB-TT" Astronomy & Astrophysics 229, 240-247
( 1990 )
[9] http://eclipse.gsfc.nasa.gov/SEcat5/deltatpoly.html
[10] Henri Arzelies - "Relativité Généralisée,
Gravitation" - Gauthier-Villars ( in French )
[11] Vincenty's
formula
[12] Paul D. Thomas - "Spheroidal Geodesics, Reference Systems
& Local Geometry" - US Naval Oceanographic Office.
[13] Marcin Ligas - Cartesian to geodetic coordinates conversion
on a triaxial ellipsoid - J Geod (2012) 86:249256
[14] If you want to visualize the geodesics on a triaxial ellipsoid,
go to this excellent
page
[15] Nachum Dershowitz & Edward M. Reingold - "Calendrical
Calculations" - Cambridge University Press - ISBN 978-0-521-70238-6
[16] Nachum Dershowitz & Edward M. Reingold - "Indian
Calendrical Calculations"