@@ -456,7 +456,12 @@ def julian_ephemeris_millennium(julian_ephemeris_century):
456
456
@jcompile (nopython = True )
457
457
def sum_mult_cos_add_mult (arr , x ):
458
458
# shared calculation used for heliocentric longitude, latitude, and radius
459
- s = (arr [:, [0 ]] * np .cos (arr [:, [1 ]] + arr [:, [2 ]] * np .expand_dims (x , axis = 0 ))).sum (axis = 0 )
459
+ if USE_NUMBA :
460
+ s = 0.
461
+ for row in range (arr .shape [0 ]):
462
+ s += arr [row , 0 ] * np .cos (arr [row , 1 ] + arr [row , 2 ] * x )
463
+ else :
464
+ s = (arr [:, [0 ]] * np .cos (arr [:, [1 ]] + arr [:, [2 ]] * np .expand_dims (x , axis = 0 ))).sum (axis = 0 )
460
465
return s
461
466
462
467
@jcompile ('float64(float64)' , nopython = True )
@@ -557,20 +562,41 @@ def moon_ascending_longitude(julian_ephemeris_century):
557
562
nopython = True )
558
563
def longitude_obliquity_nutation (julian_ephemeris_century , x0 , x1 , x2 , x3 , x4 ,
559
564
out ):
560
- a = NUTATION_ABCD_ARRAY [:, [0 ]]
561
- b = NUTATION_ABCD_ARRAY [:, [1 ]]
562
- c = NUTATION_ABCD_ARRAY [:, [2 ]]
563
- d = NUTATION_ABCD_ARRAY [:, [3 ]]
564
- arg = np .radians (
565
- NUTATION_YTERM_ARRAY [:, [0 ]]* np .expand_dims (x0 , axis = 0 ) +
566
- NUTATION_YTERM_ARRAY [:, [1 ]]* np .expand_dims (x1 , axis = 0 ) +
567
- NUTATION_YTERM_ARRAY [:, [2 ]]* np .expand_dims (x2 , axis = 0 ) +
568
- NUTATION_YTERM_ARRAY [:, [3 ]]* np .expand_dims (x3 , axis = 0 ) +
569
- NUTATION_YTERM_ARRAY [:, [4 ]]* np .expand_dims (x4 , axis = 0 )
570
- )
571
-
572
- delta_psi_sum = ((a + b * julian_ephemeris_century ) * np .sin (arg )).sum (axis = 0 )
573
- delta_eps_sum = ((c + d * julian_ephemeris_century ) * np .cos (arg )).sum (axis = 0 )
565
+
566
+ if USE_NUMBA :
567
+ delta_psi_sum = 0.0
568
+ delta_eps_sum = 0.0
569
+ for row in range (NUTATION_YTERM_ARRAY .shape [0 ]):
570
+ a = NUTATION_ABCD_ARRAY [row , 0 ]
571
+ b = NUTATION_ABCD_ARRAY [row , 1 ]
572
+ c = NUTATION_ABCD_ARRAY [row , 2 ]
573
+ d = NUTATION_ABCD_ARRAY [row , 3 ]
574
+ arg = np .radians (
575
+ NUTATION_YTERM_ARRAY [row , 0 ]* x0 +
576
+ NUTATION_YTERM_ARRAY [row , 1 ]* x1 +
577
+ NUTATION_YTERM_ARRAY [row , 2 ]* x2 +
578
+ NUTATION_YTERM_ARRAY [row , 3 ]* x3 +
579
+ NUTATION_YTERM_ARRAY [row , 4 ]* x4
580
+ )
581
+ delta_psi_sum += (a + b * julian_ephemeris_century ) * np .sin (arg )
582
+ delta_eps_sum += (c + d * julian_ephemeris_century ) * np .cos (arg )
583
+
584
+ else :
585
+ a = NUTATION_ABCD_ARRAY [:, [0 ]]
586
+ b = NUTATION_ABCD_ARRAY [:, [1 ]]
587
+ c = NUTATION_ABCD_ARRAY [:, [2 ]]
588
+ d = NUTATION_ABCD_ARRAY [:, [3 ]]
589
+ arg = np .radians (
590
+ NUTATION_YTERM_ARRAY [:, [0 ]]* np .expand_dims (x0 , axis = 0 ) +
591
+ NUTATION_YTERM_ARRAY [:, [1 ]]* np .expand_dims (x1 , axis = 0 ) +
592
+ NUTATION_YTERM_ARRAY [:, [2 ]]* np .expand_dims (x2 , axis = 0 ) +
593
+ NUTATION_YTERM_ARRAY [:, [3 ]]* np .expand_dims (x3 , axis = 0 ) +
594
+ NUTATION_YTERM_ARRAY [:, [4 ]]* np .expand_dims (x4 , axis = 0 )
595
+ )
596
+
597
+ delta_psi_sum = ((a + b * julian_ephemeris_century ) * np .sin (arg )).sum (axis = 0 )
598
+ delta_eps_sum = ((c + d * julian_ephemeris_century ) * np .cos (arg )).sum (axis = 0 )
599
+
574
600
delta_psi = delta_psi_sum * 1.0 / 36000000
575
601
delta_eps = delta_eps_sum * 1.0 / 36000000
576
602
# seems like we ought to be able to return a tuple here instead
0 commit comments