LORENE
cmp_arithm.C
1/*
2 * Arithmetical operations for class Cmp
3 *
4 */
5
6/*
7 * Copyright (c) 1999-2000 Jean-Alain Marck
8 * Copyright (c) 1999-2001 Eric Gourgoulhon
9 * Copyright (c) 1999-2001 Philippe Grandclement
10 *
11 * This file is part of LORENE.
12 *
13 * LORENE is free software; you can redistribute it and/or modify
14 * it under the terms of the GNU General Public License as published by
15 * the Free Software Foundation; either version 2 of the License, or
16 * (at your option) any later version.
17 *
18 * LORENE is distributed in the hope that it will be useful,
19 * but WITHOUT ANY WARRANTY; without even the implied warranty of
20 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
21 * GNU General Public License for more details.
22 *
23 * You should have received a copy of the GNU General Public License
24 * along with LORENE; if not, write to the Free Software
25 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
26 *
27 */
28
29
30char cmp_arithm_C[] = "$Header: /cvsroot/Lorene/C++/Source/Cmp/cmp_arithm.C,v 1.5 2014/10/13 08:52:47 j_novak Exp $" ;
31
32/*
33 * $Id: cmp_arithm.C,v 1.5 2014/10/13 08:52:47 j_novak Exp $
34 * $Log: cmp_arithm.C,v $
35 * Revision 1.5 2014/10/13 08:52:47 j_novak
36 * Lorene classes and functions now belong to the namespace Lorene.
37 *
38 * Revision 1.4 2014/10/06 15:13:03 j_novak
39 * Modified #include directives to use c++ syntax.
40 *
41 * Revision 1.3 2003/06/20 13:59:59 f_limousin
42 * L'assert pour le mapping est realise a partir du mapping lui meme et non a partir du pointeur sur le mapping.
43 *
44 * Revision 1.2 2002/10/16 14:36:34 j_novak
45 * Reorganization of #include instructions of standard C++, in order to
46 * use experimental version 3 of gcc.
47 *
48 * Revision 1.1.1.1 2001/11/20 15:19:27 e_gourgoulhon
49 * LORENE
50 *
51 * Revision 2.11 2001/05/26 15:07:47 eric
52 * Ajout de operator% : multiplication de deux Cmp avec desaliasage
53 *
54 * Revision 2.10 2000/02/15 13:39:57 phil
55 * correction -= et +=
56 *
57 * Revision 2.9 2000/01/28 16:34:57 eric
58 * Utilisation des nouvelles fonctions Cmp::check_dzpuis et Cmp::dz_nonzero
59 * dans les tests sur dzpuis.
60 *
61 * Revision 2.8 1999/12/10 16:32:52 eric
62 * Dans l'arithmetique membre (+=, -=, *=), on n'appelle desormais
63 * del_deriv() que tout a la fin.
64 *
65 * Revision 2.7 1999/11/26 14:23:38 eric
66 * Ajout du membre dzpuis.
67 *
68 * Revision 2.6 1999/11/12 17:08:35 eric
69 * Ajout de la partie manquante de l'arithmetique.
70 *
71 * Revision 2.5 1999/10/28 13:23:43 phil
72 * Deverouillage des ETATNONDEF
73 *
74 * Revision 2.4 1999/10/27 08:45:21 eric
75 * ntroduction du membre Valeur va.
76 *
77 * Revision 2.3 1999/10/22 08:14:58 eric
78 * La fonction annule() est rebaptisee annule_hard().
79 *
80 * Revision 2.2 1999/09/14 17:19:44 phil
81 * ajout de Cmp operator* (double, const Cmp&)
82 *
83 * Revision 2.1 1999/03/03 11:13:56 hyc
84 * *** empty log message ***
85 *
86 *
87 * $Header: /cvsroot/Lorene/C++/Source/Cmp/cmp_arithm.C,v 1.5 2014/10/13 08:52:47 j_novak Exp $
88 *
89 */
90
91// headers C
92#include <cassert>
93#include <cstdlib>
94
95// headers Lorene
96#include "cmp.h"
97#include "type_parite.h"
98
99 //********************//
100 // OPERATEURS UNAIRES //
101 //********************//
102
103namespace Lorene {
104Cmp operator+(const Cmp & ci) {
105 return ci ;
106}
107
108Cmp operator-(const Cmp & ci) {
109
110 // Cas particulier
111 if ((ci.get_etat() == ETATZERO) || (ci.get_etat() == ETATNONDEF)) {
112 return ci ;
113 }
114
115 // Cas general
116 assert(ci.get_etat() == ETATQCQ) ; // sinon...
117 Cmp r(ci.get_mp()) ; // Cmp resultat
118 r.set_etat_qcq() ;
119 r.va = - ci.va ;
120 r.set_dzpuis( ci.get_dzpuis() ) ;
121
122 // Termine
123 return r ;
124}
125
126 //**********//
127 // ADDITION //
128 //**********//
129// Cmp + Cmp
130// ---------
131Cmp operator+(const Cmp & c1, const Cmp & c2) {
132
133 if (c1.get_etat() == ETATNONDEF)
134 return c1 ;
135 if (c2.get_etat() == ETATNONDEF)
136 return c2 ;
137 assert(c1.get_mp() == c2.get_mp()) ;
138
139 // Cas particuliers
140 if (c1.get_etat() == ETATZERO) {
141 return c2 ;
142 }
143 if (c2.get_etat() == ETATZERO) {
144 return c1 ;
145 }
146 assert(c1.get_etat() == ETATQCQ) ; // sinon...
147 assert(c2.get_etat() == ETATQCQ) ; // sinon...
148
149 // Cas general
150
151 if ( c1.dz_nonzero() && c2.dz_nonzero() ) {
152 if ( c1.get_dzpuis() != c2.get_dzpuis() ) {
153 cout << "Operation Cmp + Cmp forbidden in the external " << endl;
154 cout << " compactified domain ! " << endl ;
155 abort() ;
156 }
157 }
158
159 Cmp r(c1) ; // Le resultat
160 r.va += c2.va ;
161
162 if (c1.dz_nonzero()) {
163 r.set_dzpuis( c1.get_dzpuis() ) ;
164 }
165 else{
166 r.set_dzpuis( c2.get_dzpuis() ) ;
167 }
168
169 // Termine
170 return r ;
171}
172
173// Cmp + double
174// ------------
175Cmp operator+(const Cmp& t1, double x)
176{
177 // Protections
178 assert(t1.get_etat() != ETATNONDEF) ;
179
180 // Cas particuliers
181 if (x == double(0)) {
182 return t1 ;
183 }
184
185 assert( t1.check_dzpuis(0) ) ;
186
187 Cmp resu(t1) ;
188
189 if (t1.get_etat() == ETATZERO) {
190 resu = x ;
191 }
192 else{
193 assert(resu.get_etat() == ETATQCQ) ; // sinon ...
194 resu.va = resu.va + x ;
195 }
196
197 resu.set_dzpuis(0) ;
198
199 return resu ;
200}
201
202// double + Cmp
203// ------------
204Cmp operator+(double x, const Cmp& t1)
205{
206 return t1 + x ;
207}
208
209// Cmp + int
210// ---------
211Cmp operator+(const Cmp& t1, int m)
212{
213 return t1 + double(m) ;
214}
215
216// int + Cmp
217// ---------
218Cmp operator+(int m, const Cmp& t1)
219{
220 return t1 + double(m) ;
221}
222
223
224
225
226
227 //**************//
228 // SOUSTRACTION //
229 //**************//
230
231// Cmp - Cmp
232// ---------
233Cmp operator-(const Cmp & c1, const Cmp & c2) {
234
235 if (c1.get_etat() == ETATNONDEF)
236 return c1 ;
237 if (c2.get_etat() == ETATNONDEF)
238 return c2 ;
239
240 assert(c1.get_mp() == c2.get_mp()) ;
241
242 // Cas particuliers
243 if (c1.get_etat() == ETATZERO) {
244 return -c2 ;
245 }
246 if (c2.get_etat() == ETATZERO) {
247 return c1 ;
248 }
249 assert(c1.get_etat() == ETATQCQ) ; // sinon...
250 assert(c2.get_etat() == ETATQCQ) ; // sinon...
251
252 // Cas general
253 if ( c1.dz_nonzero() && c2.dz_nonzero() ) {
254 if ( c1.get_dzpuis() != c2.get_dzpuis() ) {
255 cout << "Operation Cmp - Cmp forbidden in the external " << endl;
256 cout << " compactified domain ! " << endl ;
257 abort() ;
258 }
259 }
260
261 Cmp r(c1) ; // Le resultat
262 r.va -= c2.va ;
263
264 if (c1.dz_nonzero()) {
265 r.set_dzpuis( c1.get_dzpuis() ) ;
266 }
267 else{
268 r.set_dzpuis( c2.get_dzpuis() ) ;
269 }
270
271 // Termine
272 return r ;
273}
274
275// Cmp - double
276// ------------
277Cmp operator-(const Cmp& t1, double x)
278{
279 // Protections
280 assert(t1.get_etat() != ETATNONDEF) ;
281
282 // Cas particuliers
283 if (x == double(0)) {
284 return t1 ;
285 }
286
287 assert( t1.check_dzpuis(0) ) ;
288
289 Cmp resu(t1) ;
290
291 if (t1.get_etat() == ETATZERO) {
292 resu = - x ;
293 }
294 else{
295 assert(resu.get_etat() == ETATQCQ) ; // sinon ...
296 resu.va = resu.va - x ;
297 }
298
299 resu.set_dzpuis(0) ;
300
301 return resu ;
302}
303
304// double - Cmp
305// ------------
306Cmp operator-(double x, const Cmp& t1)
307{
308 return - (t1 - x) ;
309}
310
311// Cmp - int
312// ---------
313Cmp operator-(const Cmp& t1, int m)
314{
315 return t1 - double(m) ;
316}
317
318// int - Cmp
319// ---------
320Cmp operator-(int m, const Cmp& t1)
321{
322 return double(m) - t1 ;
323}
324
325
326
327
328
329
330 //****************//
331 // MULTIPLICATION //
332 //****************//
333
334// Cmp * Cmp
335// ---------
336Cmp operator*(const Cmp& c1, const Cmp& c2) {
337
338
339 // Cas particuliers
340 if ((c1.get_etat() == ETATZERO) || (c1.get_etat() == ETATNONDEF)){
341 return c1 ;
342 }
343 if ((c2.get_etat() == ETATZERO)|| (c2.get_etat() == ETATNONDEF)) {
344 return c2 ;
345 }
346 assert(c1.get_etat() == ETATQCQ) ; // sinon...
347 assert(c2.get_etat() == ETATQCQ) ; // sinon...
348
349 // Protection
350 assert(*(c1.get_mp()) == *(c2.get_mp())) ;
351
352 // Cas general
353 Cmp r(c1) ; // Le resultat
354 r.va *= c2.va ;
355
356 r.set_dzpuis( c1.get_dzpuis() + c2.get_dzpuis() ) ;
357
358 // Termine
359 return r ;
360}
361
362// Cmp % Cmp (multiplication with desaliasing)
363// -------------------------------------------
364Cmp operator%(const Cmp& c1, const Cmp& c2) {
365
366
367 // Cas particuliers
368 if ((c1.get_etat() == ETATZERO) || (c1.get_etat() == ETATNONDEF)){
369 return c1 ;
370 }
371 if ((c2.get_etat() == ETATZERO)|| (c2.get_etat() == ETATNONDEF)) {
372 return c2 ;
373 }
374 assert(c1.get_etat() == ETATQCQ) ; // sinon...
375 assert(c2.get_etat() == ETATQCQ) ; // sinon...
376
377 // Protection
378 assert(c1.get_mp() == c2.get_mp()) ;
379
380 // Cas general
381 Cmp r( c1.get_mp() ) ; // Le resultat
382 r.set_etat_qcq() ;
383 r.va = c1.va % c2.va ;
384
385 r.set_dzpuis( c1.get_dzpuis() + c2.get_dzpuis() ) ;
386
387 // Termine
388 return r ;
389}
390
391
392
393
394// double * Cmp
395// ------------
396Cmp operator*(double a, const Cmp& c1) {
397
398 // Cas particuliers
399 if ((c1.get_etat() == ETATZERO) || (c1.get_etat() == ETATNONDEF)) {
400 return c1 ;
401 }
402
403 assert(c1.get_etat() == ETATQCQ) ; // sinon...
404
405 // Cas general
406 Cmp r(c1.get_mp()) ;
407 r.set_dzpuis( c1.get_dzpuis() ) ;
408
409 if ( a == double(0) ) {
410 r.set_etat_zero() ;
411 }
412 else {
413 r.set_etat_qcq() ;
414 r.va = a * c1.va ;
415 }
416
417
418 // Termine
419 return r ;
420}
421
422
423// Cmp * double
424// ------------
425Cmp operator*(const Cmp& t1, double x)
426{
427 return x * t1 ;
428}
429
430// Cmp * int
431// ---------
432Cmp operator*(const Cmp& t1, int m)
433{
434 return t1 * double(m) ;
435}
436
437// int * Cmp
438// ---------
439Cmp operator*(int m, const Cmp& t1)
440{
441 return double(m) * t1 ;
442}
443
444
445
446
447
448
449
450 //**********//
451 // DIVISION //
452 //**********//
453
454
455// Cmp / Cmp
456// ---------
457Cmp operator/(const Cmp& c1, const Cmp& c2) {
458
459 // Protections
460 assert(c1.get_etat() != ETATNONDEF) ;
461 assert(c2.get_etat() != ETATNONDEF) ;
462 assert(c1.get_mp() == c2.get_mp()) ;
463
464 // Cas particuliers
465 if (c2.get_etat() == ETATZERO) {
466 cout << "Division by 0 in Cmp / Cmp !" << endl ;
467 abort() ;
468 }
469 if (c1.get_etat() == ETATZERO) {
470 return c1 ;
471 }
472
473 // Cas general
474
475 assert(c1.get_etat() == ETATQCQ) ; // sinon...
476 assert(c2.get_etat() == ETATQCQ) ; // sinon...
477
478 Cmp r(c1.get_mp()) ; // Le resultat
479
480 r.set_etat_qcq() ;
481 r.va = c1.va / c2.va ;
482
483 r.set_dzpuis( c1.get_dzpuis() - c2.get_dzpuis() ) ;
484
485 // Termine
486 return r ;
487}
488
489// Cmp / double
490// -------------
491Cmp operator/(const Cmp& c1, double x) {
492
493 if (c1.get_etat() == ETATNONDEF)
494 return c1 ;
495
496 // Cas particuliers
497 if ( x == double(0) ) {
498 cout << "Division by 0 in Cmp / double !" << endl ;
499 abort() ;
500 }
501 if (c1.get_etat() == ETATZERO) {
502 return c1 ;
503 }
504
505 assert(c1.get_etat() == ETATQCQ) ; // sinon...
506
507 Cmp r(c1.get_mp()) ; // Le resultat
508
509 r.set_etat_qcq() ;
510 r.va = c1.va / x ;
511
512 r.set_dzpuis( c1.get_dzpuis() ) ;
513
514 // Termine
515 return r ;
516}
517
518
519// double / Cmp
520// ------------
521Cmp operator/(double x, const Cmp& c2) {
522
523 if (c2.get_etat() == ETATNONDEF)
524 return c2 ;
525
526 if (c2.get_etat() == ETATZERO) {
527 cout << "Division by 0 in Cmp / Cmp !" << endl ;
528 abort() ;
529 }
530
531
532 assert(c2.get_etat() == ETATQCQ) ; // sinon...
533
534 Cmp r(c2.get_mp()) ; // Le resultat
535 r.set_dzpuis( - c2.get_dzpuis() ) ;
536
537 if ( x == double(0) ) {
538 r.set_etat_zero() ;
539 }
540 else {
541 r.set_etat_qcq() ;
542 r.va = x / c2.va ;
543 }
544
545 // Termine
546 return r ;
547}
548
549
550// Cmp / int
551// ---------
552Cmp operator/(const Cmp& c1, int m) {
553
554 return c1 / double(m) ;
555
556}
557
558
559// int / Cmp
560// ---------
561Cmp operator/(int m, const Cmp& c2) {
562
563 return double(m) / c2 ;
564
565}
566
567 //*******************//
568 // operateurs +=,... //
569 //*******************//
570
571//---------
572// += Cmp
573//---------
574
575void Cmp::operator+=(const Cmp & ci) {
576
577 // Protection
578 assert(mp == ci.get_mp()) ; // meme mapping
579 if (etat == ETATNONDEF)
580 return ;
581
582 // Cas particulier
583 if (ci.get_etat() == ETATZERO) {
584 return ;
585 }
586
587 if (ci.get_etat() == ETATNONDEF) {
589 return ;
590 }
591
592 // Cas general
593
594
595 if ( dz_nonzero() && ci.dz_nonzero() ) {
596 if ( dzpuis != ci.dzpuis ) {
597 cout << "Operation += Cmp forbidden in the external " << endl;
598 cout << " compactified domain ! " << endl ;
599 abort() ;
600 }
601 }
602
603 if (etat == ETATZERO) {
604 (*this) = ci ;
605 }
606 else {
607 va += ci.va ;
608
609 if( ci.dz_nonzero() ) {
610 set_dzpuis(ci.dzpuis) ;
611 }
612 }
613 // Menage (a ne faire qu'a la fin seulement)
614 del_deriv() ;
615
616
617}
618
619//---------
620// -= Cmp
621//---------
622
623void Cmp::operator-=(const Cmp & ci) {
624
625 // Protection
626 assert(mp == ci.get_mp()) ; // meme mapping
627 if (etat == ETATNONDEF)
628 return ;
629
630 // Cas particulier
631 if (ci.get_etat() == ETATZERO) {
632 return ;
633 }
634
635 if (ci.get_etat() == ETATNONDEF) {
637 return ;
638 }
639
640 // Cas general
641 if ( dz_nonzero() && ci.dz_nonzero() ) {
642 if ( dzpuis != ci.dzpuis ) {
643 cout << "Operation -= Cmp forbidden in the external " << endl;
644 cout << " compactified domain ! " << endl ;
645 abort() ;
646 }
647 }
648
649
650 if (etat == ETATZERO) {
651 (*this) = -ci ;
652 }
653 else {
654 va -= ci.va ;
655
656 if( ci.dz_nonzero() ) {
657 set_dzpuis(ci.dzpuis) ;
658 }
659 }
660 // Menage (a ne faire qu'a la fin seulement)
661 del_deriv() ;
662}
663
664//---------
665// *= Cmp
666//---------
667
668void Cmp::operator*=(const Cmp & ci) {
669
670 // Protection
671 assert(mp == ci.get_mp()) ; // meme mapping
672 if (etat == ETATNONDEF)
673 return ;
674
675 // Cas particulier
676 if (ci.get_etat() == ETATZERO) {
677 set_etat_zero() ;
678 return ;
679 }
680
681 if (etat == ETATZERO) {
682 return ;
683 }
684
685 if (ci.get_etat() == ETATNONDEF) {
687 return ;
688 }
689
690 // Cas general
691
692 assert(etat == ETATQCQ) ; // sinon....
693
694 va *= ci.va ;
695
696 dzpuis += ci.dzpuis ;
697
698 // Menage (a ne faire qu'a la fin seulement)
699 del_deriv() ;
700
701}
702}
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Definition cmp.h:446
const Map * mp
Reference mapping.
Definition cmp.h:451
int dzpuis
Power of r by which the quantity represented by this must be divided in the external compactified z...
Definition cmp.h:461
int get_etat() const
Returns the logical state.
Definition cmp.h:899
Valeur va
The numerical value of the Cmp
Definition cmp.h:464
int etat
Logical state (ETATNONDEF , ETATQCQ or ETATZERO ).
Definition cmp.h:454
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition cmp.C:304
int get_dzpuis() const
Returns dzpuis.
Definition cmp.h:903
void operator+=(const Cmp &)
+= Cmp
Definition cmp_arithm.C:575
bool dz_nonzero() const
Returns true if the last domain is compactified and *this is not zero in this domain.
Definition cmp.C:660
void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition cmp.C:289
void operator*=(const Cmp &)
*= Cmp
Definition cmp_arithm.C:668
void set_dzpuis(int)
Set a value to dzpuis.
Definition cmp.C:654
bool check_dzpuis(int dzi) const
Returns false if the last domain is compactified and *this is not zero in this domain and dzpuis is n...
Definition cmp.C:715
const Map * get_mp() const
Returns the mapping.
Definition cmp.h:901
void del_deriv()
Logical destructor of the derivatives.
Definition cmp.C:265
void operator-=(const Cmp &)
-= Cmp
Definition cmp_arithm.C:623
void set_etat_nondef()
Sets the logical state to ETATNONDEF (undefined).
Definition cmp.C:297
Base_val operator*(const Base_val &, const Base_val &)
This operator is used when calling multiplication or division of Valeur .
Cmp operator-(const Cmp &)
- Cmp
Definition cmp_arithm.C:108
Cmp operator/(const Cmp &, const Cmp &)
Cmp / Cmp.
Definition cmp_arithm.C:457
Cmp operator%(const Cmp &, const Cmp &)
Cmp * Cmp with desaliasing.
Definition cmp_arithm.C:364
Cmp operator+(const Cmp &)
Definition cmp_arithm.C:104
Lorene prototypes.
Definition app_hor.h:64