hp41programs

Template

Age of the Universe & Redshift -> Distance for the HP-41


Overview
 

1°) Empty Universes

  a)  0 < lambda < 1
  b)  General Case

2°) Tolman Universes

  a)  General Case
  b)  Cyclic Universes

3°) Our Universe ?

  a)  Via Carlson Elliptic Integrals ( 2 programs )
  b)  Via Numerical Integration
( 3 programs )

4°)  More General Cyclic Universes  ( 3 programs )

5°)  More General Universes

  a)  Cosmological Constant = 0  ( 3 programs )
  b)  Cosmological Constant # 0
 ( 4 programs )

6°)  Redshifht -> Distance

  a)  Gauss-Legendre 2-point formula
  b)  Tanh-Sinh Quadrature
  c)  Cosmological Constant = 0  ( 2 programs )
  d)  No Matter, Negative Pressure
  e)  No Pressure, Negative Matter
  f)  
Tolman Universes  ( 4 programs )
 g)  Milne Universes




-The age of the universe may be computed by the following integral:   

        t = (1/H)       § 0 1  y. [  (Omega) lambda.y 4   + ( 1-(Omega)tot ).y2 +  (Omega) mat .y + (Omega)rad  ] -1/2  dy            

  where  (Omega)tot = (Omega)mat + (Omega) lambda + (Omega)rad     ,   H = Hubble constant

-If the minimum radius of the Universe is > 0, the age of the Universe is infinite, but these programs calculate the time since the last minimum.

-In these programs below  1/H = 13.77172143  E9  years  which corresponds to H = 71 km/s/Mpc 
  or more generally,  1/H ( in GigaYears ) = 977.7922214 / H  ( if H is expressed in km/s/Mpc )


-In paragraph 6, we calculate:

          D   =  light-travel time distance ( in Giga-light-years ) = light-travel time ( in Giga-years )
          D0 = comoving distance ( in Giga-light-years )
          DL = luminosity-distance ( in Giga-light-years )

Formulae:

    D  = (c/H0)  §y(em)1  y. [  (Omega)lambda.y4 + ( 1-(Omega)tot ).y2 + (Omega)mat.y + (Omega)rad  ] -1/2 dy           y(em) = y at the instant of emission
    D0 = (c/H0)  §y(em)1      [  (Omega)lambda.y4 + ( 1-(Omega)tot ).y2 + (Omega)mat.y + (Omega)rad  ] -1/2 dy          z + 1 = 1/yem  ,  z = redshift  ,   y = R/R0

    where  (Omega)tot = (Omega)mat + (Omega)lambda + (Omega)rad  

                                                                  F(x) = Sinh(x)  if  k = -1     hyperbolic space
     DL = R0 ( z + 1 ) F(D0/R0)    where    F(x) =   x         if  k = 0      euclidean space                      R0 = current scale factor ( or radius )
                                                                 F(x) = Sin(x)    if  k = +1    spherical space


1°) Empty Universes  


 
    a)  0 < Lambda < 1
 

-Here,  (Omega)mat = (Omega)rad = 0


Data Registers:  R00: temp

Flags: /
Subroutines: /


 01 LBL "AEU"
 02 STO 00
 03 1
 04 X<>Y
 05 -
 06 SQRT
 07 1/X
 08 RCL 00
 09 SQRT
 10 STO 00
 11 RCL Y
 12 ST* Y
 13 +
 14 LN
 15 RCL 00
 16 /
 17 977.7922214
 18 R^
 19 /
 20 ST* Z
 21 *
 22 END

 
( 43 bytes / SIZE 001 )

 

          STACK          INPUTS         OUTPUTS
             Y       H0 ( km/s/Mpc )                R
             X       (Omega)lambda                T


Example:
  H0 = 71 km/s/Mpc  &  (Omega)lambda = 0.41

        71    ENTER^
      0.41   XEQ "AEU"   >>>>   16.31804618  Giga-years
                                     X<>Y   17.92925416  Giga-light-years

Notes:

-The constant k is always equal to  -1  in these cases ( hyperbolic universes ) 
-The deceleration parameter  q = - (Omega)lambda


      b)  General Case



Data Registers:   R00: temp

Flags:  F06 & F07
Subroutines: /


 01 LBL "AEU"
 02 DEG
 03 977.7922214
 04 RCL Z
 05 /
 06 STO 00
 07 X<>Y
 08 X#0?
 09 GTO 00
 10 CLST
 11 SIGN
 12 CHS
 13 RCL 00          
 14 ENTER
 15 GTO 03
 16 LBL 00
 17 ENTER
 18 CF 06
 19 X<0?
 20 SF 06
 21 ABS
 22 SQRT
 23 1/X
 24 LASTX
 25 1
 26 R^
 27 X#Y?
 28 GTO 00
 29 CLST
 30 RCL 00          
 31 90
 32 TAN
 33 R^
 34 +
 35 GTO 03
 36 LBL 00
 37 -
 38 CF 07
 39 X<0?
 40 SF 07
 41 ABS
 42 SQRT
 43 1/X
 44 STO Z
 45 *
 46 FS? 06          
 47 GTO 01
 48 ENTER
 49 X^2
 50 SIGN
 51 FS? 07
 52 CHS
 53 ST+ L
 54 X<> L
 55 SQRT
 56 +
 57 LN
 58 GTO 02
 59 LBL 01
 60 ASIN
 61 D-R
 62 LBL 02
 63 R^
 64 *
 65 RCL 00          
 66 ST* Z
 67 ST* T
 68 *
 69 R^
 70 SIGN
 71 CLX
 72 PI
 73 FC? 06
 74 CLX
 75 ST* L
 76 CLX
 77 FC?C 06
 78 FS? 07
 79 X<0?
 80 STO T
 81 ISG X
 82 CLX
 83 FC?C 07
 84 CHS
 85 X<> Z
 86 X<>Y
 87 LBL 03          
 88 END

     
         ( 133 bytes / SIZE 001 )

 

      STACK        INPUTS      OUTPUTS
           T             /     Rmin or  R max
           Z             /     k ( -1 , 0 , +1 )
           Y      H0 ( km/s/Mpc )            R
           X      (Omega)lambda            T
           L             /         0 or P

    Where  P  is the period

Example1:
     H0 = 71 km/s/Mpc  &  (Omega)lambda = -0.4

     71              ENTER^
    0.4   CHS   XEQ "AEU"    >>>>>    T = 12.27985299  Giga-years
                                               RDN      R = 11.63922896  Giga-light-years
                                               RDN       k = -1
                                               RDN  Rmax = 21.77500351  Giga-light-years
                                              LASTX    P = 68.40819107  Giga-years

Example2:      H0 = 71 km/s/Mpc  &  (Omega)lambda = 0.4

      71   ENTER^
     0.4   XEQ "AEU"    >>>>>    T = 16.23322493  Gy
                                     RDN      R = 17.77921592  Gly
                                     RDN       k = -1
                                     RDN  Rmin =  0
                                    LASTX    P =  0   ( no period )


Example3:      H0 = 71 km/s/Mpc  &  (Omega)lambda = 1.4

       71   ENTER^
      1.4   XEQ "AEU"    >>>>>     T = 14.42035714  Gy
                                      RDN      R = 21.77500351   Gly
                                      RDN       k = +1
                                      RDN  Rmin = 11.63922896  Gly
                                    LASTX    P =  0   ( no period )


Example4:      H0 = 71 km/s/Mpc  &  (Omega)lambda = 0

        71  ENTER^
         0   XEQ "AEU"    >>>>>    T = 13.77172143  Gy
                                       RDN     R = 13.77172143  Gly
                                       RDN      k = -1
                                       RDN  Rmin = 0
                                    LASTX    P =  0   ( no period )


Example5:      H0 = 71 km/s/Mpc  &  (Omega)lambda = 1

         71  ENTER^
          1   XEQ "AEU"    >>>>>   T = 9.9999999 E99  Gy  ( in fact: infinite )
                                       RDN     R = 13.77172143     Gly
                                       RDN      k = 0
                                       RDN  Rmin = 0
                                       LASTX   P =  0   ( no period )

Note:

-The deceleration parameter  q = - (Omega)lambda

2°) Tolman Universes    


      a)  General Case


-Here,  (Omega)mat = 0
-Such models have been studied by Tolman.


Data Registers:    R00 = k ( -1 , 0 , +1 )               R06 to R10:  temp

                                R01 = T            R03 = Rmin       R05 = Rmax     R11 =
H0
                                R02 = P            R04 = R
Flags:  F24
Subroutines:


 01 LBL "TOL"
 02 DEG
 03 90
 04 TAN
 05 STO 02
 06 STO 05
 07 STO 09
 08 RDN
 09 STO 06
 10 X<>Y
 11 STO 07
 12 +
 13 X<>Y
 14 STO 11
 15 CLX
 16 STO 03
 17 SIGN
 18 -
 19 STO 00
 20 CHS
 21 RCL 06
 22 X=0?
 23 SIGN
 24 ST/ 07
 25 /
 26 STO 08          
 27 RCL 06
 28 X=0?
 29 GTO 01
 30 ABS
 31 SQRT
 32 ST+ X
 33 STO 10
 34 CLX
 35 2
 36 /
 37 STO 08
 38 X^2
 39 RCL 07
 40 -
 41 X<0?
 42 GTO 02
 43 SQRT
 44 RCL 08
 45 SIGN
 46 *
 47 RCL 08
 48 +
 49 CHS
 50 RCL 07          
 51 RCL Y
 52 X=0?
 53 SIGN
 54 /
 55 X<0?
 56 CLX
 57 SQRT
 58 X<>Y
 59 X<0?
 60 CLX
 61 SQRT
 62 X<Y?
 63 X<>Y
 64 GTO 03
 65 LBL 01
 66 RCL 07
 67 CHS
 68 RCL 08
 69 X=0?
 70 GTO 02
 71 /
 72 X<0?
 73 CLX
 74 SQRT
 75 ENTER
 76 GTO 03
 77 LBL 02
 78 CLST
 79 LBL 03
 80 1
 81 X<>Y
 82 X>Y?
 83 STO 05
 84 X<> Z
 85 X>Y?
 86 STO 05          
 87 X>0?
 88 X>Y?
 89 FS? 30
 90 STO 03
 91 X<> Z
 92 X>0?
 93 X>Y?
 94 FS? 30
 95 STO 03
 96 RCL 06
 97 X=0?
 98 GTO 01
 99 X<0?
100 GTO 02
101 RCL 07
102 RCL 08
103 ST+ X
104 +
105 1
106 +
107 SQRT
108 1
109 +
110 RCL 08
111 ST+ Y
112 RCL 07
113 X<0?
114 CLX
115 SQRT
116 +
117 X=0?
118 GTO 00
119 STO 04        
120 /
121 ABS
122 LN
123 RCL 10
124 /
125 STO 01
126 RCL 09
127 RCL 05
128 X=Y?
129 GTO 03
130 X^2
131 RCL 08
132 +
133 RCL 04
134 /
135 ABS
136 LN
137 ST+ X
138 RCL 10
139 /
140 STO 02
141 GTO 03
142 LBL 01
143 1
144 RCL 07
145 X<0?
146 CLX
147 SQRT
148 STO 04
149 -
150 RCL 08        
151 X=0?
152 GTO 04
153 /
154 STO 01
155 RCL 09
156 RCL 05
157 X=Y?
158 GTO 03
159 X^2
160 RCL 08
161 *
162 RCL 07
163 +
164 SQRT
165 RCL 04
166 -
167 ST+ X
168 RCL 08
169 /
170 STO 02
171 GTO 03
172 LBL 02
173 1
174 RCL 08
175 ST+ Y
176 X^2
177 RCL 07
178 -
179 SQRT
180 /
181 ASIN
182 90
183 +
184 D-R
185 RCL 10
186 /
187 STO 01        
188 PI
189 ST+ X
190 LASTX
191 /
192 STO 02
193 GTO 03
194 LBL 04
195 .5
196 STO 01
197 GTO 03
198 LBL 00
199 RCL 09
200 STO 01
201 LBL 03
202 RCL 00
203 X#0?
204 SIGN
205 ENTER
206 X<> 00
207 ABS
208 SQRT
209 X=0?
210 SIGN
211 1/X
212 SF 24
213 977.7922214
214 RCL 11        
215 /
216 ST* 01
217 ST* 02
218 *
219 STO 04
220 ST* 03
221 ST* 05
222 RCL 01
223 RCL 02
224 RDN
225 CF 24
226 END


         ( 282 bytes / SIZE 012 )

              

           STACK           INPUTS        OUTPUTS
               T
                /
             P
               Z
               H0
             k
               Y          (Omega)rad
             R
               X        (Omega)lambda
             T


Example1:     H0 = 71 km/s/Mpc    (Omega)lambda   = 0.732      (Omega)rad  = 0. 000049

           71        ENTER^
     0.000049   ENTER^
       0.732       XEQ "TOL"   >>>>     T  = 20.19785676 Giga-years         = R10
                                               RDN    R = 26.60483306 Giga-light-years = R04
                                               RDN    k = -1                                            = R00
                                               RDN    P = 9.999999999 E99                   = R02 = period  ( infinity = no period )

-We also have   R03 = Rmin =  0   &   R05 = Rmax = 9.999999999 E99


Example2:      H0 = 71 km/s/Mpc    (Omega)lambda   = -0.1     (Omega)rad  = -0.001

         71     ENTER^
       0.001  CHS   ENTER^
        0.1     CHS   XEQ "TOL"   >>>>    T  = 13.32700845 Gy   = R10
                                                    RDN    R = 13.12485669  Gly  = R04
                                                    RDN    k = -1                           = R00
                                                    RDN    P = 136.8163821 Gy    = R02 = period

-We also have   R03 = Rmin =  0.395565882 Gly   &   R05 = Rmax = 43.54821054  Gly


Notes:

-The constant k = -1 , 0 , +1  corresponds to hyperbolic, euclidean, spherical universes respectively.
-This program also works for pulsating universes.


      b)  Cyclic Universes


-With this program, you choose the minimum, current and maximum scale factor:  Rmin R0 Rmax
  and "AUR" returns the age of the Universe T, the period P, the current Hubble "constant" H0 and the deceleration parameter q


Data Registers:       R00 = L (Gy)-2                       

                                   R01 = R2min      R04 = T
                                   R02 = R20        
                                   R03 = R2max  
Flags: /
Subroutines: /

 
 01 LBL "AUR"
 02 DEG
 03 X^2
 04 STO 00          
 05 STO 03
 06 X<>Y
 07 X^2
 08 ENTER
 09 STO 02
 10 ST- Z
 11 R^
 12 X^2
 13 ST+ 00
 14 STO 01          
 15 -
 16 STO 04
 17 R^
 18 ST* 04
 19 -
 20 ST* Y
 21 RCL 01
 22 RCL 03          
 23 -
 24 /
 25 ACOS
 26 D-R
 27 2
 28 /
 29 X<> 04
 30 ST+ Y
 31 ST/ Y
 32 SQRT
 33 RCL 02          
 34 /
 35 977.7922214
 36 *
 37 3
 38 CHS
 39 X<> 00
 40 ST/ 00
 41 SQRT
 42 ST* 04
 43 ST/ Y
 44 PI
 45 *
 46 RCL 04          
 47 END

     
         ( 76 bytes / SIZE 005 )

 

          STACK           INPUTS       OUTPUTS
               T                /              q
               Z             Rmin             H0
               Y              R0              P
               X             Rmax              T

    Where  T = age of the Universe ,  P  = the period ,  H0 = the current Hubble "constant" ( in km/s/Mpc ) and  q = the deceleration parameter

Example:
       Rmin = 1    R0 = 14   Rmax = 314   ( Gygalightyears )

      1    ENTER^
     14   ENTER^
    314  XEQ "AUR"  >>>>   T  = 13.96898881 = R04               ---Execution time = 2.8s---
                                  RDN   P  = 986.4650960
                                  RDN  H0 = 69.59427439
                                  RDN   q = -0.003136335

-We also have the cosmological constant in  R00 = L = -0.000030427 (Gy)-2

Note:

 k is always equal to -1  ( hyperbolic space )

3°) Our Universe ?
 

    a) Via Carlson Integrals


-The following program only works if the cosmological constant is positive

  and if the quartic:   (Omega)lambda.y 4 + ( 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 !




Data Registers:     R00 thru R21: temp
Flags:  F00-F01-F02-F10
Subroutines:  "RFZ"  "RJZ"  "P4"   (  cf "Elliptic Integrals for the HP41" & "Polynomials for the HP41" )


 
01 LBL "AUM"
02 STO 12
03 X<>Y
04 STO 13
05 +
06 X<>Y
07 STO 14             
08 +
09 1
10 -
11 STO 21
12 CHS
13 RCL 12
14 RCL 14
15 RCL 13             
16 ST/ Z
17 ST/ T
18 /
19 0
20 RDN
21 XEQ "P4"
22 FS?C 01
23 FS?C 02
24 SF 41
25 STO 15
26 RDN
27 STO 16
28 RDN
29 X<Y?
30 X<>Y
31 STO 12
32 STO 14             
33 X<>Y
34 ST* 14
35 -
36 STO 11
37 R^
38 RCL 15
39 R-P
40 X^2
41 STO 17
42 SIGN
43 RCL 16             
44 RCL 12
45 ST- Z
46 RCL 15
47 -
48 R-P
49 X<>Y
50 CHS
51 X<>Y
52 1/X
53 P-R
54 STO 18
55 X<>Y
56 STO 19
57 X<>Y
58 RCL 11             
59 1/X
60 R^
61 1/X
62 ST+ Y
63 ST+ Z
64 RDN
65 XEQ "RJZ"
66 STO 20
67 RCL 19
68 RCL 18
69 RCL 11              
70 1/X
71 RCL 12
72 X=0?
73 GTO 00
74 CHS
75 1/X
76 ST+ Y
77 ST+ Z
78 RDN
79 XEQ "RJZ"
80 ST- 20
81 LBL 00
82 RCL 15
83 RCL 12             
84 -
85 X^2
86 RCL 16
87 X^2
88 +
89 RCL 11
90 *
91 SQRT
92 3
93 *
94 ST/ 20
95 RCL 15               
96 X^2
97 RCL 14
98 -
99 4
100 *
101 RCL 17
102 RCL 15
103 X^2
104 -
105 *
106 SQRT
107 RCL 17
108 RCL 14
109 RCL 15            
110 ST+ X
111 STO 09
112 +
113 1
114 +
115 *
116 SQRT
117 RCL 17             
118 RCL 09
119 -
120 RCL 14
121 ST* Y
122 +
123 SQRT
124 +
125 X^2
126 RCL 15
127 X^2
128 ST+ X
129 CHS
130 RCL 17
131 -
132 RCL 14
133 -
134 X<>Y
135 ST+ Y
136 XEQ "RFZ"
137 RCL 12
138 *
139 ST+ 20
140 RCL 21             
141 X#0?
142 SIGN
143 RCL 21
144 ABS
145 SQRT
146 X=0?
147 SIGN
148 1/X
149 RCL 20
150 ST+ X
151 RCL 13            
152 SQRT
153 /
154 13.77172143 E9
155 ST* Z
156 *
157 END

       
        ( 230 bytes / SIZE 022 )

 
           STACK         INPUTS        OUTPUTS
               Z
       (Omega) rad
             k
               Y        (Omega)lamda
             R
               X        (Omega) mat
             t


  Example1:      (Omega)mat = 0.044 ,  (Omega)lambda  = 0.521 ,   (Omega)rad  = 0.000049 

        0.000049   ENTER^
        0.521         ENTER^
        0.044         XEQ "AUM"    >>>>   t = 15602214530 years
                                                   RDN   R= 20881806280 light-years
                                                   RDN   k = -1

Example2:      (Omega)mat = 0 .271 ,  (Omega)lambda  = 0.732 ,  (Omega) rad  = 0.000049 

        0.000049   ENTER^
        0.732         ENTER^
        0.271         XEQ "AUM"    >>>>   t =  13667886770 years
                                                   RDN   R= 249407504600 light-years
                                                   RDN   k = +1

Note:

-We can also use M-Code routines  RF & RJ
-The following program employs a shorter formula given in reference [1]


Data Registers:     R00 thru R15: temp
Flags:  F00-F01-F02-F10
Subroutines:  "RFZ"  "RJZ" ( M-Code ) and  "P4"   (  cf "Carlson Elliptic Integrals" & "Polynomials for the HP41" )


 
 01 LBL "AUM"
 02 STO 03
 03 X<>Y
 04 STO 09
 05 +
 06 X<>Y
 07 STO 01
 08 +
 09 977.7922214
 10 R^
 11 /
 12 STO 15          
 13 SIGN
 14 -
 15 STO 10
 16 CHS
 17 RCL 01
 18 RCL 03
 19 RCL 09
 20 ST/ Z
 21 ST/ T
 22 /
 23 0
 24 RDN
 25 XEQ "P4"
 26 FS?C 01
 27 FS?C 02
 28 SF 99
 29 STO 05
 30 RDN
 31 ABS
 32 STO 06
 33 RDN
 34 X>Y?
 35 X<>Y
 36 STO 08
 37 X<>Y
 38 STO 07          
 39 *
 40 STO 14
 41 1
 42 RCL 05
 43 -
 44 X^2
 45 RCL 06
 46 X^2
 47 STO 13
 48 +
 49 STO 04
 50 *
 51 SQRT
 52 RCL 07
 53 RCL 08
 54 1
 55 ST- Z
 56 -
 57 *
 58 STO 00
 59 RCL 05
 60 X^2
 61 RCL 13
 62 +
 63 STO 03
 64 *
 65 SQRT
 66 +
 67 X^2
 68 STO 11          
 69 RCL 00
 70 RCL 14
 71 *
 72 RCL 03
 73 *
 74 RCL 04
 75 *
 76 SQRT
 77 ST+ X
 78 RCL 14
 79 ST+ X
 80 RCL 07
 81 -
 82 RCL 08
 83 -
 84 RCL 03
 85 RCL 05
 86 -
 87 *
 88 +
 89 STO 12
 90 RCL 05
 91 RCL 07          
 92 -
 93 X^2
 94 RCL 13
 95 +
 96 CHS
 97 STO 13
 98 RCL 11
 99 +
100 STO 02
101 RCL 06
102 RCL 08
103 RCL 07
104 -
105 STO 01
106 *
107 STO 00
108 RCL 12
109 RCL 11
110 RJZ
111 RCL 13
112 *
113 RCL 01
114 *
115 1.5
116 /
117 X<> 02
118 RCL 07          
119 X=0?
120 GTO 00
121 ENTER
122 ST* Y
123 -
124 /
125 1
126 +
127 SQRT
128 1/X
129 ENTER
130 ST+ Y
131 CHS
132 1
133 +
134 /
135 LN1+X
136 ST+ 02
137 RCL 00
138 RCL 12
139 RCL 11
140 RFZ
141 ST+ X
142 RCL 07
143 *
144 ST+ 02
145 LBL 00
146 RCL 10          
147 X#0?
148 SIGN
149 RCL 10
150 ABS
151 SQRT
152 X=0?
153 SIGN
154 1/X
155 RCL 02
156 RCL 09
157 SQRT
158 /
159 RCL 15
160 ST* Z
161 *
162 END

       
        ( 203 bytes / SIZE 016 )

 
           STACK          INPUTS        OUTPUTS
               T
              H0
             /
               Z
       (Omega) mat
             k
               Y        (Omega)lamda
             R
               X        (Omega)rad
             t


  Example1:    H0 = 71 km/s/Mpc     (Omega)mat = 0.044 ,  (Omega)lambda  = 0.521 ,   (Omega)rad  = 0.000049 

           71           ENTER^
        0.044         ENTER^
        0.521         ENTER^
     0.000049      XEQ "AUM"    >>>>   t = 15.60221456 Giga-years                                         ---Execution time = 34s---
                                                   RDN   R= 20.88180628 Giga-light-years
                                                   RDN   k = -1

Example2:     H0 = 71 km/s/Mpc      (Omega)mat = 0 .271 ,  (Omega)lambda  = 0.732 ,  (Omega) rad  = 0.000049 

           71        ENTER^
        0.271      ENTER^
        0.732      ENTER^
     0.000049  XEQ "AUM"    >>>>   t = 13.66788679  Giga-years                                         ---Execution time = 34s---
                                               RDN   R= 249.4075046 Giga-light-years
                                               RDN   k = +1


   b) Via Numerical Integration  ( 3 programs )
 

