1 2c3 = <float64ne> (_c3); 3c4 = <float64ne> (_c4); 4c5 = <float64ne> (_c5); 5c6 = <float64ne> (_c6); 6c7 = <float64ne> (_c7); 7 8 9E = 0; #MAPLE 10 11zh = <float64ne> (Z); 12zl = Z - zh; #MAPLE 13 14polyHorner <float64ne>= c3 + zh * (c4 + zh * (c5 + zh * (c6 + zh * c7))); 15 16ZhSquarehl = zh * zh; #MAPLE 17zhSquareh = <float64ne> (ZhSquarehl); 18zhSquarel = <float64ne> (ZhSquarehl - zhSquareh); 19 20zhSquareHalfh = zhSquareh * (-0.5); #MAPLE 21zhSquareHalfl = zhSquarel * (-0.5); #MAPLE 22ZhSquareHalfhl = ZhSquarehl * (-0.5); #MAPLE 23 24ZhCube <float64ne>= (zh * zhSquareh); 25polyUpper <float64ne>= polyHorner * ZhCube; 26 27temp = <float64ne> (zh * zl); 28T1hl = polyUpper - temp; #MAPLE 29t1h = <float64ne> (T1hl); 30t1l = <float64ne> (T1hl - t1h); 31 32T2 = Z + ZhSquareHalfhl; #MAPLE 33t2h = <float64ne> (T2hl); 34t2l = <float64ne> (T2hl - t2h); 35 36PE = T2hl + T1hl; #MAPLE 37ph = <float64ne> (Phl); 38pl = <float64ne> (Phl - ph); 39 40 41#We can simplify the computations in the function in this case as we know that 42#all operations (add, mult) on 0 (as a double double) are exact. 43 44Loghm = Phl; #MAPLE 45 46logh = <float64ne> (Loghm); 47logm = <float64ne> (Loghm - logh); 48 49#Mathematical definition of the logarithm and the polynomial 50 51Phigher = (c3 + Z * (c4 + Z * (c5 + Z * (c6 + Z * c7)))); #MAPLE 52ZZZ = Z*Z*Z; #MAPLE 53ZZZPhigher = ZZZ * Phigher; #MAPLE 54HZZ = (-0.5*Z*Z); #MAPLE 55ZpHZZ = Z + HZZ; #MAPLE 56P = ZpHZZ + ZZZPhigher; #MAPLE 57 58#We apply the same simplification on the mathematical definition of the log 59Log = Log1pZ; #MAPLE 60 61 62#Additional useful definitions 63 64ZZ = Z*Z; #MAPLE 65ZZZPhigherPzhzl = ZZZPhigher - zh * zl; #MAPLE 66 67HZ = -0.5*Z; #MAPLE 68 69Flzhzl = temp; #MAPLE 70 71{ 72(T2hl - T2) / T2 in [-1b-103,1b-103] 73/\ (Phl - PE) / PE in [-1b-103,1b-103] 74/\ Z in [_zmin,_zmax] 75/\ (P - Log1pZ) / Log1pZ in [-_epsilonApproxQuick,_epsilonApproxQuick] 76/\ ((logh + logm) - Loghm) / Loghm in [-1b-106,1b-106] 77-> 78((logh + logm) - Log) / Log in [-5b-65,5b-65] 79} 80 81T2hl - T2 -> ((T2hl - T2) / T2) * T2; 82T2hl -> (T2hl - T2) + T2; 83 84Phl - PE -> ((Phl - PE) / PE) * PE; 85Phl -> (Phl - PE) + PE; 86 87 88(ZhSquarehl - ZZ) / ZZ -> 2 * ((zh - Z) / Z) + ((zh - Z) / Z) * ((zh - Z) / Z); 89 90(zhSquareh - ZZ) / ZZ -> ((ZhSquarehl - ZZ) / ZZ) + ((zhSquareh - ZhSquarehl) / ZZ); 91 92(zhSquareh - ZhSquarehl) / ZZ -> ((zhSquareh - ZhSquarehl) / ZhSquarehl) * (ZhSquarehl / ZZ); 93 94ZhSquarehl / ZZ -> ((ZhSquarehl - ZZ) / ZZ) + 1; 95 96(ZhCube - ZZZ) / ZZZ -> (((zh * zhSquareh) - ZZZ) / ZZZ) + ((ZhCube - (zh * zhSquareh)) / ZZZ); 97 98((zh * zhSquareh) - ZZZ) / ZZZ -> (1 + ((zh - Z) / Z)) * (1 + ((zhSquareh - ZZ) / ZZ)) - 1; 99 100((ZhCube - (zh * zhSquareh)) / ZZZ) -> ((ZhCube - (zh * zhSquareh)) / (zh * zhSquareh)) * (((zh - Z) / Z) + 1) * (((zhSquareh - ZZ) / ZZ) + 1); 101 102polyHorner / Phigher -> ((polyHorner - Phigher) / Phigher) + 1; 103 104(polyUpper - ZZZPhigher) / ZZZPhigher -> ((polyHorner - Phigher) / Phigher) + ((ZhCube - ZZZ) / ZZZ) * (polyHorner / Phigher) + 105 + ((polyUpper - (polyHorner * ZhCube)) / (polyHorner * ZhCube)) * (polyHorner / Phigher) + 106 + ((ZhCube - ZZZ) / ZZZ) * ((polyUpper - (polyHorner * ZhCube)) / (polyHorner * ZhCube)) * 107 (polyHorner / Phigher); 108 109 110((ZhSquareHalfhl - (zh * zl)) - HZZ) / HZZ -> - ((zh - Z) / Z) * ((zh - Z) / Z); 111 112(ZhSquareHalfhl - HZZ) / HZZ -> (ZhSquarehl - ZZ) / ZZ; 113 114((T2hl - (zh * zl)) - ZpHZZ) / ZpHZZ -> ((HZ * (((ZhSquareHalfhl - (zh * zl)) - HZZ) / HZZ)) + ((T2hl - T2) / T2) 115 + (HZ * ((T2hl - T2) / T2)) 116 + (HZ * ((ZhSquareHalfhl - HZZ) / HZZ) * ((T2hl - T2) / T2))) / (1 + HZ); 117 118(PE - P) / P -> (((1 + HZ) * (((T2hl - (zh * zl)) - ZpHZZ) / ZpHZZ)) + 119 ((1 + ((zh - Z) / Z)) * (Z * ((zh - Z) / Z)) * ((Flzhzl - (zh * zl)) / (zh * zl))) 120 + (ZZ * Phigher * ((polyUpper - ZZZPhigher) / ZZZPhigher))) / (1 + HZ + ZZ * Phigher); 121 122(Phl - P) / P -> ((PE - P) / P) + ((((PE - P) / P) + 1) * ((Phl - PE) / PE)); 123 124(Loghm - Log) / Log -> ((Loghm - P) / P) + ((P - Log) / Log) + ((Loghm - P) / P) * ((P - Log) / Log); 125 126(((logh + logm) - Log) / Log) -> (((logh + logm) - Loghm) / Loghm) + ((Loghm - Log) / Log) + (((logh + logm) - Loghm) / Loghm) * ((Loghm - Log) / Log); 127