C-XSC - A C++ Class Library for Extended Scientific Computing 2.5.4
cinterval.cpp
1/*
2** CXSC is a C++ library for eXtended Scientific Computing (V 2.5.4)
3**
4** Copyright (C) 1990-2000 Institut fuer Angewandte Mathematik,
5** Universitaet Karlsruhe, Germany
6** (C) 2000-2014 Wiss. Rechnen/Softwaretechnologie
7** Universitaet Wuppertal, Germany
8**
9** This library is free software; you can redistribute it and/or
10** modify it under the terms of the GNU Library General Public
11** License as published by the Free Software Foundation; either
12** version 2 of the License, or (at your option) any later version.
13**
14** This library is distributed in the hope that it will be useful,
15** but WITHOUT ANY WARRANTY; without even the implied warranty of
16** MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
17** Library General Public License for more details.
18**
19** You should have received a copy of the GNU Library General Public
20** License along with this library; if not, write to the Free
21** Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
22*/
23
24/* CVS $Id: cinterval.cpp,v 1.20 2014/01/30 17:23:44 cxsc Exp $ */
25
26#include "cinterval.hpp"
27#include "cidot.hpp"
28#include "dot.hpp"
29#include "rmath.hpp"
30#include "imath.hpp"
31
32namespace cxsc {
33
34#define CXSC_Zero 0.0
35
36cinterval::cinterval(const dotprecision &a) noexcept : re(a),im(0) {}
37cinterval::cinterval(const idotprecision &a) noexcept : re(a),im(0) {}
38cinterval::cinterval(const cdotprecision &a) noexcept : re(Re(a)),im(Im(a)) {}
40 re(rnd(InfRe(a),RND_DOWN),rnd(SupRe(a),RND_UP)),
41 im(rnd(InfIm(a),RND_DOWN),rnd(SupIm(a),RND_UP))
42{
43}
44
45
46cinterval mult_operator(const cinterval & a,const cinterval & b) noexcept
47{
48 cidotprecision akku;
49 akku=0.0;
50 accumulate(akku,a,b);
51 return rnd(akku);
52}
53
54// In complex.cpp
55void product (real, real, real, real, int&, real&, interval&);
56real quotient (real, interval, real, interval, int, int, int);
57
58// optimale komplexe Intervalldivision
59
60bool cxsc_complex_division_p[5];
61
62real cxsc_complex_division_f(real a, real b, real c, real d, int round)
63{
64 int zoverfl,noverfl;
65 real z1,n1;
66 interval z2,n2;
67
68 // f:=(a*c+b*d)/(SQR(c)+SQR(d))
69
70 product( a, c, b, d, zoverfl, z1,z2 );
71 product( c, c, d, d, noverfl, n1,n2 );
72
73 return quotient( z1,z2, n1,n2, round, zoverfl, noverfl );
74}
75
76static real minmax(int minimum, real a, real b, real y0,
77 interval x, int i, int j)
78// Calculates the inner minimum or maximum of f = (ac+bd)/(cc+dd)
79// on the interval x = [c.inf,c.sup] ( a,b,d=y0 fixated ).
80// If minimum = true the minimum will be calculated, otherwise the maximum.
81
82{
83 real q1,q2,w1,w2,t1,t2,x1,x2,ay0, minmax;
84 dotprecision akku;
85 bool scaling = false; // scaling = true <==> scaling is necessary
86
87 if (minimum)
88 minmax = MaxReal;
89 else
90 minmax = -MaxReal;
91
92 if (Inf(x) == Sup(x))
93 {
94 if (cxsc_complex_division_p[i] && cxsc_complex_division_p[j])
95 minmax = cxsc_complex_division_f( a, b, Inf(x), y0, 1-2*minimum );
96
97 cxsc_complex_division_p[i] = false;
98 cxsc_complex_division_p[j] = false;
99 } else
100 if (a == 0.0)
101 {
102 if ( b == CXSC_Zero || y0 == CXSC_Zero )
103 {
104 minmax = 0.0;
105 cxsc_complex_division_p[i] = false;
106 cxsc_complex_division_p[j] = false;
107 } else
108 {
109 if (0.0 < x) {
110 if (minimum && sign(b) != sign(y0) )
111 {
112 minmax = divd(b, y0);
113 cxsc_complex_division_p[i] = false;
114 cxsc_complex_division_p[j] = false;
115 } else
116 if (!minimum && sign(b) == sign(y0) )
117 {
118 minmax = divu(b, y0);
119 cxsc_complex_division_p[i] = false;
120 cxsc_complex_division_p[j] = false;
121 }
122 }
123 }
124 } else
125 {
126 // a != 0.0
127 if (y0 == 0.0)
128 {
129 if (minimum)
130 {
131 if (a > 0.0)
132 minmax = divd(a, Sup(x));
133 else
134 minmax = divd(a, Inf(x));
135 } else
136 {
137 if (a > 0.0)
138 minmax = divu(a, Inf(x));
139 else
140 minmax = divu(a, Sup(x));
141 }
142 cxsc_complex_division_p[i] = false;
143 cxsc_complex_division_p[j] = false;
144 } else
145 {
146 // y0 != 0.0, Calculation of extrema points and minimum|maximum
147 // values:
148 // IF NOTBEKANNT THEN
149 // Calculation of: t = sign(a) * ( |b/a| + sqrt( 1+|b/a|^2 ) )
150 // in staggered presentation: t ~ t1 + t2
151
152 // Exponent over-/undeflow in |b/a| is now considered from
153 // Blomquist/Hofschuster 06.11.02.
154
155 real invf2=1, a_skal;
156 int exf1=0, exf2=0, exinf1=0, exinf2=0;
157
158 if (sign(b)==0) { t1 = 1; t2 = 0; }
159 else
160 { // b != 0, ---> expo(b) != -2147483647;
161 if (a>0.0)
162 q2 = abs(b);
163 else
164 q2 = -abs(b);
165
166 // Skaling to avoid overflow by division b/a:
167 int expo_diff = expo(b)-expo(a), ex;
168 if (expo_diff >= 512)
169 {
170 int exdiff5 = expo_diff - 500;
171 scaling = true;
172 if (exdiff5 > 1024) // Two scaling factors necessary!
173 {
174 ex = exdiff5 / 2;
175 exf1 = ex-1;
176 exf2 = exdiff5 - ex;
177 exinf1 = -exf1;
178 invf2 = comp(0.5,1-exf2);
179 exinf2 = -exf2;
180 a_skal = a;
181 times2pown(a_skal,exf1);
182 times2pown(a_skal,exf2);
183 }
184 else // Scaling with only one factor!
185 {
186 exf2 = exdiff5-1; // exf1 = 0;
187 invf2 = comp(0.5,2-exdiff5); // invf1 = 1;
188 exinf2 = -exf2; // exinf1 = 0;
189 a_skal = a;
190 times2pown(a_skal,exf2);
191 }
192 }
193 else // Scaling not necessary!
194 { a_skal = a; }
195
196 q1 = q2/a_skal;
197 akku = q2;
198 accumulate(akku, -q1, a_skal);
199 q2 = rnd(akku) / a_skal;
200
201 akku = 0.0;
202 if (exinf1 == 0) accumulate(akku, invf2, invf2);
203 accumulate(akku, q1, q1);
204 accumulate(akku, q1, q2);
205 accumulate(akku, q1, q2);
206 accumulate(akku, q2, q2);
207
208 w1 = sqrt(rnd(akku));
209
210 accumulate(akku, -w1, w1);
211 w2 = rnd(akku) / (2.0*w1);
212
213 akku = q1;
214 akku += q2;
215 akku += w1;
216 akku += w2;
217
218 t1 = rnd(akku);
219
220 akku -= t1;
221 t2 = rnd(akku);
222 }
223 if (a<0.0) // if (a_skale<0.0)
224 {
225 t1 = -t1;
226 t2 = -t2;
227 }
228
229 // Fall differentiation for min-,max- calculation:
230 ay0 = abs(y0);
231 if (( sign(b) == sign(y0) ) == minimum)
232 {
233 // Calculation of: x1 + x2 = |y0| * ( t1 + t2 )
234 akku = 0.0;
235 accumulate(akku,ay0,t1);
236 accumulate(akku,ay0,t2);
237 x1 = rnd(akku);
238 if (expo(x1) == 2147483647) goto Ende;
239 akku -= x1;
240 x2 = rnd(akku);
241 if (scaling)
242 {
243 if (expo(x1)+exf1 > 1023) goto Ende;
244 times2pown(x1,exf1);
245 if (expo(x1)+exf2 > 1023) goto Ende;
246 times2pown(x1,exf2);
247 times2pown(x2,exf1);
248 times2pown(x2,exf2);
249 }
250 } else
251 {
252 // Calculation of: x1 + x2 = |y0| / ( t1 + t2 )
253 x1 = ay0 / t1;
254 akku = ay0;
255 accumulate(akku, -t1, x1);
256 accumulate(akku, -t2, x1);
257 x2 = rnd(akku) / t1;
258 if (scaling)
259 {
260 if (expo(x1)+exinf1 > 1023) goto Ende;
261 times2pown(x1,exinf1);
262 if (expo(x1)+exinf2 > 1023) goto Ende;
263 times2pown(x1,exinf2);
264 times2pown(x2,exinf1);
265 times2pown(x2,exinf2);
266 }
267 }
268 if (minimum)
269 {
270 x1 = -x1;
271 x2 = -x2;
272 }
273
274 if (x1 < x)
275 {
276 // Calculation of: a / ( 2*(x1+x2) )
277 q1 = a/(2*x1);
278 akku = 0.0;
279 accumulate(akku, -x1, q1); // vorher: accumulate(akku, x1, q1);
280 accumulate(akku, -x2, q1);
281
282 // exact calculation of (a + akku + akku) in new variable akku:
283 akku += akku;
284 akku += a;
285 q2 = rnd(akku) / (2.0*x1);
286 if (minimum)
287 {
288 if (sign(q2)==0 && sign(akku)!=0) minmax = pred(q1);
289 else minmax = addd(q1, q2);
290 }
291 else
292 {
293 if (sign(q2)==0 && sign(akku)!=0) minmax = succ(q1);
294 else minmax = addu(q1, q2);
295 }
296
297 cxsc_complex_division_p[i] = false;
298 cxsc_complex_division_p[j] = false;
299 }
300 Ende:;
301 } // y0 != 0.0
302 }
303 return minmax;
304} // *** minmax ***
305
306cinterval cidiv(const cinterval& A, const cinterval& B)
307{
308 real realteilINF, realteilSUP,
309 imagteilINF, imagteilSUP;
310 // da sonst eventuell zwischendurch illegale Intervalle entstehen
311 real a0,b0;
312 bool a_repeat,b_repeat;
313 int i, rep, j;
314 real AREINF, ARESUP, AIMINF, AIMSUP,
315 BREINF, BRESUP, BIMINF, BIMSUP;
316 interval ARE, AIM, BRE, BIM;
317
318 // keine Fehlerabfrage -> ist schon in CINTVAL.CPP
319 // IF ( 0.0 IN B.RE ) AND ( 0.0 IN B.IM ) THEN
320 // CIDIVISION:= COMPL( 1.0 / INTVAL(-1.0,1.0), INTVAL(0.0) );
321 // Fehlerabbruch erzwingen: Intervall enthaelt 0
322
323 // *** Berechnung des Realteils ***
324
325 AREINF = Inf(Re(A));
326 ARESUP = Sup(Re(A));
327 AIMINF = Inf(Im(A));
328 AIMSUP = Sup(Im(A));
329 BREINF = Inf(Re(B));
330 BRESUP = Sup(Re(B));
331 BIMINF = Inf(Im(B));
332 BIMSUP = Sup(Im(B));
333 ARE = Re(A);
334 AIM = Im(A);
335 BRE = Re(B);
336 BIM = Im(B);
337
338 a_repeat = ( BREINF < CXSC_Zero ) && ( CXSC_Zero < BRESUP );
339 b_repeat = ( BIMINF < CXSC_Zero ) && ( CXSC_Zero < BIMSUP );
340
341 if (a_repeat || b_repeat)
342 rep = 2;
343 else
344 rep = 1;
345
346 if (BREINF >= 0.0)
347 a0 = ARESUP;
348 else
349 a0 = AREINF;
350
351 if (BIMINF >= 0.0)
352 b0 = AIMSUP;
353 else
354 b0 = AIMINF;
355
356 realteilSUP = -MaxReal;
357
358 for (j=1; j<=rep; j++)
359 {
360 for (i=1; i<=4; cxsc_complex_division_p[i++] = true);
361
362 realteilSUP =
363 max( realteilSUP,
364 max( max( minmax( false, a0, b0, BIMINF, BRE, 1,2 ),
365 minmax( false, a0, b0, BIMSUP, BRE, 3,4 ) ),
366 max( minmax( false, b0, a0, BREINF, BIM, 1,3 ),
367 minmax( false, b0, a0, BRESUP, BIM, 2,4 ) ) )
368
369 );
370
371 if (cxsc_complex_division_p[1])
372 realteilSUP = max( realteilSUP, cxsc_complex_division_f( a0, b0, BREINF, BIMINF, +1 ) );
373 if (cxsc_complex_division_p[2])
374 realteilSUP = max( realteilSUP, cxsc_complex_division_f( a0, b0, BRESUP, BIMINF, +1 ) );
375 if (cxsc_complex_division_p[3])
376 realteilSUP = max( realteilSUP, cxsc_complex_division_f( a0, b0, BREINF, BIMSUP, +1 ) );
377 if (cxsc_complex_division_p[4])
378 realteilSUP = max( realteilSUP, cxsc_complex_division_f( a0, b0, BRESUP, BIMSUP, +1 ) );
379
380 if (a_repeat)
381 a0 = ARESUP;
382 else if (b_repeat)
383 b0 = AIMSUP;
384 }
385
386 if (BREINF >= 0.0)
387 a0 = AREINF;
388 else
389 a0 = ARESUP;
390 if (BIMINF >= 0.0)
391 b0 = AIMINF;
392 else
393 b0 = AIMSUP;
394
395 realteilINF = MaxReal;
396
397 for (j=1; j<=rep; j++)
398 {
399 for (i=1; i<=4; cxsc_complex_division_p[i++] = true);
400
401 realteilINF =
402 min( realteilINF,
403 min( min( minmax( true, a0, b0, BIMINF, BRE, 1,2 ),
404 minmax( true, a0, b0, BIMSUP, BRE, 3,4 ) ),
405 min( minmax( true, b0, a0, BREINF, BIM, 1,3 ),
406 minmax( true, b0, a0, BRESUP, BIM, 2,4 ) ) )
407 );
408 if (cxsc_complex_division_p[1])
409 realteilINF = min( realteilINF, cxsc_complex_division_f( a0, b0, BREINF, BIMINF, -1 ) );
410 if (cxsc_complex_division_p[2])
411 realteilINF = min( realteilINF, cxsc_complex_division_f( a0, b0, BRESUP, BIMINF, -1 ) );
412 if (cxsc_complex_division_p[3])
413 realteilINF = min( realteilINF, cxsc_complex_division_f( a0, b0, BREINF, BIMSUP, -1 ) );
414 if (cxsc_complex_division_p[4])
415 realteilINF = min( realteilINF, cxsc_complex_division_f( a0, b0, BRESUP, BIMSUP, -1 ) );
416
417 if (a_repeat)
418 a0 = AREINF;
419 else if (b_repeat)
420 b0 = AIMINF;
421 }
422
423
424 // Calculation of the img. part:
425 // g(a, b, c, d) = cxsc_complex_division_f(b, -a, c, d)
426
427 a_repeat = ( BIMINF < CXSC_Zero ) && ( CXSC_Zero < BIMSUP );
428 b_repeat = ( BREINF < CXSC_Zero ) && ( CXSC_Zero < BRESUP );
429
430 // IF a_repeat OR b_repeat THEN rep:= 2 ELSE rep:= 1;
431
432 if (BREINF >= 0.0)
433 b0 = AIMSUP;
434 else
435 b0 = AIMINF;
436
437 if (BIMINF >= 0.0)
438 a0 = AREINF;
439 else
440 a0 = ARESUP;
441
442 imagteilSUP = -MaxReal;
443
444 for (j=1; j<=rep; j++)
445 {
446 for (i=1; i<=4; cxsc_complex_division_p[i++] = true) ;
447
448 imagteilSUP =
449 max( imagteilSUP,
450 max( max( minmax( false, b0, -a0, BIMINF, BRE, 1,2 ),
451 minmax( false, b0, -a0, BIMSUP, BRE, 3,4 ) ),
452 max( minmax( false, -a0, b0, BREINF, BIM, 1,3 ),
453 minmax( false, -a0, b0, BRESUP, BIM, 2,4 ) ) )
454 );
455
456 if (cxsc_complex_division_p[1])
457 imagteilSUP = max( imagteilSUP, cxsc_complex_division_f( b0, -a0, BREINF, BIMINF, +1 ) );
458 if (cxsc_complex_division_p[2])
459 imagteilSUP = max( imagteilSUP, cxsc_complex_division_f( b0, -a0, BRESUP, BIMINF, +1 ) );
460 if (cxsc_complex_division_p[3])
461 imagteilSUP = max( imagteilSUP, cxsc_complex_division_f( b0, -a0, BREINF, BIMSUP, +1 ) );
462 if (cxsc_complex_division_p[4])
463 imagteilSUP = max( imagteilSUP, cxsc_complex_division_f( b0, -a0, BRESUP, BIMSUP, +1 ) );
464
465 if (b_repeat)
466 b0 = AIMSUP;
467 else if (a_repeat)
468 a0 = AREINF ;
469 }
470
471 if (BREINF >= 0.0)
472 b0 = AIMINF;
473 else
474 b0 = AIMSUP;
475
476 if (BIMINF >= 0.0)
477 a0 = ARESUP;
478 else
479 a0 = AREINF;
480
481 imagteilINF = MaxReal;
482
483 for (j=1; j<=rep; j++)
484 {
485 for (i=1; i<=4; cxsc_complex_division_p[i++] = true) ;
486
487 imagteilINF =
488 min( imagteilINF,
489 min( min( minmax( true, b0, -a0, BIMINF, BRE, 1,2 ),
490 minmax( true, b0, -a0, BIMSUP, BRE, 3,4 ) ),
491 min( minmax( true, -a0, b0, BREINF, BIM, 1,3 ),
492 minmax( true, -a0, b0, BRESUP, BIM, 2,4 ) ) )
493 );
494
495 if (cxsc_complex_division_p[1])
496 imagteilINF = min( imagteilINF, cxsc_complex_division_f( b0, -a0, BREINF, BIMINF, -1 ) );
497 if (cxsc_complex_division_p[2])
498 imagteilINF = min( imagteilINF, cxsc_complex_division_f( b0, -a0, BRESUP, BIMINF, -1 ) );
499 if (cxsc_complex_division_p[3])
500 imagteilINF = min( imagteilINF, cxsc_complex_division_f( b0, -a0, BREINF, BIMSUP, -1 ) );
501 if (cxsc_complex_division_p[4])
502 imagteilINF = min( imagteilINF, cxsc_complex_division_f( b0, -a0, BRESUP, BIMSUP, -1 ) );
503
504 if (b_repeat)
505 b0 = AIMINF;
506 else if (a_repeat)
507 a0 = ARESUP;
508 }
509
510 return cinterval(interval(realteilINF, realteilSUP),
511 interval(imagteilINF, imagteilSUP));
512} // CIDIVISION
513
514cinterval C_point_div(const cinterval& z, const cinterval& n)
515// Division of complex point intervals;
516// z,n must be point intervals!! Blomquist, 07,11.02
517// This function only for internal use!
518{
519 complex a,b,q1,q2;
520 a = complex(InfRe(z),InfIm(z));
521 b = complex(InfRe(n),InfIm(n));
522 q1 = divd(a,b);
523 q2 = divu(a,b);
524
525 interval re, im;
526 re = interval( Re(q1),Re(q2) );
527 im = interval( Im(q1),Im(q2) );
528
529 return cinterval(re,im);
530} // C_point_div
531
532cinterval div_operator (const cinterval & a, const cinterval & b)
533{
534 bool a_point, b_point;
535 a_point = InfRe(a)==SupRe(a) && InfIm(a)==SupIm(a);
536 b_point = InfRe(b)==SupRe(b) && InfIm(b)==SupIm(b);
537 if(a_point && b_point) return C_point_div(a,b); // a,b are point intervals
538 else return cidiv(a,b);
539}
540
541interval abs(const cinterval &a) noexcept
542{
543// idotakku[2]=0;
544// accumulate(idotakku[2],a.re,a.re);
545// accumulate(idotakku[2],a.im,a.im);
546// return sqrt(rnd(idotakku[2]));
547 return sqrtx2y2(a.re,a.im);
548}
549
550
551// ---- Ausgabefunkt. ---------------------------------------
552
553std::ostream & operator << (std::ostream &s, const cinterval& a) noexcept
554{
555 s << '('
556 << a.re << ','
557 << a.im
558 << ')';
559 return s;
560}
561std::string & operator << (std::string &s, const cinterval &a) noexcept
562{
563 s+='(';
564 s << a.re;
565 s+=',';
566 s << a.im;
567 s+=')';
568 return s;
569}
570
571std::istream & operator >> (std::istream &s, cinterval &a)
572{ // New version for cinterval input; Blomquist, 27.10.02;
573 char c;
574 skipeolnflag = inpdotflag = true;
575 c = skipwhitespacessinglechar (s, '(');
576 if (inpdotflag) s.putback(c);
577 c = skipwhitespacessinglechar (s, '[');
578 if (inpdotflag) s.putback(c);
579 s >> SaveOpt >> RndDown >> Inf(a.re);
580 skipeolnflag = inpdotflag = true;
581 c = skipwhitespacessinglechar (s, ',');
582 if (inpdotflag) s.putback(c);
583 s >> RndUp >> Sup(a.re);
584 c = skipwhitespacessinglechar (s, ']');
585 if (inpdotflag) s.putback(c);
586 c = skipwhitespacessinglechar (s, ',');
587 if (inpdotflag) s.putback(c);
588
589 c = skipwhitespacessinglechar (s, '[');
590 if (inpdotflag) s.putback(c);
591 s >> RndDown >> Inf(a.im);
592 skipeolnflag = inpdotflag = true;
593 c = skipwhitespacessinglechar (s, ',');
594 if (inpdotflag) s.putback(c);
595 s >> RndUp >> Sup(a.im) >> RestoreOpt;
596
597 if (!waseolnflag)
598 {
599 skipeolnflag = false, inpdotflag = true;
600 c = skipwhitespaces (s);
601 if (inpdotflag && c != ']')
602 s.putback(c);
603 }
604 if (!waseolnflag)
605 {
606 skipeolnflag = false, inpdotflag = true;
607 c = skipwhitespaces (s);
608 if (inpdotflag && c != ')')
609 s.putback(c);
610 }
611 if (Inf(a.re) > Sup(a.re) || Inf(a.im) > Sup(a.im))
612 cxscthrow(ERROR_CINTERVAL_EMPTY_INTERVAL("std::istream & operator >> (std::istream &s, cinterval &a)"));
613
614 return s;
615}
616
617std::string & operator >> (std::string &s, cinterval &a)
618{
619 s = skipwhitespacessinglechar (s, '(');
620 s = skipwhitespacessinglechar (s, '[');
621 s = s >> SaveOpt >> RndDown >> Inf(a.re);
622 s = skipwhitespacessinglechar (s, ',');
623 s = s >> RndUp >> Sup(a.re);
624 s = skipwhitespacessinglechar (s, ']');
625 s = skipwhitespacessinglechar (s, ',');
626 s = skipwhitespacessinglechar (s, '[');
627 s = s >> RndDown >> Inf(a.im);
628 s = skipwhitespacessinglechar (s, ',');
629 s = s >> RndUp >> Sup(a.im) >> RestoreOpt;
630 s = skipwhitespaces (s);
631 if (s[0] == ']')
632 s.erase(0,1);
633 s = skipwhitespaces (s);
634 if (s[0] == ')')
635 s.erase(0,1);
636
637 if (Inf(a.re) > Sup(a.re) || Inf(a.im) > Sup(a.im))
638 cxscthrow(ERROR_CINTERVAL_EMPTY_INTERVAL("std::string & operator >> (std::string &s, cinterval &a)"));
639
640 return s;
641}
642
643void operator >>(const std::string &s,cinterval &a)
644{
645 std::string r(s);
646 r>>a;
647}
648void operator >>(const char *s,cinterval &a)
649{
650 std::string r(s);
651 r>>a;
652}
653
654int in ( const cinterval& x, const cinterval& y ) // Contained-in-the-interior relation
655{ //-----------------------------------
656 return ( in(Re(x),Re(y)) && in(Im(x),Im(y)) );
657}
665cinterval Blow ( cinterval x, const real& eps ) // Epsilon inflation
666{ //------------------
667 return cinterval(Blow(Re(x),eps),Blow(Im(x),eps));
668}
669} // namespace cxsc
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
The Data Type cdotprecision.
Definition cdot.hpp:61
The Data Type cidotprecision.
Definition cidot.hpp:58
The Scalar Type cinterval.
Definition cinterval.hpp:55
cinterval(void) noexcept
Constructor of class cinterval.
Definition cinterval.hpp:64
The Data Type dotprecision.
Definition dot.hpp:112
The Data Type idotprecision.
Definition idot.hpp:48
The Scalar Type interval.
Definition interval.hpp:55
The Scalar Type real.
Definition real.hpp:114
The namespace cxsc, providing all functionality of the class library C-XSC.
Definition cdot.cpp:29
const real MaxReal
Greatest representable floating-point number.
Definition real.cpp:65
int in(const cinterval &x, const cinterval &y)
Checks if first argument is part of second argument.
ivector abs(const cimatrix_subv &mv) noexcept
Returns the absolute value of the matrix.
Definition cimatrix.inl:737
cinterval sqrt(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:1007
void times2pown(cinterval &x, int n) noexcept
Fast multiplication of reference parameter [z] with .
Definition cimath.cpp:2059
interval sqrtx2y2(const interval &x, const interval &y) noexcept
Calculates .
Definition imath.cpp:80
cinterval Blow(cinterval x, const real &eps)
Performs an epsilon inflation.