32real ln2_1067(6505485212531678.0/8796093022208.0);
34int B_lnx2y2_1[22] = {21, 31, 51, 101, 151, 201, 251, 301, 351, 401,
35 451, 501, 551, 601, 651, 701, 751, 801, 851,
37int B_lnx2y2_2[22] = {-1021,-949,-899,-849,-799,-749,-699,-649,-599,
38 -549,-499,-449,-399,-349,-299,-249,-199,-149,
41int B_lnx2y2_N1[21] ={20, 40, 61, 122, 160, 229, 259, 320, 366, 427,
42 488, 518, 549, 610, 671, 732, 763, 825, 885,
46 7804143460206699.0 / 562949953421312.0,
47 7804143460206699.0 / 281474976710656.0,
48 5950659388407608.0 / 140737488355328.0,
49 5950659388407608.0 / 70368744177664.0,
50 7804143460206699.0 / 70368744177664.0,
51 5584840163710419.0 / 35184372088832.0,
52 6316478613104797.0 / 35184372088832.0,
53 7804143460206699.0 / 35184372088832.0,
54 8925989082611412.0 / 35184372088832.0,
55 5206826964856657.0 / 17592186044416.0,
56 5950659388407608.0 / 17592186044416.0,
57 6316478613104797.0 / 17592186044416.0,
58 6694491811958559.0 / 17592186044416.0,
59 7438324235509510.0 / 17592186044416.0,
60 8182156659060461.0 / 17592186044416.0,
61 8925989082611412.0 / 17592186044416.0,
62 4652001140732587.0 / 8796093022208.0,
63 5030014339586349.0 / 8796093022208.0,
64 5395833564283538.0 / 8796093022208.0,
65 5761652788980727.0 / 8796093022208.0,
66 5950659388407608.0 / 8796093022208.0
69int uint_trail(
const unsigned int& n)
87void sqr2uv(
const real& x, real& u, real& v)
119const real expmx2_x0 = 7491658466053896.0 / 281474976710656.0;
134 if (ex<=-26) res = 1;
139 res = 1-u*( 1-v*(1-u/3) );
140 }
else if (t<=expmx2_x0) {
164 if (ex<=-26) res = 1;
170 res = 1+u*( 1+v*(1+u/3) );
191 real t(x),u,v,y,res(0);
195 if (t>=6.5) res =
expx2(t);
207 if (ex>=-8) res =
expm1(u) + v*
exp(u);
215 if(ex>=-510) res = u;
218 std::cerr <<
"expx2m1: denormalized range!"
219 << std::endl; exit(1);
233 if (expo(x) > 33)
return abs(x);
234 else return sqrt(1+x*x);
239const real sqrtp1m1_s = 9007199254740984.0 / 1125899906842624.0;
248 else if (ex>=105) y =
sqrt(x);
249 else if (ex>=53) y =
sqrt(x)-1;
250 else if (x>-0.5234375 && x<=sqrtp1m1_s ) y = x / (
sqrt(x+1) + 1);
251 else y =
sqrt(x+1)-1;
265 exa = expo(a); exb = expo(b);
269 ex = exa; exa = exb; exb = ex;
289{
const real c1 = 1.000732421875,
292 real res,ep,ep2,s1,s2,x1,x2,
arg=x;
300 s1 =
Cut26(res); s2 = res - s1;
302 arg += (ep2 - s2*(res+s1));
312 x2 = x1*x1; ep = x2-1;
314 s1 =
Cut26(res); s2 = res - s1;
315 ep2 = ep2 - s2 * (res+s1);
316 if (
arg<c3) ep -= s1*s1;
341 if (sign(t)<0) t = -t;
342 if (t>1) cxscthrow(STD_FKT_OUT_OF_DEF(
"real sqrt1mx2(const real&)"));
345 if (ex<=-26) res = 1;
366int Interval_Nr(
int* v,
const int& n,
const int& ex)
373 if (ex < v[k]) j = k-1;
391 a = sign(x)<0 ? -x : x;
392 b = sign(y)<0 ? -y : y;
393 int exa=expo(a), exb=expo(b), ex;
397 ex = exa; exa = exb; exb = ex;
401 cxscthrow(STD_FKT_OUT_OF_DEF
402 (
"real ln_sqrtx2y2(const real&, const real&)"));
405 j = Interval_Nr(B_lnx2y2_1,21,exa);
419 r =
ln(a) + B_lnx2y2_c1[j];
435 j = 20 - Interval_Nr(B_lnx2y2_2,21,exa);
454 if (s>=0.25 && s<=1.75)
455 if (s>=0.828125 && s<=1.171875)
457 if (a==1 && exb<=-28)
486typedef union {
double f;
char intern[8]; } help_real;
496 y.intern[0] = y.intern[1] = y.intern[2] = 0;
499 y.intern[7] = y.intern[6] = y.intern[5] = 0;
512 y.intern[0] = y.intern[1] = y.intern[2] = 0;
515 y.intern[7] = y.intern[6] = y.intern[5] = 0;
528 y.intern[0] = y.intern[1] = y.intern[2] = 0;
531 y.intern[7] = y.intern[6] = y.intern[5] = 0;
555 return dbl<0 ? int(dbl-0.5) : int(dbl+0.5);
563 int res = int( _double(y) );
565 if (!neg && y>0) res += 1;
575 int res = int( _double(y) );
585static real q_acoshp1[5] =
588 -6004799503160661.0 / 72057594037927936.0,
589 +5404319552844595.0 / 288230376151711744.0,
590 -6433713753386423.0 / 1152921504606846976.0,
591 +8756999275442631.0 / 4611686018427387904.0
594static const real c_ln2_B = 6243314768165359.0 / 9007199254740992.0;
605 cxscthrow(STD_FKT_OUT_OF_DEF(
"real acoshp1(const real&)"));
607 if (ex<=-50) res =
sqrt(2*x);
609 res =
sqrt(2*x)*((((q_acoshp1[4]*x+q_acoshp1[3])*x+q_acoshp1[2])
610 *x+q_acoshp1[1])*x + q_acoshp1[0]);
611 else if (ex<=0) res =
lnp1(x+
sqrt(2*x+x*x));
612 else if (ex<=50) res =
lnp1(x*(1+
sqrt(1+2/x)));
613 else if (ex<=1022) res =
ln(2*x);
614 else res =
ln(x) + c_ln2_B;
625real a_sinpix_pi[8] = {0,
626 7.450580596923828125e-9,
638 const real q_sin1_a[3] = {0.0,
639 -7408124450506663.0 / 4503599627370496.0,
640 +4595816915236340.0 / 36028797018963968.0 };
642 const real q_sin1_b[3] = {+1.0,
643 +8889749340277122.0 / 18014398509481984.0,
644 -6103596185201565.0 / 36028797018963968.0};
650 const real q_sin2_a[5] = {0.0,
651 -4829370543434630.0 / 18014398509481984.0,
652 +5229154385037560.0 / 140737488355328.0,
653 -6000119407941526.0 / 576460752303423488.0,
654 +4592842702840413.0 / 1125899906842624.0 };
656 const real q_sin2_b[5] = {-6359670668682560.0 / 576460752303423488.0,
657 -6771299865719208.0 / 1125899906842624.0,
658 +6861350950372367.0 / 1125899906842624.0,
659 -4669728680615660.0 / 2251799813685248.0,
660 +4607746187351458.0 / 2251799813685248.0 };
666 const real q_sin3_a[5] = {0.0,
667 -7954137925457403.0 / 18014398509481984.0,
668 +7534721222201852.0 / 562949953421312.0,
669 -8481413201661564.0 / 288230376151711744.0,
670 +6484759126571114.0 / 4503599627370496.0};
671 const real q_sin3_b[5] = {-8780869688499296.0 / 288230376151711744.0,
672 -7929689932517363.0 / 2251799813685248.0,
673 +8223190551048856.0 / 2251799813685248.0,
674 -5789057670884168.0 / 4503599627370496.0,
675 +5586166285064551.0 / 4503599627370496.0 };
681 const real q_sin4_a[6] = {+0.0,
682 -5469389903424399.0 / 9007199254740992.0,
683 +7704524171392577.0 / 1125899906842624.0,
684 -8513632445057187.0 / 144115188075855872.0,
685 +6447776473700967.0 / 9007199254740992.0,
686 -7362728796775514.0 / 288230376151711744.0 };
688 const real q_sin4_b[6] = {-8529339040415872.0 / 144115188075855872.0,
689 -5452397600582528.0 / 2251799813685248.0,
690 +5848879160432416.0 / 2251799813685248.0,
691 -8609222354677976.0 / 9007199254740992.0,
692 +8031029418839465.0 / 9007199254740992.0,
693 -4904826237544497.0 / 9007199254740992.0 };
699 const real q_sin5_a[6] = {0.0,
700 4*(5276398025110506.0 / 9007199254740992.0),
701 +7169490078506431.0 / 1125899906842624.0,
702 -7877917296929161.0 / 36028797018963968.0,
703 +5470058796541645.0 / 9007199254740992.0,
704 -7233501679212276.0 / 144115188075855872.0 };
706 const real q_sin5_b[6] = {+4598161403289824.0 / 144115188075855872.0,
707 +4893639363223553.0 / 2251799813685248.0,
708 -5525704702280927.0 / 2251799813685248.0,
709 +4608444394217766.0 / 4503599627370496.0,
710 -8018266172128683.0 / 9007199254740992.0,
711 +5822429161405618.0 / 9007199254740992.0 };
717 const real q_sin6_a[6] = {0.0,
718 +7461973733147728.0 / 36028797018963968.0,
719 +7979710108639590.0 / 140737488355328.0,
720 -6289496525317218.0 / 288230376151711744.0,
721 +6964694082726429.0 / 1125899906842624.0,
722 -5522978779702360.0 / 1152921504606846976.0 };
724 const real q_sin6_b[6] = {+5609829449006976.0 / 18014398509481984.0,
725 +8354019229473476.0 / 1125899906842624.0,
726 -8475200820952730.0 / 1125899906842624.0,
727 +5829377160898302.0 / 2251799813685248.0,
728 -5714036688191727.0 / 2251799813685248.0,
729 +7900583259891052.0 / 4503599627370496.0 };
759real sinpi_A1(
const real& x)
765 y = q_sin1_a[2] / ( v + q_sin1_b[2]);
766 y = q_sin1_a[1] / ( (v + q_sin1_b[1]) + y);
776real sinpi_A2(
const real& x)
784 v = 1/(x-0.08203125);
786 y = q_sin2_a[4] / ( v + q_sin2_b[4]);
787 y = q_sin2_a[3] / ( (v + q_sin2_b[3]) + y);
788 y = q_sin2_a[2] / ( (v + q_sin2_b[2]) + y);
789 y = q_sin2_a[1] / ( (v + q_sin2_b[1]) + y) + q_sin2_b[0];
800real sinpi_A3(
const real& x)
808 v = 1/(x-0.13671875);
810 y = q_sin3_a[4] / ( v + q_sin3_b[4]);
811 y = q_sin3_a[3] / ( (v + q_sin3_b[3]) + y);
812 y = q_sin3_a[2] / ( (v + q_sin3_b[2]) + y);
813 y = q_sin3_a[1] / ( (v + q_sin3_b[1]) + y) + q_sin3_b[0];
824real sinpi_A4(
const real& x)
832 v = 1/(x-0.19140625);
834 y = q_sin4_a[5] / ( v + q_sin4_b[5]);
835 y = q_sin4_a[4] / ( (v + q_sin4_b[4]) + y);
836 y = q_sin4_a[3] / ( (v + q_sin4_b[3]) + y);
837 y = q_sin4_a[2] / ( (v + q_sin4_b[2]) + y);
838 y = q_sin4_a[1] / ( (v + q_sin4_b[1]) + y) + q_sin4_b[0];
849real sinpi_A5(
const real& x)
857 v = 1/(x-0.30078125);
859 y = q_sin5_a[5] / ( v + q_sin5_b[5]);
860 y = q_sin5_a[4] / ( (v + q_sin5_b[4]) + y);
861 y = q_sin5_a[3] / ( (v + q_sin5_b[3]) + y);
862 y = q_sin5_a[2] / ( (v + q_sin5_b[2]) + y);
863 y = q_sin5_a[1] / ( (v + q_sin5_b[1]) + y) + q_sin5_b[0];
874real sinpi_A6(
const real& x)
882 v = 1/(x-0.43359375);
884 y = q_sin6_a[5] / ( v + q_sin6_b[5]);
885 y = q_sin6_a[4] / ( (v + q_sin6_b[4]) + y);
886 y = q_sin6_a[3] / ( (v + q_sin6_b[3]) + y);
887 y = q_sin6_a[2] / ( (v + q_sin6_b[2]) + y);
888 y = q_sin6_a[1] / ( (v + q_sin6_b[1]) + y) + q_sin6_b[0];
901 if (m == -2147483648)
902 cxscthrow(STD_FKT_OUT_OF_DEF(
"real sinpix_pi(const real&)"));
908 nr = int_no(a_sinpix_pi,8,xr);
912 case 0: res = xr;
break;
913 case 1: res = sinpi_A1(xr);
break;
914 case 2: res = sinpi_A2(xr);
break;
915 case 3: res = sinpi_A3(xr);
break;
916 case 4: res = sinpi_A4(xr);
break;
917 case 5: res = sinpi_A5(xr);
break;
918 case 6: res = sinpi_A6(xr);
break;
919 case 7: res = 5734161139222659.0 / 18014398509481984.0;
923 if (m%2 != 0) res = -res;
935{-0.5, 8.5, 16.5, 24.5, 35.5, 46.5, 57.5, 68.5, 79.5, 90.5,
936 101.5, 112.5, 122.5, 132.5, 142.5, 150.5, 157.5, 164.5, 171.5 };
939 const real q_gams0_a[8] =
941 -7616205496030159.0 / 18014398509481984.0,
942 6808968414757234.0 / 9007199254740992.0,
943 -7179656013820698.0 / 144115188075855872.0,
944 7646522333236288.0 / 144115188075855872.0,
945 -7301751514589296.0 / 144115188075855872.0,
946 6062274112599236.0 / 288230376151711744.0,
947 -7580146497393670.0 / 576460752303423488.0 };
949 const real q_gams0_b[8] =
951 -4965940208012024.0 / 9007199254740992.0,
952 8627036189024519.0 / 9007199254740992.0,
953 -7819065539096245.0 / 72057594037927936.0,
954 7433219784613049.0 / 36028797018963968.0,
955 7389378817674059.0 / 72057594037927936.0,
956 -6060163955665826.0 / 144115188075855872.0,
957 5769165265609039.0 / 18014398509481984.0 };
960 const real q_gams1_a[7] =
962 5262165366811426.0 / 72057594037927936.0,
963 5574236285326272.0 / 9007199254740992.0,
964 -8823729115230752.0 / 9223372036854775808.0,
965 4639548867796070.0 / 72057594037927936.0,
966 -4784760610310013.0 / 4611686018427387904.0,
967 5900435828568799.0 / 288230376151711744.0 };
969 const real q_gams1_b[7] =
970 {7108556377252296.0 / 36028797018963968.0,
971 -7219014979973550.0 / 9007199254740992.0,
972 7224885284258947.0 / 9007199254740992.0,
973 -8569030245421939.0 / 36028797018963968.0,
974 4838653997792226.0 / 18014398509481984.0,
975 -4567015707194151.0 / 36028797018963968.0,
976 5705105513288442.0 / 36028797018963968.0 };
979 const real q_gams2_a[6] =
981 7322606177885925.0 / 72057594037927936.0,
982 5856515731454725.0 / 144115188075855872.0,
983 -5420981868254561.0 / 1152921504606846976.0,
984 6209135328051357.0 / 1152921504606846976.0,
985 -8261990134869242.0 / 2305843009213693952.0 };
987 const real q_gams2_b[6] =
988 {-5267325780498998.0 / 18014398509481984.0,
989 -4688643295042637.0 / 18014398509481984.0,
990 6865435581764388.0 / 36028797018963968.0,
991 -5147410677045000.0 / 144115188075855872.0,
992 5878944809110092.0 / 144115188075855872.0,
993 5315270486582075.0 / 576460752303423488.0 };
996 const real q_gams3_a[6] =
998 4608894695736103.0 / 9007199254740992.0,
999 4622056591672246.0 / 144115188075855872.0,
1000 7936321632956658.0 / 1152921504606846976.0,
1001 7533691379187725.0 / 2305843009213693952.0,
1002 8317076627520632.0 / 4611686018427387904.0 };
1004 const real q_gams3_b[6] =
1005 {-5793090370845400.0 / 36028797018963968.0,
1006 -5993724737769479.0 / 18014398509481984.0,
1007 -7355055148590989.0 / 288230376151711744.0,
1008 -6811362569175008.0 / 288230376151711744.0,
1009 -5643865380181690.0 / 288230376151711744.0,
1010 -7257737463164532.0 / 288230376151711744.0 };
1013 const real q_gams4_a[6] =
1015 -7504135817814391.0 / 18014398509481984.0,
1016 4990516833351222.0 / 576460752303423488.0,
1017 7409944094502278.0 / 4611686018427387904.0,
1018 8239437277348520.0 / 2305843009213693952.0,
1019 -5361398870955224.0 / 4611686018427387904.0 };
1021 const real q_gams4_b[6] =
1022 {7405548010914168.0 / 18014398509481984.0,
1023 6819421126036941.0 / 36028797018963968.0,
1024 8474998846370093.0 / 288230376151711744.0,
1025 7733944242652532.0 / 144115188075855872.0,
1026 -6547033063052812.0 / 144115188075855872.0,
1027 -6636555235967309.0 / 576460752303423488.0 };
1030 const real q_gams5_a[6] =
1032 -7282435522169731.0 / 144115188075855872.0,
1033 7812862759300002.0 / 288230376151711744.0,
1034 -7591789799457117.0 / 9223372036854775808.0,
1035 7290376444377792.0 / 2305843009213693952.0,
1036 -7051827298983325.0 / 9223372036854775808.0 };
1038 const real q_gams5_b[6] =
1039 {-4692552597358020.0 / 36028797018963968.0,
1040 7065248316566075.0 / 36028797018963968.0,
1041 -5667667511562837.0 / 36028797018963968.0,
1042 7741409329322298.0 / 144115188075855872.0,
1043 -6371025742549657.0 / 144115188075855872.0,
1044 8045016330906311.0 / 288230376151711744.0 };
1047 const real q_gams6_a[6] =
1049 6719861857699617.0 / 36028797018963968.0,
1050 6146108740657636.0 / 1152921504606846976.0,
1051 -8996491371873855.0 / 4611686018427387904.0,
1052 5082150623010284.0 / 2305843009213693952.0,
1053 -4650554434211955.0 / 18446744073709551616.0 };
1055 const real q_gams6_b[6] =
1056 {6795471792542416.0 / 18014398509481984.0,
1057 -4567367130077106.0 / 36028797018963968.0,
1058 8403023343856889.0 / 288230376151711744.0,
1059 5083390457138636.0 / 144115188075855872.0,
1060 -4614092542204999.0 / 144115188075855872.0,
1061 7104477795873427.0 / 72057594037927936.0 };
1064 const real q_gams7_a[6] =
1066 5372276754422009.0 / 18014398509481984.0,
1067 4663470139182266.0 / 576460752303423488.0,
1068 8173857739398457.0 / 4611686018427387904.0,
1069 4932343261664244.0 / 4611686018427387904.0,
1070 -5822870305854791.0 / 295147905179352825856.0 };
1072 const real q_gams7_b[6] =
1073 {-4806355488125952.0 / 1152921504606846976.0,
1074 -6211400907258787.0 / 36028797018963968.0,
1075 -5429997192499743.0 / 288230376151711744.0,
1076 -5653057298064211.0 / 288230376151711744.0,
1077 -4746908442350821.0 / 1152921504606846976.0,
1078 8670733620601602.0 / 18014398509481984.0 };
1081 const real q_gams8_a[5] =
1083 -8754830070807807.0 / 72057594037927936.0,
1084 7931975738629914.0 / 2305843009213693952.0,
1085 6414800170661986.0 / 73786976294838206464.0,
1086 4558538137025832.0 / 18014398509481984.0};
1088 const real q_gams8_b[5] =
1089 {-4924952929015874.0 / 18014398509481984.0,
1090 8571263173618757.0 / 72057594037927936.0,
1091 7912939311540494.0 / 576460752303423488.0,
1092 8975599623306617.0 / 18014398509481984.0,
1093 -4569152133381960.0 / 9007199254740992.0 };
1096 const real q_gams9_a[5] =
1098 -6140046099729716.0 / 144115188075855872.0,
1099 7278997066304035.0 / 576460752303423488.0,
1100 -4759432412103962.0 / 9223372036854775808.0,
1101 6912711986126027.0 / 4611686018427387904.0 };
1103 const real q_gams9_b[5] =
1104 {-5611185402650416.0 / 72057594037927936.0,
1105 4915638659703209.0 / 36028797018963968.0,
1106 -7692469721238990.0 / 72057594037927936.0,
1107 5045754589592299.0 / 144115188075855872.0,
1108 -8290188163752846.0 / 288230376151711744.0 };
1111 const real q_gams10_a[5] =
1113 4717576521494070.0 / 72057594037927936.0,
1114 6906626643101502.0 / 1152921504606846976.0,
1115 -8147816895338065.0 / 9223372036854775808.0,
1116 7960901442281121.0 / 9223372036854775808.0 };
1118 const real q_gams10_b[5] =
1119 {7952058244493952.0 / 288230376151711744.0,
1120 -7601355771801470.0 / 72057594037927936.0,
1121 4923655207606594.0 / 72057594037927936.0,
1122 -7238499638430702.0 / 576460752303423488.0,
1123 5852764550899526.0 / 576460752303423488.0 };
1126 const real q_gams11_a[5] =
1128 5000157748740957.0 / 36028797018963968.0,
1129 6733781271808263.0 / 2305843009213693952.0,
1130 4912448110743899.0 / 18446744073709551616.0,
1131 6917678469221762.0 / 576460752303423488.0 };
1133 const real q_gams11_b[5] =
1134 {-4805173503400964.0 / 36028797018963968.0,
1135 -7686561799456867.0 / 72057594037927936.0,
1136 -6649023763933472.0 / 576460752303423488.0,
1137 -7324214555543306.0 / 72057594037927936.0,
1138 8284786467509721.0 / 72057594037927936.0 };
1141 const real q_gams12_a[5] =
1143 5226845396890227.0 / 36028797018963968.0,
1144 5602232597459763.0 / 1152921504606846976.0,
1145 4920783552645014.0 / 4611686018427387904.0,
1146 5570041179715558.0 / 9223372036854775808.0 };
1148 const real q_gams12_b[5] =
1149 {-6799658092196962.0 / 18014398509481984.0,
1150 -4810318281953418.0 / 36028797018963968.0,
1151 -8340064868499919.0 / 576460752303423488.0,
1152 -8375346434191481.0 / 576460752303423488.0,
1153 -6231299816839781.0 / 1152921504606846976.0 };
1156 const real q_gams13_a[5] =
1158 5077602289066213.0 / 18014398509481984.0,
1159 4971389438544026.0 / 576460752303423488.0,
1160 8326326060413143.0 / 4611686018427387904.0,
1161 7582214946963028.0 / 9223372036854775808.0 };
1163 const real q_gams13_b[5] =
1164 {-8336661118182000.0 / 72057594037927936.0,
1165 -6152827173580884.0 / 36028797018963968.0,
1166 -6306415240260295.0 / 576460752303423488.0,
1167 -5949724688158565.0 / 576460752303423488.0,
1168 -5625482009313801.0 / 576460752303423488.0 };
1171 const real q_gams14_a[6] =
1173 8624518348171320.0 / 18014398509481984.0,
1174 7049908374954842.0 / 576460752303423488.0,
1175 5769548496333642.0 / 2305843009213693952.0,
1176 5110826386599631.0 / 4611686018427387904.0,
1177 5902099580937297.0 / 9223372036854775808.0 };
1179 const real q_gams14_b[6] =
1180 {4591842870321392.0 / 18014398509481984.0,
1181 -7195103997153274.0 / 36028797018963968.0,
1182 -5062133254875404.0 / 576460752303423488.0,
1183 -4899776496053167.0 / 576460752303423488.0,
1184 -4703735605739871.0 / 576460752303423488.0,
1185 -8717622510435264.0 / 1152921504606846976.0 };
1188 const real q_gams15_a[6] =
1190 5047093325633351.0 / 9007199254740992.0,
1191 8838568428828895.0 / 576460752303423488.0,
1192 7169875770591066.0 / 2305843009213693952.0,
1193 6277116142443993.0 / 4611686018427387904.0,
1194 7160380168510556.0 / 9223372036854775808.0 };
1196 const real q_gams15_b[6] =
1197 {5575895624898616.0 / 18014398509481984.0,
1198 -7982725134478240.0 / 36028797018963968.0,
1199 -8683131154525203.0 / 1152921504606846976.0,
1200 -8504376590080638.0 / 1152921504606846976.0,
1201 -8272396567408614.0 / 1152921504606846976.0,
1202 -7370684458128734.0 / 1152921504606846976.0 };
1205 const real q_gams16_a[6] =
1207 8928052213955455.0 / 18014398509481984.0,
1208 5405751921217612.0 / 288230376151711744.0,
1209 8725861134087760.0 / 2305843009213693952.0,
1210 7583718370137025.0 / 4611686018427387904.0,
1211 8577302812898546.0 / 9223372036854775808.0 };
1213 const real q_gams16_b[6] =
1214 {6669445708872192.0 / 144115188075855872.0,
1215 -8769965967955460.0 / 36028797018963968.0,
1216 -7524042392635917.0 / 1152921504606846976.0,
1217 -7422547685779587.0 / 1152921504606846976.0,
1218 -7284571995868819.0 / 1152921504606846976.0,
1219 -7797018373646746.0 / 1152921504606846976.0 };
1222 const real q_gams17_a[6] =
1224 4642907317603488.0 / 9007199254740992.0,
1225 6403634366308572.0 / 288230376151711744.0,
1226 5153510478021874.0 / 1152921504606846976.0,
1227 8919203526515161.0 / 4611686018427387904.0,
1228 5015558226948384.0 / 4611686018427387904.0 };
1230 const real q_gams17_b[6] =
1231 {-6229620043098112.0 / 9223372036854775808.0,
1232 -4750296278401558.0 / 18014398509481984.0,
1233 -6639785613322219.0 / 1152921504606846976.0,
1234 -6577910457226093.0 / 1152921504606846976.0,
1235 -6491061443859352.0 / 1152921504606846976.0,
1236 -6372399542597081.0 / 1152921504606846976.0 };
1242 real gam_S0(
const real& x)
1256 y = q_gams0_a[7] / ( v + q_gams0_b[7]);
1257 y = q_gams0_a[6] / ( (v + q_gams0_b[6]) + y );
1258 y = q_gams0_a[5] / ( (v + q_gams0_b[5]) + y );
1259 y = q_gams0_a[4] / ( (v + q_gams0_b[4]) + y );
1260 y = q_gams0_a[3] / ( (v + q_gams0_b[3]) + y );
1261 y = q_gams0_a[2] / ( (v + q_gams0_b[2]) + y );
1262 y = q_gams0_a[1] / ( (v + q_gams0_b[1]) + y ) + q_gams0_b[0];
1267 real gam_S0_n0(
const real& x)
1281 y = q_gams0_a[7] / ( v + q_gams0_b[7]);
1282 y = q_gams0_a[6] / ( (v + q_gams0_b[6]) + y );
1283 y = q_gams0_a[5] / ( (v + q_gams0_b[5]) + y );
1284 y = q_gams0_a[4] / ( (v + q_gams0_b[4]) + y );
1285 y = q_gams0_a[3] / ( (v + q_gams0_b[3]) + y );
1286 y = q_gams0_a[2] / ( (v + q_gams0_b[2]) + y );
1287 y = q_gams0_a[1] / ( (v + q_gams0_b[1]) + y ) + q_gams0_b[0];
1293real gam_S0_n1(
const real& x)
1307 y = q_gams0_a[7] / ( v + q_gams0_b[7]);
1308 y = q_gams0_a[6] / ( (v + q_gams0_b[6]) + y );
1309 y = q_gams0_a[5] / ( (v + q_gams0_b[5]) + y );
1310 y = q_gams0_a[4] / ( (v + q_gams0_b[4]) + y );
1311 y = q_gams0_a[3] / ( (v + q_gams0_b[3]) + y );
1312 y = q_gams0_a[2] / ( (v + q_gams0_b[2]) + y );
1313 y = q_gams0_a[1] / ( (v + q_gams0_b[1]) + y ) + q_gams0_b[0];
1319real gammar_S0(
const real& x)
1331 case 0:
if (expo(x)<=-52) y = x;
1332 else y = x*(x+1)*gam_S0_n0(x);
break;
1333 case 1: y = x*gam_S0_n1(x);
break;
1334 case 2: y = gam_S0(x);
break;
1338 for (
int k=2; k<=p; k++) Ne *= (x-k);
1339 y = gam_S0(x-p) / Ne;
1345real gamma_S0(
const real& x)
1350 real y = 1/gammar_S0(x);
1354real gam_S1(
const real& x)
1368 y = q_gams1_a[6] / ( v + q_gams1_b[6]);
1369 y = q_gams1_a[5] / ( (v + q_gams1_b[5]) + y );
1370 y = q_gams1_a[4] / ( (v + q_gams1_b[4]) + y );
1371 y = q_gams1_a[3] / ( (v + q_gams1_b[3]) + y );
1372 y = q_gams1_a[2] / ( (v + q_gams1_b[2]) + y );
1373 y = q_gams1_a[1] / ( (v + q_gams1_b[1]) + y ) + q_gams1_b[0];
1382real gamma_S1(
const real& x)
1396 for (
int k=2; k<=p; k++)
1403 for (
int k=1; k<=p-1; k++)
1405 y = (p==0)? gam_S1(x) : gam_S1(x+p)/Pr;
1411real gam_S2(
const real& x)
1425 y = q_gams2_a[5] / ( v + q_gams2_b[5]);
1426 y = q_gams2_a[4] / ( (v + q_gams2_b[4]) + y );
1427 y = q_gams2_a[3] / ( (v + q_gams2_b[3]) + y );
1428 y = q_gams2_a[2] / ( (v + q_gams2_b[2]) + y );
1429 y = q_gams2_a[1] / ( (v + q_gams2_b[1]) + y ) + q_gams2_b[0];
1438real gamma_S2(
const real& x)
1452 for (
int k=2; k<=p; k++)
1459 for (
int k=1; k<=p-1; k++)
1461 y = (p==0)? gam_S2(x) : gam_S2(x+p)/Pr;
1467real gam_S3(
const real& x)
1481 y = q_gams3_a[5] / ( v + q_gams3_b[5]);
1482 y = q_gams3_a[4] / ( (v + q_gams3_b[4]) + y );
1483 y = q_gams3_a[3] / ( (v + q_gams3_b[3]) + y );
1484 y = q_gams3_a[2] / ( (v + q_gams3_b[2]) + y );
1485 y = q_gams3_a[1] / ( (v + q_gams3_b[1]) + y ) + q_gams3_b[0];
1494real gamma_S3(
const real& x)
1508 for (
int k=2; k<=p; k++)
1515 for (
int k=1; k<=p-1; k++)
1517 y = (p==0)? gam_S3(x) : gam_S3(x+p)/Pr;
1523real gam_S4(
const real& x)
1535 v = 1/(x-41.140625);
1537 y = q_gams4_a[5] / ( v + q_gams4_b[5]);
1538 y = q_gams4_a[4] / ( (v + q_gams4_b[4]) + y );
1539 y = q_gams4_a[3] / ( (v + q_gams4_b[3]) + y );
1540 y = q_gams4_a[2] / ( (v + q_gams4_b[2]) + y );
1541 y = q_gams4_a[1] / ( (v + q_gams4_b[1]) + y ) + q_gams4_b[0];
1550real gamma_S4(
const real& x)
1564 for (
int k=2; k<=p; k++)
1571 for (
int k=1; k<=p-1; k++)
1573 y = (p==0)? gam_S4(x) : gam_S4(x+p)/Pr;
1579real gam_S5(
const real& x)
1591 v = 1/(x-52.015625);
1593 y = q_gams5_a[5] / ( v + q_gams5_b[5]);
1594 y = q_gams5_a[4] / ( (v + q_gams5_b[4]) + y );
1595 y = q_gams5_a[3] / ( (v + q_gams5_b[3]) + y );
1596 y = q_gams5_a[2] / ( (v + q_gams5_b[2]) + y );
1597 y = q_gams5_a[1] / ( (v + q_gams5_b[1]) + y ) + q_gams5_b[0];
1606real gamma_S5(
const real& x)
1620 for (
int k=2; k<=p; k++)
1627 for (
int k=1; k<=p-1; k++)
1629 y = (p==0)? gam_S5(x) : gam_S5(x+p)/Pr;
1635real gam_S6(
const real& x)
1647 v = 1/(x-63.015625);
1649 y = q_gams6_a[5] / ( v + q_gams6_b[5]);
1650 y = q_gams6_a[4] / ( (v + q_gams6_b[4]) + y );
1651 y = q_gams6_a[3] / ( (v + q_gams6_b[3]) + y );
1652 y = q_gams6_a[2] / ( (v + q_gams6_b[2]) + y );
1653 y = q_gams6_a[1] / ( (v + q_gams6_b[1]) + y ) + q_gams6_b[0];
1662real gamma_S6(
const real& x)
1676 for (
int k=2; k<=p; k++)
1683 for (
int k=1; k<=p-1; k++)
1685 y = (p==0)? gam_S6(x) : gam_S6(x+p)/Pr;
1692real gam_S7(
const real& x)
1704 v = 1/(x-74.16015625);
1706 y = q_gams7_a[5] / ( v + q_gams7_b[5]);
1707 y = q_gams7_a[4] / ( (v + q_gams7_b[4]) + y );
1708 y = q_gams7_a[3] / ( (v + q_gams7_b[3]) + y );
1709 y = q_gams7_a[2] / ( (v + q_gams7_b[2]) + y );
1710 y = q_gams7_a[1] / ( (v + q_gams7_b[1]) + y ) + q_gams7_b[0];
1719real gamma_S7(
const real& x)
1733 for (
int k=2; k<=p; k++)
1740 for (
int k=1; k<=p-1; k++)
1742 y = (p==0)? gam_S7(x) : gam_S7(x+p)/Pr;
1748real gam_S8(
const real& x)
1760 v = 1/(x-85.1015625);
1762 y = q_gams8_a[4] / ( v + q_gams8_b[4]);
1763 y = q_gams8_a[3] / ( (v + q_gams8_b[3]) + y );
1764 y = q_gams8_a[2] / ( (v + q_gams8_b[2]) + y );
1765 y = q_gams8_a[1] / ( (v + q_gams8_b[1]) + y ) + q_gams8_b[0];
1774real gamma_S8(
const real& x)
1788 for (
int k=2; k<=p; k++)
1795 for (
int k=1; k<=p-1; k++)
1797 y = (p==0)? gam_S8(x) : gam_S8(x+p)/Pr;
1804real gam_S9(
const real& x)
1816 v = 1/(x-95.984375);
1818 y = q_gams9_a[4] / ( v + q_gams9_b[4]);
1819 y = q_gams9_a[3] / ( (v + q_gams9_b[3]) + y );
1820 y = q_gams9_a[2] / ( (v + q_gams9_b[2]) + y );
1821 y = q_gams9_a[1] / ( (v + q_gams9_b[1]) + y ) + q_gams9_b[0];
1830real gamma_S9(
const real& x)
1844 for (
int k=2; k<=p; k++)
1851 for (
int k=1; k<=p-1; k++)
1853 y = (p==0)? gam_S9(x) : gam_S9(x+p)/Pr;
1860real gam_S10(
const real& x)
1872 v = 1/(x-107.078125);
1874 y = q_gams10_a[4] / ( v + q_gams10_b[4]);
1875 y = q_gams10_a[3] / ( (v + q_gams10_b[3]) + y );
1876 y = q_gams10_a[2] / ( (v + q_gams10_b[2]) + y );
1877 y = q_gams10_a[1] / ( (v + q_gams10_b[1]) + y ) + q_gams10_b[0];
1886real gamma_S10(
const real& x)
1900 for (
int k=2; k<=p; k++)
1902 y = Pr*gam_S10(x-p);
1907 for (
int k=1; k<=p-1; k++)
1909 y = (p==0)? gam_S10(x) : gam_S10(x+p)/Pr;
1915real gam_S11(
const real& x)
1927 v = 1/(x-117.8671875);
1929 y = q_gams11_a[4] / ( v + q_gams11_b[4]);
1930 y = q_gams11_a[3] / ( (v + q_gams11_b[3]) + y );
1931 y = q_gams11_a[2] / ( (v + q_gams11_b[2]) + y );
1932 y = q_gams11_a[1] / ( (v + q_gams11_b[1]) + y ) + q_gams11_b[0];
1941real gamma_S11(
const real& x)
1955 for (
int k=2; k<=p; k++)
1957 y = Pr*gam_S11(x-p);
1962 for (
int k=1; k<=p-1; k++)
1964 y = (p==0)? gam_S11(x) : gam_S11(x+p)/Pr;
1970real gam_S12(
const real& x)
1982 v = 1/(x-126.7421875);
1984 y = q_gams12_a[4] / ( v + q_gams12_b[4]);
1985 y = q_gams12_a[3] / ( (v + q_gams12_b[3]) + y );
1986 y = q_gams12_a[2] / ( (v + q_gams12_b[2]) + y );
1987 y = q_gams12_a[1] / ( (v + q_gams12_b[1]) + y ) + q_gams12_b[0];
1996real gamma_S12(
const real& x)
2009 for (
int k=2; k<=p; k++)
2011 y = Pr*gam_S12(x-p);
2016 for (
int k=1; k<=p-1; k++)
2018 y = (p==0)? gam_S12(x) : gam_S12(x+p)/Pr;
2024real gam_S13(
const real& x)
2036 v = 1/(x-138.0390625);
2038 y = q_gams13_a[4] / ( v + q_gams13_b[4]);
2039 y = q_gams13_a[3] / ( (v + q_gams13_b[3]) + y );
2040 y = q_gams13_a[2] / ( (v + q_gams13_b[2]) + y );
2041 y = q_gams13_a[1] / ( (v + q_gams13_b[1]) + y ) + q_gams13_b[0];
2050real gamma_S13(
const real& x)
2063 for (
int k=2; k<=p; k++)
2065 y = Pr*gam_S13(x-p);
2070 for (
int k=1; k<=p-1; k++)
2072 y = (p==0)? gam_S13(x) : gam_S13(x+p)/Pr;
2078real gam_S14(
const real& x)
2086 if (x==146.94921875)
2090 v = 1/(x-146.94921875);
2092 y = q_gams14_a[5] / ( v + q_gams14_b[5]);
2093 y = q_gams14_a[4] / ( (v + q_gams14_b[4]) + y );
2094 y = q_gams14_a[3] / ( (v + q_gams14_b[3]) + y );
2095 y = q_gams14_a[2] / ( (v + q_gams14_b[2]) + y );
2096 y = q_gams14_a[1] / ( (v + q_gams14_b[1]) + y ) + q_gams14_b[0];
2105real gamma_S14(
const real& x)
2118 for (
int k=2; k<=p; k++)
2120 y = Pr*gam_S14(x-p);
2125 for (
int k=1; k<=p-1; k++)
2127 y = (p==0)? gam_S14(x) : gam_S14(x+p)/Pr;
2133real gam_S15(
const real& x)
2141 if (x==153.90234375)
2145 v = 1/(x-153.90234375);
2147 y = q_gams15_a[5] / ( v + q_gams15_b[5]);
2148 y = q_gams15_a[4] / ( (v + q_gams15_b[4]) + y );
2149 y = q_gams15_a[3] / ( (v + q_gams15_b[3]) + y );
2150 y = q_gams15_a[2] / ( (v + q_gams15_b[2]) + y );
2151 y = q_gams15_a[1] / ( (v + q_gams15_b[1]) + y ) + q_gams15_b[0];
2162real gamma_S15(
const real& x)
2175 for (
int k=2; k<=p; k++)
2177 y = Pr*gam_S15(x-p);
2182 for (
int k=1; k<=p-1; k++)
2184 y = (p==0)? gam_S15(x) : gam_S15(x+p)/Pr;
2190real gam_S16(
const real& x)
2198 if (x==161.08984375)
2202 v = 1/(x-161.08984375);
2204 y = q_gams16_a[5] / ( v + q_gams16_b[5]);
2205 y = q_gams16_a[4] / ( (v + q_gams16_b[4]) + y );
2206 y = q_gams16_a[3] / ( (v + q_gams16_b[3]) + y );
2207 y = q_gams16_a[2] / ( (v + q_gams16_b[2]) + y );
2208 y = q_gams16_a[1] / ( (v + q_gams16_b[1]) + y ) + q_gams16_b[0];
2217real gamma_S16(
const real& x)
2230 for (
int k=2; k<=p; k++)
2232 y = Pr*gam_S16(x-p);
2237 for (
int k=1; k<=p-1; k++)
2239 y = (p==0)? gam_S16(x) : gam_S16(x+p)/Pr;
2245real gam_S17(
const real& x)
2259 y = q_gams17_a[5] / ( v + q_gams17_b[5]);
2260 y = q_gams17_a[4] / ( (v + q_gams17_b[4]) + y );
2261 y = q_gams17_a[3] / ( (v + q_gams17_b[3]) + y );
2262 y = q_gams17_a[2] / ( (v + q_gams17_b[2]) + y );
2263 y = q_gams17_a[1] / ( (v + q_gams17_b[1]) + y ) + q_gams17_b[0];
2273real gamma_S17(
const real& x)
2286 for (
int k=2; k<=p; k++)
2288 y = Pr*gam_S17(x-p);
2293 for (
int k=1; k<=p-1; k++)
2295 y = (p==0)? gam_S17(x) : gam_S17(x+p)/Pr;
2301real gamma_05(
const real& x)
2309 Nr = int_no(gam_f85,19,x);
2313 case 0: y = gamma_S0(x);
break;
2314 case 1: y = gamma_S1(x);
break;
2315 case 2: y = gamma_S2(x);
break;
2316 case 3: y = gamma_S3(x);
break;
2317 case 4: y = gamma_S4(x);
break;
2318 case 5: y = gamma_S5(x);
break;
2319 case 6: y = gamma_S6(x);
break;
2320 case 7: y = gamma_S7(x);
break;
2321 case 8: y = gamma_S8(x);
break;
2322 case 9: y = gamma_S9(x);
break;
2323 case 10: y = gamma_S10(x);
break;
2324 case 11: y = gamma_S11(x);
break;
2325 case 12: y = gamma_S12(x);
break;
2326 case 13: y = gamma_S13(x);
break;
2327 case 14: y = gamma_S14(x);
break;
2328 case 15: y = gamma_S15(x);
break;
2329 case 16: y = gamma_S16(x);
break;
2330 case 17: y = gamma_S17(x);
break;
2346 if (x<-170.0 || x>171.0)
2347 cxscthrow(STD_FKT_OUT_OF_DEF(
"real gammar(const real& x)"));
2354 else y = 1/gamma_05(x);
2366 if (x>171.5 || x<-170.0)
2367 cxscthrow(STD_FKT_OUT_OF_DEF(
"real gamma(const real& x)"));
2379 void r_lfsr(
void) {;}
The Data Type dotprecision.
The namespace cxsc, providing all functionality of the class library C-XSC.
cinterval sqrtp1m1(const cinterval &z) noexcept
Calculates .
cinterval sqrt1mx2(const cinterval &z) noexcept
Calculates .
interval gamma(const interval &x)
The Gamma function.
int ifloor(const real &x) noexcept
Rounding to the greates integer smaller or equal x; -2147483649 < x <= 2147483647....
int Round(const real &x) noexcept
Rouding to the next integer; |x| < 2147483647.5.
interval sinpix_pi(const interval &x)
Calculates ;.
interval expx2(const interval &x)
Calculates .
real Cut25(const real &x)
Returns a real value, which corresponds with the first 25 mantissa bits of x.
cinterval ln(const cinterval &z) noexcept
Calculates .
int ceil(const real &x) noexcept
Rounding to the smallest integer greater or equal x; -2147483649 < x <= 2147483647....
interval expmx2(const interval &x)
Calculates .
real Cut24(const real &x)
Returns a real value, which corresponds with the first 24 mantissa bits of x.
interval acoshp1(const interval &x)
Calculates .
interval arg(const cinterval &z) noexcept
Calculates .
cinterval sqrtx2m1(const cinterval &z) noexcept
Calculates .
cinterval sqrt1px2(const cinterval &z) noexcept
Calculates .
cinterval exp(const cinterval &z) noexcept
Calculates .
interval gammar(const interval &x)
The inverse Gamma function: 1/Gamma(x)
interval expx2m1(const interval &x)
Calculates .
interval ln_sqrtx2y2(const interval &x, const interval &y) noexcept
Calculates .
cinterval expm1(const cinterval &z) noexcept
Calculates .
ivector abs(const cimatrix_subv &mv) noexcept
Returns the absolute value of the matrix.
cinterval sqrt(const cinterval &z) noexcept
Calculates .
real Cut26(const real &x)
Returns a real value, which corresponds with the first 26 mantissa bits of x.
void times2pown(cinterval &x, int n) noexcept
Fast multiplication of reference parameter [z] with .
cinterval sqr(const cinterval &z) noexcept
Calculates .
cinterval lnp1(const cinterval &z) noexcept
Calculates .
interval sqrtx2y2(const interval &x, const interval &y) noexcept
Calculates .