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 ( 4 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) §₀¹ y. [ (Omega)_lambda.y⁴ + ( 1-(Omega)_tot ).y² + (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 )
          D₀ = comoving distance ( in Giga-light-years )
          D_L = luminosity-distance ( in Giga-light-years )

Formulae:

D = (c/H₀) §_y(em)¹ y. [ (Omega)_lambda.y⁴ + ( 1-(Omega)_tot ).y² + (Omega)_mat.y + (Omega)_rad ] ^-1/2 dy y(em) = y at the instant of emission
D₀ = (c/H₀) §_y(em)¹ [ (Omega)_lambda.y⁴ + ( 1-(Omega)_tot ).y² + (Omega)_mat.y + (Omega)_rad ] ^-1/2 dy z + 1 = 1/y_em , z = redshift , y = R/R₀

where (Omega)_tot = (Omega)_mat + (Omega)_lambda + (Omega)_rad

                                                                F(x) = Sinh(x) if k = -1 hyperbolic space
     D_L = R₀ ( z + 1 ) F(D₀/R₀)    where    F(x) =   x         if k = 0 euclidean space R₀ = 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	H₀ ( km/s/Mpc )	R
X	(Omega)_lambda	T

Example:   H₀ = 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	/	R_minorR_max
Z	/	k ( -1 , 0 , +1 )
Y	H₀ ( km/s/Mpc )	R
X	(Omega)_lambda	T
L	/	0 or P

Where P is the period

Example1:    H₀ = 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 R_max = 21.77500351 Giga-light-years
LASTX    P = 68.40819107 Giga-years

Example2:    H₀ = 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 R_min = 0
LASTX    P = 0 ( no period )

Example3:    H₀ = 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 R_min = 11.63922896 Gly
LASTX    P = 0 ( no period )

Example4:    H₀ = 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 R_min = 0
LASTX    P = 0 ( no period )

Example5:    H₀ = 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 R_min = 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 = R_min R05 = R_max R11 = H₀
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	H₀	k
Y	(Omega)_rad	R
X	_{(Omega)lambda}	T

Example1: H₀ = 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 = R_min = 0 & R05 = R_max = 9.999999999 E99

Example2: H₀ = 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 = R_min = 0.395565882 Gly & R05 = R_max = 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: R_min R₀ R_max
and "AUR" returns the age of the Universe T, the period P, the current Hubble "constant" H₀ and the deceleration parameter q

Data Registers: R00 = L (Gy)^-2

R01 = R²_min R04 = T
R02 = R²₀
R03 = R²_max
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	R_min	H₀
Y	R₀	P
X	R_max	T

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

Example: R_min = 1 R₀ = 14 R_max = 314 ( Gygalightyears )

1 ENTER^
14 ENTER^
314 XEQ "AUR" >>>> T = 13.96898881 = R04 ---Execution time = 2.8s---
RDN P = 986.4650960
RDN H₀ = 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⁴ + ( 1-(Omega)_tot ).y² + (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	H₀	/
Z	(Omega)_mat	k
Y	(Omega)_lamda	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 Giga-years ---Execution time = 34s---
   RDN R= 20.88180628 Giga-light-years
   RDN k = -1

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 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 = u²

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 y² = 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¹ y. [ (Omega)_lambda.y⁴ + ( 1-(Omega)_tot ).y² + (Omega)_mat .y + (Omega)_rad ] ^-1/2 dy
-The Period P = 2 (c/H) §_ymin^ymax y. [ (Omega)_lambda.y⁴ + ( 1-(Omega)_tot ).y² + (Omega)_mat .y + (Omega)_rad ] ^-1/2 dy

-The following program assumes that the quartic equation    (Omega)_lambda.y⁴ + ( 1-(Omega)_tot ).y² + (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 y_min & 1 and 3.N subintervals between 1 & y_max ( 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 = R_min    R06 thru R18: temp
   • R02 = (Omega)_lambda < 0 R05 = R_max
   • 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 = R_min = 1.306380222 GLy
R05 = maximum radius of the Universe = R_max = 250.8400402 GLy

Notes:

-The numerical integration is performed with the change of argument: y = A + ( B - A ) Sin² u where A = R_min/R & B = R_max/R
-We must choose: 1-(Omega)_tot > 0
-The constant k = -1 ( hyperbolic universes )

-We can also choose the minimum, current & maximum scale factors R_min R₀ R_max
-The following program employs Carlson elliptic integral R_J

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

       R01 = R_min       R04 = H₀    R06 thru R15: temp
     R02 = R₀    R05 = q
     R03 = R_max
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	R_min	H₀
Y	R₀	P
X	R_max	T

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

Example:    R_min = 1 R₀ = 14 R_max = 314 ( Gigalightyears )

1 ENTER^
14 ENTER^
314 XEQ "AUC" >>>> T = 15.49057984 gigayears ---Execution time = 58s---
RDN P = 993.8651299 gigayears
RDN H₀ = 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

§_a^r x dx / sqrt [ x (b-x) (x+c) (x+d) ] = §₀^+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 R_F [ (b-r)/(b-a) , (r+c)/(a+c) , r/a ] + [(r-a)/3] R_J [ (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 = R_min       R04 = H₀ R10 = P    R06 thru R11: temp
     R02 = R₀    R05 = q
     R03 = R_max
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	R_min	H₀
Y	R₀	P
X	R_max	T

a) Cosmological Constant = 0

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

Data Registers: R00: temp R01 = R_min R02 = R₀ R03 = R_max 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	R_min	H₀
Y	R₀	P
X	R_max	T

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

Example1:    R_min = 2 R₀ = 41 R_max = 257 ( Gygalightyears )

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

Example2:    R_min = 3 R₀ = 37 R_max = 300 ( Gygalightyears )

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

Example3:    R_min = 10 R₀ = 123 R_max = 9876 ( Gygalightyears )

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

Note:

-We can also compute: (Omega)_mat & (Omega)_rad

Data Registers: R00: temp

R01 = R_min R04 = (Omega)_mat R06 = T
R02 = R₀    R05 = (Omega)_rad
R03 = R_max
Flags: /
Subroutines: /

01 LBL "AUL"
02 DEG
03 STO 00
04 STO 03
05 STO 05
06 X<>Y
07 ENTER
08 STO 02
09 ST- Z
10 R^

11 ST+ 00
12 STO 01
13 CHS
14 ST* 05
15 +
16 STO 06
17 R^
18 ST* 06
19 -
20 ST* Y

21 RCL 01
22 RCL 03
23 -
24 /
25 ACOS
26 D-R
27 RCL 00
28 STO 04
29 *
30 2

31 ST/ Z
32 /
33 X<> 06
34 ST+ Y
35 ST/ Y
36 ST/ 04
37 ST/ 05
38 SQRT
39 ST- 06
40 RCL 02

41 ST* 04
42 X^2
43 /
44 977.7922214
45 *
46 RCL 00
47 PI
48 *
49 RCL 06
50 END

( 81 bytes / SIZE 007 )

STACK	INPUTS	OUTPUTS
T	/	q
Z	R_min	H₀
Y	R₀	P
X	R_max	T

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

Example: R_min = 2 R₀ = 41 R_max = 257 ( Gygalightyears )

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

and R04 = (Omega)_mat = 1.260565052
R05 = (Omega)_rad = -0.061016144

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 = H₀

R01 = (Omega)_mat R03 = T R05 = R_min R07 = R_max R09-R10-R11: temp
R02 = (Omega)_rad R04 = P R06 = R₀ 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	H₀	R₀
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 ) R₀ = current scale factor

Example1: H₀ = 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 R₀ = 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: H₀ = 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 R₀ = 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: H₀ = 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 R₀ = 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: H₀ = 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 R₀ = 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: H₀ = 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 R₀ = 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: H₀ = 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 R₀ = 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: H₀ = 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 R₀ = 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: H₀ = 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 R₀ = 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: H₀ = 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 R₀ = 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	H₀	q
Y	(Omega)_mat	R₀
X	_(Omega)rad	T

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

Example1: H₀ = 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 R₀ = 17.18654761 Gly
RDN q = 0.2

Example2: H₀ = 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 R₀ = 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 = H₀

R01 = (Omega)_mat R04 = T R06 = R_min R08 = R₀ R10.....: temp
R02 = (Omega)_lambda R05 = P R07 = R_max 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	H₀ > 0	k
Z	(Omega)_mat	R₀
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 , R₀ = 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 R₀ = 13.10701187
RDN k = -1 ( hyperbolic space )

-We also have in R06 & R07

R_min = 1.306380222
R_max = 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 R₀ = 20.88180628
      RDN k = -1

R_min = 0
R_max = 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 R₀ = 249.4075046
   RDN k = +1 ( spherical space )

R_min = 0
R_max = 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 R₀ = 13.13082118
   RDN k = -1 ( hyperbolic space )

R_min = 0
R_max = 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 R₀ = 11.72326781
RDN k = -1 ( hyperbolic space )

R_min = 0
R_max = 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 R₀ = 435.5000702
RDN k = +1 ( spherical space )

R_min = 0
R_max = 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 , R_min = 0 and R_max = +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 N₁ for T & N₂ for P, for instance, with 30 & 100

   30 ENTER^
100 XEQ 10

Data Registers:    R00 = H₀

R01 = (Omega)_mat    R21 = T R23 = R_min    R25 = R₀
R02 = (Omega)_lambda R22 = P R24 = R_max R26 = k
R03 = (Omega)_rad

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

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

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	H₀ > 0	k
Z	(Omega)_mat	R₀
Y	(Omega)_lambda	P
X	(Omega)_rad	T

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 R₀ = 13.10701187
RDN k = -1 ( hyperbolic space )

-We also have in R23 & R24

R_min = 1.306380222
R_max = 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 R₀ = 20.88180628
      RDN k = -1

R_min = 0
R_max = 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 = H₀

R01 = (Omega)_mat    R04 = T R06 = R_min    R08 = R₀ R10.....: temp
R02 = (Omega)_lambda R05 = P R07 = R_max 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	H₀ > 0	k
Z	(Omega)_mat	R₀
Y	(Omega)_lambda	P
X	(Omega)_rad	T

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 R₀ = 13.10701187
RDN k = -1 ( hyperbolic space )

-We also have in R06 & R07

R_min = 1.306380222
R_max = 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 R₀ = 20.88180628
      RDN k = -1

R_min = 0
R_max = 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 R₀ = 249.4075046
   RDN k = +1 ( spherical space )

R_min = 0
R_max = 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 R₀ = 13.13082118
   RDN k = -1 ( hyperbolic space )

R_min = 0
R_max = 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 R₀ = 11.72326781
RDN k = -1 ( hyperbolic space )

R_min = 0
R_max = 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 R₀ = 435.5000702
RDN k = +1 ( spherical space )

R_min = 0
R_max = 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: §_a^r x dx / sqrt [(x-a) (b-x) (x+c) (x+d) ] may be re-written:

§₀^+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 = H₀

R01 = (Omega)_mat    R04 = T R06 = R_min    R08 = R₀ R10.....: temp
R02 = (Omega)_lambda R05 = P R07 = R_max 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	H₀ > 0	k
Z	(Omega)_mat	R₀
Y	(Omega)_lambda	P
X	(Omega)_rad	T

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 R₀ = 13.10701187
RDN k = -1 ( hyperbolic space )

-We also have in R06 & R07

R_min = 1.306380222
R_max = 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 R₀ = 20.88180628
      RDN k = -1

R_min = 0
R_max = 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 R₀ = 249.4075046
   RDN k = +1 ( spherical space )

R_min = 0
R_max = 9.999999999 E99 = +infinity

Example4: Hubble "constant" = 71 km/s/Mpc , (Omega)_mat = -0.6 , (Omega)_lambda = 0.4 , (Omega)_rad = 0.1

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 R₀ = 11.72326781
RDN k = -1 ( hyperbolic space )

R_min = 0
R_max = 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 R₀ = 435.5000702
RDN k = +1 ( spherical space )

R_min = 0
R_max = 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 = H₀ ( Registers R00-R01-R02-R03 are to be initialized before executing "Z-D" )

• R01 = (Omega)_mat    R04 = D R07 = R₀    R09 to R14: temp
• R02 = (Omega)_lambda R05 = D₀ R08 = k
• R03 = (Omega)_rad R06 = D_L

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	/	D_L
Y	z	D₀
X	N	D

Example1: H₀ = 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 D₀ = 35.73893637 Gly
RDN D_L = 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: R₀ = R07 = 21.80306243 Gly and k = R08 = -1 ( hyperbolic space )

Example2: H₀ = 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 D₀ = 35.94478929 Gly
RDN D_L = 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: R₀ = 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 = H₀ ( Registers R00-R01-R02-R03 are to be initialized before executing "Z-D" )

• R01 = (Omega)_mat    R04 = D R07 = R₀    R09 to R14: temp
• R02 = (Omega)_lambda R05 = D₀ R08 = k
• R03 = (Omega)_rad R06 = D_L

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	/	D_L
Y	/	D₀
X	z	D

Example1: H₀ = 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 D₀ = 35.73910110 Gly
RDN D_L = 432.2931865 Gly

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

Example2: H₀ = 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 D₀ = 35.95766592 Gly
RDN D_L = 753.5928202 Gly

-And R₀ = 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 = R_min
• R02 = R₀
• R03 = R_max
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	/	D_L
Y	/	D₀
X	z	D

Example: R_min = 2 R₀ = 41 R_max = 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 D₀ = 23.85152165 Gly
RDN D_L = 180.2301791 Gly

-In the following program, we assume that 1 - (Omega)_mat - (Omega)_rad > 0 , k = - 1 ( hyperbolic space )
-You choose H₀ , (Omega)_mat & (Omega)_rad

Data Registers: • R00 = H₀ ( 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	/	D_L
Y	/	D₀
X	z	D

Example1: H₀ = 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 D₀ = 27.30593471 Gly
RDN D_L = 383.4090767 Gly

Example2: H₀ = 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 D₀ = 20.51883416 Gly
RDN D_L = 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 = R_min    R04 = D
• R02 = R₀ R05 = D₀
• R03 = R_max R06 = D_L

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	/	D_L
Y	/	D₀
X	z	D

Example: R_min = 0.7 R₀ = 14 R_max = 314 ( Gygalightyears ) z = 7 , z = 19 = z_max

0.7 STO 01
14 STO 02
314 STO 03

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

19 XEQ "Z-D" >>>> D = 13.98718383 Gly
RDN D₀ = 51.64270084 Gly
RDN D_L = 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 = R_min
• R02 = R₀
• R03 = R_max
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: R_min = 0.7 R₀ = 14 R_max = 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 = R_min    R04 = D
• R02 = R₀ R05 = D₀
• R03 = R_max R06 = D_L

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	/	D_L
Y	/	D₀
X	z	D

Example: R_min = 1 R₀ = 14 R_max = 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 D₀ = 33.91425300 Gly
RDN D_L = 626.3433837 Gly-

f) Tolman Universes ( 4 programs )

-Here, (Omega)_mat = 0
-This program computes D , T , R₀ , & 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	H₀	k
Z	(Omega)_rad >=0	R₀
Y	(Omega)lambda > 0	T
X	_z	D

Example1: H₀ = 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 R₀ = 17.36165485 Gly
RDN k = -1

Example2: H₀ = 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 R₀ = 14.37929737 Gly
   RDN k = 0

Example3: H₀ = 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 R₀ = 32.15308638 Gly
   RDN k = +1

-A more general version:

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

                                        R01 = D            R04 = R_min             R07 to R13: temp
                                        R02 = T           R05 = R₀ R14 = H₀
                                        R03 = P             R06 = R_max
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	H₀	k
Z	(Omega)_rad	R₀
Y	(Omega)_lambda	T
X	z	D
L	/	P

Example1:    With H₀ = 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   R₀ = 13.70949896 Gly   = R05
                                          RDN    k =   -1                          = R00

                                       LASTX   P = 142.8523990 Gy = R03

And we also have:

    R04 = R_min = 0.130709543 Gly
    R06 = R_max = 45.47114300 Gly

Example2:    H₀ = 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   R₀ = 26.25288510 Gly   = R05
                                      RDN    k =   -1                          = R00

                                   LASTX   P = 9.9999999 99 Gy = R03 ( infinite )

And we also have:

    R04 = R_min = 0.000000000 Gly
    R06 = R_max = 9.9999999 99 Gly ( infinite )

Example3:   H₀ = 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   R₀ = 45.47133084 Gly   = R05
                                     RDN    k =   +1                          = R00

                                   LASTX   P = 301.6226862 Gy = R03

And we also have:

    R04 = R_min = 0.000000000 Gly
    R06 = R_max = 150.8113431 Gly

-We can use the program listed in §6°)a) to compute D₀ & D_L

-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 = H₀       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	/	D_L
Y	/	D₀
X	_z	D

Example: H₀ = 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   D₀ = 26.16037018 Gly
RDN D_L = 233.1494318 Gly

-In most cases, we cannot calculate D₀ 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	H₀	D_L
Y	0 < (Omega)lambda < 1	D₀
X	_z	D

Example: H₀ = 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   D₀ = 17.11226659 Gly
    RDN D_L = 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	/	V_R
Z	/	D_L
Y	H₀	D₀
X	_z	D
L	/	z

Example: H₀ = 71 km/s/Mpc   , z = 7

71 ENTER^
7 XEQ "MLN" >>>> D = 12.05025625 Gly ---Execution time = 1.2s---
RDN   D₀ = 28.63748965 Gly
    RDN D_L = 433.8092250 Gly
RDN V_R = 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"

hp41programs

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