-This routine employs Gauss-Legendre 3-point formula to evaluate the integral.
-The formula is applied to N subintervals


Data Registers:     R00 & R04 thru R10: temp                                  ( Registers R01 thru R03 are to be initialized before executing "AUM" )

                              •  R01 = (Omega)mat
                              •  R02 = (Omega)lambda
                              •  R03 = (Omega)rad
Flags: /
Subroutines: /


01 LBL "AUM"
02 STO 05
03 1/X
04 STO 06      
05 2
06 /
07 STO 07
08 1.6
09 STO 10
10 FRC
11 SQRT
12 *
13 STO 08
14 RCL 01
15 RCL 02
16 +
17 RCL 03      
18 +
19 1
20 -
21 STO 00
22 CLX
23 STO 04
24 GTO 01
25 LBL 00
26 STO 09
27 X^2
28 ENTER
29 ENTER
30 X^2
31 RCL 02
32 *
33 RCL 00      
34 -
35 *
36 RCL 01
37 +
38 *
39 RCL 03
40 +
41 SQRT
42 /
43 RCL 09
44 *
45 RTN
46 LBL 01
47 RCL 07
48 RCL 08      
49 -
50 XEQ 00
51 ST+ 04
52 RCL 07
53 XEQ 00
54 RCL 10
55 *
56 ST+ 04
57 RCL 07
58 RCL 08
59 +
60 XEQ 00
61 ST+ 04
62 RCL 06      
63 ST+ 07
64 DSE 05
65 GTO 01
66 RCL 00
67 X#0?
68 SIGN
69 RCL 00
70 ABS
71 SQRT
72 X=0?
73 SIGN
74 1/X
75 1.8
76 *
77 RCL 04      
78 R^
79 *
80 7650956350
81 ST* Z
82 *
83 END


      ( 118 bytes / SIZE 011 )

 
           STACK           INPUTS        OUTPUTS
               Z
                /
             k
               Y                 /
             R
               X                N
             t

    Where N = number of subintervals 

Example1:      (Omega)mat = 0 .044 ,  (Omega)lambda  = 0.521 ,   (Omega) rad  = 0.000049 

        0.000049   STO 03
        0.521         STO 02
        0.044         STO 01      and if you choose  N = 50

         50    XEQ "AUM"       >>>>   t = 15602214540 years
                                            RDN   R= 20881806270 light-years 
                                            RDN   k = -1

Example2:      (Omega)mat = 0 .271 ,  (Omega)lambda  = 0.732 ,  (Omega) rad  = 0.000049 

        0.000049   STO 03
        0.732         STO 02
        0.271         STO 01      and if you choose  N = 50       

         50         XEQ "AUM"      >>>>   t =  13667886820 years
                                                 RDN   R= 249407504600 light-years 
                                                 RDN   k = +1

Notes:

-N = 50 produces a very good precision in most cases.
-However, this program is slower with N = 50 than the version listed in paragraph 3°) a)
-On the other hand, Hubble constant is not known with an accuracy of 10 digits.
-So a smaller N-value is enough to get a 4-digit precision:   N = 10  seems quite enough !

-We can also use Gauss-Legendre 2-point formula.
-The following program makes the change of variable:     y (Omega)mat + (Omega) rad  =  u2


Data Registers:     R00 & R04 thru R11: temp         

                                R01 = (Omega)mat
                                R02 = (Omega)lambda
                                R03 = (Omega)rad
Flags: /
Subroutines: /


 01 LBL "AUM"
 02 STO 03
 03 X<>Y
 04 STO 02
 05 +
 06 X<>Y
 07 STO 01          
 08 +
 09 977.7922214
 10 R^
 11 /
 12 STO 00
 13 SIGN
 14 -
 15 STO 04
 16 RCL 01
 17 RCL 03
 18 +
 19 SQRT
 20 RCL 03          
 21 SQRT
 22 STO 07
 23 -
 24 50
 25 STO 09
 26 ST+ 09
 27 /
 28 STO Y
 29 3
 30 SQRT
 31 STO 11
 32 /
 33 STO 08
 34 -
 35 STO 05
 36 2
 37 /
 38 ST- 07
 39 CLX
 40 LBL 01
 41 RCL 07          
 42 RCL 08
 43 X<> 05
 44 STO 08
 45 +
 46 STO 07
 47 X^2
 48 STO 10
 49 RCL 03
 50 -
 51 RCL 01
 52 /
 53 ENTER
 54 X^2
 55 STO 06
 56 RCL 02          
 57 *
 58 RCL 04
 59 -
 60 RCL 06
 61 *
 62 RCL 10
 63 +
 64 SQRT
 65 /
 66 RCL 07
 67 *
 68 +
 69 DSE 09
 70 GTO 01
 71 RCL 11          
 72 *
 73 RCL 08
 74 *
 75 RCL 00
 76 *
 77 RCL 01
 78 /
 79 RCL 04
 80 X#0?
 81 SIGN
 82 RCL 00          
 83 RCL 04
 84 ABS
 85 SQRT
 86 X=0?
 87 SIGN
 88 /
 89 R^
 90 END


      ( 115 bytes / SIZE 012 )

 
           STACK           INPUTS        OUTPUTS
               T
       H ( km/s/Mpc )
             /
               Z
      (Omega)mat > 0
             k
               Y        (Omega)lambda
             R
               X         (Omega)rad
             T


Example1:
   H = 71 km/s/Mpc  (Omega)mat = 0.044 ,  (Omega)lambda  = 0.521 ,   (Omega) rad  = 0.000049 

        71       ENTER^
     0.044     ENTER^
     0.521     ENTER^
  0.000049   XEQ "AUM"   >>>>   T = 15.60221456   ( Gy )                                         ---Execution time = 108s---
                                            RDN    R = 20.88180628  ( Gly )
                                            RDN    k = - 1  ( hyperbolic space )

Example2:    H = 71 km/s/Mpc   (Omega)mat = 0.271 ,  (Omega)lambda  = 0.732 ,  (Omega) rad  = 0.000049 

        71       ENTER^
     0.271     ENTER^
     0.732     ENTER^
  0.000049   XEQ "AUM"   >>>>   T = 13.66788679   ( Gy )
                                            RDN    R = 249.4075046  ( Gly )
                                            RDN    k = + 1  ( spherical space )


Example3:    H = 71 km/s/Mpc   (Omega)mat = 0.271 ,  (Omega)lambda  = 0.729 ,  (Omega) rad  = 0

        71       ENTER^
     0.271     ENTER^
     0.729     ENTER^
        0        XEQ "AUM"   >>>>   T = 13.65710890   ( Gy )
                                          RDN    R = 13.77172143  ( Gly )
                                          RDN    k =  0  ( Euclidean space )

Notes:

 50 subintervals seem enough for a very good precision
 Change line 24 if you prefer another value.

-We may simplify the program if  (Omega) rad  = 0
-The following routine employs again Gauss-Legendre 2-point formula ( and the change of variable  y2 = u )


Data Registers:   R00 to R07: temp
Flags: /
Subroutines: /


 01 LBL "AUM"
 02 STO 02
 03 2
 04 /
 05 X<>Y
 06 STO T
 07 ST+ T
 08 -
 09 STO 04          
 10 3
 11 *
 12 X<>Y
 13 STO 00
 14 SIGN
 15 -
 16 +
 17 STO 05
 18 CLX
 19 SIGN
 20 50
 21 STO 03
 22 ST+ 03
 23 /
 24 STO Y
 25 3
 26 SQRT
 27 /
 28 STO 01          
 29 -
 30 STO 06
 31 2
 32 /
 33 CHS
 34 STO 07
 35 CLX
 36 LBL 01
 37 RCL 07
 38 RCL 01
 39 X<> 06
 40 STO 01
 41 +
 42 STO 07
 43 X^2
 44 ENTER
 45 X^2
 46 RCL 04          
 47 *
 48 RCL 05
 49 -
 50 RCL Y
 51 *
 52 RCL 02
 53 +
 54 SQRT
 55 /
 56 +
 57 DSE 03
 58 GTO 01
 59 RCL 01          
 60 *
 61 3
 62 SQRT
 63 *
 64 977.7922214
 65 RCL 00
 66 /
 67 STO 06
 68 *
 69 RCL 05
 70 X#0?
 71 SIGN
 72 RCL 06          
 73 RCL 05
 74 ABS
 75 SQRT
 76 X=0?
 77 SIGN
 78 /
 79 R^
 80 END

     
         ( 107 bytes / SIZE 008 )

 
           STACK           INPUTS        OUTPUTS
               Z
       H ( km/s/Mpc )
             k
               Y                q
             R
               X         (Omega)mat
             T


Example1:
   H = 71 km/s/Mpc   q = -0.5    (Omega)mat = 0.044

        71       ENTER^
      -0.5       ENTER^
     0.044     XEQ "AUM"   >>>>    T  = 15.62753949   ( Gy )                                         ---Execution time = 82s---
                                            RDN    R = 20.90467233  ( Gly )
                                            RDN    k = - 1  ( hyperbolic space )

Example2:    H = 71 km/s/Mpc   q = -0.6   (Omega)mat = 0.271 

        71       ENTER^
      -0.6      ENTER^
     0.271    XEQ "AUM"   >>>>   T = 13.68753141   ( Gy )                                         ---Execution time = 82s---
                                         RDN    R = 170.8171812  ( Gly )
                                         RDN    k = + 1  ( spherical space )

Example3:    H = 71 km/s/Mpc    q = -0.595    (Omega)mat = 0.27

        71      ENTER^
    -0.595    ENTER^
      0.27     XEQ "AUM"   >>>>   T = 13.67100708   ( Gy )                                         ---Execution time = 82s---
                                          RDN    R = 13.77172143  ( Gly )
                                          RDN    k =  0  ( Euclidean space )

Example4:    H = 71 km/s/Mpc    q = -0.5    (Omega)mat = 0

        71     ENTER^
      -0.5    ENTER^
        0      XEQ "AUM"   >>>>   T = 17.16576879   ( Gy )                                         ---Execution time = 82s---
                                       RDN    R = 19.47615522  ( Gly )
                                       RDN    k = - 1  ( hyperbolic space )

Notes:

-N = 50 subintervals seem enough for a very good accuracy ( line 20 )

-Here we have  q = (1/2) (Omega)mat -  (Omega)lambda

-Though this program works well if  (Omega)mat = 0 , it's much slower than the routine listed in paragraph 1°)


4°)  More General Cyclic Universes


-The age of the Universe  T = (c/H)   § ymin 1  y. [  (Omega) lambda.y 4   + ( 1-(Omega)tot ).y2 +  (Omega) mat .y + (Omega)rad  ] -1/2  dy            
-The Period  P = 2 (c/H) § ymin ymax  y. [  (Omega) lambda.y 4   + ( 1-(Omega)tot ).y2 +  (Omega) mat .y + (Omega)rad  ] -1/2  dy            

-The following program assumes that the quartic equation    (Omega) lambda.y 4   + ( 1-(Omega)tot ).y2 +  (Omega) mat .y + (Omega)rad  = 0  has (at least ) 2 positive real roots.
-They are computed by Newton's method.

-The integrals are then calculated by Gauss-Legendre 3-point formula, applied to N subintervals between ymin & 1  and 3.N subintervals between 1 & ymax ( lines 266-267 )


Data Registers:     R00 = H ( km/s/Mpc )                                                         ( Registers R00 thru R03 are to be initialized before executing "ACU" )

                               •  R01 = (Omega)mat                     R04 = Rmin                        R06 thru R18: temp
                               •  R02 = (Omega)lambda < 0          R05 = Rmax
                               •  R03 = (Omega)rad
Flags: /
Subroutines: /


 01 LBL "ACU"
 02 STO 16
 03 1
 04 RCL 01          
 05 -
 06 RCL 02
 07 -
 08 RCL 03
 09 -
 10 STO 06
 11 STO 17
 12 2 E-8
 13 STO 11
 14 CLX
 15 STO 09
 16 STO 10
 17 LBL 01
 18 RCL 10
 19 ENTER
 20 ENTER
 21 X^2
 22 RCL 02
 23 *
 24 RCL 06
 25 +
 26 *
 27 RCL 01
 28 +
 29 *
 30 RCL 03          
 31 +
 32 X<>Y
 33 X^2
 34 RCL 02
 35 *
 36 ST+ X
 37 RCL 06
 38 +
 39 R^
 40 *
 41 ST+ X
 42 RCL 01
 43 +
 44 /
 45 -
 46 STO 10          
 47 ST- Y
 48 X=0?
 49 SIGN
 50 /
 51 ABS
 52 RCL 11
 53 X<Y?
 54 GTO 01
 55 RCL 10
 56 STO 04
 57 ENTER
 58 X^2
 59 RCL 02
 60 *
 61 RCL 06
 62 +
 63 STO 07
 64 *
 65 RCL 01
 66 +
 67 STO 08
 68 LBL 02
 69 RCL 10
 70 RCL 10
 71 RCL 10
 72 RCL 04
 73 +
 74 *
 75 RCL 02
 76 *
 77 RCL 07
 78 +
 79 *
 80 RCL 08
 81 +
 82 X<>Y
 83 3
 84 *
 85 RCL 04
 86 ST+ X
 87 +
 88 RCL 02          
 89 *
 90 R^
 91 *
 92 RCL 07
 93 +
 94 /
 95 -
 96 STO 10
 97 ST- Y
 98 X=0?
 99 SIGN
100 /
101 ABS
102 RCL 11
103 X<Y?
104 GTO 02
105 RCL 10
106 STO 05
107 RCL 04
108 +
109 CHS
110 2
111 /
112 ENTER
113 X^2
114 RCL 04
115 RCL 05
116 +
117 RCL 04
118 *
119 RCL 05
120 X^2
121 +
122 STO 07
123 -
124 RCL 06
125 RCL 02
126 /
127 ST+ 07
128 -
129 X<0?
130 GTO 04
131 SQRT
132 ST+ Z
133 -
134 X>Y?
135 X<>Y
136 RCL 04          
137 X>Y?
138 X<>Y
139 RCL 05
140 X>Y?
141 X<>Y
142 STO 07
143 RDN
144 X<Y?
145 X<>Y
146 RDN
147 X<Y?
148 X<>Y
149 R^
150 X<Y?
151 X<>Y
152 STO 05
153 RDN
154 STO 04
155 X<>Y
156 STO 06
157 1
158 RCL 04
159 X<Y?
160 GTO 03
161 RCL 06
162 X<> 04
163 X<> 05
164 STO 06
165 LBL 03
166 RCL 06
167 RCL 07
168 +
169 CHS
170 X<> 06
171 ST* 07
172 GTO 05
173 LBL 04
174 RCL 04
175 RCL 05
176 X>Y?
177 X<>Y
178 STO 04          
179 X<>Y
180 STO 05
181 +
182 STO 06
183 LBL 05
184 1
185 RCL 05
186 RCL 04
187 X<=0?
188 LN
189 ST- Z
190 -
191 STO 08
192 /
193 SQRT
194 RAD
195 ASIN
196 STO 11
197 1.6
198 STO 15
199 GTO 14
200 LBL 00
201 SIN
202 X^2
203 RCL 08
204 *
205 RCL 04
206 +
207 ENTER
208 STO Z
209 RCL 06
210 +
211 *
212 RCL 07
213 +
214 SQRT
215 /
216 RTN
217 LBL 10
218 RCL 11
219 RCL 09          
220 -
221 X<>Y
222 STO 13
223 /
224 STO 14
225 2
226 /
227 ST+ 09
228 .6
229 SQRT
230 *
231 STO 10
232 CLX
233 STO 12
234 LBL 12
235 RCL 09
236 RCL 10
237 -
238 XEQ 00
239 ST+ 12
240 RCL 09
241 XEQ 00
242 RCL 15
243 *
244 ST+ 12
245 RCL 09
246 RCL 10
247 +
248 XEQ 00
249 ST+ 12
250 RCL 14
251 ST+ 09
252 DSE 13
253 GTO 12
254 RCL 12
255 *
256 RTN
257 LBL 14
258 RCL 16
259 XEQ 10
260 STO 18          
261 SIGN
262 ASIN
263 X<> 11
264 STO 09
265 RCL 16
266 3
267 *
268 XEQ 10
269 RCL 18
270 ST+ Y
271 X<>Y
272 ST+ X
273 1.8
274 RCL 02
275 CHS
276 SQRT
277 *
278 ST/ Z
279 /
280 977.7922214
281 RCL 00
282 /
283 ST* Y
284 ST* Z
285 RCL 17
286 SQRT
287 /
288 ST* 04
289 ST* 05
290 X<> Z
291 DEG
292 END


         ( 370 bytes / SIZE 019 )

 
           STACK           INPUTS        OUTPUTS
               Z
                /
             R
               Y                 /
             P
               X                N
             T

   Where  N = number of subintervals   ,  T = Age of the Universe in Gigayears  ,
                P = Period of the Universe in Gigayears   ,   R = radius of the Universe in Gigalightyears

Example:      Hubble "constant" = 71 km/s/Mpc  ,  (Omega)mat = -0.1 ,  (Omega)lambda  = -0.003 ,   (Omega) rad  = -0.001 

            71      STO 00

          -0.1     STO 01
        -0.003   STO 02
        -0.001   STO 03      and if you choose  N = 8   subintervals

             8       XEQ "ACU"    >>>>   T = 14.61256009   Gy                                           ---Execution time = 2mn10s---
                                              RDN    P = 796.1455569   Gy
                                              RDN    R =  13.10701187  GLy

