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