26#include "lx_cinterval.hpp"
34l_cinterval & l_cinterval::operator = (
const lx_cinterval &a)
noexcept
75std::string & operator >> (std::string &s, lx_cinterval &a)
noexcept
83 str = skipwhitespacessinglechar (str,
'(');
86 str >> SaveOpt >> Lar;
90 s = skipwhitespacessinglechar (s,
',');
91 s >> Lai >> RestoreOpt;
94 a = lx_cinterval(Lar,Lai);
99void operator >> (
const std::string &s, lx_cinterval &a)
noexcept
106void operator >> (
const char *s, lx_cinterval &a)
noexcept
112std::istream & operator >> (std::istream &s, lx_cinterval &a)
noexcept
114 lx_interval Lar, Lai;
117 std::cerr <<
"Real part: {Exponent to base 10, [a,b]} = ?"
120 std::cerr <<
"Img. part: {Exponent to base 10, [a,b]} = ?"
122 s >> Lai >> RestoreOpt;
123 a = lx_cinterval(Lar,Lai);
127 skipeolnflag =
false; inpdotflag =
true;
128 c = skipwhitespaces (s);
129 if (inpdotflag && c !=
')')
144 int max(
int a,
int b)
156lx_cinterval
sqr(
const lx_cinterval &z)
noexcept
158 lx_interval rez(Re(z)), reza(
abs(rez)),
159 imz(Im(z)), imza(
abs(imz));
166 hxl(irez), hxu(srez), hyl(iimz), hyu(simz);
167 lx_real resxl, resxu;
168 resxl = Inf( (hxl-hyu)*(hxl+hyu) );
169 resxu = Sup( (hxu-hyl)*(hxu+hyl) );
174 return lx_cinterval( lx_interval(resxl,resxu), hxl );
181 lx_interval Sqrt_zpx_m2(
const lx_interval &x,
const lx_interval &y )
191 lx_interval Sqrt_zpx_d2(
const lx_interval &x,
const lx_interval &y )
201lx_interval Re_Sqrt_Point(
const lx_interval &rez,
const lx_interval &imz)
208 lx_real irez = Inf( rez ), iimz = Inf( imz );
211 res = (ge_zero(irez))?
sqrt(rez) : lx_interval(0);
213 res = (ge_zero(irez))? Sqrt_zpx_d2(rez,imz) :
214 abs(iimz)/Sqrt_zpx_m2(rez,imz);
218lx_interval Im_Sqrt_Point(
const lx_interval &rez,
const lx_interval &imz)
225 lx_real irez = Inf( rez ), iimz = Inf( imz );
228 res = (ge_zero(irez))? lx_interval(0) :
sqrt(-rez);
230 if (ge_zero(irez)) res = iimz/Sqrt_zpx_m2(rez,imz);
233 res = Sqrt_zpx_d2(rez,imz);
234 if (se_zero(iimz)) res = -res;
239lx_cinterval
sqrt(
const lx_cinterval& z)
noexcept
248 hxl( irez ), hxu( srez ), hyl( iimz ), hyu( simz );
250 resxl, resxu, resyl, resyu;
252 if ( sm_zero(irez) && sm_zero(iimz) && ge_zero(simz) )
253 cxscthrow(STD_FKT_OUT_OF_DEF(
254 "lx_cinterval sqrt(const lx_cinterval& z); z not in principal branch."));
258 resxl = Inf( Re_Sqrt_Point( hxl,hyl ) );
259 resxu = Sup( Re_Sqrt_Point( hxu,hyu ) );
261 resyl = Inf( Im_Sqrt_Point( hxu,hyl ) );
262 resyu = Sup( Im_Sqrt_Point( hxl,hyu ) );
267 resxl = Inf( Re_Sqrt_Point( hxl,hyu ) );
268 resxu = Sup( Re_Sqrt_Point( hxu,hyl ) );
270 resyl = Inf( Im_Sqrt_Point( hxl,hyl ) );
271 resyu = Sup( Im_Sqrt_Point( hxu,hyu ) );
275 resxl = Inf(
sqrt(hxl) );
276 resxu = ( -iimz>simz ? Sup( Re_Sqrt_Point( hxu,hyl ) )
277 : Sup( Re_Sqrt_Point( hxu,hyu ) ) );
279 resyl = Inf( Im_Sqrt_Point( hxl,hyl ) );
280 resyu = Sup( Im_Sqrt_Point( hxl,hyu ) );
282 y = lx_cinterval( lx_interval(resxl,resxu), lx_interval(resyl,resyu) );
289std::list<lx_cinterval>
sqrt_all(
const lx_cinterval& z )
noexcept
296 lx_interval hxl( irez ), hxu( srez ), hyl( iimz ), hyu( simz );
297 lx_real resxl, resxu, resyl, resyu;
300 if( irez < 0.0 && iimz <= 0.0 && simz >= 0.0 )
309 resxl = Inf( Re_Sqrt_Point( hxl, hyl ) );
310 resxu = Sup( Re_Sqrt_Point( hxu, hyu ) );
313 resyl = Inf( Im_Sqrt_Point( hxu, hyl ) );
314 resyu = Sup( Im_Sqrt_Point( hxl, hyu ) );
326 resxu = Sup( Re_Sqrt_Point( hxu, hyl ) );
329 resyl = Inf( Im_Sqrt_Point( hxl, hyl ) );
330 resyu = ( srez > 0.0 ? lx_real(0.0,
l_real(0)) : -Inf(
sqrt(-hxu) ) );
342 resxl = Inf( Im_Sqrt_Point(-hxu, hyl ) );
343 resxu = Sup( Re_Sqrt_Point( hxu, hyu ) );
346 resyl = Inf(
sqrt( -hxu ) );
347 resyu = ( - iimz > simz ? Sup( Re_Sqrt_Point( -hxl, hyl ) ) :
348 Sup( Im_Sqrt_Point( hxl, hyu ) ) );
356 resxu = ( - iimz > simz ? Sup( Re_Sqrt_Point( hxu, hyl ) ) :
357 Sup( Re_Sqrt_Point( hxu, hyu ) ) );
360 resyl = Inf( Im_Sqrt_Point( hxl, hyl ) );
361 resyu = Sup( Im_Sqrt_Point( hxl, hyu ) );
365 w = lx_cinterval( lx_interval(resxl,resxu), lx_interval(resyl,resyu ) );
371 std::list<lx_cinterval> res;
418 lx_cinterval
exp(
const lx_cinterval& z)
noexcept
420 int stagsave = stagprec,
422 if (stagprec > stagmax)
425 lx_interval lreal(Re(z)), limg(Im(z));
430 y = lx_cinterval( A*
cos( B ) , A*
sin( B ) );
438 lx_cinterval
exp2(
const lx_cinterval& z)
noexcept
443 lx_cinterval
exp10(
const lx_cinterval& z)
noexcept
450 lx_cinterval
cos(
const lx_cinterval& z)
noexcept
452 int stagsave = stagprec,
454 if (stagprec > stagmax)
457 lx_interval lreal(Re(z)),limg(Im(z));
472 lx_cinterval
sin(
const lx_cinterval& z)
noexcept
474 int stagsave = stagprec,
476 if (stagprec > stagmax)
479 lx_interval lreal(Re(z)),limg(Im(z));
494 lx_cinterval
cosh(
const lx_cinterval& z)
noexcept
496 int stagsave = stagprec,
498 if (stagprec > stagmax)
501 lx_interval lreal(Re(z)),limg(Im(z));
516 lx_cinterval
sinh(
const lx_cinterval& z)
noexcept
518 int stagsave = stagprec,
520 if (stagprec > stagmax)
523 lx_interval lreal(Re(z)),limg(Im(z));
544lx_interval Atan(
const lx_interval& y,
const lx_interval& x)
550 lx_interval res(0),u;
551 l_interval yl(li_part(y));
553 signy(sign(Inf(y))), signx(sign(Inf(x)));
554 real ex_y(expo(Inf(y))), ex_x(expo(Inf(x)));
557 if (ex_yl > -1000000)
559 neg = signy * signx == -1;
560 if (ex_y > ex_x+4197)
568 if (ex_x >= 9007199254738946.0)
573 res = lx_interval( lx_real(0,l_real(0)),
574 lx_real(-Max_Int_R,l_real(
minreal)) );
587 if (ex_x <= -9007199254738891.0)
603lx_interval Atan(
const lx_interval& y,
const lx_real& x)
605 return Atan(y,lx_interval(x));
623lx_interval
Arg(
const lx_cinterval& z)
noexcept
632 hxl(irez), hxu(srez), hyl(iimz), hyu(simz), Pid2;
641 resl = ( srez > 0.0 ? Inf( Atan( hyl,hxu ) ) :
644 resu = ( irez > 0.0 ? Sup( Atan( hyu,hxl ) ) :
647 return lx_interval( resl, resu );
655 ( irez > 0.0 ? Inf( Atan( hyl,hxl ) ) : -Sup( Pid2 ) ) );
657 ( srez > 0.0 ? Sup( Atan( hyu,hxu ) ) : -Inf(Pid2) ) );
658 return lx_interval( resl, resu );
667 resl = iimz < 0.0 ? Inf( Atan( hyl,hxl ) ) : lx_real(0.0);
668 return lx_interval( resl, Sup( Atan( hyu,hxl ) ) );
676 cxscthrow (STD_FKT_OUT_OF_DEF(
"lx_interval Arg(const lx_cinterval& z); z contains negative real numbers"));
677 return lx_interval(0.0);
685 resl = iimz < 0.0 ? -Sup(Pid2) : lx_real(0.0);
686 resu = simz > 0.0 ? Sup(Pid2) : lx_real(0.0);
687 return lx_interval( resl, resu );
692 if( eq_zero(iimz) && eq_zero(simz) )
694 return lx_interval(0.0);
697 resl = ( iimz < 0.0 ? - Sup(Pid2) : Inf(Pid2) );
698 resu = ( simz > 0.0 ? Sup(Pid2) : -Inf(Pid2) );
699 return lx_interval( resl, resu );
718lx_interval
arg(
const lx_cinterval& z )
noexcept
730 if( sm_zero(irez) && se_zero(iimz) && ge_zero(simz) )
736 resl = ( sm_zero(iimz)? - Sup(
Pi_lx_interval() ) : lx_real(0.0) );
737 resu = ( ( sm_zero(iimz) && eq_zero(simz) ) ? lx_real(0.0) :
739 return lx_interval( resl, resu );
751 resl = ( gr_zero(simz)? Inf(Pid2) :
753 resu = ( sm_zero(iimz)?
754 ( gr_zero(simz)? Sup(3*Pid2) : -Inf(Pid2) ) :
756 return lx_interval( resl, resu );
761 lx_interval hyl(iimz), hyu(simz);
762 resl = ( simz > 0.0 ? Inf( Atan( hyu,srez ) +
Pi_lx_interval() ) :
764 resu = ( iimz < 0.0 ? ( simz > 0.0 ? Sup( Atan( hyl,srez ) +
767 return lx_interval( resl, resu );
781 lx_interval Re_Sqrt_point(
const lx_interval& rez,
const lx_interval& imz,
793 if( eq_zero(Sup(a)) )
795 return lx_interval(0);
798 cos(
Arg( lx_cinterval( rez, imz ) ) / n );
801 lx_interval Im_Sqrt_point(
const lx_interval& rez,
const lx_interval& imz,
812 if( eq_zero(Sup(a)) )
814 return lx_interval(0);
817 sin(
Arg( lx_cinterval( rez, imz ) ) / n );
820lx_cinterval
sqrt(
const lx_cinterval& z,
int n )
noexcept
826 if( n == 0 )
return lx_cinterval(lx_interval(1));
827 if( n == 1 )
return z;
828 if( n == 2 )
return sqrt(z);
836 lx_interval hxl( irez ), hxu( srez ), hyl( iimz ), hyu( simz );
837 lx_real resxl, resxu, resyl, resyu;
839 if( sm_zero(irez) && se_zero(iimz) && ge_zero(simz) )
842 cxscthrow(STD_FKT_OUT_OF_DEF(
"lx_cinterval sqrt(const lx_cinterval& z, int n ); z contains negative real values."));
850 lx_cinterval hres =
sqrt( lx_cinterval( Re(z), -Im(z) ), n );
851 return lx_cinterval( Re(hres), -Im(hres) );
860 lx_real tanglel = Inf( tangle ), tangleu = Sup( tangle );
864 if ( ge_zero(irez) || Sup( hyl / irez ) <= tanglel )
867 resxl = Inf( Re_Sqrt_point( hxl, hyl, n ) );
870 if( sm_zero(srez) && Inf( hyl / srez ) >= tangleu )
873 resxl = Inf( Re_Sqrt_point( hxu, hyl, n ) );
877 resxl = Inf( Re_Sqrt_point( iimz / tangle , hyl, n ) );
882 if ( ge_zero(irez) || Sup( hyu / irez ) <= tanglel )
885 resxu = Sup( Re_Sqrt_point( lx_interval(srez),
886 lx_interval(simz), n ) );
889 if ( sm_zero(srez) && Inf( hyu / srez ) >= tangleu )
892 resxu = Sup( Re_Sqrt_point( hxl, hyu, n ) );
896 resxu = max( Sup( Re_Sqrt_point( hxl, hyu, n ) ),
897 Sup( Re_Sqrt_point( hxu, hyu, n ) ) );
902 if ( ge_zero(srez) || Sup( hyl / srez ) <= tanglel )
905 resyl = Inf( Im_Sqrt_point( hxu, hyl, n ) );
908 if( Inf( hyu / srez ) >= tangleu )
911 resyl = Inf( Im_Sqrt_point( hxu, hyu, n ) );
915 resyl = Inf(Im_Sqrt_point( hxu, srez * tangle, n ));
920 if( ge_zero(irez) || Sup( hyl / irez ) <= tanglel )
923 resyu = Sup( Im_Sqrt_point( hxl, hyu, n ) );
926 if( Inf( hyu / irez ) >= tangleu )
929 resyu = Sup( Im_Sqrt_point( hxl, hyl, n ) );
933 resyu = max( Sup( Im_Sqrt_point( hxl, hyl, n ) ),
934 Sup( Im_Sqrt_point( hxl, hyu, n ) ) );
942 resxl = ( eq_zero(irez)? lx_real(0.0) : Inf(
sqrt( hxl, n ) ) );
943 resxu = ( - iimz > simz ? Sup( Re_Sqrt_point( hxu, hyl, n ) ) :
944 Sup( Re_Sqrt_point( hxu, hyu, n ) ) );
947 resyl = Inf( Im_Sqrt_point( hxl, hyl, n ) );
948 resyu = Sup( Im_Sqrt_point( hxl, hyu, n ) );
950 return lx_cinterval( lx_interval( resxl, resxu ),
951 lx_interval( resyl, resyu ) );
961 std::list<lx_cinterval>
sqrt_all(
const lx_cinterval& z,
int n )
noexcept
980 std::list<lx_cinterval> res;
984 res.push_back( lx_cinterval( lx_interval(0,
l_interval(1)),
993 else if( n == 2 )
return sqrt_all( z );
997 arg_z =
arg( z ), root_abs_z =
sqrt(
abs(z), n );
999 for(
int k = 0; k < n; k++)
1001 lx_interval arg_k = ( arg_z + 2 * k *
Pi_l_interval() ) / n;
1003 res.push_back( lx_cinterval( root_abs_z *
cos( arg_k ),
1004 root_abs_z *
sin( arg_k ) ) );
1024lx_cinterval
Ln(
const lx_cinterval& z)
noexcept
1026 int stagsave = stagprec,
1028 if (stagprec > stagmax) stagprec = stagmax;
1032 srez = Sup( Re(z) ),
1033 simz = Sup( Im(z) ),
1034 iimz = Inf( Im(z) );
1035 lx_interval a1(
abs(Re(z)) ), a2(
abs(Im(z)) );
1036 if ( eq_zero(Inf(a1)) && eq_zero(Inf(a2)) )
1038 cxscthrow(STD_FKT_OUT_OF_DEF(
1039 "lx_cinterval Ln(const lx_cinterval& z); z contains 0"));
1040 if ( sm_zero(srez) && sm_zero(iimz) && ge_zero(simz) )
1041 cxscthrow(STD_FKT_OUT_OF_DEF(
1042 "lx_cinterval Ln(const lx_cinterval& z); z not allowed"));
1046 stagprec = stagsave;
1058lx_cinterval
ln(
const lx_cinterval& z)
noexcept
1060 int stagsave = stagprec,
1062 if (stagprec > stagmax) stagprec = stagmax;
1065 lx_interval a1(
abs(Re(z)) ), a2(
abs(Im(z)) );
1066 if ( eq_zero(Inf(a1)) && eq_zero(Inf(a2)) )
1068 cxscthrow(STD_FKT_OUT_OF_DEF
1069 (
"lx_cinterval ln(const lx_cinterval& z); z contains 0"));
1072 stagprec = stagsave;
1085lx_cinterval
log2(
const lx_cinterval& z)
noexcept
1090lx_cinterval
log10(
const lx_cinterval& z)
noexcept
1109 if( n == 0 )
return lx_cinterval(lx_interval(1));
1110 else if( n == 1 )
return z;
1111 else if( n == -1 )
return 1 / z;
1112 else if( n == 2 )
return sqr(z);
1130 lx_interval abs_z =
abs(z);
1132 if( ((n < 0) && eq_zero(Inf(abs_z))) || !(
Is_Integer(n)))
1134 cxscthrow (STD_FKT_OUT_OF_DEF(
"lx_cinterval power_fast(const lx_cinterval& z, const real& n); z contains 0 or n is not integer."));
1136 if( eq_zero(Sup(abs_z)) )
return lx_cinterval( lx_interval(0));
1139 lx_interval arg_z =
arg(z);
1140 lx_interval abs_z_n =
power(abs_z,n);
1142 return lx_cinterval( abs_z_n *
cos( n * arg_z ),
1143 abs_z_n *
sin( n * arg_z ) );
1156 lx_cinterval
power(
const lx_cinterval& x,
const real& n )
noexcept
1160 cxscthrow(STD_FKT_OUT_OF_DEF(
"lx_cinterval power(const lx_cinterval& z, const real& n); n is not integer."));
1162 real one(1.0), zhi(2.0), N(n), r;
1164 lx_cinterval y,neu,X(x);
1166 if (x == one) y = x;
1168 if (N == 0.0) y = one;
1173 if (N == 2) y =
sqr(x);
1187 dbl = _double(N/zhi);
1208 lx_cinterval
pow(
const lx_cinterval& z,
const lx_interval& p )
noexcept
1210 return exp( p*
Ln(z) );
1221 lx_cinterval
pow(
const lx_cinterval& z,
const lx_cinterval& p)
noexcept
1223 return exp( p*
Ln(z) );
1246void horizontal_check(
const lx_interval& hy,
const lx_interval& cosh_2y,
1247 lx_real irez, lx_real srez,
const lx_interval& hxl,
1248 const lx_interval& hxu, lx_real& resxl, lx_real& resxu )
1262 bool both =
false, left =
false, right =
false;
1263 lx_interval hlp, hlp1;
1266 hlp1 = lx_interval(srez) - lx_interval(irez);
1267 if (Inf(hlp1) > Inf(hlp))
1273 res_l =
cos( 2 * hxl ) - cosh_2y,
1274 res_u =
cos( 2 * hxu ) - cosh_2y;
1276 if( Inf( res_l * res_u ) > 0.0 )
1281 if( Sup( res_l * res_u ) < 0.0 )
1285 if( Inf( res_l ) > 0.0 )
1300 sin_2xl =
sin( 2 * hxl ),
1301 sin_2xu =
sin( 2 * hxu );
1303 if( !
Disjoint( lx_interval(0), res_l ) )
1306 if( ge_zero(Inf(sin_2xl)) )
1311 res_l = -lx_interval(1);
1315 if( se_zero(Sup(sin_2xl)) )
1320 res_l = lx_interval(1);
1330 if( !
Disjoint( lx_interval(0), res_u ) )
1333 if( ge_zero(Inf(sin_2xu)) )
1338 res_u = lx_interval(1);
1342 if( se_zero(Sup(sin_2xu)) )
1347 res_u = -lx_interval(1);
1359 if( Inf( res_l * res_u ) < 0.0 )
1367 lx_interval re_tan = 1 /
sinh( 2 *
abs( hy ) );
1372 resxl = min( resxl, Inf( re_tan ) );
1373 resxu = max( resxu, Sup( re_tan ) );
1379 resxl = min( resxl, -Sup( re_tan ) );
1380 resxu = max( resxu, -Inf( re_tan ) );
1385lx_cinterval Tan(
const lx_cinterval& z)
noexcept
1398 lx_interval hxl(irez), hxu(srez), hyl(iimz), hyu(simz);
1400 lx_real resxl, resxu, resyl, resyu;
1404 if( ( !
Disjoint( lx_interval(0), imz ) ) &&
1406 cxscthrow (STD_FKT_OUT_OF_DEF(
1407 "lx_cinterval tan(const lx_cinterval& z); Pole(s) in z"));
1414 cos_2rez =
cos( 2 * rez ), sinh_imz_2 =
sqr(
sinh( imz ) );
1417 re_tan_l=
sin( 2*hxl ) / ( 2*(
sqr(
cos( hxl ) ) + sinh_imz_2 ) ),
1418 re_tan_u=
sin( 2*hxu ) / ( 2*(
sqr(
cos( hxu ) ) + sinh_imz_2 ) );
1420 resxl = min( Inf( re_tan_l ), Inf( re_tan_u ) );
1421 resxu = max( Sup( re_tan_l ), Sup( re_tan_u ) );
1433 cosh_2yl = - 1 /
cosh( 2 * hyl ),
1434 cosh_2yu = - 1 /
cosh( 2 * hyu );
1436 if( !
Disjoint( cos_2rez, cosh_2yl ) && iimz != 0.0 )
1438 horizontal_check(hyl,cosh_2yl,irez,srez,hxl,hxu,resxl,resxu);
1440 if( !
Disjoint( cos_2rez, cosh_2yu ) && simz != 0.0 )
1442 horizontal_check(hyu,cosh_2yu,irez,srez,hxl,hxu,resxl,resxu);
1449 lx_interval cos_rez_2 =
sqr(
cos( rez ) );
1452 im_tan_l =
sinh( 2*hyl ) / (2*(cos_rez_2 +
sqr(
sinh( hyl ) ))),
1453 im_tan_u =
sinh( 2*hyu ) / (2*(cos_rez_2 +
sqr(
sinh( hyu ) )));
1455 resyl = min( Inf( im_tan_l ), Inf( im_tan_u ) );
1456 resyu = max( Sup( im_tan_l ), Sup( im_tan_u ) );
1465 cos_2xl =
cos( 2 * hxl ), cos_2xu =
cos( 2 * hxu );
1472 imz_h = lx_interval( iimz, min( simz, lx_real(0.0) ) ),
1473 cosh_2imz = - 1 /
cosh( 2 * imz_h );
1475 if( ( !
Disjoint( cosh_2imz, cos_2xl ) ) )
1479 im_tan = - 1 /
abs(
sin( 2 * hxl ) );
1480 resyl = min( resyl, Inf( im_tan ) );
1481 resyu = max( resyu, Sup( im_tan ) );
1483 if( ( !
Disjoint( cosh_2imz, cos_2xu ) ) )
1487 im_tan = - 1 /
abs(
sin( 2 * hxu ) );
1488 resyl = min( resyl, Inf( im_tan ) );
1489 resyu = max( resyu, Sup( im_tan ) );
1496 imz_h = lx_interval( max( iimz, lx_real(0.0) ), simz ),
1497 cosh_2imz = - 1 /
cosh( 2 * imz_h );
1499 if( ( !
Disjoint( cosh_2imz, cos_2xl ) ) )
1503 im_tan = + 1 /
abs(
sin( 2 * hxl ) );
1504 resyl = min( resyl, Inf( im_tan ) );
1505 resyu = max( resyu, Sup( im_tan ) );
1507 if( ( !
Disjoint( cosh_2imz, cos_2xu ) ) )
1511 im_tan = + 1 /
abs(
sin( 2 * hxu ) );
1512 resyl = min( resyl, Inf( im_tan ) );
1513 resyu = max( resyu, Sup( im_tan ) );
1517 y = lx_cinterval( lx_interval(resxl,resxu ),lx_interval(resyl,resyu ) );
1522lx_cinterval
tan(
const lx_cinterval& z)
noexcept
1537 const real S = 1e-15;
1538 int stagsave = stagprec,
1540 if (stagprec > stagmax)
1551 real x(
mid(re)),k, pi(3.14159265358979323);
1558 if (k == 9007199254740992.0)
1559 cxscthrow (STD_FKT_OUT_OF_DEF(
1560 "lx_cinterval tan(const lx_cinterval& z); z out of range"));
1567 u = Re(eps); u =
abs(u);
1568 v = Im(eps); v =
abs(v);
1571 y = (Sup(u_)<S && Sup(v_)<S)? -lx_cinterval(1) / Tan(eps) : Tan(z);
1573 stagprec = stagsave;
1579 lx_cinterval
cot(
const lx_cinterval& z)
noexcept
1594 const real S = 1e-15;
1595 int stagsave = stagprec,
1597 if (stagprec > stagmax)
1603 lx_interval rez = Re(z);
1607 real x(
mid(re)), k, pi(3.14159265358979323);
1614 if (k == 9007199254740992.0)
1615 cxscthrow (STD_FKT_OUT_OF_DEF(
1616 "lx_cinterval cot(const lx_cinterval& z); z out of range"));
1618 u = Re(eps); u =
abs(u);
1619 v = Im(eps); v =
abs(v);
1622 if (Sup(u)<S && Sup(v)<S)
1631 stagprec = stagsave;
1641lx_cinterval
tanh(
const lx_cinterval& z)
noexcept
1643 int stagsave = stagprec,
1645 if (stagprec > stagmax) stagprec = stagmax;
1647 lx_cinterval res =
tan( lx_cinterval( Im(z), Re(z) ) ),y;
1648 y = lx_cinterval( Im(res), Re(res) );
1650 stagprec = stagsave;
1663lx_cinterval
coth(
const lx_cinterval& z)
noexcept
1665 int stagsave = stagprec,
1667 if (stagprec > stagmax) stagprec = stagmax;
1669 lx_cinterval zh = lx_cinterval( -Im(z), Re(z) );
1670 lx_cinterval res =
cot(zh);
1671 zh = lx_cinterval( -Im(res), Re(res) );
1673 stagprec = stagsave;
1700lx_interval f_aux_asin(
const lx_interval& x,
const lx_interval& y)
1712 if (y != 0.0 || Inf(res) < 1.0)
1722 lx_real hlb = max( lx_real(1.0),
abs( Sup(x) ) );
1723 if( Inf(res) < hlb )
1724 res = lx_interval( hlb, Sup(res) );
1728 lx_interval ACOSH_f_aux(
const lx_interval& x,
const lx_interval& y)
1742 lx_interval res, xa(
abs(x)), ya(
abs(y)), S1, S2, S3;
1743 lx_real rx((Inf(xa))), ry(Inf(ya));
1745 if (rx>2.0 || ry>2.0) res =
acosh( f_aux_asin(x,y) );
1750 res = ya/2; S2 = res/2;
1751 S1 = ya*(0.5 + S2/(
sqrt1px2(res)+1));
1771 if (eq_zero(Inf(S1)))
1772 S1 = lx_interval(Inf(lx_interval(-2097,l_interval(1))),
1778 if (eq_zero(Inf(S1)))
1779 S1 = lx_interval(Inf(lx_interval(-2097,l_interval(1))),
1783 if (Inf(S1)>Inf(S2))
1792 S1 = (S3 + S2) + S1;
1802 lx_interval Asin_beta(
const lx_interval& x,
const lx_interval& y )
1810 const real c1 = 0.75;
1811 bool neg_x, fertig(
false);
1812 lx_real Infxa(Inf(x)), L_y(Inf(y));
1813 int ex_xl(
expo_gr(lr_part(Infxa))),
1815 real ex_x(expo(Infxa)), ex_y(expo(L_y));
1816 lx_interval res,beta,abs_beta,delta,w_delta,xa,Ne,Kla;
1822 if ( ex_yl>-1000000 && ex_y>=9007199254739972.0-ex_yl &&
1823 ex_x <= 2097-ex_xl )
1829 if (neg_x) xa = -xa;
1830 res = lx_interval( lx_real(0), lx_real(-9007199254737870.0,l_real(1)) );
1838 beta = x / ( Ne/2 );
1839 if (Inf(beta)<-1) Inf(beta) = -1;
1840 if (Sup(beta)> 1) Sup(beta) = 1;
1841 abs_beta =
abs(beta);
1842 if (Inf(abs_beta) < c1)
1850 if (neg_x) xa = -xa;
1860 delta = lx_interval(lx_real(-2100,7.4699)) * delta;
1862 if (Sup(
abs(y)) < Inf(delta))
1866 delta =
sqr(y/beta);
1895 res -=
asin(
sqrt(delta*(2-delta)) );
1896 if (neg_x) res = -res;
1904lx_cinterval
asin(
const lx_cinterval& z)
noexcept
1906 int stagsave = stagprec,
1908 if (stagprec>stagmax) stagprec = stagmax;
1911 lx_interval rez = Re(z), imz = Im(z);
1919 lx_interval hxl(irez), hxu(srez), hyl(iimz), hyu(simz);
1921 lx_real resxl, resxu, resyl, resyu;
1923 bool bl = (sm_zero(iimz)) && (gr_zero(simz)),
1924 raxis = (eq_zero(iimz)) && (eq_zero(simz));
1929 if ( (irez<-1 && (bl || (sm_zero(iimz) && eq_zero(simz)))) ||
1930 (srez >1 && (bl || (eq_zero(iimz) && gr_zero(simz)))) )
1931 cxscthrow(STD_FKT_OUT_OF_DEF(
1932 "lx_cinterval asin(const lx_cinterval& z); z contains singularities."));
1937 if (sm_zero(iimz) && gr_zero(simz))
1941 resxl = Inf(
asin(hxl) );
1943 resxl = Inf(Asin_beta(hxl,lx_interval( max(-iimz,simz) )) );
1945 resxu = Sup(Asin_beta(hxu,lx_interval( max(-iimz,simz) )) );
1947 resxu = Sup(
asin(hxu) );
1951 if ( (ge_zero(iimz) && ge_zero(irez)) ||
1952 (se_zero(simz) && se_zero(irez)) )
1955 resxl = Inf( Asin_beta(hxl,hyu) );
1959 resxl = Inf( Asin_beta(hxl,hyl) );
1960 if ( (ge_zero(iimz) && ge_zero(srez)) || (se_zero(simz) && se_zero(srez) ) )
1963 resxu = Sup( Asin_beta(hxu,hyl) );
1967 resxu = Sup( Asin_beta(hxu,hyu) );
1976 if (sm_zero(srez)) resyl = Inf( ACOSH_f_aux( hxu, hyu ));
1977 else resyl = -Sup( ACOSH_f_aux( hxu, hyu ));
1978 if (gr_zero(irez)) resyu = -Inf( ACOSH_f_aux( hxl, hyu ));
1979 else resyu = Sup( ACOSH_f_aux( hxl, hyu ));
1990 resyl = -Sup( ACOSH_f_aux( hxl, hyl ) );
1992 resyu = -Inf( ACOSH_f_aux( hxu, hyu ) );
1994 resyu = -Inf( ACOSH_f_aux( lx_interval(0),hyu ) );
1999 resyl = -Sup( ACOSH_f_aux( hxu, hyl ) );
2001 resyu = -Inf( ACOSH_f_aux( hxl, hyu ) );
2003 resyu = -Inf( ACOSH_f_aux( lx_interval(0),hyu ) );
2015 resyu = Sup( ACOSH_f_aux( hxl, hyu ) );
2017 resyl = Inf( ACOSH_f_aux( hxu, hyl ) );
2019 resyl = Inf( ACOSH_f_aux( lx_interval(0), hyl ) );
2024 resyu = Sup( ACOSH_f_aux( hxu, hyu ) );
2026 resyl = Inf( ACOSH_f_aux( hxl, hyl ) );
2028 resyl = Inf( ACOSH_f_aux( lx_interval(0), hyl ) );
2039 resyl = - Sup( ACOSH_f_aux( hxl, hyl ) );
2040 resyu = Sup( ACOSH_f_aux( hxl, hyu ) );
2044 resyl = - Sup( ACOSH_f_aux( hxu, hyl ) );
2045 resyu = Sup( ACOSH_f_aux( hxu, hyu ) );
2049 res = lx_cinterval( lx_interval(resxl,resxu),lx_interval(resyl,resyu) );
2050 stagprec = stagsave;
2060lx_interval BETA_xy(
const lx_interval& x,
const lx_interval& y)
2072 xl = li_part(xa); yl = li_part(y);
2075 exx = expo(xa); exxl =
expo_gr(Inf(xl));
2076 exy = expo(y); exyl =
expo_gr(Inf(yl));
2079 if (Inf(xa)<1) res = xa;
2086 z = interval(exy) + (exyl-exxl-2103);
2087 if (exx<Inf(z)) res = 0;
2089 res = xa / f_aux_asin(xa,y);
2091 if (Sup(res)>1) res = 1.0;
2093 if (neg) res = -res;
2097lx_interval Asin_arg(
const lx_interval& x,
const lx_interval& y )
2105 const real c2 = -9007199254739968.0;
2106 lx_interval res,xa,F_xa,S1,S2,xm1,xp1,w_delta,T;
2111 if (neg_x) xa = -xa;
2123 if (Inf(
abs(y)) > Inf(F_xa))
2131 xm1 = lx_interval(Inf( lx_interval(-2097,l_interval(1)) ),
2136 w_delta =
sqrt(2 + S1 + S2);
2137 w_delta = w_delta *
sqrt(S1+xm1);
2140 T = w_delta *
sqrt(T);
2148 real exx,exxl,exy,exyl;
2152 exx = expo(xa); exy = expo(y);
2156 z = interval(exy)-interval(exx) + (exyl-exxl+1079);
2163 if (
diam(z)==0) exyl = Sup(z);
2166 dbl = floor(_double(exx));
2167 exyl = real(dbl) + 1;
2171 res = lx_interval(lx_real(0.0),
2172 lx_real(exyl,l_real(
minreal)));
2176 T =
sqrt(xa*(xa-1.0));
2192 res = lx_interval(lx_real(0.0),Sup(res));
2201 res = 2 + T +
sqrtx2y2(lx_interval(0,l_interval(2)),T);
2222 T =
sqrt( T*(2-T) );
2230lx_interval Acos_beta(
const lx_interval& x,
const lx_interval& y)
2234 const real c1 = 0.75;
2235 lx_interval res(0),beta;
2236 beta = BETA_xy(x,y);
2242 else res =
acos(beta);
2244 res = Asin_arg(x,y);
2249lx_cinterval
acos(
const lx_cinterval& z)
noexcept
2262 hxl(irez), hxu(srez), hyl(iimz), hyu(simz);
2264 bool bl = iimz< 0.0 && simz>0.0,
2265 raxis = iimz==0.0 && simz==0;
2266 lx_real resxl, resxu, resyl, resyu;
2270 if( (irez<-1 && (bl || (iimz<0 && simz==0))) ||
2271 (srez >1 && (bl || (iimz==0 && simz>0))) )
2272 cxscthrow(STD_FKT_OUT_OF_DEF(
"lx_cinterval acos(const lx_cinterval& z); z contains singularities."));
2278 if( iimz < 0.0 && simz > 0.0 )
2281 if( irez <= 0.0 ) resxu = Sup(
acos( hxl ) );
2282 else resxu = Sup( Acos_beta(hxl,lx_interval( max(-iimz,simz) )) );
2285 resxl = Inf( Acos_beta(hxu,lx_interval(max(-iimz,simz))) );
2286 else resxl = Inf(
acos( hxu ) );
2290 if (irez<0 && srez>0)
2294 resxl = Inf( Acos_beta(hxu,hyl) );
2295 resxu = Sup( Acos_beta(hxl,hyl) );
2299 resxl = Inf( Acos_beta(hxu,hyu) );
2300 resxu = Sup( Acos_beta(hxl,hyu) );
2304 if ( (ge_zero(iimz) && ge_zero(irez)) ||
2305 (se_zero(simz) && sm_zero(irez)) )
2308 resxl = Inf( Acos_beta(hxu,hyl) );
2312 resxl = Inf( Acos_beta(hxu,hyu) );
2314 if ( (ge_zero(iimz) && gr_zero(srez)) ||
2315 (se_zero(simz) && se_zero(srez)) )
2318 resxu = Sup( Acos_beta(hxl,hyu) );
2322 resxu = Sup( Acos_beta(hxl,hyl) );
2332 if (srez<0.0) resyl = Inf( ACOSH_f_aux( hxu, hyu ));
2333 else resyl = -Sup( ACOSH_f_aux( hxu, hyu ));
2334 if (irez>0.0) resyu = -Inf( ACOSH_f_aux( hxl, hyu ));
2335 else resyu = Sup( ACOSH_f_aux( hxl, hyu ));
2346 resyl = -Sup( ACOSH_f_aux( hxl, hyl ) );
2348 resyu = -Inf( ACOSH_f_aux( hxu, hyu ) );
2350 resyu = -Inf( ACOSH_f_aux( lx_interval(0), hyu ) );
2355 resyl = -Sup( ACOSH_f_aux( hxu, hyl ) );
2357 resyu = -Inf( ACOSH_f_aux( hxl, hyu ) );
2359 resyu = -Inf( ACOSH_f_aux( lx_interval(0), hyu ) );
2371 resyu = Sup( ACOSH_f_aux( hxl, hyu ) );
2373 resyl = Inf( ACOSH_f_aux( hxu, hyl ) );
2375 resyl = Inf( ACOSH_f_aux( lx_interval(0), hyl ) );
2380 resyu = Sup( ACOSH_f_aux( hxu, hyu ) );
2382 resyl = Inf( ACOSH_f_aux( hxl, hyl ) );
2384 resyl = Inf( ACOSH_f_aux( lx_interval(0), hyl ) );
2395 resyl = -Sup( ACOSH_f_aux( hxl, hyl ) );
2396 resyu = Sup( ACOSH_f_aux( hxl, hyu ) );
2400 resyl = -Sup( ACOSH_f_aux( hxu, hyl ) );
2401 resyu = Sup( ACOSH_f_aux( hxu, hyu ) );
2405 return lx_cinterval( lx_interval( resxl, resxu ),-lx_interval( resyl, resyu ) );
2413lx_cinterval
asinh(
const lx_cinterval& z )
noexcept
2418 lx_cinterval res =
asin( lx_cinterval( Im(z), -Re(z) ) );
2419 return lx_cinterval( -Im(res), Re(res) );
2426lx_cinterval
acosh(
const lx_cinterval& z )
noexcept
2441 lx_interval hxl(irez), hxu(srez), hyl(iimz), hyu(simz);
2443 lx_real resxl, resxu, resyl, resyu;
2448 if( ( se_zero(iimz) && ge_zero(simz) ) && ( irez < 1.0 ) )
2449 cxscthrow(STD_FKT_OUT_OF_DEF(
2450 "lx_cinterval acosh( const lx_cinterval& z ); z contains singularities."));
2463 lx_cinterval res =
acos(z);
2464 return lx_cinterval( -Im(res),Re(res) );
2473 lx_cinterval res =
acos(z);
2474 return lx_cinterval( Im(res), -Re(res) );
2482 resxl = Inf(
acosh( hxl ) );
2483 lx_interval ytilde( max( -iimz, simz ) );
2485 resxu = Sup( ACOSH_f_aux(hxu,ytilde) );
2491 resyl = -Sup( Acos_beta(hxl,hyl) );
2493 resyu = Sup( Acos_beta(hxl,hyu) );
2495 return lx_cinterval(lx_interval( resxl, resxu ),lx_interval( resyl, resyu ));
2504bool sign_test(
const lx_interval& x,
int s_org)
2508 if (
diam(li_part(x))>0)
2510 bl1 = eq_zero(Sup(x)) && (s_org==1);
2511 bl2 = eq_zero(Inf(x)) && (s_org==-1);
2515 alter = sign(Sup(li_part(x))) != s_org;
2519void re_atan(
const lx_interval& y, lx_interval& x, lx_interval& res )
2524 const real c1 = 4503599627369982.0;
2525 lx_interval ya(
abs(y)), One(lx_interval(0,l_interval(1))), xa;
2527 real ex_x,ex_xl,ex_y,ex_yl,s;
2532 ex_xl =
expo_gr(Inf(li_part(x)));
2534 ex_yl =
expo_gr(Sup(li_part(ya)));
2536 signx = sign(Inf(li_part(x)));
2538 if (ex_xl < -1000000)
2541 if (ex_yl < -1000000)
2551 res = lx_interval(Inf(res),lx_real(1.0));
2557 s = c1 - max(ex_x,ex_y);
2558 times2pown_neg(One,s);
2559 times2pown_neg(x,s);
2561 res = (One-x)*(One+x) -
sqr(ya);
2562 times2pown_neg(x,s);
2566 res =
sqr(y); xa =
abs(x);
2575 ex_xl = ex_xl - 1051;
2576 if (ex_y+ex_yl < 2*(ex_x+ex_xl))
2583 x = (Sup(x)>0)? -1.0 : 1.0;
2590 R = Sup(
abs(res) );
2591 times2pown_neg(x,s);
2598 times2pown_neg(res,s);
2599 times2pown_neg(res,s);
2603 times2pown_neg(x,s);
2617 res = (x-1)*(x+1) +
sqr(ya);
2637 res = (1-x)*(1+x) -
sqr(y);
2638 alter = sign_test(x,signx);
2646void re_vert(
const lx_real& x,
const lx_interval& hx,
2647 const lx_real& rew_inf,
const lx_real& rew_sup,
2648 lx_real& resxl, lx_real& resxu )
2662 lx_interval hx2(hx),Pid2,Pid4;
2669 resxl = gr_zero(rew_sup)? Inf( Atan( hx2,rew_sup )/2.0 )
2670 : ( sm_zero(rew_sup)? Inf( (Atan( hx2,rew_sup ) +
Pi_lx_interval() )/2.0 )
2673 resxu = gr_zero(rew_inf)? Sup( Atan( hx2,rew_inf )/2.0 )
2674 : ( sm_zero(rew_inf)? Sup( (Atan( hx2,rew_inf ) +
Pi_lx_interval())/2.0 )
2680 resxl = sm_zero(rew_inf)? Inf( (Atan( hx2,rew_inf ) -
Pi_lx_interval())/2.0 )
2681 : ( gr_zero(rew_inf)? Inf( Atan( hx2,rew_inf )/2.0 )
2683 resxu = sm_zero(rew_sup)? Sup( (Atan( hx2,rew_sup ) -
Pi_lx_interval())/2.0 )
2684 : ( gr_zero(rew_sup)? Sup( Atan( hx2,rew_sup )/2.0 )
2690lx_interval T_atan(
const lx_real& x)
2696 const real c1 = 4503599627367859.0;
2699 real ex_x(expo(ix));
2718lx_interval Q_atan(
const lx_interval& x,
const lx_interval& y)
2725 const real c1 = 4503599627367859.0,
2726 c2 = 9007199254740990.0,
2727 c3 = 4503599627370495.0;
2728 const lx_real S = lx_real(-c2,
MinReal);
2733 lx_interval res(0.0);
2734 real ex_y(expo(y)), ex_yl(
expo_gr(li_part(y))),
2735 ex_x(expo(Sup(x))),ex_xl,exx,exxl,dbl,s,up,exy,exyl;
2738 ex_xl =
expo_gr(lr_part(Sup(x)));
2739 if (ex_xl<-100000) ex_x = 0;
2748 res =
ln(4+
sqr(x)) - res;
2756 if (514-exxl < -c3+exx)
2757 res = lx_interval(lx_real(0),S);
2759 if (3 - exxl >= -c3 + exx)
2761 dbl = -2*(exx+exxl-3);
2762 res = lx_interval(lx_real(0),lx_real(dbl,l_real(1)));
2766 s = c2-exx; s = (s - exx) -2*exxl + 6;
2767 r = (int) _double(s);
2768 res = lx_interval(lx_real(0),lx_real(-c2,comp(0.5,r+1)));
2776 if (ex_x>c1 || ex_y>c1 )
2779 exx = expo(R); exxl =
expo_gr(lr_part(R));
2781 exy = expo(R); exyl =
expo_gr(lr_part(R));
2784 dbli = 2*( interval(exy) + (exyl-2) );
2787 z = 2*( interval(exx) + (exxl-2) );
2788 dbli = 2*( interval(exy) + (exyl-2) );
2789 if (Sup(z) > Sup(dbli)) dbli = z;
2791 dbli = (interval(exy) - dbli) + (exyl + 3);
2796 up = Inf( interval(-c2) - 1022);
2798 res = lx_interval(lx_real(0),S);
2804 Dbl = floor(_double(dbl)) + 1;
2807 res = lx_interval(lx_real(0),lx_real(dbl,l_real(1)));
2811 s = Sup(interval(c2) + dbl);
2814 Dbl = floor(_double(s)) + 1;
2817 else r = (int) _double(s);
2818 res = lx_interval(lx_real(0),lx_real(-c2,comp(0.5,r+1)));
2829 res = res/(
sqr(x) +
sqr(1-y));
2837lx_cinterval
atan(
const lx_cinterval& z )
noexcept
2839 int stagsave = stagprec,
2841 if (stagprec>stagmax) stagprec = stagmax;
2855 hxl(irez), hxu(srez), hyl(iimz), hyu(simz);
2858 resxl, resxu, resyl, resyu;
2862 if( (se_zero(irez) && ge_zero(srez)) && ( iimz <= -1.0 || simz >= 1.0 ) )
2863 cxscthrow(STD_FKT_OUT_OF_DEF(
"lx_cinterval atan( const lx_cinterval& z ); points of the branch cuts are not allowed in z."));
2879 bool sqrImz_1 = (iimz==simz) && (iimz==1.0 || iimz==-1.0);
2884 rew_l = -
abs(hxl); hxl = lx_interval(0,
l_interval(sign(irez)));
2885 rew_u = -
abs(hxu); hxu = lx_interval(0,
l_interval(sign(srez)));
2890 re_atan(imz,hxl,rew_l);
2892 re_atan(imz,hxu,rew_u);
2899 lx_real rew_inf = Inf( rew_l );
2900 lx_real rew_sup = Sup( rew_l );
2901 re_vert( irez, hxl, rew_inf, rew_sup, resxl, resxu );
2906 rew_inf = Inf( rew_u );
2907 rew_sup = Sup( rew_u );
2908 lx_real res_l, res_u;
2909 re_vert( srez, hxu, rew_inf, rew_sup, res_l, res_u );
2911 if (res_l<resxl) resxl = res_l;
2912 if (res_u>resxu) resxu = res_u;
2919 lx_real abs_y_min = Inf(
abs( imz ) );
2920 if( abs_y_min > 1.0 )
2923 abs_hyl = lx_interval( abs_y_min ),
2927 if( Sup( abs_hxl ) > irez && Inf( abs_hxl ) < srez )
2933 if ( -Inf( abs_hxl ) > irez && -Sup( abs_hxl ) < srez )
2944 im_atan_l, im_atan_u;
2946 if ( sm_zero(iimz) )
2949 im_atan_l = -Q_atan(abs_rez,-hyl);
2952 im_atan_l = Q_atan(abs_rez,hyl);
2954 if ( sm_zero(simz) )
2957 im_atan_u = -Q_atan(abs_rez,-hyu);
2960 im_atan_u = Q_atan(abs_rez,hyu);
2962 resyl = min( Inf( im_atan_l ), Inf( im_atan_u ) );
2963 resyu = max( Sup( im_atan_l ), Sup( im_atan_u ) );
2968 lx_real abs_x_min = Inf(
abs( rez ) );
2970 x_extr = lx_interval( abs_x_min ),
2974 if( Inf( y_extr ) < simz && Sup( y_extr ) > iimz )
2980 rez = T_atan(abs_x_min);
2985 if( -Sup(y_extr) < simz && -Inf(y_extr) > iimz )
2991 rez = T_atan(abs_x_min);
2996 res = lx_cinterval( lx_interval(resxl,resxu),lx_interval(resyl,resyu) );
2997 stagprec = stagsave;
3004lx_cinterval
acot(
const lx_cinterval& z )
noexcept
3006 int stagsave = stagprec,
3008 if (stagprec>stagmax) stagprec = stagmax;
3022 hxl(irez), hxu(srez), hyl(iimz), hyu(simz);
3025 resxl, resxu, resyl, resyu;
3029 if( (se_zero(irez) && ge_zero(srez)) && ( iimz <= 1.0 || simz >= -1.0 ) )
3030 cxscthrow(STD_FKT_OUT_OF_DEF(
"lx_cinterval acot( const lx_cinterval& z ); points of the branch cuts are not allowed in z."));
3045 bool sqrImz_1 = (iimz==simz) && (iimz==1.0 || iimz==-1.0);
3050 rew_l =
abs(hxl); hxl = lx_interval(0,
l_interval(sign(irez)));
3051 rew_u =
abs(hxu); hxu = lx_interval(0,
l_interval(sign(srez)));
3056 re_atan(imz,hxl,rew_l);
3059 re_atan(imz,hxu,rew_u);
3067 lx_real rew_inf = Inf( rew_l );
3068 lx_real rew_sup = Sup( rew_l );
3069 re_vert( irez, hxl, rew_inf, rew_sup, resxl, resxu );
3074 rew_inf = Inf( rew_u );
3075 rew_sup = Sup( rew_u );
3076 lx_real res_l, res_u;
3077 re_vert( srez, hxu, rew_inf, rew_sup, res_l, res_u );
3079 if (res_l<resxl) resxl = res_l;
3080 if (res_u>resxu) resxu = res_u;
3087 lx_real abs_y_min = Inf(
abs( imz ) );
3089 if( abs_y_min > 1.0 )
3092 abs_hyl = lx_interval( abs_y_min ),
3096 if( Sup( abs_hxl ) > irez && Inf( abs_hxl ) < srez )
3100 resxu = Sup(
atan( 1.0 / abs_hxl ) / 2.0 );
3101 if ( -Inf( abs_hxl ) > irez && -Sup( abs_hxl ) < srez )
3104 resxl = -Sup(
atan( 1.0 / abs_hxl ) /2.0 );
3113 im_atan_l, im_atan_u;
3115 if ( sm_zero(iimz) )
3118 im_atan_l = -Q_atan(abs_rez,-hyl);
3121 im_atan_l = Q_atan(abs_rez,hyl);
3124 if ( sm_zero(simz) )
3127 im_atan_u = -Q_atan(abs_rez,-hyu);
3130 im_atan_u = Q_atan(abs_rez,hyu);
3133 resyl = min( Inf( im_atan_l ), Inf( im_atan_u ) );
3134 resyu = max( Sup( im_atan_l ), Sup( im_atan_u ) );
3139 lx_real abs_x_min = Inf(
abs( rez ) );
3141 x_extr = lx_interval( abs_x_min ),
3145 if( Inf( y_extr ) < simz && Sup( y_extr ) > iimz )
3151 rez = T_atan(abs_x_min);
3156 if( -Sup(y_extr) < simz && -Inf(y_extr) > iimz )
3162 rez = T_atan(abs_x_min);
3167 res = lx_cinterval( lx_interval(resxl,resxu),lx_interval(-resyu,-resyl) );
3168 stagprec = stagsave;
3177lx_cinterval
atanh(
const lx_cinterval& z )
noexcept
3182 lx_cinterval res =
atan( lx_cinterval( -Im(z), Re(z) ) );
3183 return lx_cinterval( Im(res), -Re(res) );
3190lx_cinterval
acoth(
const lx_cinterval& z )
noexcept
3195 lx_cinterval res =
acot( lx_cinterval( -Im(z), Re(z) ) );
3196 return lx_cinterval( -Im(res), Re(res) );
3207 const lx_real c = lx_real(1600,
l_real(1));
3208 int stagsave = stagprec,
3210 if (stagprec>stagmax) stagprec = stagmax;
3213 lx_interval absz(
abs(z));
3214 lx_real Inf_absz(Inf(absz));
3218 absz = 1 / lx_interval(Inf_absz);
3219 Inf_absz = Sup(absz);
3220 res = lx_cinterval(lx_interval(-Inf_absz,Inf_absz),
3221 lx_interval(-Inf_absz,Inf_absz));
3224 res = (Inf(Re(z))>=0)? z + res : -z + res;
3228 res = lx_cinterval(lx_interval(0,
l_interval(0)),
3230 if ( Sup(
abs(z-res))<0.5 || Sup(
abs(z+res))<0.5 )
3232 res = lx_cinterval(-Im(z),Re(z));
3234 res =
sqrt( (1-res)*(1+res) );
3239 if (sm_zero(Inf(Re(res))))
3240 res = lx_cinterval(lx_interval(lx_real(0),Sup(Re(res))),Im(res));
3241 stagprec = stagsave;
3251 const lx_real c = lx_real(1600,
l_real(1));
3252 int stagsave = stagprec,
3254 if (stagprec>stagmax) stagprec = stagmax;
3257 lx_interval absz(
abs(z));
3258 lx_real Inf_absz(Inf(absz));
3262 absz = 1 / lx_interval(Inf_absz);
3263 Inf_absz = Sup(absz);
3264 res = lx_cinterval(lx_interval(-Inf_absz,Inf_absz),
3265 lx_interval(-Inf_absz,Inf_absz));
3266 u = lx_cinterval(-Im(z),Re(z));
3269 res = (Inf(Im(z))>=0)? -u + res : u + res;
3276 if (sm_zero(Inf(Re(res))))
3277 res = lx_cinterval(lx_interval(lx_real(0),Sup(Re(res))),Im(res));
3278 stagprec = stagsave;
3287 const lx_real c = lx_real(1600,
l_real(1));
3288 int stagsave = stagprec,
3290 if (stagprec>stagmax) stagprec = stagmax;
3293 lx_interval absz(
abs(z));
3294 lx_real Inf_absz(Inf(absz));
3298 absz = 1 / lx_interval(Inf_absz);
3299 Inf_absz = Sup(absz);
3300 res = lx_cinterval(lx_interval(-Inf_absz,Inf_absz),
3301 lx_interval(-Inf_absz,Inf_absz));
3303 res = (Inf(Re(z))>=0)? z + res : -z + res;
3311 if (sm_zero(Inf(Re(res))))
3312 res = lx_cinterval(lx_interval(lx_real(0),Sup(Re(res))),Im(res));
3314 stagprec = stagsave;
3323 const real c = 0.125;
3324 int stagsave = stagprec,
3326 if (stagprec>stagmax) stagprec = stagmax;
3329 lx_interval absz(
abs(z));
3330 lx_real Sup_absz(Sup(absz));
3333 res = z / (
sqrt(z+1) + 1);
3335 res =
sqrt(z+1) - 1;
3337 stagprec = stagsave;
3342lx_cinterval
expm1(
const lx_cinterval& z)
noexcept
3346 int stagsave = stagprec,
3348 if (stagprec>stagmax) stagprec = stagmax;
3350 const lx_interval cancl_test = lx_interval(0,
l_interval(0.995,1.005));
3351 lx_interval rez(Re(z)), imz(Im(z));
3352 lx_interval exp_x, sin_y, cos_y, h, Xt;
3360 if (h < cancl_test && cos_y < cancl_test)
3369 imz = exp_x * sin_y;
3370 res = lx_cinterval(h,imz);
3372 stagprec = stagsave;
3378lx_cinterval
lnp1(
const lx_cinterval& z)
noexcept
3381 int stagsave = stagprec,
3383 if (stagprec > stagmax) stagprec = stagmax;
3384 const real c1 = 1e-128;
3386 lx_interval abs_z(
abs(z));
3388 srez = Sup( Re(z) ),
3389 simz = Sup( Im(z) ),
3390 iimz = Inf( Im(z) );
3393 cxscthrow(STD_FKT_OUT_OF_DEF(
3394 "lx_cinterval lnp1(const lx_cinterval& z); z contains -1"));
3395 if ( srez<-1 && iimz<0 && simz>=0 )
3396 cxscthrow(STD_FKT_OUT_OF_DEF(
3397 "lx_cinterval lnp1(const lx_cinterval& z); z not allowed"));
3398 if (Sup(abs_z) < c1)
3401 abs_z =
lnp1( abs_z*(2+abs_z) +
sqr(Im(z)) );
3403 y = lx_cinterval( abs_z,
arg(1+z) );
3407 stagprec = stagsave;
3423std::list<lx_cinterval>
pow_all(
const lx_cinterval& z,
3424 const lx_interval& p )
noexcept
3426 lx_interval abs_z =
abs(z);
3428 if( 0.0 < Inf( abs_z ) )
3430 lx_interval abs_z_p =
exp( p *
ln( abs_z ) );
3435 lx_interval rad_2 = lx_interval(Sup( abs_z_p ));
3437 std::list<lx_cinterval> res;
3440 res.push_back(lx_cinterval(lx_interval( Inf( rad_1 ), Sup( rad_2 ) ),
3441 lx_interval( -Sup( rad_1 ), Sup( rad_2 ) )));
3442 res.push_back(lx_cinterval(lx_interval( -Sup( rad_2 ), Sup( rad_1 ) ),
3443 lx_interval( Inf( rad_1 ), Sup( rad_2 ) ) ));
3444 res.push_back(lx_cinterval(lx_interval( -Sup( rad_2 ), -Inf( rad_1 ) ),
3445 lx_interval( -Sup( rad_2 ), Sup(rad_1) ) ) );
3446 res.push_back(lx_cinterval(lx_interval( -Sup( rad_1 ), Sup( rad_2 ) ),
3447 lx_interval( -Sup( rad_2 ), -Inf(rad_1) ) ) );
3454 if( Inf( p ) > 0.0 )
3459 lx_interval abs_z_p =
exp( p *
ln( lx_interval( Sup( abs_z ) ) ) );
3460 lx_real rad_p = Sup( abs_z_p );
3462 std::list<lx_cinterval> res;
3464 res.push_back( lx_cinterval( lx_interval( -rad_p, rad_p ),
3465 lx_interval( -rad_p, rad_p ) ) );
3474 cxscthrow(STD_FKT_OUT_OF_DEF(
3475 "pow_all(lx_cinterval, lx_interval); 0^p is undefined for p <= 0."));
3476 std::list<lx_cinterval> res;
The Scalar Type cinterval.
cinterval & operator=(const real &) noexcept
Implementation of standard assigning operator.
The Scalar Type interval.
The Multiple-Precision Data Type l_cinterval.
The Multiple-Precision Data Type l_interval.
The Multiple-Precision Data Type l_real.
The namespace cxsc, providing all functionality of the class library C-XSC.
cinterval sqrtp1m1(const cinterval &z) noexcept
Calculates .
cinterval exp2(const cinterval &z) noexcept
Calculates .
cinterval sqrt1mx2(const cinterval &z) noexcept
Calculates .
const real MinReal
Smallest normalized representable floating-point number.
cinterval asinh(const cinterval &z) noexcept
Calculates .
cinterval coth(const cinterval &z) noexcept
Calculates .
cinterval log2(const cinterval &z) noexcept
Calculates .
l_interval Pi_l_interval() noexcept
Enclosure-Interval for .
cinterval power(const cinterval &z, int n) noexcept
Calculates .
cinterval log10(const cinterval &z) noexcept
Calculates .
cinterval Ln(const cinterval &z) noexcept
Calculates .
int Disjoint(const interval &a, const interval &b)
Checks arguments for disjointness.
cinterval ln(const cinterval &z) noexcept
Calculates .
lx_interval Sqrt2_lx_interval() noexcept
Enclosure-Interval for .
bool Is_Integer(const real &x)
Returns 1 if x is an integer value and if .
int expo_gr(const l_interval &x)
lx_interval Ln2_lx_interval() noexcept
Enclosure-Interval for .
lx_interval Pid4_lx_interval() noexcept
Enclosure-Interval for .
cvector diam(const cimatrix_subv &mv) noexcept
Returns the diameter of the matrix.
const real minreal
Smallest positive denormalized representable floating-point number.
cinterval pow(const cinterval &z, const interval &p) noexcept
Calculates .
lx_interval Ln10_lx_interval() noexcept
Enclosure-Interval for .
cinterval sinh(const cinterval &z) noexcept
Calculates .
cinterval asin(const cinterval &z) noexcept
Calculates .
cinterval tan(const cinterval &z) noexcept
Calculates .
cinterval exp10(const cinterval &z) noexcept
Calculates .
interval arg(const cinterval &z) noexcept
Calculates .
std::list< cinterval > sqrt_all(const cinterval &z)
Calculates and returns all possible solutions.
cinterval acos(const cinterval &z) noexcept
Calculates .
l_interval Sqrt2r_l_interval() noexcept
Enclosure-Interval for .
cinterval sqrtx2m1(const cinterval &z) noexcept
Calculates .
real cutint(const real &x) noexcept
Returns the truncated integer part of x.
cinterval acosh(const cinterval &z) noexcept
Calculates .
cinterval cosh(const cinterval &z) noexcept
Calculates .
cinterval cos(const cinterval &z) noexcept
Calculates .
lx_interval One_p_lx_interval() noexcept
Enclosure-Interval for .
cinterval sqrt1px2(const cinterval &z) noexcept
Calculates .
cinterval exp(const cinterval &z) noexcept
Calculates .
lx_interval Pi_lx_interval() noexcept
Enclosure-Interval for .
cinterval tanh(const cinterval &z) noexcept
Calculates .
interval ln_sqrtx2y2(const interval &x, const interval &y) noexcept
Calculates .
lx_interval Pid2_lx_interval() noexcept
Enclosure-Interval for .
std::list< cinterval > pow_all(const cinterval &z, const interval &p) noexcept
Calculates and returns all possible solutions.
cinterval expm1(const cinterval &z) noexcept
Calculates .
cinterval cot(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 .
cinterval power_fast(const cinterval &z, int n) noexcept
Calculates .
cinterval acot(const cinterval &z) noexcept
Calculates .
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 .
cvector mid(const cimatrix_subv &mv) noexcept
Returns the middle of the matrix.
cinterval atan(const cinterval &z) noexcept
Calculates .
cinterval atanh(const cinterval &z) noexcept
Calculates .
interval Arg(const cinterval &z) noexcept
Calculates .
cinterval acoth(const cinterval &z) noexcept
Calculates .
interval sqrtx2y2(const interval &x, const interval &y) noexcept
Calculates .
lx_interval One_m_lx_interval() noexcept
Enclosure-Interval for .
cinterval sin(const cinterval &z) noexcept
Calculates .