GCC Code Coverage Report


Directory: ./
File: lib/geogram/numerics/expansion_nt.h
Date: 2026-09-07 02:28:19
Exec Total Coverage
Lines: 51 66 77.3%
Functions: 6 6 100.0%
Branches: 211 486 43.4%

Line Branch Exec Source
1 /*
2 * Copyright (c) 2000-2022 Inria
3 * All rights reserved.
4 *
5 * Redistribution and use in source and binary forms, with or without
6 * modification, are permitted provided that the following conditions are met:
7 *
8 * * Redistributions of source code must retain the above copyright notice,
9 * this list of conditions and the following disclaimer.
10 * * Redistributions in binary form must reproduce the above copyright notice,
11 * this list of conditions and the following disclaimer in the documentation
12 * and/or other materials provided with the distribution.
13 * * Neither the name of the ALICE Project-Team nor the names of its
14 * contributors may be used to endorse or promote products derived from this
15 * software without specific prior written permission.
16 *
17 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
18 * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
19 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
20 * ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE
21 * LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
22 * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
23 * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
24 * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
25 * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
26 * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
27 * POSSIBILITY OF SUCH DAMAGE.
28 *
29 * Contact: Bruno Levy
30 *
31 * https://www.inria.fr/fr/bruno-levy
32 *
33 * Inria,
34 * Domaine de Voluceau,
35 * 78150 Le Chesnay - Rocquencourt
36 * FRANCE
37 *
38 */
39
40 #ifndef GEOGRAM_NUMERICS_EXPANSION_NT
41 #define GEOGRAM_NUMERICS_EXPANSION_NT
42
43 #include <geogram/basic/common.h>
44 #include <geogram/numerics/multi_precision.h>
45 #include <geogram/basic/matrix.h>
46 #include <geogram/basic/rationalg.h>
47
48 /**
49 * \file geogram/numerics/expansion_nt.h
50 * \brief High-level interface to multi-precision arithmetics
51 * \details
52 * This file provides a "number-type" that encapsulates a (low-level)
53 * GEO::expansion object.
54 */
55
56 namespace GEO {
57
58 class expansion_nt;
59
60 /**
61 * \brief Expansion_nt (expansion Number Type) is used to compute the
62 * sign of polynoms exactly.
63 * \details Expansion_nt can be used like float and double. It supports
64 * three arithmetic operations (+,-,*), comparisons (>,>=,<,<=,==,!=)
65 * and exact sign computation. expansion_nt is a wrapper around
66 * an \ref expansion allocated on the heap. When
67 * performance is a concern, the lower-level expansion class may be
68 * used instead.
69 */
70 class GEOGRAM_API expansion_nt {
71 public:
72
73 /**
74 * \brief This type is used by the constructor that
75 * takes two expansion%s.
76 */
77 enum Operation {
78 SUM, DIFF, PRODUCT
79 };
80
81 /**
82 * \brief Constructs an uninitialized expansion_nt.
83 */
84
4/10
✓ Branch 1 taken 5684 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✓ Branch 4 taken 5684 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 240137 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 46515 times.
✗ Branch 11 not taken.
298020 expansion_nt() : rep_(nullptr) {
85 }
86
87 /**
88 * \brief Constructs a new expansion_nt from a double.
89 * \param[in] x the value to initialize this expansion.
90 */
91
0/8
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
3551 explicit expansion_nt(double x) {
92 908671 rep_ = expansion::new_expansion_on_heap(1);
93
0/8
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
3551 rep()[0] = x;
94 rep().set_length(1);
95 }
96
97 /**
98 * \brief Constructs a new expansion_nt from an expansion
99 * \details A new expansion is created on the heap and its
100 * content is initialized from \p rhs
101 * \param[in] rhs the expansion to be copied
102 */
103 222801 explicit expansion_nt(const expansion& rhs) {
104 222801 rep_ = expansion::new_expansion_on_heap(rhs.length());
105 rep().assign(rhs);
106 222801 }
107
108 /**
109 * \brief Constructs a new expansion_nt from two expansions
110 * \details A new expansion is created on the heap and its
111 * content is initialized from \p x \p op \p y. This function
112 * is used by code that combines the low-level API (expansion)
113 * with the high-level number type (expansion_nt). When returning
114 * the result of operations that combine expansion as an expansion_nt,
115 * it makes it possible to avoid copying the result of the
116 * last operation by directly assigning it to an expansion_nt.
117 * \param[in] x , y the two operands
118 * \param[in] op one of
119 * expansion_nt::SUM, expansion_nt::DIFF, expansion_nt::PRODUCT
120 */
121 668403 explicit expansion_nt(
122 Operation op, const expansion& x, const expansion& y
123 668403 ) {
124
1/4
✓ Branch 0 taken 668403 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
668403 switch(op) {
125 case SUM:
126 668403 rep_ = expansion::new_expansion_on_heap(
127 expansion::sum_capacity(x,y)
128 );
129 668403 rep_->assign_sum(x,y);
130 668403 break;
131 case DIFF:
132 rep_ = expansion::new_expansion_on_heap(
133 expansion::diff_capacity(x,y)
134 );
135 rep_->assign_diff(x,y);
136 break;
137 case PRODUCT:
138 rep_ = expansion::new_expansion_on_heap(
139 expansion::product_capacity(x,y)
140 );
141 rep_->assign_product(x,y);
142 break;
143 }
144 668403 }
145
146 /**
147 * \brief Constructs a new expansion_nt from three expansions
148 * \details A new expansion is created on the heap and its
149 * content is initialized from \p x \p op \p y \p op \p z.
150 * This function is used by code that combines the low-level
151 * API (expansion) with the high-level number type (expansion_nt).
152 * When returning the result of operations that combine expansion
153 * as an expansion_nt, it makes it possible to avoid copying
154 * the result of the last operation by directly assigning it
155 * to an expansion_nt.
156 * \param[in] x , y , z the three operands
157 * \param[in] op one of expansion_nt::SUM, expansion_nt::PRODUCT
158 */
159 197602 explicit expansion_nt(
160 Operation op,
161 const expansion& x, const expansion& y, const expansion& z
162 197602 ) {
163
1/4
✓ Branch 0 taken 197602 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
197602 switch(op) {
164 case SUM:
165 197602 rep_ = expansion::new_expansion_on_heap(
166 expansion::sum_capacity(x,y,z)
167 );
168 197602 rep_->assign_sum(x,y,z);
169 197602 break;
170 case DIFF:
171 geo_assert_not_reached;
172 break;
173 case PRODUCT:
174 rep_ = expansion::new_expansion_on_heap(
175 expansion::product_capacity(x,y,z)
176 );
177 rep_->assign_product(x,y,z);
178 break;
179 }
180 197602 }
181
182 /**
183 * \brief Constructs a new expansion_nt from four expansions
184 * \details A new expansion is created on the heap and its
185 * content is initialized from \p x \p op \p y \p op \p z \p op \p t.
186 * This function is used by code that combines the
187 * low-level API (expansion) with the high-level number type
188 * (expansion_nt). When returning
189 * the result of operations that combine expansion as an expansion_nt,
190 * it makes it possible to avoid copying the result of the
191 * last operation by directly assigning it to an expansion_nt.
192 * \param[in] x , y , z , t the four operands
193 * \param[in] op one of expansion_nt::SUM, expansion_nt::PRODUCT
194 */
195 explicit expansion_nt(
196 Operation op,
197 const expansion& x, const expansion& y,
198 const expansion& z, const expansion& t
199 ) {
200 switch(op) {
201 case SUM:
202 rep_ = expansion::new_expansion_on_heap(
203 expansion::sum_capacity(x,y,z,t)
204 );
205 rep_->assign_sum(x,y,z,t);
206 break;
207 case DIFF:
208 geo_assert_not_reached;
209 break;
210 case PRODUCT:
211 // HERE: TODO CHECK SIZE
212 const expansion& p1 = expansion_product(x,y);
213 const expansion& p2 = expansion_product(z,t);
214 rep_ = expansion::new_expansion_on_heap(
215 expansion::product_capacity(p1,p2)
216 );
217 rep_->assign_sum(p1,p2);
218 break;
219 }
220 }
221
222 /**
223 * \brief Constructs a new expansion_nt from two doubles
224 * \details A new expansion is created on the heap and its
225 * content is initialized from \p x \p op \p y. This function
226 * is used by code that combines the low-level API (expansion)
227 * with the high-level number type (expansion_nt). When returning
228 * the result of operations that combine expansion as an expansion_nt,
229 * it makes it possible to avoid copying the result of the
230 * last operation by directly assigning it to an expansion_nt.
231 * \param[in] x , y the two operands
232 * \param[in] op one of
233 * expansion_nt::SUM, expansion_nt::DIFF, expansion_nt::PRODUCT
234 */
235 3889407 explicit expansion_nt(Operation op, double x, double y) {
236
1/4
✗ Branch 0 not taken.
✓ Branch 1 taken 3889407 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
3889407 switch(op) {
237 case SUM:
238 rep_ = expansion::new_expansion_on_heap(
239 expansion::sum_capacity(x,y)
240 );
241 rep_->assign_sum(x,y);
242 break;
243 case DIFF:
244 3889407 rep_ = expansion::new_expansion_on_heap(
245 expansion::diff_capacity(x,y)
246 );
247 rep_->assign_diff(x,y);
248 break;
249 case PRODUCT:
250 rep_ = expansion::new_expansion_on_heap(
251 expansion::product_capacity(x,y)
252 );
253 rep_->assign_product(x,y);
254 break;
255 }
256 3889407 }
257
258 /**
259 * \brief Copy-constructor.
260 * \param[in] rhs the expansion to be copied
261 */
262 14260385 expansion_nt(const expansion_nt& rhs) {
263 5677870 copy(rhs);
264 }
265
266 /**
267 * \brief Move-constructor.
268 * \details Steals the expansion from \p rhs
269 * \param[in] rhs the victim expansion_nt
270 */
271 24896 expansion_nt(expansion_nt&& rhs) {
272 rep_ = nullptr;
273 std::swap(rep_, rhs.rep_);
274 }
275
276 /**
277 * \brief Assignment operator.
278 * \param[in] rhs the expansion to be copied
279 * \return the new value of this expansion (rhs)
280 */
281 1169344 expansion_nt& operator= (const expansion_nt& rhs) {
282
1/2
✓ Branch 0 taken 1169344 times.
✗ Branch 1 not taken.
1169344 if(&rhs != this) {
283 cleanup();
284 1169344 copy(rhs);
285 }
286 1169344 return *this;
287 }
288
289 /**
290 * \brief Assignment operator with move semantics
291 * \param[in] rhs the expansion to be copied
292 * \return the new value of this expansion (rhs)
293 */
294 expansion_nt& operator= (expansion_nt&& rhs) {
295
4/8
✓ Branch 0 taken 8139 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 8139 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 8139 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 4319 times.
✗ Branch 7 not taken.
28736 if(&rhs != this) {
296 cleanup();
297 std::swap(rep_, rhs.rep_);
298 }
299 return *this;
300 }
301
302 /**
303 * \brief Expansion_nt destructor.
304 * \details The stored expansion is deallocated whenever
305 * reference counting reaches 0.
306 */
307 ~expansion_nt() {
308 cleanup();
309
37/68
✓ Branch 0 taken 32061 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 108006 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 32568 times.
✓ Branch 5 taken 949920 times.
✓ Branch 6 taken 32568 times.
✓ Branch 7 taken 949920 times.
✓ Branch 8 taken 81 times.
✓ Branch 9 taken 427318 times.
✗ Branch 10 not taken.
✓ Branch 11 taken 427318 times.
✗ Branch 12 not taken.
✓ Branch 13 taken 89003 times.
✗ Branch 14 not taken.
✓ Branch 15 taken 618357 times.
✗ Branch 16 not taken.
✓ Branch 17 taken 618357 times.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✗ Branch 21 not taken.
✗ Branch 22 not taken.
✗ Branch 23 not taken.
✓ Branch 24 taken 222801 times.
✗ Branch 25 not taken.
✓ Branch 26 taken 222801 times.
✗ Branch 27 not taken.
✓ Branch 28 taken 222801 times.
✗ Branch 29 not taken.
✓ Branch 30 taken 8439 times.
✓ Branch 31 taken 238479 times.
✓ Branch 32 taken 8520 times.
✗ Branch 33 not taken.
✓ Branch 34 taken 8439 times.
✓ Branch 35 taken 81 times.
✗ Branch 36 not taken.
✓ Branch 37 taken 36144 times.
✗ Branch 38 not taken.
✓ Branch 39 taken 36144 times.
✗ Branch 40 not taken.
✓ Branch 41 taken 40401 times.
✓ Branch 42 taken 4 times.
✓ Branch 43 taken 79295 times.
✓ Branch 44 taken 305071 times.
✓ Branch 45 taken 102623 times.
✓ Branch 46 taken 305071 times.
✗ Branch 47 not taken.
✓ Branch 48 taken 12735 times.
✓ Branch 49 taken 150982 times.
✓ Branch 50 taken 12735 times.
✗ Branch 51 not taken.
✓ Branch 53 taken 67558 times.
✗ Branch 54 not taken.
✓ Branch 56 taken 67558 times.
✗ Branch 57 not taken.
✓ Branch 62 taken 5684 times.
✗ Branch 63 not taken.
✓ Branch 65 taken 1365 times.
✗ Branch 66 not taken.
✓ Branch 68 taken 1365 times.
✗ Branch 69 not taken.
✗ Branch 70 not taken.
✗ Branch 71 not taken.
✗ Branch 72 not taken.
✗ Branch 73 not taken.
✓ Branch 74 taken 465729 times.
✗ Branch 75 not taken.
9054557 }
310
311 /**
312 * \brief Optimizes the internal representation without changing the
313 * represented value
314 * \details this function can reduce the length of an expansion
315 */
316 void optimize() {
317 1844677 rep().optimize();
318 254583 }
319
320 /**
321 * \brief Flips the sign of an expansion.
322 */
323 void negate() {
324 rep().negate();
325 }
326
327 /********************************************************************/
328
329 /**
330 * \brief Adds an expansion_nt to this expansion_nt
331 * \param[in] rhs the expansion_nt to be added to this expansion_nt
332 * \return the new value of this expansion_nt
333 */
334 expansion_nt& operator+= (const expansion_nt& rhs);
335
336 /**
337 * \brief Subtracts an expansion_nt to this expansion_nt
338 * \param[in] rhs the expansion_nt to be subtracted
339 * \return the new value of this expansion_nt
340 */
341 expansion_nt& operator-= (const expansion_nt& rhs);
342
343 /**
344 * \brief Multiplies this expansion_nt by an expansion_nt
345 * \param[in] rhs the expansion_nt to multiply this expansion_nt by
346 * \return the new value of this expansion_nt
347 */
348 expansion_nt& operator*= (const expansion_nt& rhs);
349
350 /**
351 * \brief Adds a double to this expansion_nt
352 * \param[in] rhs the double to be added to this expansion_nt
353 * \return the new value of this expansion_nt
354 */
355 expansion_nt& operator+= (double rhs);
356
357 /**
358 * \brief Subtracts a double from this expansion_nt
359 * \param[in] rhs the double to be subtracted from this expansion_nt
360 * \return the new value of this expansion_nt
361 */
362 expansion_nt& operator-= (double rhs);
363
364 /**
365 * \brief Multiplies this expansion_nt by a double
366 * \details If the double is a constant (possibly negative) power of
367 * two (e.g. 0.125, 0.5, 2.0, 4.0 ...), one may use
368 * scale_fast() instead.
369 * \param[in] rhs the double to multiply this expansion_nt with
370 * \return the new value of this expansion_nt
371 */
372 expansion_nt& operator*= (double rhs);
373
374 /********************************************************************/
375
376 /**
377 * \brief Computes the sum of two expansion_nt%s
378 * \param[in] rhs the expansion_nt to be added to this expansion_nt
379 * \return the sum of this expansion_nt and \p rhs
380 */
381 expansion_nt operator+ (const expansion_nt& rhs) const;
382
383 /**
384 * \brief Computes the difference between two expansion_nt%s
385 * \param[in] rhs the expansion_nt to be subtracted from
386 * this expansion_nt
387 * \return the difference between this expansion_nt and \p rhs
388 */
389 expansion_nt operator- (const expansion_nt& rhs) const;
390
391 /**
392 * \brief Computes the product between two expansion_nt%s
393 * \param[in] rhs the expansion_nt to be multiplied by
394 * this expansion_nt
395 * \return the product between this expansion_nt and \p rhs
396 */
397 expansion_nt operator* (const expansion_nt& rhs) const;
398
399 /**
400 * \brief Computes the sum of an expansion_nt and a double.
401 * \param[in] rhs the double to be added to this expansion_nt
402 * \return the sum of this expansion_nt and \p rhs
403 */
404 expansion_nt operator+ (double rhs) const;
405
406 /**
407 * \brief Computes the difference between an expansion_nt and a double.
408 * \param[in] rhs the double to be subtracted from this expansion_nt
409 * \return the difference between this expansion_nt and \p rhs
410 */
411 expansion_nt operator- (double rhs) const;
412
413 /**
414 * \brief Computes the product between an expansion_nt and a double.
415 * \param[in] rhs the double to be multiplied by this expansion_nt
416 * \return the product between this expansion_nt and \p rhs
417 */
418 expansion_nt operator* (double rhs) const;
419
420 /********************************************************************/
421
422 /**
423 * \brief Computes the opposite of this expansion_nt.
424 * \return the opposite of this expansion_nt
425 */
426 expansion_nt operator- () const;
427
428 /********************************************************************/
429
430 /**
431 * \brief Compares two expansion_nt
432 * \return the sign of this expansion minus rhs
433 */
434 Sign compare(const expansion_nt& rhs) const {
435 return rep().compare(rhs.rep());
436 }
437
438 /**
439 * \brief Compares an expansion_nt with a double
440 * \return the sign of this expansion minus rhs
441 */
442 Sign compare(double rhs) const {
443 return rep().compare(rhs);
444 }
445
446 /**
447 * \brief Compares this expansion_nt with another one.
448 * \details Internally computes the sign of the difference
449 * between this expansion_nt and \p rhs.
450 * \return true if this expansion_nt is greater than \p rhs,
451 * false otherwise
452 */
453 bool operator> (const expansion_nt& rhs) const {
454 return (int(compare(rhs))>0);
455 }
456
457 /**
458 * \brief Compares this expansion_nt with another one.
459 * \details Internally computes the sign of the difference
460 * between this expansion_nt and \p rhs.
461 * \return true if this expansion_nt is greater or equal than \p rhs,
462 * false otherwise
463 */
464 bool operator>= (const expansion_nt& rhs) const {
465 return (int(compare(rhs))>=0);
466 }
467
468 /**
469 * \brief Compares this expansion_nt with another one.
470 * \details Internally computes the sign of the difference
471 * between this expansion_nt and \p rhs.
472 * \return true if this expansion_nt is smaller than \p rhs,
473 * false otherwise
474 */
475 bool operator< (const expansion_nt& rhs) const {
476 return (int(compare(rhs))<0);
477 }
478
479 /**
480 * \brief Compares this expansion_nt with another one.
481 * \details Internally computes the sign of the difference
482 * between this expansion_nt and \p rhs.
483 * \return true if this expansion_nt is smaller or equal than \p rhs,
484 * false otherwise
485 */
486 bool operator<= (const expansion_nt& rhs) const {
487 return (int(compare(rhs))<=0);
488 }
489
490 /**
491 * \brief Compares this expansion_nt with another one.
492 * \details Internally computes the sign of the difference
493 * between this expansion_nt and \p rhs.
494 * \return true if this expansion_nt is greater than \p rhs,
495 * false otherwise
496 */
497 bool operator> (double rhs) const {
498 return (int(compare(rhs))>0);
499 }
500
501 /**
502 * \brief Compares this expansion_nt with another one.
503 * \details Internally computes the sign of the difference
504 * between this expansion_nt and \p rhs.
505 * \return true if this expansion_nt is greater or equal than \p rhs,
506 * false otherwise
507 */
508 bool operator>= (double rhs) const {
509 return (int(compare(rhs))>=0);
510 }
511
512 /**
513 * \brief Compares this expansion_nt with another one.
514 * \details Internally computes the sign of the difference
515 * between this expansion_nt and \p rhs.
516 * \return true if this expansion_nt is smaller than \p rhs,
517 * false otherwise
518 */
519 bool operator< (double rhs) const {
520 return (int(compare(rhs))<0);
521 }
522
523 /**
524 * \brief Compares this expansion_nt with another one.
525 * \details Internally computes the sign of the difference
526 * between this expansion_nt and \p rhs.
527 * \return true if this expansion_nt is smaller or equal than \p rhs,
528 * false otherwise
529 */
530 bool operator<= (double rhs) const {
531 return (int(compare(rhs))<=0);
532 }
533
534 /********************************************************************/
535
536 /**
537 * \brief Computes an approximation of the stored
538 * value in this expansion.
539 * \return an approximation of the stored value.
540 */
541 double estimate() const {
542 return rep().estimate();
543 }
544
545 /**
546 * \brief Gets the sign of this expansion_nt.
547 * \return the sign of this expansion_nt, computed exactly.
548 */
549 Sign sign() const {
550 return rep().sign();
551 }
552
553 /**
554 * \brief Gets the length of this expansion.
555 * \return the number of components used internally
556 * to represend this expansion
557 * \note most client code will not need to use this
558 * (advanced use only).
559 */
560 index_t length() const {
561 return rep().length();
562 }
563
564 /**
565 * \brief Gets the i-th component of this expansion.
566 * \param i index of the component
567 * \return the i-th component of this expansion
568 * \pre i < length()
569 * \note most client code will not need to use this
570 * (advanced use only).
571 */
572 double component(index_t i) const {
573 geo_debug_assert(i < length());
574
2/2
✓ Branch 0 taken 17483882 times.
✓ Branch 1 taken 17601269 times.
82975314 return rep()[i];
575 }
576
577 /**
578 * \brief Constructs a new expansion_nt from an expansion.
579 * \details Used internally
580 * \param[in] rep should be a reference-counted expansion, created
581 * by new_expansion_on_heap(). Its reference counter is incremented.
582 * \note most client code will not need to use this
583 * (advanced use only).
584 */
585 17052 expansion_nt(expansion* rep) :
586 5895422 rep_(rep) {
587 }
588
589 /**
590 * \brief Gets the internal expansion that represents this
591 * expansion_nt.
592 * \return a reference to the expansion that represents
593 * this expansion_nt
594 * \note most client code will not need to use this
595 * (advanced use only).
596 */
597 expansion& rep() {
598
12/28
✓ Branch 4 taken 36063 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 36063 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 36063 times.
✗ Branch 11 not taken.
✓ Branch 13 taken 36063 times.
✗ Branch 14 not taken.
✓ Branch 16 taken 36063 times.
✗ Branch 17 not taken.
✓ Branch 19 taken 36063 times.
✗ Branch 20 not taken.
✓ Branch 22 taken 12735 times.
✗ Branch 23 not taken.
✓ Branch 25 taken 12735 times.
✗ Branch 26 not taken.
✓ Branch 28 taken 12735 times.
✗ Branch 29 not taken.
✗ Branch 31 not taken.
✗ Branch 32 not taken.
✗ Branch 34 not taken.
✗ Branch 35 not taken.
✓ Branch 36 taken 7495 times.
✗ Branch 37 not taken.
✓ Branch 38 taken 7495 times.
✗ Branch 39 not taken.
✓ Branch 40 taken 7495 times.
✗ Branch 41 not taken.
777212 return *rep_;
599 }
600
601 /**
602 * \brief Gets the internal expansion that represents
603 * this expansion_nt.
604 * \return a const reference to the expansion that represents
605 * this expansion_nt
606 * \note most client code will not need to use this
607 * (advanced use only).
608 */
609 const expansion& rep() const {
610
42/90
✓ Branch 0 taken 2235047 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2235047 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2235047 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2263207 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 239043 times.
✓ Branch 10 taken 153932 times.
✗ Branch 11 not taken.
✓ Branch 12 taken 35851 times.
✓ Branch 13 taken 1704036 times.
✓ Branch 14 taken 30418 times.
✓ Branch 15 taken 1704036 times.
✓ Branch 16 taken 153932 times.
✗ Branch 17 not taken.
✓ Branch 18 taken 36556 times.
✗ Branch 19 not taken.
✓ Branch 20 taken 31129 times.
✗ Branch 21 not taken.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
✓ Branch 25 taken 153932 times.
✗ Branch 26 not taken.
✓ Branch 27 taken 285798 times.
✗ Branch 28 not taken.
✓ Branch 29 taken 52927 times.
✗ Branch 30 not taken.
✓ Branch 31 taken 5684 times.
✓ Branch 32 taken 11583 times.
✓ Branch 33 taken 8439 times.
✗ Branch 34 not taken.
✓ Branch 35 taken 8439 times.
✗ Branch 36 not taken.
✓ Branch 37 taken 8439 times.
✗ Branch 38 not taken.
✗ Branch 39 not taken.
✗ Branch 40 not taken.
✗ Branch 41 not taken.
✗ Branch 42 not taken.
✗ Branch 43 not taken.
✗ Branch 44 not taken.
✓ Branch 45 taken 36063 times.
✗ Branch 46 not taken.
✓ Branch 47 taken 36063 times.
✗ Branch 48 not taken.
✓ Branch 49 taken 36063 times.
✗ Branch 50 not taken.
✓ Branch 51 taken 36063 times.
✗ Branch 52 not taken.
✓ Branch 53 taken 12735 times.
✗ Branch 54 not taken.
✓ Branch 55 taken 12735 times.
✗ Branch 56 not taken.
✓ Branch 57 taken 12735 times.
✗ Branch 58 not taken.
✓ Branch 59 taken 12735 times.
✗ Branch 60 not taken.
✗ Branch 61 not taken.
✗ Branch 62 not taken.
✗ Branch 63 not taken.
✗ Branch 64 not taken.
✓ Branch 66 taken 3319597 times.
✗ Branch 67 not taken.
✓ Branch 68 taken 3319597 times.
✗ Branch 69 not taken.
✓ Branch 70 taken 3319597 times.
✗ Branch 71 not taken.
✓ Branch 72 taken 861190 times.
✗ Branch 73 not taken.
✓ Branch 74 taken 861190 times.
✗ Branch 75 not taken.
✓ Branch 76 taken 861190 times.
✗ Branch 77 not taken.
✓ Branch 78 taken 7495 times.
✗ Branch 79 not taken.
✓ Branch 81 taken 632111 times.
✗ Branch 82 not taken.
✓ Branch 83 taken 632111 times.
✗ Branch 84 not taken.
✓ Branch 85 taken 632111 times.
✗ Branch 86 not taken.
✓ Branch 87 taken 7953 times.
✗ Branch 88 not taken.
✓ Branch 89 taken 7953 times.
✗ Branch 90 not taken.
✓ Branch 91 taken 7953 times.
✗ Branch 92 not taken.
87054765 return *rep_;
611 }
612
613 /**
614 * \brief Gets a string representation of this expansion
615 * \return a string with the length and components or "null"
616 * if this expansion_nt was explicitely set to uninitialized.
617 */
618 std::string to_string() const {
619 return (rep_ == nullptr) ?
620 std::string("null") :
621 rep_->to_string() ;
622 }
623
624 /**
625 * \brief Tests whether an expansion_nt is equal to one.
626 * \details Optimized using the low-level API
627 * \retval true if \p x is equal to one
628 * \retval false otherwise
629 */
630 bool is_one() const {
631 return rep().equals(1.0);
632 }
633
634 /**
635 * \brief Tests whether an expansion_nt is equal to zero
636 * \details Optimized using the low-level API
637 * \retval true if \p x is equal to zero
638 * \retval false otherwise
639 */
640 bool is_zero() const {
641 return sign() == ZERO;
642 }
643
644 protected:
645
646 /**
647 * \brief Copies an expansion into this one.
648 * \details current rep_ pointer is supposed to be
649 * uninitialized or freed before calling this function.
650 * \param[in] rhs a const reference to the expansion to be copied
651 */
652 16428745 void copy(const expansion_nt& rhs) {
653
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 16428745 times.
16428745 if(rhs.rep_ == nullptr) {
654 rep_ = nullptr;
655 } else {
656 16428745 rep_ = expansion::new_expansion_on_heap(rhs.rep().capacity());
657 rep_->set_length(rhs.rep().length());
658
2/2
✓ Branch 0 taken 37996282 times.
✓ Branch 1 taken 16428745 times.
54425027 for(index_t i=0; i<rep_->length(); ++i) {
659 37996282 (*rep_)[i] = rhs.rep()[i];
660 }
661 }
662 16428745 }
663
664 /**
665 * \brief Cleanups the memory associated with this expansion_nt.
666 */
667 void cleanup() {
668
82/202
✓ Branch 0 taken 32061 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 32061 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 32061 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 81 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 76026 times.
✗ Branch 9 not taken.
✓ Branch 10 taken 76026 times.
✗ Branch 11 not taken.
✓ Branch 12 taken 32568 times.
✗ Branch 13 not taken.
✓ Branch 14 taken 32568 times.
✗ Branch 15 not taken.
✓ Branch 16 taken 32568 times.
✗ Branch 17 not taken.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✗ Branch 21 not taken.
✗ Branch 22 not taken.
✗ Branch 23 not taken.
✓ Branch 24 taken 427318 times.
✗ Branch 25 not taken.
✓ Branch 26 taken 427318 times.
✗ Branch 27 not taken.
✓ Branch 28 taken 427318 times.
✗ Branch 29 not taken.
✓ Branch 30 taken 89003 times.
✗ Branch 31 not taken.
✓ Branch 32 taken 89003 times.
✗ Branch 33 not taken.
✓ Branch 34 taken 92823 times.
✗ Branch 35 not taken.
✓ Branch 36 taken 92823 times.
✗ Branch 37 not taken.
✓ Branch 38 taken 3820 times.
✗ Branch 39 not taken.
✓ Branch 40 taken 1479274 times.
✗ Branch 41 not taken.
✓ Branch 42 taken 1479274 times.
✗ Branch 43 not taken.
✓ Branch 44 taken 1479274 times.
✗ Branch 45 not taken.
✗ Branch 46 not taken.
✗ Branch 47 not taken.
✗ Branch 48 not taken.
✗ Branch 49 not taken.
✗ Branch 50 not taken.
✗ Branch 51 not taken.
✗ Branch 52 not taken.
✗ Branch 53 not taken.
✗ Branch 54 not taken.
✗ Branch 55 not taken.
✗ Branch 56 not taken.
✗ Branch 57 not taken.
✓ Branch 58 taken 222801 times.
✗ Branch 59 not taken.
✓ Branch 60 taken 222801 times.
✗ Branch 61 not taken.
✓ Branch 62 taken 222801 times.
✗ Branch 63 not taken.
✓ Branch 64 taken 222801 times.
✗ Branch 65 not taken.
✗ Branch 66 not taken.
✗ Branch 67 not taken.
✗ Branch 68 not taken.
✗ Branch 69 not taken.
✗ Branch 70 not taken.
✗ Branch 71 not taken.
✓ Branch 72 taken 19962 times.
✗ Branch 73 not taken.
✓ Branch 74 taken 246918 times.
✗ Branch 75 not taken.
✓ Branch 76 taken 8439 times.
✗ Branch 77 not taken.
✓ Branch 78 taken 8439 times.
✗ Branch 79 not taken.
✓ Branch 80 taken 8439 times.
✗ Branch 81 not taken.
✓ Branch 82 taken 8439 times.
✗ Branch 83 not taken.
✗ Branch 84 not taken.
✗ Branch 85 not taken.
✗ Branch 86 not taken.
✗ Branch 87 not taken.
✗ Branch 88 not taken.
✗ Branch 89 not taken.
✗ Branch 90 not taken.
✗ Branch 91 not taken.
✓ Branch 92 taken 81 times.
✗ Branch 93 not taken.
✓ Branch 94 taken 81 times.
✗ Branch 95 not taken.
✓ Branch 96 taken 81 times.
✗ Branch 97 not taken.
✓ Branch 98 taken 81 times.
✗ Branch 99 not taken.
✓ Branch 100 taken 81 times.
✗ Branch 101 not taken.
✓ Branch 102 taken 81 times.
✗ Branch 103 not taken.
✓ Branch 104 taken 36144 times.
✗ Branch 105 not taken.
✓ Branch 106 taken 36144 times.
✗ Branch 107 not taken.
✓ Branch 108 taken 36144 times.
✗ Branch 109 not taken.
✓ Branch 110 taken 40405 times.
✗ Branch 111 not taken.
✓ Branch 112 taken 40405 times.
✗ Branch 113 not taken.
✓ Branch 114 taken 36063 times.
✗ Branch 115 not taken.
✓ Branch 116 taken 36063 times.
✓ Branch 117 taken 4319 times.
✓ Branch 118 taken 36063 times.
✓ Branch 119 taken 4319 times.
✓ Branch 120 taken 36063 times.
✓ Branch 121 taken 4319 times.
✓ Branch 122 taken 40382 times.
✗ Branch 123 not taken.
✓ Branch 124 taken 328399 times.
✗ Branch 125 not taken.
✓ Branch 126 taken 328399 times.
✗ Branch 127 not taken.
✓ Branch 128 taken 328399 times.
✗ Branch 129 not taken.
✓ Branch 130 taken 292336 times.
✗ Branch 131 not taken.
✓ Branch 132 taken 292336 times.
✗ Branch 133 not taken.
✓ Branch 134 taken 371631 times.
✗ Branch 135 not taken.
✓ Branch 136 taken 79295 times.
✗ Branch 137 not taken.
✓ Branch 138 taken 79295 times.
✗ Branch 139 not taken.
✓ Branch 140 taken 79295 times.
✗ Branch 141 not taken.
✓ Branch 142 taken 79295 times.
✗ Branch 143 not taken.
✓ Branch 144 taken 79295 times.
✗ Branch 145 not taken.
✓ Branch 146 taken 79295 times.
✗ Branch 147 not taken.
✓ Branch 148 taken 150982 times.
✗ Branch 149 not taken.
✓ Branch 150 taken 150982 times.
✗ Branch 151 not taken.
✗ Branch 152 not taken.
✗ Branch 153 not taken.
✓ Branch 154 taken 67558 times.
✗ Branch 155 not taken.
✓ Branch 156 taken 80293 times.
✗ Branch 157 not taken.
✓ Branch 158 taken 80293 times.
✗ Branch 159 not taken.
✓ Branch 160 taken 12735 times.
✗ Branch 161 not taken.
✓ Branch 162 taken 12735 times.
✗ Branch 163 not taken.
✓ Branch 164 taken 18419 times.
✗ Branch 165 not taken.
✓ Branch 166 taken 18419 times.
✗ Branch 167 not taken.
✓ Branch 168 taken 18419 times.
✗ Branch 169 not taken.
✗ Branch 170 not taken.
✓ Branch 171 taken 5684 times.
✗ Branch 172 not taken.
✓ Branch 173 taken 5684 times.
✗ Branch 174 not taken.
✓ Branch 175 taken 1365 times.
✗ Branch 176 not taken.
✓ Branch 177 taken 1365 times.
✗ Branch 178 not taken.
✓ Branch 179 taken 1365 times.
✗ Branch 180 not taken.
✓ Branch 181 taken 1365 times.
✗ Branch 182 not taken.
✓ Branch 183 taken 1365 times.
✗ Branch 184 not taken.
✓ Branch 185 taken 1365 times.
✗ Branch 186 not taken.
✗ Branch 187 not taken.
✗ Branch 188 not taken.
✗ Branch 189 not taken.
✗ Branch 190 not taken.
✗ Branch 191 not taken.
✗ Branch 192 not taken.
✗ Branch 193 not taken.
✗ Branch 194 not taken.
✗ Branch 195 not taken.
✓ Branch 196 taken 465729 times.
✗ Branch 197 not taken.
✓ Branch 198 taken 465729 times.
✗ Branch 199 not taken.
✗ Branch 200 not taken.
✓ Branch 201 taken 1169344 times.
13151113 if(rep_ != nullptr) {
669 expansion::delete_expansion_on_heap(rep_);
670 15779 rep_ = nullptr;
671 }
672 }
673
674 private:
675 expansion* rep_;
676 friend expansion_nt operator- (double a, const expansion_nt& b);
677
678 friend expansion_nt expansion_nt_sq_dist(
679 const double* a, const double* b, coord_index_t dim
680 );
681
682 friend expansion_nt expansion_nt_dot_at(
683 const double* a, const double* b, const double* c,
684 coord_index_t dim
685 );
686 // friend class rational_nt;
687 };
688
689 /**
690 * \brief Computes the sum of a double and an expansion_nt
691 * \param[in] a the double to be added
692 * \param[in] b the expansion_nt to be added
693 * \return an expansion_nt that represents \p a + \p b
694 * \relates expansion_nt
695 */
696 inline expansion_nt operator+ (double a, const expansion_nt& b) {
697 return b + a;
698 }
699
700 /**
701 * \brief Computes the difference between a double and an expansion_nt
702 * \param[in] a the double
703 * \param[in] b the expansion_nt to be subtracted
704 * \return an expansion_nt that represents \p a - \p b
705 * \relates expansion_nt
706 */
707 inline expansion_nt operator- (double a, const expansion_nt& b) {
708 expansion_nt result = b - a;
709 result.rep().negate();
710 return result;
711 }
712
713 /**
714 * \brief Computes the product of a double and an expansion_nt
715 * \param[in] a the double
716 * \param[in] b the expansion_nt to be multiplied
717 * \return an expansion_nt that represents \p a * \p b
718 * \relates expansion_nt
719 */
720 inline expansion_nt operator* (double a, const expansion_nt& b) {
721 return b * a;
722 }
723
724 /**
725 * \brief Tests equality between two expansion_nt%s.
726 * \details Implemented by testing whether the difference between
727 * \p a and \p b is 0.
728 * \return true if \p a and \p b represent exactly the same value, false
729 * otherwise
730 * \relates expansion_nt
731 */
732 inline bool operator== (const expansion_nt& a, const expansion_nt& b) {
733 return a.rep().equals(b.rep());
734 }
735
736 /**
737 * \brief Tests equality between an expansion_nt and a double.
738 * \details Implemented by testing whether the difference between
739 * \p a and \p b is 0.
740 * \return true if \p a and \p b represent exactly the same value, false
741 * otherwise
742 * \relates expansion_nt
743 */
744 inline bool operator== (const expansion_nt& a, double b) {
745 return a.rep().equals(b);
746 }
747
748 /**
749 * \brief Tests equality between a double and an expansion_nt.
750 * \details Implemented by testing whether the difference between
751 * \p a and \p b is 0.
752 * \return true if \p a and \p b represent exactly the same value, false
753 * otherwise
754 * \relates expansion_nt
755 */
756 inline bool operator== (double a, const expansion_nt& b) {
757 return b.rep().equals(a);
758 }
759
760 /**
761 * \brief Tests whether two expansion_nt%s differ.
762 * \details Implemented by testing whether the difference between
763 * \p a and \p b is different from 0.
764 * \return true if \p a and \p b do not represent the same exact value,
765 * false otherwise
766 * \relates expansion_nt
767 */
768 inline bool operator!= (const expansion_nt& a, const expansion_nt& b) {
769 return !a.rep().equals(b.rep());
770 }
771
772 /**
773 * \brief Tests whether an expansion_nt differs from a double.
774 * \details Implemented by testing whether the difference between
775 * \p a and \p b is different from 0.
776 * \return true if \p a and \p b do not represent the same exact value,
777 * false otherwise
778 * \relates expansion_nt
779 */
780 inline bool operator!= (const expansion_nt& a, double b) {
781 return !a.rep().equals(b);
782 }
783
784 /**
785 * \brief Tests whether a double differs from an expansion_nt.
786 * \details Implemented by testing whether the difference between
787 * \p a and \p b is different from 0.
788 * \return true if \p a and \p b do not represent the same exact value,
789 * false otherwise
790 * \relates expansion_nt
791 */
792 inline bool operator!= (double a, const expansion_nt& b) {
793 return !b.rep().equals(a);
794 }
795
796 /**
797 * \brief Computes an expansion that represents the square distance between
798 * two points.
799 * \param[in] a an array of \p dim doubles
800 * \param[in] b an array of \p dim doubles
801 * \param[in] dim the dimension of the points
802 * \return an expansion_nt that represent the exact squared
803 * distance between \p a and \p b
804 * \relates expansion_nt
805 */
806 inline expansion_nt expansion_nt_sq_dist(
807 const double* a, const double* b, coord_index_t dim
808 ) {
809 expansion* result = expansion::new_expansion_on_heap(
810 expansion::sq_dist_capacity(dim)
811 );
812 result->assign_sq_dist(a, b, dim);
813 return expansion_nt(result);
814 }
815
816 /**
817 * \brief Computes an expansion that represents the dot product of
818 * two vectors determined by three points.
819 * \param[in] a an array of \p dim doubles
820 * \param[in] b an array of \p dim doubles
821 * \param[in] c an array of \p dim doubles
822 * \param[in] dim the dimension of the points
823 * \return an expansion_nt that represents the exact
824 * dot product (\p a - \p c).(\p b - \p c)
825 * \relates expansion_nt
826 */
827 inline expansion_nt expansion_nt_dot_at(
828 const double* a, const double* b, const double* c, coord_index_t dim
829 ) {
830 expansion* result = expansion::new_expansion_on_heap(
831 expansion::dot_at_capacity(dim)
832 );
833 result->assign_dot_at(a, b, c, dim);
834 return expansion_nt(result);
835 }
836
837 /************************************************************************/
838
839 /**
840 * \brief Specialization of geo_sgn() for expansion_nt.
841 * \param x a const reference to an expansion_nt
842 * \return the (exact) sign of x (one of POSITIVE, ZERO, NEGATIVE)
843 */
844 template <> inline Sign geo_sgn(const expansion_nt& x) {
845 return x.sign();
846 }
847
848 /**
849 * \brief Specialization of geo_cmp() for expansion_nt.
850 * \param x , y const references to two expansion_nt
851 * \retval POSITIVE if x > y
852 * \retval ZERO if x == y
853 * \retval NEGATIVE if x < y
854 */
855 template <> inline Sign geo_cmp(
856 const expansion_nt& x, const expansion_nt& y
857 ) {
858 return x.compare(y);
859 }
860
861 /************************************************************************/
862
863 /**
864 * \brief Tests whether an expansion_nt is zero.
865 * \details Optimized using the low-level API
866 * \param[in] x a const reference to the expansion_nt to be tested
867 * \retval true if \p x is equal to zero
868 * \retval false otherwise
869 */
870 inline bool expansion_nt_is_zero(const expansion_nt& x) {
871 return (x.sign() == GEO::ZERO);
872 }
873
874 /**
875 * \brief Tests whether an expansion_nt is equal to one.
876 * \details Optimized using the low-level API
877 * \param[in] x a const reference to the expansion_nt to be tested
878 * \retval true if \p x is equal to one
879 * \retval false otherwise
880 */
881 inline bool expansion_nt_is_one(const expansion_nt& x) {
882 return x.rep().equals(1.0);
883 }
884
885
886 /**
887 * \brief Compares two expansion_nt
888 * \details Optimized using the low-level API
889 * \param [in] x , y the two expansion_nt to compare
890 * \retval POSITIVE if \p x is greater than \p y
891 * \retval ZERO if \p x equals \p y
892 * \retval NEGATIVE if \p x is smaller than \p y
893 */
894 inline Sign expansion_nt_compare(
895 const expansion_nt& x, const expansion_nt& y
896 ) {
897 const expansion& diff = expansion_diff(x.rep(), y.rep());
898 return diff.sign();
899 }
900
901 /**
902 * \brief Computes the square of an expansion_nt
903 * \details Optimized using the low-level API
904 * \param[in] x the expansion_nt to be squared
905 * \return \p x * \p x
906 */
907 inline expansion_nt expansion_nt_square(const expansion_nt& x) {
908 expansion_nt result(
909 expansion::new_expansion_on_heap(
910 expansion::square_capacity(x.rep()
911 ))
912 );
913 result.rep().assign_square(x.rep());
914 return result;
915 }
916
917
918 /**
919 * \brief Computes a 2x2 determinant
920 * \details Specialization using the low-evel API for expansions.
921 * This gains some performance as compared to using CGAL's
922 * determinant template with expansion_nt.
923 */
924 expansion_nt GEOGRAM_API expansion_nt_determinant(
925 const expansion_nt& a00,const expansion_nt& a01,
926 const expansion_nt& a10,const expansion_nt& a11
927 );
928
929 /**
930 * \brief Computes a 3x3 determinant
931 * \details Specialization using the low-evel API for expansions.
932 * This gains some performance as compared to using CGAL's determinant
933 * template with expansion_nt.
934 */
935 expansion_nt GEOGRAM_API expansion_nt_determinant(
936 const expansion_nt& a00,const expansion_nt& a01,const expansion_nt& a02,
937 const expansion_nt& a10,const expansion_nt& a11,const expansion_nt& a12,
938 const expansion_nt& a20,const expansion_nt& a21,const expansion_nt& a22
939 );
940
941 /**
942 * \brief Computes a 4x4 determinant
943 * \details Specialization using the low-evel API for expansions.
944 * This gains some performance as compared to using CGAL's determinant
945 * template with expansion_nt.
946 */
947 expansion_nt GEOGRAM_API expansion_nt_determinant(
948 const expansion_nt& a00,const expansion_nt& a01,
949 const expansion_nt& a02,const expansion_nt& a03,
950 const expansion_nt& a10,const expansion_nt& a11,
951 const expansion_nt& a12,const expansion_nt& a13,
952 const expansion_nt& a20,const expansion_nt& a21,
953 const expansion_nt& a22,const expansion_nt& a23,
954 const expansion_nt& a30,const expansion_nt& a31,
955 const expansion_nt& a32,const expansion_nt& a33
956 );
957
958 // Make things a bit faster if target OS has large stack size
959 #ifdef GEO_HAS_BIG_STACK
960
961 /**
962 * \brief Specialization of det2x2
963 * \details Calls the optimized implementation for expansion_nt
964 */
965
966 template <> inline expansion_nt det2x2(
967 const expansion_nt& a11, const expansion_nt& a12,
968 const expansion_nt& a21, const expansion_nt& a22
969 ) {
970 return expansion_nt_determinant(
971 a11,a12,
972 a21,a22
973
17/36
✓ Branch 1 taken 32487 times.
✓ Branch 2 taken 31980 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 32487 times.
✓ Branch 5 taken 31980 times.
✗ Branch 6 not taken.
✓ Branch 7 taken 89084 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 89084 times.
✓ Branch 11 taken 529354 times.
✗ Branch 12 not taken.
✓ Branch 13 taken 89003 times.
✓ Branch 14 taken 529354 times.
✗ Branch 15 not taken.
✓ Branch 16 taken 246918 times.
✗ Branch 17 not taken.
✓ Branch 19 taken 19962 times.
✗ Branch 20 not taken.
✓ Branch 22 taken 11583 times.
✗ Branch 23 not taken.
✗ Branch 25 not taken.
✗ Branch 26 not taken.
✗ Branch 28 not taken.
✗ Branch 29 not taken.
✗ Branch 31 not taken.
✗ Branch 32 not taken.
✓ Branch 34 taken 36063 times.
✗ Branch 35 not taken.
✓ Branch 37 taken 36063 times.
✗ Branch 38 not taken.
✓ Branch 40 taken 36063 times.
✗ Branch 41 not taken.
✓ Branch 44 taken 12735 times.
✗ Branch 45 not taken.
✓ Branch 47 taken 12735 times.
✗ Branch 48 not taken.
1010247 );
974 }
975
976 /**
977 * \brief Specialization of det3x3
978 * \details Calls the optimized implementation for expansion_nt
979 */
980
981 template <> inline expansion_nt det3x3(
982 const expansion_nt& a11, const expansion_nt& a12,
983 const expansion_nt& a13,
984 const expansion_nt& a21, const expansion_nt& a22,
985 const expansion_nt& a23,
986 const expansion_nt& a31, const expansion_nt& a32,
987 const expansion_nt& a33
988 ) {
989 return expansion_nt_determinant(
990 a11,a12,a13,
991 a21,a22,a23,
992 a31,a32,a33
993
4/8
✓ Branch 1 taken 13179 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 1365 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 1365 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 1365 times.
✗ Branch 11 not taken.
17274 );
994 }
995
996 /**
997 * \brief Specialization of det4x4
998 * \details Calls the optimized implementation for expansion_nt
999 */
1000
1001 template <> inline expansion_nt det4x4(
1002 const expansion_nt& a11, const expansion_nt& a12,
1003 const expansion_nt& a13, const expansion_nt& a14,
1004 const expansion_nt& a21, const expansion_nt& a22,
1005 const expansion_nt& a23, const expansion_nt& a24,
1006 const expansion_nt& a31, const expansion_nt& a32,
1007 const expansion_nt& a33, const expansion_nt& a34,
1008 const expansion_nt& a41, const expansion_nt& a42,
1009 const expansion_nt& a43, const expansion_nt& a44
1010 ) {
1011 return expansion_nt_determinant(
1012 a11,a12,a13,a14,
1013 a21,a22,a23,a24,
1014 a31,a32,a33,a34,
1015 a41,a42,a43,a44
1016 );
1017 }
1018
1019 #endif
1020
1021 /************************************************************************/
1022 }
1023
1024 /**
1025 * \brief Displays the approximated value of an expansion_nt to a stream.
1026 * \param[out] os the stream
1027 * \param[in] a the expansion_nt to be sent to the stream
1028 * \return a reference to the stream
1029 */
1030 inline std::ostream& operator<< (
1031 std::ostream& os, const GEO::expansion_nt& a
1032 ) {
1033 return os << a.estimate();
1034 }
1035
1036 /**
1037 * \brief Reads a double precision number from a stream and converts it to
1038 * an approximation.
1039 * \param[in] is the stream
1040 * \param[out] a the read expansion_nt
1041 * \return a reference to the stream
1042 */
1043 inline std::istream& operator>> ( std::istream& is, GEO::expansion_nt& a) {
1044 double d;
1045 is >> d;
1046 if (is) {
1047 a = GEO::expansion_nt(d);
1048 }
1049 return is;
1050 }
1051
1052 /*****************************************************************************/
1053
1054 namespace GEO {
1055
1056 /**************************************************************************/
1057
1058 namespace Numeric {
1059
1060 template<> inline void optimize_number_representation(expansion_nt& x) {
1061 x.optimize();
1062 }
1063
1064 /**
1065 * \brief Compares two rational numbers given as separate
1066 * numerators and denominators. Specialization for exact_nt.
1067 * \param[in] a_num , a_denom defines a = \p a_num / \p a_denom
1068 * \param[in] b_num , b_denom defines b = \p b_num / \p b_denom
1069 * \return the sign of a - b
1070 */
1071 template<> Sign GEOGRAM_API ratio_compare(
1072 const expansion_nt& a_num, const expansion_nt& a_denom,
1073 const expansion_nt& b_num, const expansion_nt& b_denom
1074 );
1075 }
1076
1077 /**************************************************************************/
1078
1079 typedef rationalg<expansion_nt> rational_nt;
1080
1081 /**************************************************************************/
1082
1083 /** \brief Specialization of GEO::is_scalar */
1084 template <> struct is_scalar<expansion_nt> {
1085 typedef expansion_nt type;
1086 static constexpr bool value = true;
1087 };
1088
1089 /**************************************************************************/
1090 }
1091
1092 #endif
1093