-We also have:   R04 = minimum radius of the Universe  = Rmin = 1.306380222  GLy
                          R05 = maximum radius of the Universe = Rmax = 250.8400402 GLy


Notes:

-The numerical integration is performed with the change of argument:  y = A + ( B - A ) Sin2 u   where  A = Rmin/R  &  B = Rmax/R
-We must choose:   1-(Omega)tot > 0
-The constant k = -1  ( hyperbolic universes )

-We can also choose the minimum, current & maximum scale factors  Rmin R0 Rmax
-The following program employs Carlson elliptic integral  RJ


Data Registers:     R00 =  L = cosmological constant ( Gy-2 )

                                 R01 = Rmin            R04 = H0                        R06 thru R15: temp
                                 R02 = R0               R05 =  q
                                 R03 = Rmax
Flag:  F07
Subroutines: /


 01 LBL "ACU"
 02 DEG
 03 STO 03
 04 RDN
 05 STO 02
 06 X<>Y
 07 STO 01
 08 R^
 09 +
 10 STO 04
 11 LASTX
 12 SF 07
 13 XEQ 00
 14 ST+ X
 15 STO 09
 16 RCL 02
 17 LBL 00          
 18 STO 10          
 19 RCL 03
 20 RCL 01
 21 -
 22 STO 11
 23 *
 24 RCL 01
 25 ST* 11
 26 STO 12
 27 RCL 04
 28 +
 29 ST* 12
 30 *
 31 X<> 11
 32 RCL 04
 33 RCL 10
 34 +
 35 *
 36 X<> 12
 37 RCL 03
 38 RCL 10
 39 -
 40 *
 41 RCL 10
 42 RCL 01
 43 -
 44 ST/ 11
 45 ST/ 12
 46 /
 47 STO 13
 48 RCL 11
 49 RCL 03
 50 RCL 04          
 51 *
 52 -
 53 STO 14
 54 STO 15
 55 CLX
 56 STO 00
 57 STO 08
 58 SIGN
 59 STO 07
 60 LBL 01
 61 RCL 11
 62 SQRT
 63 ENTER
 64 STO 05
 65 RCL 12
 66 SQRT
 67 ST* Z
 68 STO 06
 69 RCL 13
 70 SQRT
 71 ST* T
 72 ST* 06
 73 +
 74 ST* 05
 75 +
 76 RCL 14
 77 *
 78 +
 79 X^2
 80 RCL 14
 81 RCL 05
 82 RCL 06          
 83 +
 84 ST+ 11
 85 ST+ 12
 86 ST+ 13
 87 RCL 14
 88 +
 89 STO 14
 90 X^2
 91 *
 92 X=Y?
 93 GTO 02
 94 -
 95 STO 05
 96 X<>Y
 97 /
 98 X<0?
 99 GTO 03
100 SQRT
101 ENTER
102 ST+ Y
103 SIGN
104 LASTX
105 -
106 /
107 LN1+X
108 2
109 /
110 GTO 04
111 LBL 02
112 SQRT
113 1/X
114 GTO 05
115 LBL 03
116 CHS
117 SQRT
118 ATAN
119 D-R
120 LBL 04
121 RCL 05          
122 ABS
123 SQRT
124 /
125 LBL 05
126 RCL 07
127 *
128 ST+ 00
129 4
130 ST/ 07
131 ST/ 11
132 ST/ 12
133 ST/ 13
134 ST/ 14
135 RCL 08
136 RCL 07
137 RCL 11
138 RCL 12
139 +
140 RCL 13
141 +
142 RCL 14
143 ST+ X
144 +
145 5
146 /
147 ENTER
148 SQRT
149 *
150 /
151 RCL 00          
152 3
153 *
154 +
155 STO 08
156 X>Y?
157 GTO 01
158 RCL 01
159 *
160 RCL 03
161 *
162 RCL 04
163 *
164 3
165 /
166 CHS
167 RCL 15
168 RCL 01
169 RCL 10
170 *
171 /
172 SQRT
173 1/X
174 FS? 07
175 SIGN
176 ASIN
177 D-R
178 +
179 ST+ X
180 FS?C 07
181 RTN
182 RCL 03
183 RCL 04          
184 *
185 RCL 01
186 X^2
187 +
188 STO 00
189 SQRT
190 ST* 09
191 *
192 3
193 RCL 00
194 /
195 CHS
196 STO 00
197 RCL 02
198 RCL 01
199 -
200 *
201 RCL 02
202 RCL 03
203 -
204 *
205 RCL 02
206 RCL 04
207 +
208 *
209 3
210 RCL 02
211 ST/ Z
212 X^2
213 *
214 /
215 STO 04          
216 RCL 02
217 X^2
218 1/X
219 -
220 RCL 00
221 -
222 RCL 04
223 ST+ X
224 /
225 STO 05
226 RCL 04
227 SQRT
228 977.7922214
229 *
230 STO 04
231 RCL 09
232 R^
233 END


         ( 288 bytes / SIZE 016 )

 
          STACK           INPUTS       OUTPUTS
               T                /              q
               Z             Rmin             H0
               Y              R0              P
               X            Rmax              T

    Where  T = age of the Universe ,  P  = the period ,  H0 = the current Hubble "constant" and  q = the deceleration parameter

Example:
       Rmin = 1    R0 = 14   Rmax = 314   ( Gigalightyears )

   1    ENTER^
  14   ENTER^
 314  XEQ "AUC"   >>>>   T = 15.49057984   gigayears                                 ---Execution time = 58s---
                                 RDN   P = 993.8651299   gigayears
                                 RDN  H0 = 67.22990237  km/s/Mpc
                                 RDN  q  =  -0.036404801

  and  R00 = L = -0.000030330 Gy-2


-We can use M-Code routines  RF  &  RJ

-I've also remarked that the change of variable:    x = (a.y+r) / (y+1)   yields

     §ar x dx / sqrt [ x (b-x) (x+c) (x+d) ]  =  §0+infinity  sqrt(r-a)  sqrt (a.y+r)  /  (y+1) / sqrt [ ((b-a).y + (b-r)) ((a+c).y + (r+c)) ]  dy  

  =  ...  =  2 SQRT [ (r-a) / a / (b-a) / (a+c) ]  {  a RF [ (b-r)/(b-a) , (r+c)/(a+c) , r/a ]  + [(r-a)/3]  RJ [ (b-r)/(b-a) , (r+c)/(a+c) , r/a ; 1 ] }

  and with this formula, the program becomes slightly shorter.


Data Registers:     R00 =  L = cosmological constant ( Gy-2 )

                                 R01 = Rmin            R04 = H0           R10 = P             R06 thru R11: temp
                                 R02 = R0               R05 =  q
                                 R03 = Rmax
Flag:  F07
Subroutines:  M-Code routines:  RF & RJ  ( cf "Carlson Elliptic Integrals for the HP41" )


 01 LBL "AUM"
 02 STO 03          
 03 STO 04
 04 STO 05
 05 RDN
 06 STO 02
 07 X<>Y
 08 STO 01
 09 ST- 04
 10 ST+ 05
 11 RCL 05
 12 +
 13 STO 09
 14 RCL 03
 15 SF 07
 16 XEQ 01
 17 ST+ X
 18 STO 10
 19 RCL 02
 20 LBL 01
 21 STO 00          
 22 RCL 01
 23 /
 24 STO 11
 25 RCL 03
 26 RCL 00
 27 -
 28 RCL 04
 29 /
 30 STO 06
 31 RCL 00
 32 RCL 05
 33 +
 34 RCL 09
 35 /
 36 STO 07
 37 RF
 38 RCL 01
 39 *
 40 STO 08          
 41 SIGN
 42 RCL 06
 43 RCL 07
 44 RCL 11
 45 RJ
 46 RCL 00
 47 RCL 01
 48 -
 49 STO 06
 50 *
 51 3
 52 /
 53 RCL 08
 54 +
 55 RCL 06
 56 RCL 01
 57 /
 58 RCL 04
 59 /
 60 RCL 09          
 61 /
 62 SQRT
 63 ST+ X
 64 *
 65 FS?C 07
 66 RTN
 67 RCL 03
 68 RCL 05
 69 *
 70 RCL 01
 71 X^2
 72 +
 73 STO 00
 74 SQRT
 75 ST* 10
 76 *
 77 3
 78 RCL 00          
 79 /
 80 CHS
 81 STO 00
 82 RCL 02
 83 ST+ 05
 84 RCL 01
 85 -
 86 *
 87 RCL 02
 88 RCL 03
 89 -
 90 *
 91 RCL 05
 92 *
 93 RCL 02
 94 3
 95 ST/ Z
 96 Y^X
 97 /
 98 STO 04          
 99 RCL 02
100 X^2
101 1/X
102 -
103 RCL 00
104 -
105 RCL 04
106 ST+ X
107 /
108 STO 05
109 RCL 04
110 SQRT
111 977.7922214
112 *
113 STO 04
114 RCL 10
115 R^
116 END


      ( 148 bytes / SIZE 012 )

 
          STACK           INPUTS       OUTPUTS
               T                /              q
               Z             Rmin             H0
               Y              R0              P
               X             Rmax              T

    Where  T = age of the Universe ,  P  = the period ,  H0 = the current Hubble "constant" and  q = the deceleration parameter

Example:
       Rmin = 1    R0 = 14   Rmax = 314   ( Gigalightyears )

   1    ENTER^
  14   ENTER^
 314  XEQ "AUC"   >>>>   T = 15.49057983   gigayears                                 ---Execution time = 42s---
                                 RDN   P = 993.8651299   gigayears
                                 RDN  H0 = 67.22990236  km/s/Mpc
                                 RDN  q  =  -0.036404801

  and  R00 = L = -0.000030330 Gy-2

Note:

-Cf also "Cyclic Universes for the HP41"


5°) More General Universes

      a)  Cosmological Constant = 0


-This program uses Einstein's equation with  Lambda = 0  &  k = +1 ( spherical space )
-You choose Rmin R0 Rmax and "AUL" returns the age of the Universe T, the period P, the current Hubble "constant" H0 and the deceleration parameter q



Data Registers:       R00: temp  R01 = Rmin  R02 = R0  R03 = Rmax   R04 = T
Flags: /
Subroutines: /


 01 LBL "AUL"
 02 DEG
 03 STO 00
 04 STO 03          
 05 X<>Y
 06 ENTER
 07 STO 02
 08 ST- Z
 09 R^
 10 ST+ 00
 11 STO 01
 12 -
 13 STO 04          
 14 R^
 15 ST* 04
 16 -
 17 ST* Y
 18 RCL 01
 19 RCL 03
 20 -
 21 /
 22 ACOS
 23 D-R
 24 RCL 00          
 25 *
 26 2
 27 ST/ Z
 28 /
 29 X<> 04          
 30 ST+ Y
 31 ST/ Y
 32 SQRT
 33 ST- 04
 34 977.7922214
 35 *
 36 RCL 02
 37 X^2
 38 /
 39 RCL 00
 40 PI
 41 *
 42 RCL 04          
 43 END

     
         ( 70 bytes / SIZE 005 )

 

          STACK           INPUTS       OUTPUTS
               T                /              q
               Z             Rmin             H0
               Y              R0              P
               X             Rmax              T

    Where  T = age of the Universe ,  P  = the period ,  H0 = the current Hubble "constant" and  q = the deceleration parameter

Example1:
       Rmin = 2    R0 = 41   Rmax = 257   ( Gygalightyears )

   2    ENTER^
  41   ENTER^
 257  XEQ "AUL"   >>>>   T  = 12.28422136   gigayears                                 ---Execution time = 2.5s---
                                 RDN   P = 813.6724974   gigayears
                                 RDN  H0 = 53.38731059  km/s/Mpc
                                 RDN  q  =  0.569266382


Example2:        Rmin = 3    R0 = 37   Rmax = 300   ( Gygalightyears )

   3    ENTER^
  37   ENTER^
 300  XEQ "AUL"   >>>>   T  = 10.02097732   gigayears                              
                                 RDN   P = 951.9025742   gigayears
                                 RDN  H0 = 67.53990768  km/s/Mpc
                                 RDN  q  =  0.526224558


Example3:        Rmin = 10    R0 = 123   Rmax = 9876   ( Gygalightyears )

   10    ENTER^
  123   ENTER^
 9876  XEQ "AUL"   >>>>  T  = 10.23464300   gigayears                              
                                  RDN   P  = 31057.78498   gigayears
                                  RDN  H0 = 67.84919334  km/s/Mpc
                                  RDN  q  =  0.462057965

Note:

-The following program is more general:
-You choose the parameters Hubble "constant" ,  (Omega)mat &  (Omega) rad 

-Since  (Omega) lambda = 0  ,  we can calculate the integrals with elementary functions ( cf reference [2] )



Data Registers:    R00 = H0    

                                R01 = 
(Omega)mat      R03 = T          R05 = Rmin      R07 = Rmax      R09-R10-R11:  temp
                                R02 = 
(Omega)rad       R04 = P          R06 = R0         R08 = k
Flags:  F24
Subroutines:  /


-Lines 27-78-102  are three-byte GTOs


 01 LBL "AUL"
 02 DEG
 03 STO 02
 04 X<>Y
 05 STO 01
 06 +
 07 X<>Y
 08 STO 00
 09 SIGN
 10 X<>Y
 11 -
 12 STO 06
 13 90
 14 TAN
 15 STO 07
 16 CLX
 17 STO 04
 18 STO 05
 19 X<>Y
 20 X#0?
 21 GTO 03
 22 RCL 01
 23 X#0?
 24 GTO 01
 25 2
 26 1/X
 27 GTO 06
 28 LBL 01
 29 CHS
 30 X<=0?
 31 GTO 02
 32 RCL 02          
 33 X<>Y
 34 /
 35 STO 07
 36 RCL 02
 37 ENTER
 38 SQRT
 39 *
 40 8
 41 *
 42 RCL 01
 43 X^2
 44 3
 45 *
 46 /
 47 STO 04
 48 LBL 02
 49 1
 50 RCL 02
 51 RCL 01
 52 CHS
 53 /
 54 X<=Y?
 55 X<0?
 56 CLX
 57 STO 05
 58 RCL 01
 59 *
 60 RCL 02
 61 +
 62 SQRT
 63 RCL 02
 64 LASTX
 65 3
 66 /
 67 -
 68 *
 69 RCL 02          
 70 -
 71 3
 72 1/X
 73 +
 74 ST+ X
 75 RCL 01
 76 X^2
 77 /
 78 GTO 06
 79 LBL 03
 80 ABS
 81 SQRT
 82 ST+ X
 83 STO 09
 84 RCL 02
 85 ABS
 86 SQRT
 87 STO 08
 88 RCL 01
 89 CHS
 90 RCL 06
 91 ST+ X
 92 STO 10
 93 /
 94 ENTER
 95 X^2
 96 RCL 02
 97 RCL 06          
 98 /
 99 STO T
100 -
101 X<0?
102 GTO 05
103 SQRT
104 RCL Y
105 SIGN
106 *
107 +
108 X#0?
109 ST/ Y
110 X>Y?
111 X<>Y
112 STO 11
113 X<0?
114 CLX
115 STO 05
116 1
117 X<Y?
118 GTO 04
119 X<>Y
120 RDN
121 X<Y?
122 GTO 00
123 X<>Y
124 X>0?
125 STO 05
126 GTO 03
127 LBL 00
128 X<>Y
129 STO 07
130 LBL 03
131 RCL 05
132 0
133 X#Y?
134 STO 08
135 RCL 06        
136 X>0?
137 GTO 05
138 RCL 05
139 RCL 11
140 -
141 RCL 05
142 +
143 RCL 07
144 ST- Y
145 RCL 11
146 -
147 STO 10
148 /
149 ASIN
150 STO 03
151 ST+ X
152 D-R
153 PI
154 -
155 RCL 01
156 *
157 RCL 09
158 /
159 RCL 08
160 ST+ X
161 -
162 STO 04
163 RCL 07
164 RCL 11
165 +
166 2
167 -
168 RCL 10
169 /
170 ASIN
171 RCL 03
172 +
173 D-R
174 RCL 01
175 *
176 RCL 09        
177 /
178 1
179 +
180 RCL 08
181 -
182 RCL 06
183 ST/ 04
184 /
185 GTO 06
186 LBL 04
187 RCL 08
188 RCL 09
189 *
190 RCL 10
191 RCL 11
192 STO 07
193 *
194 RCL 01
195 ST+ Z
196 +
197 /
198 LN
199 RCL 01
200 *
201 RCL 09
202 /
203 RCL 08
204 -
205 RCL 06
206 /
207 ST+ X
208 STO 04
209 CLX
210 STO 05
211 LBL 05
212 1
213 RCL 05
214 RCL 10        
215 *
216 RCL 01
217 +
218 RCL 08
219 ST- Z
220 RCL 09
221 *
222 +
223 RCL 01
224 RCL 10
225 +
226 RCL 09
227 +
228 /
229 LN
230 RCL 01
231 *
232 RCL 09
233 /
234 +
235 RCL 06
236 /
237 LBL 06
238 STO 03
239 RCL 06
240 CHS
241 X#0?
242 SIGN
243 STO 08
244 977.7922214
245 RCL 00
246 /
247 ST* 03
248 ST* 04
249 RCL 06        
250 X=0?
251 SIGN
252 ABS
253 SQRT
254 /
255 SF 24
256 ST* 05
257 STO 06
258 ST* 07
259 RCL 04
260 RCL 03
261 CF 24
262 END


         ( 314 bytes / SIZE 012 )

              

           STACK          INPUTS        OUTPUTS
               T
               /
             k
               Z
              H0
             R0
               Y        (Omega)mat
             P
               X        (Omega)rad
             T

  Where    T = age of the Universe ( or the time since the last minimum scale factor )   P = period  ( or 0 if no period )   R0 = current scale factor

Example1:
   H0 = 68 km/s/Mpc     (Omega)mat = 1.4       (Omega)rad = -0.2

    68  ENTER^
   1.4  ENTER^
  -0.2  XEQ "AUL"   >>>>    T  = 10.22400044  Gy                                       ---Execution time = 5s---
                                 RDN     P  = 707.0833001  Gy
                                 RDN    R0 =  32.15308639  Gly
                                 RDN    k  =  +1  ( spherical space )

   and    R05 = 4.691072085 Gly  = minimum radius of this Universe
             R07 = 220.3805326 Gly = maximum radius of this Universe.


Example2:    H0 = 68 km/s/Mpc     (Omega)mat = 0.2       (Omega)rad = 0.1

    68  ENTER^
   0.2  ENTER^
   0.1  XEQ "AUL"   >>>>    T  = 10.35899203  Gy                                    
                                 RDN     P  = 0.000000000  Gy   ( no periodic )
                                 RDN    R0 =  17.18654760  Gly
                                 RDN    k  =  -1  ( hyperbolic space )

   and    R05 = 0.000000000 Gly     = minimum radius of this Universe
             R07 = 9.9999999 E99 Gly = maximum radius of this Universe = infinity


Example3:    H0 = 68 km/s/Mpc     (Omega)mat = 1.2       (Omega)rad = -2

    68  ENTER^
   1.2  ENTER^
   -2   XEQ "AUL"   >>>>    T  =   5.732772984  Gy                                    
                                 RDN     P  = 0.000000000  Gy   ( no periodic )
                                 RDN    R0 =  10.71769547  Gly
                                 RDN    k  =  -1  ( hyperbolic space )

   and    R05 = 8.276293007 Gly     = minimum radius of this Universe
             R07 = 9.9999999 E99 Gly = maximum radius of this Universe = infinity


Example4:    H0 = 68 km/s/Mpc     (Omega)mat = 1.2       (Omega)rad = -0.2

    68  ENTER^
   1.2  ENTER^
  -0.2  XEQ "AUL"   >>>>    T  = 10.65133139  Gy                                    
                                 RDN     P  = 0.000000000  Gy   ( no periodic )
                                 RDN    R0 =  14.37929737  Gly
                                 RDN    k  =  0  ( Euclidean space )

   and    R05 = 2.396549562 Gly     = minimum radius of this Universe
             R07 = 9.9999999 E99 Gly = maximum radius of this Universe = infinity


Example5:    H0 = 68 km/s/Mpc     (Omega)mat = 1       (Omega)rad = 0

    68  ENTER^
     1   ENTER^
     0  XEQ "AUL"   >>>>    T  = 9.586198246  Gy                                    
                               RDN     P  = 0.000000000  Gy   ( no periodic )
                               RDN    R0 =  14.37929737  Gly
                               RDN    k  =  0  ( Euclidean space )

   and    R05 = 0.000000000 Gly     = minimum radius of this Universe
             R07 = 9.9999999 E99 Gly = maximum radius of this Universe = infinity


Example6:    H0 = 68 km/s/Mpc     (Omega)mat = 0.4       (Omega)rad = 0.6

    68   ENTER^
   0.4   ENTER^
   0.6  XEQ "AUL"   >>>>   T  = 7.759788005  Gy                                    
                               RDN     P  = 0.000000000  Gy   ( no periodic )
                               RDN    R0 =  14.37929737  Gly
                               RDN    k  =  0  ( Euclidean space )

   and    R05 = 0.000000000 Gly     = minimum radius of this Universe
             R07 = 9.9999999 E99 Gly = maximum radius of this Universe = infinity


Example7:    H0 = 68 km/s/Mpc     (Omega)mat = 0       (Omega)rad = 1

    68   ENTER^
     0   ENTER^
     1  XEQ "AUL"   >>>>   T  = 7.189648685  Gy                                    
                               RDN     P  = 0.000000000  Gy   ( no periodic )
                               RDN    R0 =  14.37929737  Gly
                               RDN    k  =  0  ( Euclidean space )

   and    R05 = 0.000000000 Gly     = minimum radius of this Universe
             R07 = 9.9999999 E99 Gly = maximum radius of this Universe = infinity


Example8:    H0 = 68 km/s/Mpc     (Omega)mat = -0.4       (Omega)rad = 1.4

    68   ENTER^
  -0.4   ENTER^
   1.4  XEQ "AUL"   >>>>   T  = 6.770532843  Gy                                    
                               RDN     P  = 396.9889956  Gy  
                               RDN    R0 =  14.37929737  Gly
                               RDN    k  =  0  ( Euclidean space )

   and    R05 = 0.000000000 Gly     = minimum radius of this Universe
             R07 = 50.32754080 Gly = maximum radius of this Universe


Example9:    H0 = 68 km/s/Mpc     (Omega)mat = -0.4       (Omega)rad = 1.3999

    68     ENTER^
  -0.4     ENTER^
1.3999  XEQ "AUL"   >>>>    T  = 6.770851649  Gy                            
                                   RDN     P  = 397.7828431  Gy
                                   RDN    R0 =  1437.929737  Gly
                                   RDN    k  =  -1  ( hyperbolic space )

   and    R05 =  0                   Gly  = minimum radius of this Universe
             R07 = 5036.805349 Gly = maximum radius of this Universe.


Notes:

-In examples 1 to 8, the  results have a good precision.
-But in example 9, free42 gives   T = 6.770673328 & P = 397.7820792
-When  1 - (Omega)mat - (Omega)rad  is very small but # 0, this program does not give accurate results with an HP41 ( the outputs T & P may even be meaningless... )

-We can make a much smaller program if we only use:  

      (Omega)mat > 0     (Omega)rad > 0   with    1 - (Omega)mat - (Omega)rad  > 0


Data Registers:       R00 to R03: temp
Flags: /
Subroutines: /


 01 LBL "AUL"
 02 STO 02
 03 X<>Y
 04 STO 01          
 05 2
 06 /
 07 +
 08 X<>Y
 09 977.7922214
 10 /
 11 STO 00
 12 SIGN
 13 RCL 01          
 14 -
 15 RCL 02
 16 -
 17 1/X
 18 STO 03
 19 ST* 02
 20 RCL 01
 21 *
 22 2
 23 /
 24 STO 01          
 25 RCL 02
 26 SQRT
 27 STO 02
 28 X<>Y
 29 ST+ Y
 30 1
 31 +
 32 RCL 03
 33 SQRT
 34 STO 03          
 35 +
 36 /
 37 LN
 38 RCL 01
 39 *
 40 RCL 02
 41 -
 42 RCL 03
 43 ST+ Y
 44 STO Z
 45 *
 46 RCL 00          
 47 ST/ Z
 48 /
 49 END

     
         ( 72 bytes / SIZE 004 )

              

           STACK          INPUTS        OUTPUTS
               Z
              H0
             q
               Y        (Omega)mat
             R0
               X        (Omega)rad
             T

  Where    T = age of the Universe     R0 = current scale factor     q = deceleration parameter

Example1:
   H0 = 68 km/s/Mpc     (Omega)mat = 0.2       (Omega)rad = 0.1

    68  ENTER^
   0.2  ENTER^
   0.1  XEQ "AUL"   >>>>    T  = 10.35899205  Gy                                    ---Execution time = 2s---
                                 RDN    R0 = 17.18654761  Gly
                                 RDN     q  = 0.2


Example2:    H0 = 71 km/s/Mpc     (Omega)mat = 0.044       (Omega)rad = 0.000049

     71        ENTER^
   0.044     ENTER^
0.000049  XEQ "AUL"   >>>>    T  = 12.93861786  Gy                                    
                                       RDN    R0 = 14.08543985  Gly
                                       RDN     q  =  0.022049

Note:

-Always k  =  -1  ( hyperbolic space )

      b)  Cosmological Constant # 0


-This program employs Romberg integration.
-The precision is controlled by the display format ( lines 395...398 ).


Data Registers:    R00 = H0    

                                R01 = 
(Omega)mat            R04 = T          R06 = Rmin         R08 = R0      R10.....:  temp
                                R02 = 
(Omega)lambda       R05 = P          R07 = Rmax         R09 = k
                                R03 = (Omega)rad

Flags:  F08-F09-F10-F24
Subroutines: /

-Lines 152-230-239-248 are three-byte GTOs


 01 LBL "AUM"
 02 CF 08
 03 CF 09
 04 CF 10
 05 R^
 06 STO 00
 07 SIGN
 08 X<>Y
 09 STO 03
 10 STO 05
 11 -
 12 X<>Y
 13 STO 02
 14 -
 15 X<>Y
 16 STO 01
 17 STO 04
 18 -
 19 STO 10
 20 RCL 02
 21 ST/ 04
 22 ST/ 05
 23 ST/ 10
 24 CLX
 25 STO 06
 26 STO 07
 27 LBL 01
 28 RCL 06
 29 X^2
 30 STO 11
 31 RCL 10
 32 +
 33 STO 13
 34 RCL 07
 35 ST+ X
 36 STO 12
 37 -
 38 STO 08
 39 RCL 04
 40 RCL 06
 41 RCL 08
 42 *
 43 -
 44 STO 09
 45 *
 46 RCL 07
 47 RCL 13
 48 -
 49 RCL 07
 50 *
 51 RCL 05
 52 +
 53 STO 14
 54 RCL 06
 55 ST+ X
 56 STO 15
 57 *
 58 +
 59 RCL 10
 60 RCL 11
 61 3
 62 *
 63 +
 64 RCL 12
 65 -
 66 ST* 08
 67 RCL 14          
 68 *
 69 RCL 06
 70 RCL 09
 71 *
 72 RCL 12
 73 *
 74 -
 75 RCL 07
 76 RCL 15
 77 X^2
 78 *
 79 RCL 08
 80 +
 81 ST/ Z
 82 /
 83 ST+ 07
 84 X<>Y
 85 ST+ 06
 86 R-P
 87 VIEW X
 88  E-8
 89 X<Y?
 90 GTO 01
 91 RCL 06
 92 CHS
 93 STO 08
 94 X^2
 95 RCL 10
 96 RCL 07
 97 -
 98 +
 99 STO 09
100 RCL 07
101 RCL 06
102 2
103 /
104 CHS
105 ENTER
106 X^2
107 RCL 07
108 -
109 CLD
110 X<0?
111 SF 09
112 X<0?
113 GTO 03
114 SQRT
115 RCL Y
116 SIGN
117 *
118 +
119 X#0?
120 ST/ Y
121 STO 06        
122 X<>Y
123 STO 07
124 LBL 03
125 RCL 09
126 RCL 08
127 2
128 /
129 CHS
130 ENTER
131 X^2
132 RCL 09
133 -
134 X<0?
135 SF 10
136 X<0?
137 GTO 07
138 SQRT
139 RCL Y
140 SIGN
141 *
142 +
143 X#0?
144 ST/ Y
145 STO 08
146 X<>Y
147 STO 09
148 LBL 07
149 FS? 09
150 FC? 10
151 FS? 30
152 GTO 13
153 FC? 09
154 GTO 03
155 RCL 08
156 X<> 06
157 STO 08
158 RCL 09
159 X<> 07
160 STO 09
161 LBL 03
162 FC?C 10
163 FS? 09
164 GTO 04
165 RCL 06
166 RCL 07
167 X>Y?
168 X<>Y
169 RCL 08
170 X>Y?
171 X<>Y
172 RCL 09
173 X>Y?
174 X<>Y
175 STO 06
176 RDN
177 X<Y?
178 X<>Y
179 RDN
180 X>Y?
181 X<>Y
182 STO 07
183 X<> T
184 X<Y?
185 X<>Y
186 STO 09
187 X<>Y
188 STO 08
189 1
190 X>Y?
191 GTO 03
192 RCL 07        
193 X>Y?
194 GTO 11
195 X<> 06
196 X<> 08
197 STO 07
198 GTO 11
199 LBL 03
200 RCL 08
201 X<> 06
202 STO 08
203 RCL 09
204 X<> 07
205 STO 09
206 LBL 11
207 RCL 08
208 RCL 09
209 +
210 CHS
211 X<> 08
212 ST* 09
213 GTO 07
214 LBL 04
215 RCL 06
216 RCL 07
217 X>Y?
218 X<>Y
219 STO 06
220 X<>Y
221 STO 07
222 CF 09
223 LBL 07
224 RCL 06
225 1
226 X<Y?
227 GTO 03
228 RCL 07
229 X<Y?
230 GTO 04
231 RCL 06
232 ST- Z
233 -
234 STO 05
235 /
236 SQRT
237 RAD
238 ASIN
239 GTO 14
240 LBL 03
241 SF 08
242 RCL 06
243 1
244 -
245 SQRT
246 RCL 06
247 SQRT
248 GTO 12
249 LBL 00
250 FS? 08
251 GTO 11
252 FS? 09
253 GTO 09
254 FS? 10
255 GTO 08
256 SIN
257 X^2
258 RCL 05
259 *
260 RCL 06
261 +
262 ENTER
263 STO Z
264 RCL 08
265 +
266 *
267 RCL 09
268 +
269 SQRT
270 /
271 RTN
272 LBL 08
273 X^2
274 RCL 07
275 +
276 ENTER
277 STO Z
278 RCL 08
279 +
280 *
281 RCL 09        
282 +
283 X<>Y
284 RCL 06
285 -
286 *
287 SQRT
288 /
289 RTN
290 LBL 09
291 ENTER
292 ENTER
293 X^2
294 RCL 10
295 +
296 *
297 RCL 04
298 +
299 *
300 RCL 05
301 +
302 SQRT
303 /
304 RTN
305 LBL 11
306 X^2
307 RCL 06
308 X<>Y
309 -
310 ENTER
311 STO Z
312 RCL 08
313 +
314 *
315 RCL 09
316 +
317 RCL 07
318 R^
319 -
320 *
321 SQRT
322 /
323 CHS
324 RTN
325 LBL 10
326 STO 11
327 X<>Y
328 STO 12
329 STO 13
330 XEQ 00
331 STO 15
332 RCL 11
333 ST- 13
334 XEQ 00
335 RCL 15
336 +
337 2
338 /
339 STO 15
340 RCL 13
341 *
342 STO 20
343 20
344 STO 19
345 SIGN
346 STO 17
347 LBL 02
348 RCL 17
349 STO 16
350 RCL 11
351 RCL 13
352 2
353 /
354 +
355 LBL 05
356 STO 14        
357 XEQ 00
358 ST+ 15
359 RCL 13
360 RCL 14
361 +
362 DSE 16
363 GTO 05
364 2
365 ST/ 13
366 ST* 17
367 X^2
368 STO 18
369 RCL 19
370 INT
371  E3
372 /
373 20
374 +
375 STO 19
376 RCL 13
377 RCL 15
378 *
379 LBL 06
380 ENTER
381 ENTER
382 X<> IND 19
383 STO 14
384 -
385 RCL 18
386 4
387 ST* 18
388 SIGN
389 -
390 /
391 +
392 ISG 19
393 GTO 06
394 STO IND 19
395 VIEW X
396 RND
397 RCL 14
398 RND
399 X#Y?
400 GTO 02
401 RCL IND 19
402 RTN
403 LBL 04
404 -
405 SQRT
406 RCL 07
407 CHS
408 X<0?
409 CLX
410 SF 10
411 SQRT
412 XEQ 10
413 ST+ X
414 STO 04
415 RCL 07
416 STO 06
417 CLX
418 STO 07
419 GTO 03
420 LBL 13
421 CLX
422 SIGN
423 LASTX
424 XEQ 10
425 STO 04
426 CLX
427 STO 07        
428 GTO 03
429 LBL 12
430 XEQ 10
431 STO 04
432 CLST
433 RCL 12
434 XEQ 10
435 RCL 06
436 STO 07
437 CLX
438 STO 06
439 X<>Y
440 GTO 12
441 LBL 14
442 RCL 06
443 RCL 06
444 RCL 07
445 -
446 /
447 X<0?
448 CLX
449 SQRT
450 ASIN
451 XEQ 10
452 STO 04
453 SIGN
454 ASIN
455 RCL 12
456 XEQ 10
457 LBL 12
458 RCL 04
459 ST+ 04
460 +
461 4
462 *
463 LBL 03
464 STO 05
465 RCL 06
466 X<0?
467 CLX
468 STO 06
469 DEG
470 90
471 TAN
472 RCL 07
473 X<=0?
474 X<>Y
475 STO 07
476 SF 24
477 RCL 02
478 ABS
479 SQRT
480 ST/ 04
481 ST/ 05
482 RCL 02
483 RCL 10
484 *
485 STO 14
486 CHS
487 X#0?
488 SIGN
489 STO 09
490 977.7922214
491 RCL 00
492 /
493 ST* 04
494 ST* 05
495 RCL 14        
496 ABS
497 SQRT
498 X=0?
499 SIGN
500 /
501 ST* 06
502 ST* 07
503 STO 08
504 RCL 05
505 RCL 04
506 CLD
507 END


         ( 649 bytes / SIZE 020+??? )

 
           STACK           INPUTS        OUTPUTS
               T
            H0 > 0
            k
               Z
        (Omega) mat
            R0
               Y         (Omega)lambda
            P
               X         (Omega)rad
            T

     T = Age of the Universe in Gigayears  ,       k = +1 , 0 , -1  =  spherical, euclidean , hyperbolic Universes
     P = Period of the Universe in Gigayears   ,   R0 = current scale factor of the Universe in Gigalightyears

Example1:      Hubble "constant" = 71 km/s/Mpc  ,  (Omega)mat = -0.1 ,  (Omega)lambda  = -0.003 ,   (Omega) rad  = -0.001 

                      FIX 7

             71    ENTER^
           -0.1   ENTER^
        -0.003  ENTER^
        -0.001  XEQ "AUM"   >>>>   T  = 14.61256011                            ---Execution time = 2m52s---
                                            RDN   P =  796.1455562
                                            RDN   R0 = 13.10701187
                                            RDN    k = -1 ( hyperbolic space )

-We also have in R06 & R07

    Rmin = 1.306380222
    Rmax = 250.8400402

Example2:      Hubble "constant" = 71 km/s/Mpc  ,  (Omega)mat = 0.044 ,  (Omega)lambda  = 0.521 ,   (Omega) rad  = 0.000049   ( our Universe ? ) 

                        FIX 7

              71      ENTER^
           0.044    ENTER^
           0.521    ENTER^
       0.000049  XEQ "AUM"   >>>>   T  = 15.60221461                            ---Execution time = 2m37s---
                                                RDN   P =  0   ( no periodic universe )
                                                RDN   R0 = 20.88180628
                                                RDN    k = -1

    Rmin = 0
    Rmax = 9.999999999 E99 = +infinity

Example3:      Hubble "constant" = 71 km/s/Mpc  ,  (Omega)mat = 0.271 ,  (Omega)lambda  = 0.732 ,   (Omega) rad  = 0.000049   ( our Universe ? ) 

                         FIX 7

              71       ENTER^
           0.271     ENTER^
           0.732     ENTER^
       0.000049   XEQ "AUM"   >>>>   T  = 13.66788682                            ---Execution time = 2m31s---
                                                 RDN   P =  0   ( no periodic universe )
                                                 RDN   R0 = 249.4075046
                                                 RDN    k = +1   ( spherical space )

    Rmin = 0
    Rmax = 9.999999999 E99 = +infinity

Example4:      Hubble "constant" = 71 km/s/Mpc  ,  (Omega)mat = -0.6 ,  (Omega)lambda  = 0.4 ,   (Omega) rad  = 0.1

                       FIX 7

              71     ENTER^
            -0.6    ENTER^
              0.4    ENTER^
              0.1    XEQ "AUM"   >>>>   T  = 17.27498133                             ---Execution time = 4m13s---
                                               RDN   P =  0   ( no periodic universe )
                                               RDN   R0 = 13.13082118
                                               RDN    k = -1   ( hyperbolic space )

    Rmin = 0
    Rmax = 9.999999999 E99 = +infinity

Example5:      Hubble "constant" = 71 km/s/Mpc  ,  (Omega)mat = 0.01 ,  (Omega)lambda  = -0.4 ,   (Omega) rad  = 0.01 

                    FIX 7

            71    ENTER^
          0.01   ENTER^
          -0.4   ENTER^
          0.01   XEQ "AUM"   >>>>   T  = 11.30457454                            ---Execution time = 7m18s---
                                            RDN   P =  66.21313843
                                            RDN   R0 = 11.72326781
                                            RDN    k = -1 ( hyperbolic space )

    Rmin = 0
    Rmax = 21.83998798

Example6:      Hubble "constant" = 71 km/s/Mpc  ,  (Omega)mat = -0.4 ,  (Omega)lambda  = +0.001 ,   (Omega) rad  = 1.4 

                    FIX 7

            71    ENTER^
          -0.4   ENTER^
         0.001  ENTER^
           1.4    XEQ "AUM"   >>>>   T  = 6.484937796                            ---Execution time = 2m32s---
                                            RDN   P =  900.0920456
                                            RDN   R0 = 435.5000702
                                            RDN    k = +1 ( spherical space )

    Rmin = 0
    Rmax = 1889.773264




Notes:

-We use Bairstow method ( LBL 01 ) to factorize the quartic polynomial before finding its roots.
-Romberg method = LBL 10
-FIX 9 would often produce infinite loops

-This program also works if (Omega)mat = 0 , but in this case, "AUR" listed in paragraph 2°) is much faster !
-It doesn't work if (Omega)lambda  = 1 and  (Omega)mat = 0 = (Omega) rad but in this case,  T = +infinity ,  Rmin = 0  and   Rmax = +infinity  ( de Sitter model )
-We could add a few lines to take this case into account...

-We can save several bytes if we use Gauss-Legendre 2-point formula.

-The first estimations are computed with 15 subintervals for T and 55 subintervals for P ( lines 414 & 415 )
-After these results, you can use other values  N1 for T & N2  for P, for instance, with  30 & 100

   30   ENTER^
  100  XEQ 10


Data Registers:    R00 = H0    

                                R01 = 
(Omega)mat            R21 = T          R23 = Rmin         R25 = R0      
                                R02 = 
(Omega)lambda       R22 = P          R24 = Rmax         R26 = k
                                R03 = (Omega)rad

Flags:  F08-F09-F10-F24
Subroutines: /

-Lines 152-230-239-250 are three-byte GTOs


 01 LBL "AUM"
 02 CF 08
 03 CF 09
 04 CF 10
 05 R^
 06 STO 00
 07 SIGN
 08 X<>Y
 09 STO 03
 10 STO 05
 11 -
 12 X<>Y
 13 STO 02
 14 -
 15 X<>Y
 16 STO 01
 17 STO 04
 18 -
 19 STO 10
 20 RCL 02
 21 ST/ 04
 22 ST/ 05
 23 ST/ 10
 24 CLX
 25 STO 06
 26 STO 07
 27 LBL 01
 28 RCL 06
 29 X^2
 30 STO 11
 31 RCL 10
 32 +
 33 STO 13
 34 RCL 07
 35 ST+ X
 36 STO 12
 37 -
 38 STO 08
 39 RCL 04
 40 RCL 06
 41 RCL 08
 42 *
 43 -
 44 STO 09
 45 *
 46 RCL 07
 47 RCL 13
 48 -
 49 RCL 07
 50 *
 51 RCL 05
 52 +
 53 STO 14
 54 RCL 06
 55 ST+ X
 56 STO 15
 57 *
 58 +
 59 RCL 10
 60 RCL 11
 61 3
 62 *
 63 +
 64 RCL 12          
 65 -
 66 ST* 08
 67 RCL 14
 68 *
 69 RCL 06
 70 RCL 09
 71 *
 72 RCL 12
 73 *
 74 -
 75 RCL 07
 76 RCL 15
 77 X^2
 78 *
 79 RCL 08
 80 +
 81 ST/ Z
 82 /
 83 ST+ 07
 84 X<>Y
 85 ST+ 06
 86 R-P
 87 VIEW X
 88  E-8
 89 X<Y?
 90 GTO 01
 91 RCL 06
 92 CHS
 93 STO 08
 94 X^2
 95 RCL 10
 96 RCL 07
 97 -
 98 +
 99 STO 09
100 RCL 07
101 RCL 06
102 2
103 /
104 CHS
105 ENTER
106 X^2
107 RCL 07
108 -
109 CLD
110 X<0?
111 SF 09
112 X<0?
113 GTO 03
114 SQRT
115 RCL Y
116 SIGN
117 *
118 +
119 X#0?
120 ST/ Y
121 STO 06
122 X<>Y
123 STO 07
124 LBL 03
125 RCL 09
126 RCL 08
127 2
128 /
129 CHS
130 ENTER
131 X^2
132 RCL 09        
133 -
134 X<0?
135 SF 10
136 X<0?
137 GTO 07
138 SQRT
139 RCL Y
140 SIGN
141 *
142 +
143 X#0?
144 ST/ Y
145 STO 08
146 X<>Y
147 STO 09
148 LBL 07
149 FS? 09
150 FC? 10
151 FS? 30
152 GTO 13
153 FC? 09
154 GTO 03
155 RCL 08
156 X<> 06
157 STO 08
158 RCL 09
159 X<> 07
160 STO 09
161 LBL 03
162 FC?C 10
163 FS? 09
164 GTO 04
165 RCL 06
166 RCL 07
167 X>Y?
168 X<>Y
169 RCL 08
170 X>Y?
171 X<>Y
172 RCL 09
173 X>Y?
174 X<>Y
175 STO 06
176 RDN
177 X<Y?
178 X<>Y
179 RDN
180 X>Y?
181 X<>Y
182 STO 07
183 X<> T
184 X<Y?
185 X<>Y
186 STO 09
187 X<>Y
188 STO 08
189 1
190 X>Y?
191 GTO 03
192 RCL 07
193 X>Y?
194 GTO 11
195 X<> 06
196 X<> 08
197 STO 07        
198 GTO 11
199 LBL 03
200 RCL 08
201 X<> 06
202 STO 08
203 RCL 09
204 X<> 07
205 STO 09
206 LBL 11
207 RCL 08
208 RCL 09
209 +
210 CHS
211 X<> 08
212 ST* 09
213 GTO 07
214 LBL 04
215 RCL 06
216 RCL 07
217 X>Y?
218 X<>Y
219 STO 06
220 X<>Y
221 STO 07
222 CF 09
223 LBL 07
224 RCL 06
225 1
226 X<Y?
227 GTO 03
228 RCL 07
229 X<Y?
230 GTO 04
231 RCL 06
232 ST- Z
233 -
234 STO 05
235 /
236 SQRT
237 RAD
238 ASIN
239 GTO 14
240 LBL 03
241 SF 08
242 RCL 06
243 SQRT
244 STO 16
245 LASTX
246 1
247 -
248 SQRT
249 STO 17
250 GTO 12
251 LBL 06
252 STO 12
253 -
254 X<>Y
255 STO 11
256 ST+ 11
257 /
258 STO Y
259 3
260 SQRT
261 /
262 STO 13
263 -
264 STO 14        
265 2
266 /
267 ST- 12
268 CLX
269 STO 15
270 LBL 05
271 RCL 12
272 RCL 13
273 X<> 14
274 STO 13
275 +
276 STO 12
277 XEQ 00
278 ST+ 15
279 DSE 11
280 GTO 05
281 RCL 13
282 RCL 15
283 *
284 3
285 SQRT
286 *
287 RTN
288 LBL 00
289 FS? 08
290 GTO 11
291 FS? 09
292 GTO 09
293 FS? 10
294 GTO 08
295 SIN
296 X^2
297 RCL 05
298 *
299 RCL 06
300 +
301 ENTER
302 STO Z
303 RCL 08
304 +
305 *
306 RCL 09
307 +
308 SQRT
309 /
310 RTN
311 LBL 08
312 X^2
313 RCL 07
314 +
315 ENTER
316 STO Z
317 RCL 08
318 +
319 *
320 RCL 09
321 +
322 X<>Y
323 RCL 06
324 -
325 *
326 SQRT
327 /
328 RTN
329 LBL 09
330 ENTER
331 ENTER
332 X^2
333 RCL 10        
334 +
335 *
336 RCL 04
337 +
338 *
339 RCL 05
340 +
341 SQRT
342 ST+ X
343 /
344 RTN
345 LBL 11
346 X^2
347 RCL 06
348 X<>Y
349 -
350 ENTER
351 STO Z
352 RCL 08
353 +
354 *
355 RCL 09
356 +
357 RCL 07
358 R^
359 -
360 *
361 SQRT
362 /
363 CHS
364 RTN
365 LBL 04
366 -
367 SQRT
368 STO 17
369 RCL 07
370 STO 23
371 CHS
372 X<0?
373 CLX
374 SF 10
375 SQRT
376 STO 16
377 CLX
378 STO 24
379 GTO 04
380 LBL 13
381 CLX
382 STO 16
383 STO 24
384 SIGN
385 STO 17
386 RCL 06
387 STO 23
388 GTO 04
389 LBL 12
390 CLX
391 STO 18
392 STO 23
393 RCL 06
394 STO 24
395 GTO 04
396 LBL 14
397 STO 17
398 RCL 06
399 RCL 06
400 STO 23        
401 RCL 07
402 STO 24
403 -
404 /
405 X<0?
406 CLX
407 SQRT
408 ASIN
409 STO 16
410 SIGN
411 ASIN
412 STO 18
413 LBL 04
414 15
415 55
416 LBL 10
417 RAD
418 STO 20
419 X<>Y
420 STO 19
421 RCL 17
422 RCL 16
423 XEQ 06
424 STO 21
425 CLX
426 FC? 09
427 FS? 10
428 GTO 00
429 RCL 20
430 RCL 18
431 RCL 16
432 XEQ 06
433 ST+ X
434 LBL 00
435 STO 22
436 RCL 23
437 X<0?
438 CLX
439 STO 23
440 DEG
441 90
442 TAN
443 RCL 24
444 X<=0?
445 X<>Y
446 STO 24
447 SF 24
448 RCL 02
449 ABS
450 SQRT
451 ST/ 21
452 ST/ 22
453 RCL 02
454 RCL 10
455 *
456 STO 14
457 CHS
458 X#0?
459 SIGN
460 STO 26
461 977.7922214
462 RCL 00
463 /
464 ST* 21
465 ST* 22
466 RCL 14        
467 ABS
468 SQRT
469 X=0?
470 SIGN
471 /
472 ST* 23
473 ST* 24
474 STO 25
475 RCL 22
476 RCL 21
477 CF 24
478 END


         ( 631 bytes / SIZE 027 )

 
           STACK           INPUTS        OUTPUTS
               T
            H0 > 0
            k
               Z
        (Omega) mat
            R0
               Y         (Omega)lambda
            P
               X         (Omega)rad
            T

     T = Age of the Universe in Gigayears  ,       k = +1 , 0 , -1  =  spherical, euclidean , hyperbolic Universes
     P = Period of the Universe in Gigayears   ,   R0 = current scale factor of the Universe in Gigalightyears

Example1:      Hubble "constant" = 71 km/s/Mpc  ,  (Omega)mat = -0.1 ,  (Omega)lambda  = -0.003 ,   (Omega) rad  = -0.001 


             71    ENTER^
           -0.1   ENTER^
        -0.003  ENTER^
        -0.001  XEQ "AUM"   >>>>   T  = 14.61256012                            ---Execution time = 3m33s---
                                            RDN   P =  796.1455560
                                            RDN   R0 = 13.10701187
                                            RDN    k = -1 ( hyperbolic space )

-We also have in R23 & R24

    Rmin = 1.306380222
    Rmax = 250.8400402

Example2:      Hubble "constant" = 71 km/s/Mpc  ,  (Omega)mat = 0.044 ,  (Omega)lambda  = 0.521 ,   (Omega) rad  = 0.000049   ( our Universe ? ) 


              71      ENTER^
           0.044    ENTER^
           0.521    ENTER^
       0.000049  XEQ "AUM"   >>>>   T  = 15.60221156                            ---Execution time = 48s---
                                                RDN   P =  0   ( no periodic universe )
                                                RDN   R0 = 20.88180628
                                                RDN    k = -1

    Rmin = 0
    Rmax = 9.999999999 E99 = +infinity

-With another number of subintervals to compute T - for instance 30 -

   30   ENTER^  ( any other number in X-register since P = 0 )   XEQ 10  >>>>   T= 15.60221435

-With N = 60 , it yields:   T = 15.60221460

 ... and so on ...

Note:

-Of course, we can use Carlson elliptic integrals to obtain very accurate results:


Data Registers:    R00 = H0    

                                R01 = 
(Omega)mat            R04 = T          R06 = Rmin         R08 = R0      R10.....:  temp
                                R02 = 
(Omega)lambda       R05 = P          R07 = Rmax         R09 = k
                                R03 = (Omega)rad

Flags:  F08-F09-F10-F24
Subroutines:  ( M-Code routines )  RF  RFZ  RJ  RJZ  ( cf "Carlson Elliptic Integrals for the HP41" )

-Lines 156-334 are three-byte GTOs


 01 LBL "AUM"
 02 CF 08
 03 CF 09
 04 CF 10
 05 R^
 06 STO 00
 07 SIGN
 08 X<>Y
 09 STO 03
 10 STO 05
 11 -
 12 X<>Y
 13 STO 02
 14 -
 15 X<>Y
 16 STO 01
 17 STO 04
 18 -
 19 STO 10
 20 STO 20
 21 RCL 02
 22 ST/ 04
 23 ST/ 05
 24 ST/ 10
 25 CLX
 26 STO 06
 27 STO 07
 28 LBL 01
 29 RCL 06
 30 X^2
 31 STO 11
 32 RCL 10
 33 +
 34 STO 13
 35 RCL 07
 36 ST+ X
 37 STO 12
 38 -
 39 STO 08
 40 RCL 04
 41 RCL 06
 42 RCL 08
 43 *
 44 -
 45 STO 09
 46 *
 47 RCL 07
 48 RCL 13
 49 -
 50 RCL 07
 51 *
 52 RCL 05
 53 +
 54 STO 14
 55 RCL 06
 56 ST+ X
 57 STO 15
 58 *
 59 +
 60 RCL 10
 61 RCL 11
 62 3
 63 *
 64 +
 65 RCL 12
 66 -
 67 ST* 08
 68 RCL 14
 69 *
 70 RCL 06
 71 RCL 09
 72 *
 73 RCL 12
 74 *
 75 -
 76 RCL 07
 77 RCL 15
 78 X^2
 79 *
 80 RCL 08
 81 +
 82 ST/ Z
 83 /
 84 ST+ 07
 85 X<>Y
 86 ST+ 06
 87 R-P
 88 VIEW X
 89  E-8
 90 X<Y?
 91 GTO 01
 92 RCL 06
 93 CHS
 94 STO 08
 95 X^2
 96 RCL 10
 97 RCL 07
 98 -
 99 +
100 STO 09
101 RCL 07
102 RCL 06          
103 2
104 /
105 CHS
106 STO 11
107 ENTER
108 X^2
109 RCL 07
110 -
111 STO 12
112 CLD
113 X<0?
114 SF 09
115 X<0?
116 GTO 02
117 SQRT
118 RCL Y
119 SIGN
120 *
121 +
122 X#0?
123 ST/ Y
124 STO 06
125 X<>Y
126 STO 07
127 LBL 02
128 RCL 09
129 RCL 08
130 2
131 /
132 CHS
133 STO 13
134 ENTER
135 X^2
136 RCL 09
137 -
138 STO 14
139 X<0?
140 SF 10
141 X<0?
142 GTO 02
143 SQRT
144 RCL Y
145 SIGN
146 *
147 +
148 X#0?
149 ST/ Y
150 STO 08
151 X<>Y
152 STO 09
153 LBL 02
154 FS? 09
155 FC? 10
156 GTO 03
157 1
158 RCL 07
159 +
160 RCL 06
161 +
162 SQRT
163 STO 11
164 SIGN
165 RCL 06
166 -
167 RCL 09
168 +
169 SQRT
170 STO 12
171 RCL 09
172 SQRT
173 +
174 X^2
175 1
176 -
177 SQRT
178 RCL 11
179 RCL 07
180 SQRT
181 +
182 X^2
183 1
184 -
185 SQRT
186 *
187 X^2
188 STO 15
189 RCL 07
190 RCL 09
191 +
192 ST+ X
193 RCL 06
194 X^2
195 STO 14
196 +
197 STO 16
198 X^2
199 RCL 07
200 4
201 *
202 RCL 14
203 -
204 STO 08
205 RCL 09
206 4
207 *
208 RCL 14
209 -
210 *
211 -
212 SQRT
213 ENTER
214 CHS
215 RCL 15
216 RCL 16
217 +
218 STO 18
219 ST+ Z
220 +
221 STO 14          
222 X<>Y
223 STO 13
224 RCL 15
225 RCL 08
226 +
227 STO 17
228 X<> Z
229 RCL 15
230 RJ
231 3
232 /
233 X<> 18
234 2
235 /
236 RCL 11
237 RCL 09
238 SQRT
239 *
240 RCL 12
241 RCL 07
242 SQRT
243 *
244 +
245 STO 19
246 X^2
247 -
248 STO 16
249 RCL 17
250 *
251 RCL 19
252 /
253 RCL 08
254 ST+ X
255 RCL 19
256 *
257 +
258 X^2
259 RCL 16
260 RCL 19
261 /
262 X^2
263 RCL 08
264 +
265 RCL 17
266 X^2
267 *
268 ENTER
269 RF
270 2
271 /
272 RCL 18
273 +
274 RCL 08
275 *
276 RCL 06
277 ST+ X
278 *
279 STO 18
280 RCL 16
281 RCL 11
282 ST/ 12
283 RCL 07
284 SQRT
285 *
286 STO 11
287 /
288 2
289 +
290 X^2
291 RCL 16
292 X^2
293 RCL 19
294 X^2
295 RCL 08
296 *
297 +
298 RCL 11
299 X^2
300 /
301 ENTER
302 RF
303 RCL 06
304 2
305 +
306 RCL 12
307 *
308 RCL 09
309 RCL 07
310 /
311 SQRT
312 RCL 06
313 *
314 +
315 *
316 2
317 /
318 ST+ 18
319 RCL 15
320 RCL 14
321 RCL 13
322 RF
323 RCL 06
324 *
325 CHS
326 RCL 18        
327 +
328 ST+ X
329 STO 04
330 CLX
331 STO 05
332 STO 06
333 STO 07
334 GTO 07
335 LBL 05
336 STO 11
337 RCL 06
338 -
339 RCL 11
340 RCL 07
341 -
342 FS? 09
343 CHS
344 STO 13
345 *
346 X<>Y
347 STO 12
348 RCL 08
349 -
350 FS? 10
351 GTO 00
352 *
353 RCL 12
354 RCL 09
355 -
356 *
357 SQRT
358 RCL 12
359 RCL 06
360 -
361 RCL 12
362 RCL 07
363 -
364 FS? 09
365 CHS
366 STO 14
367 *
368 RCL 11
369 RCL 08
370 -
371 *
372 RCL 11
373 RCL 09
374 -
375 *
376 SQRT
377 +
378 X^2
379 STO 15
380 STO 18
381 RCL 11
382 RCL 06
383 -
384 RCL 11
385 RCL 08
386 -
387 *
388 RCL 14
389 *
390 RCL 12
391 RCL 09
392 -
393 *
394 SQRT
395 RCL 12
396 RCL 06
397 -
398 RCL 12
399 RCL 08
400 -
401 *
402 RCL 13
403 *
404 RCL 11
405 RCL 09
406 -
407 *
408 SQRT
409 +
410 X^2
411 STO 16
412 RCL 11
413 RCL 06
414 -
415 RCL 11
416 RCL 09
417 -
418 *
419 RCL 14
420 *
421 RCL 12
422 RCL 08
423 -
424 *
425 SQRT
426 RCL 12
427 RCL 06
428 -
429 RCL 12
430 RCL 09
431 -
432 *
433 RCL 13
434 *
435 RCL 11        
436 RCL 08
437 -
438 *
439 SQRT
440 +
441 X^2
442 STO 17
443 GTO 02
444 LBL 00
445 X^2
446 RCL 09
447 X^2
448 +
449 *
450 SQRT
451 RCL 12
452 RCL 06
453 -
454 RCL 12
455 RCL 07
456 -
457 FS? 09
458 CHS
459 STO 14
460 *
461 RCL 11
462 RCL 08
463 -
464 X^2
465 RCL 09
466 X^2
467 +
468 *
469 SQRT
470 +
471 X^2
472 STO 15
473 STO 18
474 RCL 11
475 RCL 12
476 -
477 RCL 09
478 *
479 STO 05
480 RCL 11
481 RCL 08
482 -
483 RCL 12
484 LASTX
485 -
486 *
487 RCL 09
488 X^2
489 +
490 STO 19
491 RCL 11
492 RCL 06
493 -
494 RCL 14
495 *
496 ST* Z
497 *
498 R-P
499 X<>Y
500 2
501 /
502 X<>Y
503 SQRT
504 P-R
505 STO 16
506 X<>Y
507 STO 17
508 RCL 05
509 CHS
510 RCL 19
511 RCL 12
512 RCL 06
513 -
514 RCL 13
515 *
516 ST* Z
517 *
518 R-P
519 X<>Y
520 2
521 /
522 X<>Y
523 SQRT
524 P-R
525 X<>Y
526 RCL 17
527 +
528 X<>Y
529 RCL 16
530 +
531 R-P
532 X<>Y
533 ST+ X
534 X<>Y
535 X^2
536 P-R
537 STO 16
538 X<>Y
539 STO 17
540 LBL 02
541 RCL 11
542 RCL 12
543 -
544 X^2
545 ST/ 15
546 ST/ 16
547 ST/ 17
548 ST/ 18
549 RCL 09
550 RCL 06        
551 -
552 RCL 08
553 LASTX
554 -
555 ST* Y
556 X^2
557 RCL 09
558 X^2
559 +
560 FC? 10
561 X<>Y
562 STO 19
563 FC? 09
564 CHS
565 ST+ 18
566 RCL 18
567 RCL 17
568 RCL 16
569 RCL 15
570 FS? 10
571 RJZ
572 FC? 10
573 RJ
574 RCL 06
575 RCL 07
576 -
577 FC? 08
578 FS? 09
579 CHS
580 *
581 3
582 /
583 ST* 19
584 RCL 18
585 RCL 11
586 RCL 06
587 -
588 RCL 12
589 LASTX
590 -
591 *
592 X=0?
593 GTO 00
594 /
595 ENTER
596 STO Z
597 1
598 FS? 09
599 CHS
600 +
601 XEQ "RF"
602 LBL 00
603 ST+ 19
604 RCL 17
605 RCL 16
606 RCL 15
607 FS? 10
608 RFZ
609 FC? 10
610 RF
611 RCL 06
612 *
613 RCL 19
614 FS? 08
615 CHS
616 +
617 ST+ X
618 RTN
619 LBL 03
620 FC? 09
621 GTO 03
622 RCL 08
623 STO 06
624 RCL 09
625 STO 07
626 RCL 11
627 STO 08
628 RCL 12
629 CHS
630 SQRT
631 STO 09
632 LBL 03
633 FC? 10
634 GTO 03
635 RCL 13
636 STO 08
637 RCL 14
638 CHS
639 SQRT
640 STO 09
641 LBL 03
642 FC? 10
643 FS? 09
644 GTO 04
645 RCL 06
646 RCL 07
647 X>Y?
648 X<>Y
649 RCL 08
650 X>Y?
651 X<>Y
652 RCL 09
653 X>Y?
654 X<>Y
655 STO 06
656 RDN
657 X<Y?
658 X<>Y
659 RDN
660 X>Y?
661 X<>Y
662 STO 07
663 X<> T
664 X<Y?
665 X<>Y
666 STO 09        
667 X<>Y
668 STO 08
669 1
670 X>Y?
671 GTO 02
672 RCL 07
673 X>Y?
674 GTO 03
675 X<> 06
676 X<> 08
677 STO 07
678 GTO 03
679 LBL 02
680 RCL 08
681 X<> 06
682 STO 08
683 RCL 09
684 X<> 07
685 STO 09
686 GTO 03
687 LBL 04
688 RCL 06
689 RCL 07
690 X>Y?
691 X<>Y
692 STO 06
693 X<>Y
694 STO 07
695 LBL 03
696 FS?C 09
697 SF 10
698 RCL 06
699 1
700 X<Y?
701 GTO 02
702 RCL 07
703 X>Y?
704 GTO 06
705 CF 09
706 RCL 06
707 X<> 07
708 STO 06
709 X<0?
710 0
711 ENTER
712 SIGN
713 XEQ 05
714 STO 04
715 CLX
716 STO 05
717 STO 07
718 GTO 07
719 LBL 02
720 CF 09
721 SF 08
722 CLST
723 SIGN
724 XEQ 05
725 STO 04
726 CLST
727 RCL 06
728 XEQ 05
729 ST+ X
730 STO 05
731 CLX
732 X<> 06
733 STO 07
734 GTO 07
735 LBL 06
736 SF 09
737 RCL 06
738 X<0?
739 0
740 ENTER
741 SIGN
742 XEQ 05
743 STO 04
744 RCL 06
745 X<0?
746 0
747 RCL 07
748 XEQ 05
749 ST+ X
750 STO 05
751 LBL 07
752 RCL 06
753 X<0?
754 CLX
755 STO 06
756 DEG
757 90
758 TAN
759 RCL 07
760 X<=0?
761 X<>Y
762 STO 07
763 SF 24
764 RCL 02
765 ABS
766 SQRT
767 ST/ 04
768 ST/ 05
769 RCL 20
770 CHS
771 X#0?
772 SIGN
773 STO 09
774 977.7922214
775 RCL 00        
776 /
777 ST* 04
778 ST* 05
779 RCL 20
780 ABS
781 SQRT
782 X=0?
783 SIGN
784 /
785 ST* 06
786 ST* 07
787 STO 08
788 RCL 05
789 RCL 04
790 CF 24
791 END


         ( 977 bytes / SIZE 021 )

 
           STACK           INPUTS        OUTPUTS
               T
            H0 > 0
            k
               Z
        (Omega) mat
            R0
               Y         (Omega)lambda
            P
               X         (Omega)rad
            T

     T = Age of the Universe in Gigayears  ,       k = +1 , 0 , -1  =  spherical, euclidean , hyperbolic Universes
     P = Period of the Universe in Gigayears   ,   R0 = current scale factor of the Universe in Gigalightyears

Example1:      Hubble "constant" = 71 km/s/Mpc  ,  (Omega)mat = -0.1 ,  (Omega)lambda  = -0.003 ,   (Omega) rad  = -0.001 


             71    ENTER^
           -0.1   ENTER^
        -0.003  ENTER^
        -0.001  XEQ "AUM"   >>>>   T  = 14.61256009                           ---Execution time = 59s---
                                            RDN   P =  796.1455565
                                            RDN   R0 = 13.10701187
                                            RDN    k = -1 ( hyperbolic space )

-We also have in R06 & R07

    Rmin = 1.306380222
    Rmax = 250.8400402

Example2:      Hubble "constant" = 71 km/s/Mpc  ,  (Omega)mat = 0.044 ,  (Omega)lambda  = 0.521 ,   (Omega) rad  = 0.000049   ( our Universe ? ) 


              71      ENTER^
           0.044    ENTER^
           0.521    ENTER^
       0.000049  XEQ "AUM"   >>>>   T  = 15.60221452                            ---Execution time = 48s---
                                                RDN   P =  0   ( no periodic universe )
                                                RDN   R0 = 20.88180628
                                                RDN    k = -1

    Rmin = 0
    Rmax = 9.999999999 E99 = +infinity

Example3:      Hubble "constant" = 71 km/s/Mpc  ,  (Omega)mat = 0.271 ,  (Omega)lambda  = 0.732 ,   (Omega) rad  = 0.000049   ( our Universe ? ) 


              71       ENTER^
           0.271     ENTER^
           0.732     ENTER^
       0.000049   XEQ "AUM"   >>>>   T  = 13.66788681                            ---Execution time = 102s---
                                                 RDN   P =  0   ( no periodic universe )
                                                 RDN   R0 = 249.4075046
                                                 RDN    k = +1   ( spherical space )

    Rmin = 0
    Rmax = 9.999999999 E99 = +infinity

Example4:      Hubble "constant" = 71 km/s/Mpc  ,  (Omega)mat = -0.6 ,  (Omega)lambda  = 0.4 ,   (Omega) rad  = 0.1


              71     ENTER^
            -0.6    ENTER^
              0.4    ENTER^
              0.1    XEQ "AUM"   >>>>   T  = 17.27498133                             ---Execution time = 44s---
                                               RDN   P =  0   ( no periodic universe )
                                               RDN   R0 = 13.13082118
                                               RDN    k = -1   ( hyperbolic space )

    Rmin = 0
    Rmax = 9.999999999 E99 = +infinity

Example5:      Hubble "constant" = 71 km/s/Mpc  ,  (Omega)mat = 0.01 ,  (Omega)lambda  = -0.4 ,   (Omega) rad  = 0.01 


            71    ENTER^
          0.01   ENTER^
          -0.4   ENTER^
          0.01   XEQ "AUM"   >>>>   T  = 11.30457447                            ---Execution time = 90s---
                                            RDN   P =  66.21313811
                                            RDN   R0 = 11.72326781
                                            RDN    k = -1 ( hyperbolic space )

    Rmin = 0
    Rmax = 21.83998798

Example6:      Hubble "constant" = 71 km/s/Mpc  ,  (Omega)mat = -0.4 ,  (Omega)lambda  = +0.001 ,   (Omega) rad  = 1.4 


            71    ENTER^
          -0.4   ENTER^
         0.001  ENTER^
           1.4    XEQ "AUM"   >>>>   T  = 6.484937782                           ---Execution time = 114s---
                                            RDN   P =  900.0920456
                                            RDN   R0 = 435.5000702
                                            RDN    k = +1 ( spherical space )

    Rmin = 0
    Rmax = 1889.773264


Notes:

-Lines 157 to 330 employ the formula given in reference [3]
-Lines 335 to 618 use the formula given in reference [1]         ( for elliptic integrals [-1,-1,-1;-1,+2] )

-But the integral:     §ar x dx / sqrt [(x-a) (b-x) (x+c) (x+d) ]      may be re-written:

     §0+infinity  sqrt(r-a)  (a.y+r)  /  (y+1) / sqrt [ ((b-a).y + (b-r)) ((a+c).y + (r+c)) ((a+d).y + (r+d)) ]  dy    by the change of variable:    x = (a.y+r) / (y+1)

-This reduces the number of bytes ( 758 instead of 976 ) , but the execution times will be a little increased... ( except in examples 1 & 4 )



Data Registers:    R00 = H0    

                                R01 = 
(Omega)mat            R04 = T          R06 = Rmin         R08 = R0      R10.....:  temp
                                R02 = 
(Omega)lambda       R05 = P          R07 = Rmax         R09 = k
                                R03 = (Omega)rad

Flags:  F09-F10-F24
Subroutines:  ( M-Code routines )  RF  RFZ  RJ  RJZ  ( cf "Carlson Elliptic Integrals for the HP41" )

-Lines 155-333 are three-byte GTOs


 01 LBL "AUM"
 02 CF 09
 03 CF 10
 04 R^
 05 STO 00
 06 SIGN
 07 X<>Y
 08 STO 03
 09 STO 05
 10 -
 11 X<>Y
 12 STO 02
 13 -
 14 X<>Y
 15 STO 01
 16 STO 04
 17 -
 18 STO 10
 19 STO 20
 20 RCL 02
 21 ST/ 04
 22 ST/ 05
 23 ST/ 10
 24 CLX
 25 STO 06
 26 STO 07
 27 LBL 01
 28 RCL 06
 29 X^2
 30 STO 11
 31 RCL 10
 32 +
 33 STO 13
 34 RCL 07
 35 ST+ X
 36 STO 12
 37 -
 38 STO 08
 39 RCL 04
 40 RCL 06
 41 RCL 08
 42 *
 43 -
 44 STO 09
 45 *
 46 RCL 07
 47 RCL 13
 48 -
 49 RCL 07
 50 *
 51 RCL 05
 52 +
 53 STO 14
 54 RCL 06
 55 ST+ X
 56 STO 15
 57 *
 58 +
 59 RCL 10
 60 RCL 11
 61 3
 62 *
 63 +
 64 RCL 12
 65 -
 66 ST* 08
 67 RCL 14
 68 *
 69 RCL 06
 70 RCL 09          
 71 *
 72 RCL 12
 73 *
 74 -
 75 RCL 07
 76 RCL 15
 77 X^2
 78 *
 79 RCL 08
 80 +
 81 ST/ Z
 82 /
 83 ST+ 07
 84 X<>Y
 85 ST+ 06
 86 R-P
 87 VIEW X
 88  E-8
 89 X<Y?
 90 GTO 01
 91 RCL 06
 92 CHS
 93 STO 08
 94 X^2
 95 RCL 10
 96 RCL 07
 97 -
 98 +
 99 STO 09
100 RCL 07
101 RCL 06
102 2
103 /
104 CHS
105 STO 11
106 ENTER
107 X^2
108 RCL 07
109 -
110 STO 12
111 CLD
112 X<0?
113 SF 09
114 X<0?
115 GTO 02
116 SQRT
117 RCL Y
118 SIGN
119 *
120 +
121 X#0?
122 ST/ Y
123 STO 06
124 X<>Y
125 STO 07
126 LBL 02
127 RCL 09
128 RCL 08
129 2
130 /
131 CHS
132 STO 13
133 ENTER
134 X^2
135 RCL 09
136 -
137 STO 14
138 X<0?
139 SF 10
140 X<0?
141 GTO 02
142 SQRT
143 RCL Y
144 SIGN
145 *
146 +
147 X#0?
148 ST/ Y
149 STO 08
150 X<>Y
151 STO 09
152 LBL 02
153 FS? 09
154 FC? 10
155 GTO 03
156 1
157 RCL 07
158 +
159 RCL 06
160 +
161 SQRT
162 STO 11
163 SIGN
164 RCL 06
165 -
166 RCL 09
167 +
168 SQRT
169 STO 12
170 RCL 09        
171 SQRT
172 +
173 X^2
174 1
175 -
176 SQRT
177 RCL 11
178 RCL 07
179 SQRT
180 +
181 X^2
182 1
183 -
184 SQRT
185 *
186 X^2
187 STO 15
188 RCL 07
189 RCL 09
190 +
191 ST+ X
192 RCL 06
193 X^2
194 STO 14
195 +
196 STO 16
197 X^2
198 RCL 07
199 4
200 *
201 RCL 14
202 -
203 STO 08
204 RCL 09
205 4
206 *
207 RCL 14
208 -
209 *
210 -
211 SQRT
212 ENTER
213 CHS
214 RCL 15
215 RCL 16
216 +
217 STO 18
218 ST+ Z
219 +
220 STO 14
221 X<>Y
222 STO 13
223 RCL 15
224 RCL 08
225 +
226 STO 17
227 X<> Z
228 RCL 15
229 RJ
230 3
231 /
232 X<> 18
233 2
234 /
235 RCL 11
236 RCL 09
237 SQRT
238 *
239 RCL 12
240 RCL 07
241 SQRT
242 *
243 +
244 STO 19
245 X^2
246 -
247 STO 16
248 RCL 17
249 *
250 RCL 19
251 /
252 RCL 08
253 ST+ X
254 RCL 19
255 *
256 +
257 X^2
258 RCL 16         
259 RCL 19
260 /
261 X^2
262 RCL 08
263 +
264 RCL 17
265 X^2
266 *
267 ENTER
268 RF
269 2
270 /
271 RCL 18
272 +
273 RCL 08
274 *
275 RCL 06
276 ST+ X
277 *
278 STO 18
279 RCL 16
280 RCL 11
281 ST/ 12
282 RCL 07
283 SQRT
284 *
285 STO 11
286 /
287 2
288 +
289 X^2
290 RCL 16
291 X^2
292 RCL 19
293 X^2
294 RCL 08
295 *
296 +
297 RCL 11
298 X^2
299 /
300 ENTER
301 RF
302 RCL 06
303 2
304 +
305 RCL 12
306 *
307 RCL 09
308 RCL 07
309 /
310 SQRT
311 RCL 06
312 *
313 +
314 *
315 2
316 /
317 ST+ 18
318 RCL 15
319 RCL 14
320 RCL 13
321 RF
322 RCL 06
323 *
324 CHS
325 RCL 18
326 +
327 ST+ X
328 STO 04
329 CLX
330 STO 05
331 STO 06        
332 STO 07
333 GTO 07
334 LBL 05
335 STO 11
336 RCL 06
337 RCL 08
338 -
339 RCL 06
340 RCL 09
341 -
342 *
343 STO 15
344 RCL 07
345 RCL 11
346 -
347 RCL 07
348 RCL 06
349 -
350 /
351 STO 12
352 RCL 11
353 RCL 08
354 -
355 RCL 06
356 LASTX
357 -
358 FS? 10
359 GTO 00
360 /
361 STO 13
362 RCL 11
363 RCL 09
364 -
365 RCL 06
366 LASTX
367 -
368 /
369 STO 14
370 GTO 02
371 LBL 00
372 *
373 RCL 09
374 X^2
375 +
376 STO 13
377 RCL 11
378 RCL 06
379 ST- Y
380 RCL 08
381 -
382 X^2
383 RCL 09
384 ST* Z
385 X^2
386 +
387 STO 15
388 ST/ 13
389 /
390 STO 14
391 LBL 02
392 1
393 X<>Y
394 RCL 13
395 RCL 12
396 FS? 10
397 RJZ
398 FC? 10
399 RJ
400 RCL 11
401 RCL 06
402 -
403 *
404 3
405 /
406 X<> 14
407 RCL 13
408 RCL 12
409 FS? 10
410 RFZ
411 FC? 10
412 RF
413 RCL 06
414 *
415 RCL 14        
416 +
417 RCL 15
418 SQRT
419 /
420 RCL 11
421 RCL 06
422 -
423 RCL 07
424 LASTX
425 -
426 /
427 ABS
428 SQRT
429 ST+ X
430 *
431 RTN
432 LBL 03
433 FC? 09
434 GTO 03
435 RCL 08
436 STO 06
437 RCL 09
438 STO 07
439 RCL 11
440 STO 08
441 RCL 12
442 CHS
443 SQRT
444 STO 09
445 LBL 03
446 FC? 10
447 GTO 03
448 RCL 13
449 STO 08
450 RCL 14
451 CHS
452 SQRT
453 STO 09
454 LBL 03
455 FC? 10
456 FS? 09
457 GTO 04
458 RCL 06
459 RCL 07
460 X>Y?
461 X<>Y
462 RCL 08
463 X>Y?
464 X<>Y
465 RCL 09
466 X>Y?
467 X<>Y
468 STO 06
469 RDN
470 X<Y?
471 X<>Y
472 RDN
473 X>Y?
474 X<>Y
475 STO 07
476 X<> T
477 X<Y?
478 X<>Y
479 STO 09
480 X<>Y
481 STO 08
482 1
483 X>Y?
484 GTO 02
485 RCL 07
486 X>Y?
487 GTO 03
488 X<> 06
489 X<> 08
490 STO 07
491 GTO 03
492 LBL 02
493 RCL 08
494 X<> 06
495 STO 08
496 RCL 09
497 X<> 07
498 STO 09
499 GTO 03
500 LBL 04
501 RCL 06
502 RCL 07        
503 X>Y?
504 X<>Y
505 STO 06
506 X<>Y
507 STO 07
508 LBL 03
509 FS?C 09
510 SF 10
511 RCL 06
512 1
513 X<Y?
514 GTO 02
515 RCL 07
516 X>Y?
517 GTO 06
518 X<> 06
519 STO 07
520 1
521 XEQ 05
522 STO 04
523 RCL 06
524 CHS
525 X<=0?
526 GTO 00
527 CLX
528 XEQ 05
529 ST- 04
530 LBL 00
531 CLX
532 STO 05
533 STO 07
534 GTO 07
535 LBL 02
536 CLX
537 XEQ 05
538 STO 04
539 ST+ X
540 STO 05
541 1
542 XEQ 05
543 ST- 04
544 CLX
545 X<> 06
546 STO 07
547 GTO 07
548 LBL 06
549 1
550 XEQ 05
551 STO 04
552 RCL 07
553 XEQ 05
554 ST+ X
555 STO 05
556 RCL 06
557 CHS
558 X<=0?
559 GTO 07
560 CLX
561 XEQ 05
562 ST- 04
563 ST+ X
564 ST- 05
565 LBL 07
566 RCL 06
567 X<0?
568 CLX
569 STO 06
570 DEG
571 90
572 TAN
573 RCL 07
574 X<=0?
575 X<>Y
576 STO 07
577 SF 24
578 RCL 02
579 ABS
580 SQRT
581 ST/ 04
582 ST/ 05
583 RCL 20
584 CHS
585 X#0?
586 SIGN
587 STO 09
588 977.7922214
589 RCL 00        
590 /
591 ST* 04
592 ST* 05
593 RCL 20
594 ABS
595 SQRT
596 X=0?
597 SIGN
598 /
599 ST* 06
600 ST* 07
601 STO 08
602 RCL 05
603 RCL 04
604 CF 24
605 END


         ( 758 bytes / SIZE 021 )

 
           STACK           INPUTS        OUTPUTS
               T
            H0 > 0
            k
               Z
        (Omega) mat
            R0
               Y         (Omega)lambda
            P
               X         (Omega)rad
            T

     T = Age of the Universe in Gigayears  ,       k = +1 , 0 , -1  =  spherical, euclidean , hyperbolic Universes
     P = Period of the Universe in Gigayears   ,   R0 = current scale factor of the Universe in Gigalightyears

Example1:      Hubble "constant" = 71 km/s/Mpc  ,  (Omega)mat = -0.1 ,  (Omega)lambda  = -0.003 ,   (Omega) rad  = -0.001 


             71    ENTER^
           -0.1   ENTER^
        -0.003  ENTER^
        -0.001  XEQ "AUM"   >>>>   T  = 14.61256009                           ---Execution time = 53s---
                                            RDN   P =  796.1455560
                                            RDN   R0 = 13.10701187
                                            RDN    k = -1 ( hyperbolic space )

-We also have in R06 & R07

    Rmin = 1.306380222
    Rmax = 250.8400402

Example2:      Hubble "constant" = 71 km/s/Mpc  ,  (Omega)mat = 0.044 ,  (Omega)lambda  = 0.521 ,   (Omega) rad  = 0.000049   ( our Universe ? ) 


              71      ENTER^
           0.044    ENTER^
           0.521    ENTER^
       0.000049  XEQ "AUM"   >>>>   T  = 15.60221457                            ---Execution time = 54s---
                                                RDN   P =  0   ( no periodic universe )
                                                RDN   R0 = 20.88180628
                                                RDN    k = -1

    Rmin = 0
    Rmax = 9.999999999 E99 = +infinity

Example3:      Hubble "constant" = 71 km/s/Mpc  ,  (Omega)mat = 0.271 ,  (Omega)lambda  = 0.732 ,   (Omega) rad  = 0.000049   ( our Universe ? ) 


              71       ENTER^
           0.271     ENTER^
           0.732     ENTER^
       0.000049   XEQ "AUM"   >>>>   T  = 13.66788679                            ---Execution time = 106s---
                                                 RDN   P =  0   ( no periodic universe )
                                                 RDN   R0 = 249.4075046
                                                 RDN    k = +1   ( spherical space )

    Rmin = 0
    Rmax = 9.999999999 E99 = +infinity

Example4:      Hubble "constant" = 71 km/s/Mpc  ,  (Omega)mat = -0.6 ,  (Omega)lambda  = 0.4 ,   (Omega) rad  = 0.1


              71     ENTER^
            -0.6    ENTER^
              0.4    ENTER^
              0.1    XEQ "AUM"   >>>>   T  = 17.27498133                             ---Execution time = 44s---
                                               RDN   P =  0   ( no periodic universe )
                                               RDN   R0 = 13.13082118
                                               RDN    k = -1   ( hyperbolic space )

    Rmin = 0
    Rmax = 9.999999999 E99 = +infinity

Example5:      Hubble "constant" = 71 km/s/Mpc  ,  (Omega)mat = 0.01 ,  (Omega)lambda  = -0.4 ,   (Omega) rad  = 0.01 


            71    ENTER^
          0.01   ENTER^
          -0.4   ENTER^
          0.01   XEQ "AUM"   >>>>   T  = 11.30457458                                ---Execution time = 107s---
                                            RDN   P =  66.21313789
                                            RDN   R0 = 11.72326781
                                            RDN    k = -1 ( hyperbolic space )

    Rmin = 0
    Rmax = 21.83998798

Example6:      Hubble "constant" = 71 km/s/Mpc  ,  (Omega)mat = -0.4 ,  (Omega)lambda  = +0.001 ,   (Omega) rad  = 1.4 


            71    ENTER^
          -0.4   ENTER^
         0.001  ENTER^
           1.4    XEQ "AUM"   >>>>   T  = 6.484938348                           ---Execution time = 109s---
                                            RDN   P =  900.0920468
                                            RDN   R0 = 435.5000702
                                            RDN    k = +1 ( spherical space )

    Rmin = 0
    Rmax = 1889.773264

-In example 6, T is calculated by a difference between 2 numbers close to P
-So, the precision may be only 6 decimals with an HP41:
-With Free42 , it yields:  T = 6.4849377713


6°)
Redshift -> Distance  


 
    a)  Gauss-Legendre 2-Point Formula




Data Registers:     R00 = H0                                                 ( Registers R00-R01-R02-R03 are to be initialized before executing "Z-D" )

                               
 R01 = (Omega)mat            R04 = D          R07 = R0       R09 to R14:  temp
                               
 R02 = (Omega)lambda       R05 = D0         R08 = k    
                                 R03 = (Omega)rad            R06 = DL

Flags:  /
Subroutines:  /


 01 LBL "Z-D"
 02 DEG
 03 STO 14
 04 SIGN
 05 +
 06 STO 10          
 07 RCL 01
 08 LASTX
 09 -
 10 RCL 02
 11 +
 12 RCL 03
 13 +
 14 STO 11
 15 X#0?
 16 SIGN
 17 STO 08
 18 977.7922214
 19 RCL 00
 20 /
 21 STO 09
 22 RCL 11          
 23 ABS
 24 SQRT
 25 X=0?
 26 SIGN
 27 /
 28 STO 07
 29 1
 30 RCL 10
 31 1/X
 32 STO 12
 33 -
 34 RCL 14
 35 ST+ 14
 36 /
 37 STO Y
 38 3
 39 SQRT
 40 /
 41 STO 13          
 42 -
 43 STO 06
 44 2
 45 /
 46 ST- 12
 47 CLX
 48 STO 04
 49 STO 05
 50 LBL 01
 51 RCL 12
 52 RCL 13
 53 X<> 06
 54 STO 13
 55 +
 56 STO 12          
 57 ENTER
 58 ENTER
 59 X^2
 60 RCL 02
 61 *
 62 RCL 11
 63 -
 64 *
 65 RCL 01
 66 +
 67 *
 68 RCL 03
 69 +
 70 SQRT
 71 1/X
 72 ST+ 05
 73 *
 74 ST+ 04
 75 DSE 14
 76 GTO 01
 77 RCL 09          
 78 RCL 13
 79 *
 80 60
 81 SIN
 82 *
 83 ST* 04
 84 ST* 05
 85 RCL 05
 86 RCL 07          
 87 /
 88 ENTER
 89 E^X-1
 90 LASTX
 91 CHS
 92 E^X-1
 93 -
 94 2
 95 /
 96 RCL Y
 97 R-D
 98 SIN
 99 RCL 08          
100 X#0?
101 ENTER
102 X>0?
103 ENTER
104 R^
105 RCL 10
106 *
107 RCL 07
108 *
109 STO 06
110 RCL 05
111 RCL 04
112 END


         ( 142 bytes / SIZE 015 )

 
           STACK           INPUTS        OUTPUTS
               Z
               /
            DL
               Y                z
            D0
               X               N
            D
    

Example1:      H0 = 68 km/s/Mpc  ,  (Omega)mat = 0.044 ,  (Omega)lambda  = 0.521 ,   (Omega) rad  = 0.000049        z = 7

       68       STO 00
    0.044     STO 01
    0.521     STO 02
 0.000049  STO 03          and if  you choose   N = 20

      7   ENTER^
     20  XEQ "Z-D"     >>>>     D = 14.75577029  Gly                             ---Execution time = 42s---
                                  RDN     D0 = 35.73893637  Gly
                                  RDN     DL = 432.2896646  Gly

-With N = 40, it yields:   14.75576636     35.73909040     432.2929577  
-With N = 100 -------    14.75576611     35.73910080     432.2931801

-And we have:  R0 = R07 = 21.80306243 Gly  and  k = R08 = -1 ( hyperbolic space )


Example2:      H0 = 68 km/s/Mpc  ,  (Omega)mat = -0.1 ,  (Omega)lambda  = -0.003 ,   (Omega) rad  = -0.001         z = 7

       68        STO 00
     -0.1       STO 01
   -0.003     STO 02
   -0.001     STO 03          and if  we choose   N = 20

      7   ENTER^
     20  XEQ "Z-D"     >>>>     D = 13.88126253  Gly                             ---Execution time = 42s---
                                  RDN     D0 = 35.94478929  Gly
                                  RDN     DL = 752.8766458  Gly

-With N = 40, it yields:   13.88233261     35.95641236     753.5230703
-With N = 100 -------    13.88244594     35.95762620     753.5906104

-And we have:  R0 = R07 = 13.68526238 Gly  and  k = R08 = -1 ( hyperbolic space )

Note:

-The precision is relatively good, provided z is not too close to the maximum redshift.


      b)  Tanh-Sinh Quadrature  


-The precision is controlled by the display format.


Data Registers:     R00 = H0                                                 ( Registers R00-R01-R02-R03 are to be initialized before executing "Z-D" )

                               
 R01 = (Omega)mat            R04 = D          R07 = R0       R09 to R14:  temp
                               
 R02 = (Omega)lambda       R05 = D0         R08 = k    
                                 R03 = (Omega)rad            R06 = DL

Flags:  /
Subroutines:  /


 01 LBL "Z-D"
 02 DEG
 03 1
 04 STO 06
 05 STO 13
 06 STO 15
 07 +
 08 STO 16
 09 1/X
 10 RCL 01
 11 RCL 02
 12 +
 13 RCL 03
 14 +
 15 RCL 06
 16 -
 17 STO 09
 18 CLX
 19 STO 07          
 20 STO 08
 21 SIGN
 22 X<>Y
 23 ST- 15
 24 +
 25 2
 26 ST/ 15
 27 /
 28 STO 14
 29 XEQ 00
 30 STO 05
 31 X<>Y
 32 STO 04
 33 GTO 02
 34 LBL 00
 35 ENTER
 36 ENTER
 37 X^2
 38 RCL 02
 39 *
 40 RCL 09
 41 -
 42 *
 43 RCL 01
 44 +
 45 *
 46 RCL 03          
 47 +
 48 SQRT
 49 RCL 15
 50 X<>Y
 51 /
 52 ST* Y
 53 RTN
 54 LBL 01
 55 2
 56 ST/ 04
 57 ST/ 05
 58 ST/ 06
 59 STO 13
 60 SIGN
 61 CHS
 62 STO 08
 63 LBL 02
 64 RCL 13
 65 ST+ 08
 66 RCL 06
 67 RCL 08          
 68 *
 69 E^X-1
 70 LASTX
 71 CHS
 72 E^X-1
 73 -
 74 STO 12
 75 E^X-1
 76 STO Y
 77 2
 78 ST/ 12
 79 +
 80 /
 81 RCL 15
 82 *
 83 STO 11
 84 RCL 14
 85 +
 86 XEQ 00
 87 STO 10
 88 X<>Y
 89 X<> 11
 90 CHS
 91 RCL 14
 92 +
 93 XEQ 00
 94 RCL 10          
 95 +
 96 X<>Y
 97 RCL 11
 98 +
 99 RCL 06
100 ST* Z
101 ST* Y
102 RCL 08
103 *
104 E^X
105 ENTER
106 1/X
107 +
108 ST+ X
109 ST* Z
110 *
111 RCL 12
112 E^X
113 ENTER
114 1/X
115 +
116 X^2
117 ST/ Z
118 /
119 RCL 04
120 +
121 STO 04        
122 LASTX
123 -
124 ABS
125 X<>Y
126 RCL 05
127 +
128 STO 05
129 LASTX
130 -
131 ABS
132 +
133 X#0?
134 GTO 02
135 VIEW 05
136 RCL 05
137 RND
138 LASTX
139 X<> 07
140 RND
141 X#Y?
142 GTO 01
143 977.7922214
144 RCL 00
145 /
146 ST* 04
147 ST* 05
148 RCL 05        
149 X<>Y
150 RCL 09
151 ABS
152 SQRT
153 X=0?
154 SIGN
155 /
156 STO 07
157 /
158 ENTER
159 E^X-1
160 LASTX
161 CHS
162 E^X-1
163 -
164 2
165 /
166 RCL Y
167 R-D
168 SIN
169 RCL 09
170 X#0?
171 SIGN
172 STO 08        
173 X#0?
174 ENTER
175 X>0?
176 ENTER
177 R^
178 RCL 07
179 *
180 RCL 16
181 *
182 STO 06
183 RCL 05
184 RCL 04
185 CLD
186 END


         ( 235 bytes / SIZE 017 )

 
           STACK           INPUTS        OUTPUTS
               Z
               /
            DL
               Y                /
            D0
               X                z
            D
    

Example1:      H0 = 68 km/s/Mpc  ,  (Omega)mat = 0.044 ,  (Omega)lambda  = 0.521 ,   (Omega) rad  = 0.000049        z = 7

       68       STO 00
    0.044     STO 01
    0.521     STO 02
 0.000049  STO 03      and with  FIX 7

      7   XEQ "Z-D"     >>>>     D = 14.75576607  Gly                             ---Execution time = 142s---
                                  RDN     D0 = 35.73910110  Gly
                                  RDN     DL = 432.2931865  Gly


-And we have:  R0 = R07 = 21.80306243 Gly  and  k = R08 = -1 ( hyperbolic space )


Example2:      H0 = 68 km/s/Mpc  ,  (Omega)mat = -0.1 ,  (Omega)lambda  = -0.003 ,   (Omega) rad  = -0.001         z = 7

       68        STO 00
     -0.1       STO 01
   -0.003     STO 02
   -0.001     STO 03       FIX 7

      7   XEQ "Z-D"     >>>>     D = 13.88244965  Gly                             ---Execution time = 145s---
                                  RDN     D0 = 35.95766592  Gly
                                  RDN     DL = 753.5928202  Gly


-And   R0 = R07 = 13.68526238 Gly  and  k = R08 = -1 ( hyperbolic space )


      c)  Cosmological Constant = 0  ( 2 programs )


-With positive matter & negative pressure.

-These universes have spherical spaces ( k = +1 )
-All the distances are computed with elementary functions.

-You choose the minimum, current & maximum scale factors in R01 to R03


Data Registers:       R00-R04-R05-R06: temp                              ( Registers R01-R02-R03 are to be initialized before executing "Z-D" )

                               
 R01 = Rmin      
                               
 R02 = R0           
                                 R03 = Rmax      
Flags:  /
Subroutines:  /


 01 LBL "Z-D"
 02 RAD
 03 1
 04 +
 05 STO 00
 06 RCL 03
 07 RCL 02          
 08 -
 09 STO 04
 10 LASTX
 11 RCL 01
 12 -
 13 ST* 04
 14 -
 15 RCL 01
 16 RCL 03
 17 -
 18 STO 06          
 19 /
 20 ASIN
 21 RCL 01
 22 RCL 02
 23 R^
 24 /
 25 -
 26 STO 05
 27 LASTX
 28 RCL 03
 29 -
 30 ST* 05
 31 -
 32 RCL 06          
 33 /
 34 ASIN
 35 -
 36 ENTER
 37 SIN
 38 RCL 00
 39 *
 40 RCL 01
 41 RCL 03          
 42 +
 43 2
 44 /
 45 R^
 46 *
 47 RCL 04
 48 SQRT
 49 -
 50 RCL 05          
 51 SQRT
 52 +
 53 R^
 54 RCL 02
 55 ST* T
 56 *
 57 X<>Y
 58 DEG
 59 END

     
         ( 70 bytes / SIZE 007 )

 
           STACK           INPUTS        OUTPUTS
               Z
               /
            DL
               Y                /
            D0
               X                z
            D
   

Example:         Rmin = 2   R0 = 41   Rmax = 257  ( Gygalightyears )       z = 7    

    2    STO 01
   41   STO 02
  257  STO 03

    7   XEQ "Z-D"   >>>>     D = 11.60902212   Gly                         ---Execution time = 3s---
                              RDN     D0 = 23.85152165  Gly
                              RDN     DL = 180.2301791  Gly


-In the following program, we assume that   1 - (Omega)mat - (Omega)rad  > 0  ,  k = - 1 ( hyperbolic space )
-You choose  H0 ,  (Omega)mat & (Omega)rad


Data Registers:     R00 = H0                                                 ( Registers R00-R01-R02 are to be initialized before executing "Z-D" )

                               
 R01 = (Omega)mat       R03 = z+1    R04-R06:  temp
                               
 R02 = (Omega)rad       R05 = D    
                               
Flags:  /
Subroutines:  /


 01 LBL "Z-D"
 02 1
 03 +
 04 1/X
 05 STO 03
 06 RCL 01          
 07 1
 08 RCL 01
 09 -
 10 RCL 02
 11 -
 12 1/X
 13 STO 04
 14 *
 15 STO 05
 16 +
 17 *
 18 RCL 02          
 19 RCL 04
 20 *
 21 +
 22 SQRT
 23 STO 06
 24 +
 25 RCL 05
 26 2
 27 /
 28 STO Z
 29 ST+ Y
 30 1
 31 +
 32 RCL 04          
 33 SQRT
 34 STO 04
 35 ST- 06
 36 +
 37 /
 38 LN
 39 STO 05
 40 *
 41 RCL 06
 42 -
 43 X<> 05
 44 CHS
 45 ENTER
 46 E^X-1
 47 LASTX
 48 CHS
 49 E^X-1
 50 -
 51 RCL 03          
 52 ST+ X
 53 /
 54 X<>Y
 55 RCL 04
 56 RCL 00          
 57 /
 58 977.7922214
 59 *
 60 ST* 05
 61 ST* Z
 62 *
 63 RCL 05
 64 END


         ( 89 bytes / SIZE 007 )

 
           STACK           INPUTS        OUTPUTS
               Z
               /
            DL
               Y                /
            D0
               X                z
            D
    

Example1:      H0 = 71 km/s/Mpc  ,  (Omega)mat = 0.044 ,   (Omega) rad  = 0.000049        z = 7

       71       STO 00
    0.044     STO 01
 0.000049  STO 02  

      7   XEQ "Z-D"     >>>>     D =   11.71594217  Gly                             ---Execution time = 3.5s---
                                  RDN     D0 =  27.30593471  Gly
                                  RDN     DL = 383.4090767  Gly


Example2:    H0 = 68 km/s/Mpc     (Omega)mat = 0.2       (Omega)rad = 0.1      z = 7

      68  STO 00
     0.2  STO 01
     0.1  STO 02

      7   XEQ "Z-D"     >>>>     D =   10.03687104  Gly                      
                                  RDN     D0 =  20.51883416  Gly
                                  RDN     DL = 206.0220177  Gly


      d)  No Matter, Negative Pressure


-These universes have hyperbolic spaces ( k = -1 )
-The light-travel time distance is computed with elementary functions.
-The comoving distance employs Carlson elliptic integrals RF

-You choose the minimum, current & maximum scale factors in R01 to R03


Data Registers:       R00 thru R11: temp                                               ( Registers R01-R02-R03 are to be initialized before executing "Z-D" )

                               
 R01 = Rmin            R04 = D       
                               
 R02 = R0               R05 = D0          
                                 R03 = Rmax            R06 = DL

Flags:  /
Subroutine:  ( M-Code )  RF  ( cf "Carlson Elliptic Integrals for the HP41" )


 01 LBL "Z-D"
 02 DEG
 03 1
 04 +
 05 STO 10
 06 RCL 02
 07 X<>Y
 08 /
 09 STO 00
 10 RCL 02
 11 X^2
 12 STO 06          
 13 ST+ X
 14 RCL 01
 15 X^2
 16 STO 04
 17 STO 07
 18 RCL 03
 19 X^2
 20 STO 08
 21 ST- 04
 22 +
 23 STO 11
 24 -
 25 RCL 04
 26 /
 27 ACOS
 28 RCL 00
 29 X^2
 30 STO 09
 31 ST+ X
 32 RCL 11          
 33 -
 34 RCL 04
 35 /
 36 ACOS
 37 -
 38 D-R
 39 2
 40 /
 41 STO 04
 42 RCL 06
 43 RCL 07
 44 -
 45 RCL 08
 46 RCL 09
 47 -
 48 *
 49 SQRT
 50 RCL 09
 51 RCL 07          
 52 -
 53 RCL 08
 54 RCL 06
 55 -
 56 *
 57 SQRT
 58 +
 59 X^2
 60 STO 08
 61 RCL 01
 62 RCL 02
 63 +
 64 STO 09
 65 RCL 03
 66 RCL 02
 67 -
 68 STO 06
 69 *
 70 RCL 00          
 71 RCL 01
 72 -
 73 ST* 09
 74 *
 75 RCL 00
 76 RCL 03
 77 +
 78 ST* 06
 79 *
 80 SQRT
 81 RCL 00
 82 RCL 01
 83 +
 84 ST* 06
 85 RCL 03
 86 RCL 00
 87 -
 88 ST* 09
 89 *
 90 RCL 02          
 91 RCL 01
 92 -
 93 ST* 06
 94 *
 95 RCL 02
 96 RCL 03
 97 +
 98 ST* 09
 99 *
100 SQRT
101 +
102 X^2
103 RCL 09
104 SQRT
105 RCL 06
106 SQRT
107 +
108 X^2
109 RCL 02
110 RCL 00        
111 -
112 X^2
113 ST/ 08
114 ST/ Z
115 /
116 RCL 08
117 RF
118 ST+ X
119 RCL 11
120 SQRT
121 ST* 04
122 *
123 STO 05
124 E^X-1
125 LASTX
126 CHS
127 E^X-1
128 -
129 2
130 /
131 RCL 10        
132 *
133 RCL 02
134 ST* 05
135 *
136 STO 06
137 RCL 05
138 RCL 04
139 END


         ( 162 bytes / SIZE 012 )

 
           STACK           INPUTS        OUTPUTS
               Z
               /
            DL
               Y                /
            D0
               X                z
            D
    

Example:         Rmin = 0.7   R0 = 14   Rmax = 314  ( Gygalightyears )       z = 7    ,    z = 19 = zmax

   0.7  STO 01
   14   STO 02
  314  STO 03

    7   XEQ "Z-D"   >>>>     D = 12.38326745   Gly                                            ---Execution time = 11s---
                              RDN     D0 = 29.70733010  Gly
                              RDN     DL = 460.7466900  Gly

   19   XEQ "Z-D"   >>>>     D = 13.98718383   Gly       
                              RDN     D0 = 51.64270084  Gly
                              RDN     DL = 5595.856028  Gly

-With z = 19,  D has the same value as the age of this universe ( more exactly the time since the last minimum scale factor...

-If we only want to compute D, the program becomes smaller:


Data Registers:       R00-R04: temp                                    ( Registers R01-R02-R03 are to be initialized before executing "Z-D" )

                               
 R01 = Rmin         
                               
 R02 = R0                     
                                 R03 = Rmax        
Flags:  /
Subroutines:  /



 01 LBL "Z-D"
 02 DEG
 03 1
 04 +
 05 RCL 02          
 06 X^2
 07 ST+ X
 08 RCL 01
 09 X^2
 10 STO 04          
 11 RCL 03
 12 X^2
 13 ST- 04
 14 +
 15 STO 00          
 16 -
 17 RCL 04
 18 /
 19 ACOS
 20 RCL 02
 21 R^
 22 /
 23 X^2
 24 ST+ X
 25 RCL 00          
 26 -
 27 RCL 04
 28 /
 29 ACOS
 30 -
 31 D-R
 32 RCL 00          
 33 SQRT
 34 *
 35 2
 36 /
 37 END

     
         ( 48 bytes / SIZE 005 )

 
           STACK           INPUTS        OUTPUTS
               Y                /
           z+1
               X                z
            D
   

Example:             Rmin = 0.7   R0 = 14   Rmax = 314  ( Gygalightyears )       z = 7  

   0.7  STO 01
   14   STO 02
  314  STO 03

    7   XEQ "Z-D"   >>>>     D = 12.38326745   Gly                                            ---Execution time = 2.5s---


      e)  No Pressure, Negative Matter


-These universes have hyperbolic spaces ( k = -1 )
-The light-travel time distance and the comoving distance employ Carlson elliptic integrals RF & RJ

-You choose the minimum, current & maximum scale factors in R01 to R03


Data Registers:       R00 thru R12: temp                                               ( Registers R01-R02-R03 are to be initialized before executing "Z-D" )

                               
 R01 = Rmin            R04 = D       
                               
 R02 = R0               R05 = D0          
                                 R03 = Rmax            R06 = DL

Flags:  /
Subroutines:  ( M-Code )  RF & RJ  ( cf "Carlson Elliptic Integrals for the HP41" )


 01 LBL "Z-D"
 02 DEG
 03 1
 04 +
 05 STO 00
 06 RCL 02
 07 X<>Y
 08 /
 09 ENTER
 10 STO 04
 11 RCL 01          
 12 -
 13 STO 12
 14 *
 15 RCL 03
 16 RCL 02
 17 -
 18 STO 11
 19 *
 20 RCL 01
 21 RCL 03
 22 +
 23 STO 06
 24 STO 10
 25 RCL 02
 26 +
 27 STO 08
 28 *
 29 SQRT
 30 RCL 02          
 31 RCL 02
 32 RCL 01
 33 -
 34 STO 05
 35 *
 36 RCL 03
 37 RCL 04
 38 ST+ 06
 39 -
 40 ST* 08
 41 *
 42 RCL 06
 43 ST* 11
 44 *
 45 SQRT
 46 +
 47 RCL 02
 48 RCL 04
 49 -
 50 STO 07          
 51 /
 52 X^2
 53 STO 09
 54 RCL 04
 55 RCL 11
 56 *
 57 RCL 05
 58 *
 59 SQRT
 60 RCL 02
 61 RCL 08
 62 *
 63 RCL 12
 64 ST* 11
 65 *
 66 SQRT
 67 +
 68 RCL 07
 69 /
 70 X^2
 71 X<> 08
 72 RCL 04          
 73 *
 74 RCL 05
 75 *
 76 SQRT
 77 RCL 02
 78 ST* 04
 79 RCL 11
 80 *
 81 SQRT
 82 +
 83 RCL 07
 84 /
 85 X^2
 86 STO 07
 87 RCL 01
 88 RCL 10
 89 *
 90 -
 91 STO 05
 92 RCL 07          
 93 RCL 08
 94 RCL 09
 95 RJ
 96 RCL 01
 97 RCL 03
 98 *
 99 RCL 10
100 *
101 *
102 3
103 /
104 X<> 05
105 RCL 04
106 STO Z
107 -
108 /
109 SQRT
110 ATAN
111 D-R
112 RCL 05          
113 -
114 STO 04
115 RCL 07
116 RCL 08
117 RCL 09
118 RF
119 STO 05
120 RCL 03
121 RCL 10
122 *
123 RCL 01
124 X^2
125 +
126 SQRT
127 ST+ X
128 ST* 04
129 ST* 05
130 RCL 05          
131 E^X-1
132 LASTX
133 CHS
134 E^X-1
135 -
136 2
137 /
138 RCL 00
139 *
140 RCL 02
141 ST* 05
142 *
143 STO 06
144 RCL 05
145 RCL 04
146 END


         ( 168 bytes / SIZE 0123)

 
           STACK           INPUTS        OUTPUTS
               Z
               /
            DL
               Y                /
            D0
               X                z
            D
    

Example:         Rmin = 1   R0 = 14   Rmax = 314  ( Gygalightyears )   ,    z = 7  

    1    STO 01
   14   STO 02
  314  STO 03

    7   XEQ "Z-D"   >>>>     D =  13.56149794   Gly                                                ---Execution time = 21s---
                              RDN     D0 = 33.91425300  Gly
                              RDN     DL = 626.3433837  Gly-


      f)  Tolman Universes  ( 4 programs )


-Here,  (Omega)mat = 0
-This program computes D , T , R0 , & k


Data Registers:       R00-R06: temp
Flags:  /
Subroutines:  /


 01 LBL "Z-D"
 02 R^
 03 1/X
 04 STO 06
 05 SIGN
 06 +
 07 X^2
 08 1/X
 09 STO 00          
 10 SIGN
 11 X<>Y
 12 STO 01
 13 -
 14 X<>Y
 15 STO 02
 16 -
 17 STO 03
 18 X#0?
 19 SIGN
 20 CHS
 21 RCL 03
 22 RCL 01
 23 ST+ X
 24 /
 25 STO 04          
 26 977.7922214
 27 ST* 06
 28 SIGN
 29 +
 30 RCL 01
 31 SQRT
 32 1/X
 33 +
 34 STO 05
 35 RCL 00
 36 RCL 03
 37 *
 38 RCL 02
 39 +
 40 RCL 01          
 41 /
 42 RCL 00
 43 X^2
 44 +
 45 SQRT
 46 RCL 04
 47 +
 48 RCL 00
 49 +
 50 /
 51 SQRT
 52 LN
 53 X<> 05
 54 RCL 02
 55 RCL 01          
 56 /
 57 SQRT
 58 RCL 04
 59 +
 60 /
 61 SQRT
 62 LN
 63 RCL 01
 64 SQRT 
 65 ST/ 05
 66 /
 67 RCL 03
 68 ABS
 69 X=0?
 70 SIGN
 71 SQRT
 72 1/X
 73 X<>Y
 74 RCL 06          
 75 ST* 05
 76 ST* Z
 77 *
 78 RCL 05
 79 END

     
         ( 103 bytes / SIZE 007 )

              

           STACK          INPUTS        OUTPUTS
               T
              H0
             k
               Z
     (Omega)rad >=0
             R0
               Y     (Omega)lambda > 0
             T
               X                z
             D


Example1:    H0 = 68 km/s/Mpc    (Omega)rad  = 0.000049      (Omega)lambda   = 0.314      z = 7


       68        ENTER^
 0.000049   ENTER^
    0.314      ENTER^
        7         XEQ "Z-D"   >>>>     D = 14.08131707 Gly                             ---Execution time = 3.7s---
                                        RDN      T =  16.10720928 Gy
                                        RDN     R0 =  17.36165485 Gly
                                        RDN      k  = -1

Example2:    H0 = 68 km/s/Mpc    (Omega)rad  = 0.3      (Omega)lambda   = 0.7        z = 7


       68   ENTER^
      0.3   ENTER^
      0.7   ENTER^
        7    XEQ "Z-D"   >>>>     D = 10.19222258 Gly
                                   RDN      T =  10.39730378 Gy
                                   RDN     R0 =  14.37929737 Gly
                                   RDN      k  =  0

Example3:    H0 = 68 km/s/Mpc    (Omega)rad  = 0.4      (Omega)lambda   = 0.8        z = 7


       68   ENTER^
      0.4   ENTER^
      0.8   ENTER^
        7    XEQ "Z-D"   >>>>     D =  9.715222804 Gly
                                   RDN      T =  9.893178897 Gy
                                   RDN     R0 =  32.15308638 Gly
                                   RDN      k  =  +1

-A more general version:


Data Registers:             R00 = k  ( -1 , 0 , +1 )

                                        R01 = D             R04 = Rmin             R07 to R13: temp
                                        R02 = T             R05 = R0                 R14 = H0
                                        R03 = P             R06 = Rmax
Flag:  F24
Subroutines: /
 
 
 

 01 LBL "TOL"
 02 DEG
 03 R^
 04 STO 14
 05 SIGN
 06 +
 07 X^2
 08 1/X
 09 STO 12
 10 RDN
 11 STO 13
 12 X<>Y
 13 STO 07
 14 +
 15 90
 16 TAN
 17 STO 03
 18 STO 06
 19 STO 09
 20 CLX
 21 STO 04
 22 SIGN
 23 -
 24 STO 00
 25 CHS
 26 RCL 13
 27 X=0?
 28 SIGN
 29 ST/ 07
 30 /
 31 STO 08
 32 RCL 13          
 33 X=0?
 34 GTO 01
 35 ABS
 36 SQRT
 37 ST+ X
 38 STO 10
 39 CLX
 40 2
 41 /
 42 STO 08
 43 X^2
 44 RCL 07
 45 -
 46 X<0?
 47 GTO 02
 48 SQRT
 49 RCL 08
 50 SIGN
 51 *
 52 RCL 08
 53 +
 54 CHS
 55 RCL 07
 56 RCL Y
 57 X=0?
 58 SIGN
 59 /
 60 X<0?
 61 CLX
 62 SQRT
 63 X<>Y
 64 X<0?
 65 CLX
 66 SQRT
 67 X<Y?
 68 X<>Y
 69 GTO 03
 70 LBL 01
 71 RCL 07
 72 CHS
 73 RCL 08          
 74 X=0?
 75 GTO 02
 76 /
 77 X<0?
 78 CLX
 79 SQRT
 80 ENTER
 81 GTO 03
 82 LBL 02
 83 CLST
 84 LBL 03
 85 1
 86 X<>Y
 87 X>Y?
 88 STO 06
 89 X<> Z
 90 X>Y?
 91 STO 06
 92 X>0?
 93 X>Y?
 94 FS? 30
 95 STO 04
 96 X<> Z
 97 X>0?
 98 X>Y?
 99 FS? 30
100 STO 04
101 RCL 13
102 X=0?
103 GTO 01
104 X<0?
105 GTO 02
106 RCL 07
107 RCL 08
108 ST+ X
109 +
110 1
111 +
112 SQRT
113 1
114 +
115 RCL 08
116 +
117 STO 11        
118 LASTX
119 RCL 07
120 X<0?
121 CLX
122 SQRT
123 +
124 X=0?
125 GTO 00
126 STO 05
127 /
128 ABS
129 LN
130 RCL 10
131 /
132 STO 02
133 RCL 11
134 RCL 08
135 ST+ X
136 RCL 12
137 ST* Y
138 X^2
139 +
140 RCL 07
141 +
142 SQRT
143 RCL 08
144 +
145 RCL 12
146 +
147 /
148 ABS
149 LN
150 RCL 10
151 /
152 STO 01
153 RCL 09
154 RCL 06
155 X=Y?
156 GTO 03
157 X^2
158 RCL 08
159 +
160 RCL 05        
161 /
162 ABS
163 LN
164 ST+ X
165 RCL 10
166 /
167 STO 03
168 GTO 03
169 LBL 00
170 RCL 09
171 STO 02
172 GTO 03
173 LBL 01
174 1
175 RCL 07
176 X<0?
177 CLX
178 SQRT
179 STO 05
180 -
181 RCL 08
182 X=0?
183 GTO 04
184 /
185 STO 02
186 1
187 RCL 08
188 RCL 12
189 *
190 RCL 07
191 +
192 SQRT
193 -
194 RCL 08
195 /
196 STO 01
197 RCL 09
198 RCL 06
199 X=Y?
200 GTO 03
201 X^2
202 RCL 08        
203 *
204 RCL 07
205 +
206 SQRT
207 RCL 05
208 -
209 ST+ X
210 RCL 08
211 /
212 STO 03
213 GTO 03
214 LBL 02
215 1
216 RCL 08
217 ST+ Y
218 X^2
219 RCL 07
220 -
221 SQRT
222 STO 01
223 /
224 ASIN
225 STO 09
226 90
227 +
228 D-R
229 RCL 10
230 /
231 STO 02
232 RCL 09
233 RCL 08
234 RCL 12
235 +
236 RCL 01
237 /
238 ASIN
239 -
240 D-R
241 RCL 10
242 /
243 STO 01        
244 PI
245 ST+ X
246 RCL 10
247 /
248 STO 03
249 GTO 03
250 LBL 04
251 .5
252 ENTER
253 STO 02
254 RCL 12
255 *
256 -
257 STO 01
258 LBL 03
259 RCL 00
260 X#0?
261 SIGN
262 ENTER
263 X<> 00
264 ABS
265 SQRT
266 X=0?
267 SIGN
268 1/X
269 SF 24
270 977.7922214
271 RCL 14
272 /
273 ST* 01
274 ST* 02
275 ST* 03
276 *
277 STO 05
278 ST* 04
279 ST* 06
280 RCL 02        
281 RCL 03
282 SIGN
283 CLX
284 RCL 01
285 CF 24
286 END

 
        ( 344 bytes / SIZE 015 )
 
 

          STACK           INPUTS          OUTPUTS
              T               H0                k
              Z         (Omega)rad               R0
              Y        (Omega)lambda               T
              X                z               D
              L                /               P



Example1:    With  H0 = 68 km/s/Mpc   (Omega)rad = -0.0001  ,    (Omega)lambda  = -0.1  ,  z = 7

           68      ENTER^
      -0.0001  ENTER^
         -0.1     ENTER^
            7     XEQ "TOL"   >>>>   D  = 12.21653956  Gly    = R01                             ---Execution time = 6.7s---
                                          RDN    T =  13.92565130  Gy    = R02
                                          RDN   R0 =  13.70949896 Gly   = R05
                                          RDN    k  =   -1                          = R00

                                       LASTX   P  =  142.8523990  Gy  = R03

 And we also have:

    R04 =  Rmin =  0.130709543  Gly
    R06 =  Rmax =  45.47114300  Gly

Example2:    H0 = 68 km/s/Mpc  (Omega)rad = 0.1  ,   (Omega)lambda = 0.6  ,   z  = 7

         68    ENTER^
         0.1   ENTER^
         0.6   ENTER^
           7    XEQ "TOL"  >>>>  D  = 12.18609378  Gly    = R01
                                      RDN    T =  12.53718802  Gy    = R02
                                      RDN   R0 =  26.25288510 Gly   = R05
                                      RDN    k  =   -1                          = R00

                                   LASTX   P  =  9.9999999 99  Gy  = R03   ( infinite )

 And we also have:

    R04 =  Rmin =  0.000000000  Gly
    R06 =  Rmax =  9.9999999 99  Gly   ( infinite )
 
Example3:   H0 = 68 km/s/Mpc   (Omega)rad = 1.1  ,   (Omega)lambda  = 0  ,   z  = 7

          68   ENTER^
         1.1   ENTER^
          0     ENTER^
          7     XEQ "TOL"  >>>>  D  =  6.911221056 Gly    = R01
                                      RDN   T =  7.018369397  Gy    = R02
                                     RDN   R0 =  45.47133084 Gly   = R05
                                     RDN    k  =   +1                          = R00

                                   LASTX   P  =  301.6226862 Gy   = R03

 And we also have:

    R04 =  Rmin =  0.000000000  Gly
    R06 =  Rmax =  150.8113431  Gly


-We can use the program listed in §6°)a) to compute D0 & DL

-Here is another program with Gauss-Legendre 3-point formula:


Data Registers:     R00 = N                                                 ( Registers R00-R01-R02-R03 are to be initialized before executing "Z-D" )

                               
 R01 = H0                                                R04 to R10:  temp
                               
 R02 = (Omega)lambda  > 0  
                                 R03 = (Omega)rad       > 0                      0 <  1 -  (Omega)lambda - (Omega)rad  < 1

Flags:  /
Subroutines:  /


 01 LBL "Z-D"
 02 1
 03 +
 04 1/X
 05 STO 04
 06 1
 07 RCL 02
 08 -
 09 RCL 03
 10 -
 11 STO 05
 12 1
 13 RCL 04          
 14 STO 07
 15 -
 16 RCL 00
 17 STO 08
 18 /
 19 STO 06
 20 2
 21 /
 22 ST+ 07
 23 .6
 24 SQRT
 25 *
 26 STO 09
 27 CLX
 28 STO 10
 29 GTO 01
 30 LBL 00
 31 X^2
 32 STO Y
 33 RCL 02
 34 *
 35 RCL 05
 36 +
 37 *
 38 RCL 03
 39 +
 40 SQRT
 41 1/X
 42 RTN
 43 LBL 01
 44 RCL 07          
 45 RCL 09
 46 -
 47 XEQ 00
 48 ST+ 10
 49 RCL 07
 50 XEQ 00
 51 1.6
 52 *
 53 ST+ 10
 54 RCL 07
 55 RCL 09
 56 +
 57 XEQ 00
 58 ST+ 10
 59 RCL 06
 60 ST+ 07
 61 DSE 08
 62 GTO 01
 63 RCL 10
 64 *
 65 3.6
 66 /
 67 STO 07
 68 RCL 05
 69 RCL 02          
 70 ST+ X
 71 /
 72 STO 06
 73 1
 74 +
 75 X^2
 76 RCL 03
 77 RCL 02
 78 /
 79 RCL 06
 80 X^2
 81 -
 82 STO 08
 83 +
 84 SQRT
 85 1
 86 +
 87 RCL 04
 88 X^2
 89 STO 09
 90 RCL 06
 91 +
 92 X^2
 93 RCL 08
 94 +
 95 SQRT
 96 RCL 06          
 97 ST+ Z
 98 +
 99 RCL 09
100 +
101 /
102 LN
103 RCL 02
104 SQRT
105 ST+ X
106 /
107 STO 06
108 RCL 07
109 RCL 05
110 SQRT
111 STO 05
112 *
113 E^X-1
114 LASTX
115 CHS
116 E^X-1
117 -
118 RCL 04
119 ST+ X
120 /
121 RCL 05
122 /
123 977.7922214
124 RCL 01        
125 /
126 ST* 06
127 ST* 07
128 *
129 RCL 07
130 RCL 06
131 END

     
         ( 175 bytes / SIZE 011 )

              

           STACK           INPUTS        OUTPUTS
               Z
                /
             DL
               Y                 /
             D0
               X                 z
             D


Example:    H0 = 68 km/s/Mpc     (Omega)lambda   = 0.7    (Omega)rad  = 0.1     z = 7

-With N = 10

      10  STO 00
 
     68   STO 01
     0.7  STO 02
     0.1  STO 03

        7         XEQ "Z-D"   >>>>     D  =  12.55193503 Gly                             ---Execution time = 26s---
                                        RDN     D0 =  26.16037018 Gly
                                        RDN     DL =  233.1494318 Gly


-In most cases, we cannot calculate D0 with elementary functions.

-But if  (Omega)rad  = 1 + (Omega)lambda - 2 sqrt (Omega)lambda  , it becomes possible !

-Here is such a program:


Data Registers:       R00-R07: temp
Flags:  /
Subroutines:  /


 01 LBL "Z-D"
 02 DEG
 03 1
 04 +
 05 1/X
 06 STO 00
 07 RDN
 08 STO 01
 09 ENTER
 10 SQRT
 11 STO 04          
 12 ST+ X
 13 STO 05
 14 -
 15 1
 16 +
 17 STO 02
 18 X<>Y
 19 STO 06
 20 RCL 04
 21 1/X
 22 1
 23 -
 24 STO 03
 25 2
 26 +
 27 RCL 00          
 28 X^2
 29 RCL 03
 30 ST+ Z
 31 +
 32 ST+ X
 33 /
 34 LN
 35 RCL 05
 36 /
 37 STO 05
 38 RCL 03          
 39 SQRT
 40 1/X
 41 ATAN
 42 RCL 00
 43 LASTX
 44 *
 45 ATAN
 46 -
 47 D-R
 48 RCL 01
 49 ST- 04
 50 RCL 03
 51 *
 52 SQRT
 53 /
 54 STO 07
 55 RCL 04
 56 ST+ X
 57 SQRT
 58 STO 04          
 59 *
 60 E^X-1
 61 LASTX
 62 CHS
 63 E^X-1
 64 -
 65 RCL 00
 66 ST+ X
 67 RCL 04
 68 *
 69 /
 70 977.7922214
 71 RCL 06          
 72 /
 73 ST* 05
 74 ST* 07
 75 *
 76 RCL 02
 77 X<>Y
 78 RCL 07
 79 RCL 05
 80 END

     
         ( 106 bytes / SIZE 008 )


             STACK              INPUTS           OUTPUTS
                  T
                  /
          (Omega)rad
                  Z
                H0
                DL
                  Y    0 < (Omega)lambda < 1
                D0
                  X                  z
                D


Example:    H0 = 68 km/s/Mpc    (Omega)lambda = 0.16      z = 7

    68   ENTER^
  0.16  ENTER^
     7    XEQ "Z-D"   >>>>    D  =   8.995379940 Gly                             ---Execution time = 4.6s---
                                RDN     D0 =  17.11226659 Gly
                                RDN     DL =  152.9445025 Gly
                                RDN    (Omega)rad  = 0.36


      f)  Milne Universes


-Here,  (Omega)mat = (Omega)lambda = (Omega)rad  = 0  &  k = -1


Data Registers:       R00: temp
Flags:  /
Subroutines:  /


 01 LBL "MLN"
 02 RCL X
 03 977.7922214
 04 R^
 05 /
 06 STO 00
 07 SIGN
 08 +
 09 ST/ Y
 10 LASTX
 11 +
 12 2
 13 /
 14 R^
 15 ST* Y
 16 LN1+X
 17 RCL 00
 18 ST* Z
 19 ST* T
 20 X<>Y
 21 ST* Y
 22 X<> T
 23 END

     
( 48 bytes / SIZE 001 )

            

           STACK          INPUTS        OUTPUTS
               T
               /
             VR
               Z
               /
             DL
               Y               H0
             D0
               X                z
             D
               L
               /
              z


Example:    H0 = 71 km/s/Mpc   ,   z = 7

    71  ENTER^
     7   XEQ "MLN"   >>>>   D  =  12.05025625 Gly                             ---Execution time = 1.2s---
                                RDN     D0 =  28.63748965 Gly
                                RDN     DL =  433.8092250 Gly
                                RDN     VR =  2.079441542 = recessional velocity ( speed of light = 1 )



References:


[1]  B. C. Carlson - "A Table of Elliptic Integrals of the Third Kind"
[2]  Abramowitz and Stegun - "Handbook of Mathematical Functions" -  Dover Publications -  ISBN  0-486-61272-4
[3]  B. C. Carlson - "A Table of Elliptic Integrals: two quadratic factors"