LORENE
valeur.C
1 /*
2  * Methods of class Valeur
3  *
4  * (see file valeur.h for documentation)
5  *
6  */
7 
8 /*
9  * Copyright (c) 1999-2000 Jean-Alain Marck
10  * Copyright (c) 1999-2003 Eric Gourgoulhon
11  * Copyright (c) 1999-2001 Philippe Grandclement
12  *
13  * This file is part of LORENE.
14  *
15  * LORENE is free software; you can redistribute it and/or modify
16  * it under the terms of the GNU General Public License as published by
17  * the Free Software Foundation; either version 2 of the License, or
18  * (at your option) any later version.
19  *
20  * LORENE is distributed in the hope that it will be useful,
21  * but WITHOUT ANY WARRANTY; without even the implied warranty of
22  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
23  * GNU General Public License for more details.
24  *
25  * You should have received a copy of the GNU General Public License
26  * along with LORENE; if not, write to the Free Software
27  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
28  *
29  */
30 
31 
32 char valeur_C[] = "$Header: /cvsroot/Lorene/C++/Source/Valeur/valeur.C,v 1.8 2014/10/13 08:53:49 j_novak Exp $" ;
33 
34 /*
35  * $Id: valeur.C,v 1.8 2014/10/13 08:53:49 j_novak Exp $
36  * $Log: valeur.C,v $
37  * Revision 1.8 2014/10/13 08:53:49 j_novak
38  * Lorene classes and functions now belong to the namespace Lorene.
39  *
40  * Revision 1.7 2014/10/06 15:13:22 j_novak
41  * Modified #include directives to use c++ syntax.
42  *
43  * Revision 1.6 2005/10/25 08:56:40 p_grandclement
44  * addition of std_spectral_base in the case of odd functions near the origin
45  *
46  * Revision 1.5 2004/11/23 12:45:00 f_limousin
47  * Add function filtre_tp(int nn, int nz1, int nz2).
48  *
49  * Revision 1.4 2003/10/19 19:52:56 e_gourgoulhon
50  * Added new method display_coef.
51  *
52  * Revision 1.3 2002/10/16 14:37:15 j_novak
53  * Reorganization of #include instructions of standard C++, in order to
54  * use experimental version 3 of gcc.
55  *
56  * Revision 1.2 2001/12/04 21:27:54 e_gourgoulhon
57  *
58  * All writing/reading to a binary file are now performed according to
59  * the big endian convention, whatever the system is big endian or
60  * small endian, thanks to the functions fwrite_be and fread_be
61  *
62  * Revision 1.1.1.1 2001/11/20 15:19:27 e_gourgoulhon
63  * LORENE
64  *
65  * Revision 2.37 2000/09/11 13:53:07 eric
66  * Ajout des membres p_mult_cp et p_mult_sp.
67  *
68  * Revision 2.36 2000/09/08 10:07:13 eric
69  * Ajout des methodes set_base_r, etc...
70  *
71  * Revision 2.35 1999/12/29 13:11:23 eric
72  * Ajout de la fonction val_point_jk.
73  *
74  * Revision 2.34 1999/12/20 16:35:23 eric
75  * Ajout de la fonction set_base.
76  *
77  * Revision 2.33 1999/12/10 16:09:47 eric
78  * Annulation des membres derives dans la fonction annule(int,int).
79  *
80  * Revision 2.32 1999/12/07 14:53:00 eric
81  * Changement ordre des arguments (phi,theta,xi) --> (xi,theta,phi)
82  * dans la routine val_point.
83  *
84  * Revision 2.31 1999/12/06 16:47:48 eric
85  * Ajout de la fonction val_point.
86  *
87  * Revision 2.30 1999/11/30 12:41:53 eric
88  * Le membre base est desormais un objet de type Base_val et non plus
89  * un pointeur vers une Base_val.
90  *
91  * Revision 2.29 1999/11/23 16:19:33 eric
92  * Suppression du constructeur par defaut.
93  *
94  * Revision 2.28 1999/11/23 14:32:16 novak
95  * Ajout des membres mult_ct et mult_st
96  *
97  * Revision 2.27 1999/11/22 15:41:33 eric
98  * Ajout de la fonction annule(int l).
99  *
100  * Revision 2.26 1999/11/19 11:51:01 eric
101  * Ajout du membre p_stdsdp.
102  *
103  * Revision 2.25 1999/11/19 10:53:17 eric
104  * *** empty log message ***
105  *
106  * Revision 2.24 1999/11/19 10:34:52 eric
107  * Corrections de diverses erreurs dans l'affectation (operator=):
108  * 1/ l'auto-affectation est desormais interdite
109  * 2/ l'affectation a des elements derives est desormais possible
110  *
111  * Revision 2.23 1999/11/19 09:32:15 eric
112  * Ajout du membre p_lapang.
113  *
114  * Revision 2.22 1999/11/16 13:28:32 novak
115  * Ajout de la gestion des pointeurs sur mult_x et scost
116  *
117  * Revision 2.21 1999/10/29 15:14:59 eric
118  * *** empty log message ***
119  *
120  * Revision 2.20 1999/10/28 13:24:15 phil
121  * copie passe avec ETATNONDEF
122  *
123  * Revision 2.19 1999/10/27 09:54:34 eric
124  * *** empty log message ***
125  *
126  * Revision 2.18 1999/10/22 08:42:53 eric
127  * *** empty log message ***
128  *
129  * Revision 2.17 1999/10/22 08:32:55 eric
130  * Anglicisation de operator<<
131  *
132  * Revision 2.16 1999/10/21 14:33:05 eric
133  * *** empty log message ***
134  *
135  * Revision 2.15 1999/10/21 14:21:23 eric
136  * Constructeur par lecture de fichier.
137  * Fonction sauve().
138  *
139  * Revision 2.14 1999/10/20 15:32:20 eric
140  * Ajout de la routine Valeur::std_base_scal().
141  * Ajout de operator=(const Mtbl_cf&).
142  *
143  * Revision 2.13 1999/10/19 15:30:33 eric
144  * Ajout de la fonction affiche_seuil.
145  *
146  * Revision 2.12 1999/10/18 15:16:17 eric
147  * *** empty log message ***
148  *
149  * Revision 2.11 1999/10/18 15:09:01 eric
150  * La fonction membre annule() est rebaptisee annule_hard().
151  * Introduction de la fonction membre annule(int, int).
152  *
153  * Revision 2.10 1999/10/18 13:43:29 eric
154  * Suppression de sxdsdx() et p_sxdsdx (car non-implementes).
155  *
156  * Revision 2.9 1999/10/13 15:52:35 eric
157  * Depoussierage.
158  *
159  * Revision 2.8 1999/04/12 13:05:40 phil
160  * *** empty log message ***
161  *
162  * Revision 2.7 1999/04/12 12:47:56 phil
163  * Correction du constructeur par recopie.
164  *
165  * Revision 2.6 1999/04/12 12:13:26 phil
166  * Correction de l'operateur Valeur = Valeur
167  *
168  * Revision 2.5 1999/03/01 15:10:42 eric
169  * *** empty log message ***
170  *
171  * Revision 2.4 1999/02/26 11:43:07 hyc
172  * *** empty log message ***
173  *
174  * Revision 2.3 1999/02/24 15:26:31 hyc
175  * *** empty log message ***
176  *
177  * Revision 2.2 1999/02/23 11:46:40 hyc
178  * *** empty log message ***
179  *
180  * $Header: /cvsroot/Lorene/C++/Source/Valeur/valeur.C,v 1.8 2014/10/13 08:53:49 j_novak Exp $
181  *
182  */
183 
184 // headers C
185 #include <cassert>
186 #include <cstdlib>
187 
188 // headers Lorene
189 #include "valeur.h"
190 #include "utilitaires.h"
191 #include "proto.h"
192 
193 
194 
195  //---------------//
196  // Constructeurs //
197  //---------------//
198 
199 namespace Lorene {
200 Valeur::Valeur(const Mg3d& mgi) : mg(&mgi), base(mgi.get_nzone()) {
201 
202  // C'est nouveau
203  nouveau() ;
204 }
205 
206 Valeur::Valeur(const Mg3d* mgi) : mg(mgi), base(mgi->get_nzone()) {
207 
208  // C'est nouveau
209  nouveau() ;
210 }
211 
212 // Copie
213 Valeur::Valeur(const Valeur& vc) : base(vc.base) {
214 
215  // Tout par default
216  mg = vc.get_mg() ;
217  nouveau() ;
218 
219  // Que faire ?
220  switch(vc.get_etat()) {
221  case ETATZERO:
222  set_etat_zero() ;
223  break ;
224 
225  case ETATNONDEF:
226  set_etat_nondef() ;
227  break ;
228 
229  case ETATQCQ:
230  assert((vc.c != 0x0) || (vc.c_cf != 0x0) ) ;
231  del_deriv() ;
232 
233  if (vc.c != 0x0) {
234  if (c == 0x0) {
235  c = new Mtbl( *(vc.c) ) ;
236  }
237  else{
238  *c = *(vc.c) ;
239  }
240  }
241 
242  if (vc.c_cf != 0x0) {
243  if (c_cf == 0x0) {
244  c_cf = new Mtbl_cf( *(vc.c_cf) ) ;
245  }
246  else{
247  *c_cf = *(vc.c_cf) ;
248  }
249  }
250 
251  etat = ETATQCQ ;
252  break ;
253 
254  default:
255  cout << "Etat pas possible" << endl ;
256  abort() ;
257  break ;
258  }
259 }
260 
261 // Constructeur a partir d'une grille et d'un fichier
262 // --------------------------------------------------
263 Valeur::Valeur(const Mg3d& g, FILE* fd) : mg(&g), base(fd) {
264 
265  // La multi-grille
266  Mg3d* mg_tmp = new Mg3d(fd) ; // la multi-grille d'origine
267  if (*mg != *mg_tmp) {
268  cout <<
269  "Valeur::Valeur(const Mg3d& g, FILE* fd): grid not consistent !"
270  << endl ;
271  abort() ;
272  }
273  delete mg_tmp ;
274 
275  fread_be(&etat, sizeof(int), 1, fd) ; // L'etat
276 
277  if (etat == ETATQCQ) {
278  c = new Mtbl(g, fd) ; // Les valeurs ds l'espace des conf.
279 
280  // Tous les autres pointeurs sont mis a zero :
281  c_cf = 0x0 ;
282  set_der_0x0() ;
283  }
284  else {
285  c = 0x0 ;
286  c_cf = 0x0 ;
287  set_der_0x0() ;
288  }
289 
290 }
291 
292  //--------------//
293  // Destructeurs //
294  //--------------//
295 
296 // Destructeur
298  del_t() ;
299 }
300 
301  //-------------//
302  // Affectation //
303  //-------------//
304 
305 // From double
306 // -----------
307 void Valeur::operator=(double x) {
308 
309  // Cas x = 0
310  if (x == 0) {
311  set_etat_zero() ;
312  return ;
313  }
314 
315  // Cas general
316  // les dependances
317  set_etat_c_qcq() ;
318 
319  // Les bases
321 
322  *c = x ;
323 }
324 
325 // From Valeur
326 // -----------
327 void Valeur::operator=(const Valeur& v) {
328 
329  // Protections:
330  // -----------
331  assert(&v != this) ; // pour eviter l'auto-affectation
332  assert(mg == v.mg) ; // meme grille
333 
334  // Copie de v :
335  // ------------
336  etat = v.etat ; // l'etat
337 
338  base = v.base ; // Les bases
339 
340  delete c ;
341  c = 0x0 ; // eventuellement provisoire...
342 
343  delete c_cf ;
344  c_cf = 0x0 ; // eventuellement provisoire...
345 
346  // La valeur eventuelle
347  switch(v.etat) {
348  case ETATNONDEF: {
349  set_etat_nondef() ;
350  break ; // valeur par defaut
351  }
352 
353  case ETATZERO: {
354  set_etat_zero() ;
355  break ;
356  }
357 
358  case ETATQCQ: {
359  if (v.c != 0x0) {
360  c = new Mtbl( *(v.c) ) ;
361  }
362 
363  if (v.c_cf != 0x0) {
364  c_cf = new Mtbl_cf( *(v.c_cf) ) ;
365  }
366 
367  etat = ETATQCQ ;
368 
369  // On detruit les dependances (seulement lorsque tout est fini !)
370  del_deriv() ;
371 
372  break ;
373  }
374 
375  default: {
376  cout << "Unkwon state in Valeur::operator=(const Valeur&) !"
377  << endl ;
378  abort() ;
379  break ;
380  }
381  }
382 
383 }
384 
385 // From Mtbl
386 // ---------
387 void Valeur::operator=(const Mtbl& mt) {
388 
389  // Protections
390  // -----------
391  assert(mt.get_etat() != ETATNONDEF) ;
392  assert(mg == mt.get_mg()) ; // meme grille
393  assert(&mt != c) ; // pour eviter l'autoaffectation
394 
395  // Les bases
397 
398  delete c_cf ;
399  c_cf = 0x0 ;
400 
401  // Suivant le cas...
402  switch(mt.get_etat()) {
403  case ETATZERO: {
404  set_etat_zero() ;
405  break ;
406  }
407 
408  case ETATQCQ: {
409  delete c ;
410  c = new Mtbl(mt) ;
411 
412  del_deriv() ; // On detruit les dependances...
413 
414  etat = ETATQCQ ;
415  break ;
416  }
417 
418  default: {
419  cout << "Unkwon state in Valeur::operator=(const Mtbl&) !"
420  << endl ;
421  abort() ;
422  break ;
423  }
424  }
425 
426 }
427 
428 // From Mtbl_cf
429 // ------------
430 void Valeur::operator=(const Mtbl_cf& mt) {
431 
432  // Protections
433  // -----------
434  assert(mt.get_etat() != ETATNONDEF) ;
435  assert(mg == mt.get_mg()) ; // meme grille
436  assert(&mt != c_cf) ; // pour eviter l'autoaffectation
437 
438  // Les bases
439  base = mt.base ;
440 
441  delete c ;
442  c = 0x0 ;
443 
444  // Suivant le cas...
445  switch(mt.get_etat()) {
446  case ETATZERO: {
447  set_etat_zero() ;
448  break ;
449  }
450 
451  case ETATQCQ: {
452  delete c_cf ;
453  c_cf = new Mtbl_cf(mt) ;
454 
455  del_deriv() ; // On detruit les dependances...
456 
457  etat = ETATQCQ ;
458  break ;
459  }
460 
461  default: {
462  cout << "Unkwon state in Valeur::operator=(const Mtbl_cf&) !"
463  << endl ;
464  abort() ;
465  break ;
466  }
467  }
468 
469 }
470 
471  //------------//
472  // Sauvegarde //
473  //------------//
474 
475 void Valeur::sauve(FILE* fd) const {
476 
477  base.sauve(fd) ; // la base
478  mg->sauve(fd) ; // la multi-grille
479  fwrite_be(&etat, sizeof(int), 1, fd) ; // l'etat
480 
481  if (etat == ETATQCQ) {
482  if (c == 0x0) { // La sauvegarde s'effectue dans l'espace
483  coef_i() ; // des configurations
484  }
485  c->sauve(fd) ;
486  }
487 
488 }
489 
490  //------------//
491  // Impression //
492  //------------//
493 
494 // Operator <<
495 // -----------
496 ostream& operator<<(ostream& o, const Valeur & vi) {
497 
498  switch(vi.etat) {
499  case ETATNONDEF: {
500  o << "*** Valeur in UNDEFINED STATE" ;
501  break ;
502  }
503 
504  case ETATZERO: {
505  o << "*** Valeur IDENTICALLY ZERO" ;
506  break ;
507  }
508 
509  case ETATQCQ: {
510  if (vi.c != 0x0) {
511  o << "*** Valeur (configuration space) :" << endl ;
512  o << *(vi.c) << endl ;
513  }
514  if (vi.c_cf != 0x0) {
515  o << "*** Valeur (coefficients) :" << endl ;
516  o << *(vi.c_cf) << endl ;
517  }
518  break ;
519  }
520 
521  default: {
522  cout << "operator<<(ostream&, const Valeur&) : unknown state !"
523  << endl ;
524  abort() ;
525  break ;
526  }
527 
528  }
529 
530  // Termine
531  return o ;
532 }
533 
534 // display_coef
535 // ------------
536 void Valeur::display_coef(double thres, int precis, ostream& ost) const {
537 
538  if (etat == ETATNONDEF) {
539  ost << " state: UNDEFINED" << endl ;
540  return ;
541  }
542 
543  if (etat == ETATZERO) {
544  ost << " state: ZERO" << endl ;
545  return ;
546  }
547 
548  coef() ; // the coefficients are required
549 
550  c_cf->display(thres, precis, ost) ;
551 
552 }
553 
554 
555 // affiche_seuil
556 //---------------
557 
558 void Valeur::affiche_seuil(ostream& ost, int type, int precis,
559  double seuil) const {
560  ost << "*** Valeur " << endl ;
561 
562 
563  // Cas particuliers
564  //-----------------
565 
566  if (etat == ETATNONDEF) {
567  ost << " state: UNDEFINED" << endl ;
568  return ;
569  }
570 
571  if (etat == ETATZERO) {
572  ost << " state: ZERO" << endl ;
573  return ;
574  }
575 
576  switch(type) {
577  case 0 : {
578  coef() ; // Mise a jour eventuelle des coefficients
579  ost << " Coefficients of the Valeur : " << endl ;
580  c_cf->affiche_seuil(ost, precis, seuil) ;
581  break ;
582  }
583 
584  case 1 : {
585  coef_i() ; // Mise a jour eventuelle dans l'espace des conf.
586  ost << " Values of the Valeur at the collocation points: " << endl ;
587  c->affiche_seuil(ost, precis, seuil) ;
588  break ;
589  }
590 
591  case 2 : {
592  coef() ; // Mise a jour eventuelle des coefficients
593  coef_i() ; // Mise a jour eventuelle dans l'espace des conf.
594  ost << " Coefficients of the Valeur : " << endl ;
595  c_cf->affiche_seuil(ost, precis, seuil) ;
596  ost << " Values of the Valeur at the collocation points: " << endl ;
597  c->affiche_seuil(ost, precis, seuil) ;
598  break ;
599  }
600 
601  default : {
602  cout << "Unknown type in Valeur::affiche_seuil !" << endl ;
603  abort() ;
604  break ;
605  }
606  }
607 
608 
609 }
610 
611 
612 
613  //-----------------------//
614  // Gestion de la memoire //
615  //-----------------------//
616 
618  // Les donnees
619  c = 0x0 ;
620  c_cf = 0x0 ;
621  set_der_0x0() ;
622  // L'etat
623  etat = ETATNONDEF ;
624 }
625 
627 
628  delete c ;
629  c = 0x0 ;
630 
631  delete c_cf ;
632  c_cf = 0x0 ;
633 
634  del_deriv() ;
635 
636 }
637 
639  // Destruction
640  delete p_dsdx ;
641  delete p_d2sdx2 ;
642  delete p_sx ;
643  delete p_sx2 ;
644  delete p_mult_x ;
645 
646  delete p_dsdt ;
647  delete p_d2sdt2 ;
648  delete p_ssint ;
649  delete p_scost ;
650  delete p_mult_ct ;
651  delete p_mult_st ;
652 
653  delete p_dsdp ;
654  delete p_stdsdp ;
655  delete p_d2sdp2 ;
656  delete p_mult_cp ;
657  delete p_mult_sp ;
658 
659  delete p_lapang ;
660 
661  // Pointeurs a 0x0
662  set_der_0x0() ;
663 }
664 
666  p_dsdx = 0x0 ;
667  p_d2sdx2 = 0x0 ;
668  p_sx = 0x0 ;
669  p_sx2 = 0x0 ;
670  p_mult_x = 0x0 ;
671 
672  p_dsdt = 0x0 ;
673  p_d2sdt2 = 0x0 ;
674  p_ssint = 0x0 ;
675  p_scost = 0x0 ;
676  p_mult_ct = 0x0 ;
677  p_mult_st = 0x0 ;
678 
679  p_dsdp = 0x0 ;
680  p_stdsdp = 0x0 ;
681  p_d2sdp2 = 0x0 ;
682  p_mult_cp = 0x0 ;
683  p_mult_sp = 0x0 ;
684 
685  p_lapang = 0x0 ;
686 }
687 
688 // ETATZERO
690  if (etat == ETATZERO) return ;
691  del_t() ;
692  etat = ETATZERO ;
693 }
694 // ETATNONDEF
696  if (etat == ETATNONDEF) return ;
697  del_t() ;
698  etat = ETATNONDEF ;
699 }
700 // ETATQCQ physique
702  // Detruit les dependances
703  del_deriv() ;
704  delete c_cf ; c_cf = 0x0 ;
705 
706  if (c == 0x0) {
707  c = new Mtbl(mg) ;
708  }
709  etat = ETATQCQ ;
710 }
711 // ETATQCQ coefficients
713  // Detruit les dependances
714  del_deriv() ;
715  delete c ; c = 0x0 ;
716 
717  if (c_cf == 0x0) {
718  c_cf = new Mtbl_cf(mg, base) ;
719  }
720  etat = ETATQCQ ;
721 }
722 // ZERO hard
724  // Detruit les dependances
725  del_deriv() ;
726 
727  if (c == 0x0) {
728  c = new Mtbl(mg) ;
729  }
730  if (c_cf == 0x0) {
731  c_cf = new Mtbl_cf(mg, base) ;
732  }
733 
734  c->annule_hard() ;
735  c_cf->annule_hard() ;
736 
737  etat = ETATQCQ ;
738 }
739 
740 
741 // Sets the Valeur to zero in a given domain
742 // -----------------------------------------
743 
744 void Valeur::annule(int l) {
745 
746  annule(l, l) ;
747 }
748 
749 
750 
751 // Sets the Valeur to zero in several domains
752 // ------------------------------------------
753 
754 void Valeur::annule(int l_min, int l_max) {
755 
756  // Cas particulier: annulation globale :
757  if ( (l_min == 0) && (l_max == mg->get_nzone()-1) ) {
758  set_etat_zero() ;
759  return ;
760  }
761 
762  assert( etat != ETATNONDEF ) ;
763 
764  if ( etat == ETATZERO ) {
765  return ; // rien n'a faire si c'est deja zero
766  }
767  else {
768  assert( etat == ETATQCQ ) ; // sinon...
769 
770  if (c != 0x0) {
771  c->annule(l_min, l_max) ; // Annule le Mtbl
772  }
773 
774  if (c_cf != 0x0) {
775  c_cf->annule(l_min, l_max) ; // Annule le Mtbl_cf
776  }
777 
778  // Annulation des membres derives
779 
780  if (p_dsdx != 0x0) p_dsdx->annule(l_min, l_max) ;
781  if (p_d2sdx2 != 0x0) p_d2sdx2->annule(l_min, l_max) ;
782  if (p_sx != 0x0) p_sx->annule(l_min, l_max) ;
783  if (p_sx2 != 0x0) p_sx2->annule(l_min, l_max) ;
784  if (p_mult_x != 0x0) p_mult_x->annule(l_min, l_max) ;
785 
786  if (p_dsdt != 0x0) p_dsdt->annule(l_min, l_max) ;
787  if (p_d2sdt2 != 0x0) p_d2sdt2->annule(l_min, l_max) ;
788  if (p_ssint != 0x0) p_ssint->annule(l_min, l_max) ;
789  if (p_scost != 0x0) p_scost->annule(l_min, l_max) ;
790  if (p_mult_ct != 0x0) p_mult_ct->annule(l_min, l_max) ;
791  if (p_mult_st != 0x0) p_mult_st->annule(l_min, l_max) ;
792 
793  if (p_dsdp != 0x0) p_dsdp->annule(l_min, l_max) ;
794  if (p_stdsdp != 0x0) p_stdsdp->annule(l_min, l_max) ;
795  if (p_d2sdp2 != 0x0) p_d2sdp2->annule(l_min, l_max) ;
796  if (p_mult_cp != 0x0) p_mult_cp->annule(l_min, l_max) ;
797  if (p_mult_sp != 0x0) p_mult_sp->annule(l_min, l_max) ;
798 
799  if (p_lapang != 0x0) p_lapang->annule(l_min, l_max) ;
800 
801  }
802 
803 }
804 
805 
806  //--------------------------------------//
807  // Spectral bases manipulation //
808  //--------------------------------------//
809 
810 void Valeur::set_base(const Base_val& base_i) {
811 
812  base = base_i ;
813 
814  if (c_cf != 0x0) {
815  if ( !(c_cf->base == base_i) ) {
816  delete c_cf ;
817  c_cf = 0x0 ;
818  }
819  }
820 
821 }
822 
823 
825 
826  set_base( mg->std_base_scal() ) ;
827 
828 }
829 
831 
832  set_base( mg->std_base_scal_odd() ) ;
833 
834 }
835 
836 void Valeur::set_base_r(int l, int base_r) {
837 
838  base.set_base_r(l, base_r) ;
839 
840  if (c_cf != 0x0) {
841  if ( !(c_cf->base == base) ) {
842  delete c_cf ;
843  c_cf = 0x0 ;
844  }
845  }
846 
847 }
848 
849 void Valeur::set_base_t(int base_t) {
850 
851  base.set_base_t(base_t) ;
852 
853  if (c_cf != 0x0) {
854  if ( !(c_cf->base == base) ) {
855  delete c_cf ;
856  c_cf = 0x0 ;
857  }
858  }
859 
860 }
861 
862 void Valeur::set_base_p(int base_p) {
863 
864  base.set_base_p(base_p) ;
865 
866  if (c_cf != 0x0) {
867  if ( !(c_cf->base == base) ) {
868  delete c_cf ;
869  c_cf = 0x0 ;
870  }
871  }
872 
873 }
874 
875 
876 
877 
878  //-----------------------------------------------//
879  // Value at an arbitrary point //
880  //-----------------------------------------------//
881 
882 double Valeur::val_point(int l, double x, double theta, double phi) const {
883 
884  assert(etat != ETATNONDEF) ;
885 
886  double resu ;
887 
888  if (etat == ETATZERO) {
889  resu = 0 ;
890  }
891  else{
892  assert(etat == ETATQCQ) ;
893  coef() ; // The coefficients are required
894  resu = c_cf->val_point(l, x, theta, phi) ; // Call to the Mtbl_cf version
895  }
896 
897  return resu ;
898 }
899 
900 double Valeur::val_point_jk(int l, double x, int j, int k) const {
901 
902  assert(etat != ETATNONDEF) ;
903 
904  double resu ;
905 
906  if (etat == ETATZERO) {
907  resu = 0 ;
908  }
909  else{
910  assert(etat == ETATQCQ) ;
911  coef() ; // The coefficients are required
912  resu = c_cf->val_point_jk(l, x, j, k) ; // Call to the Mtbl_cf version
913  }
914 
915  return resu ;
916 }
917 
918 
919  //-----------------------------------------------//
920  // Filtres //
921  //-----------------------------------------------//
922 
923 
924 void Valeur::filtre_tp(int nn, int nz1, int nz2) {
925 
926  int nz = mg->get_nzone() ;
927  int nr = mg->get_nr(0) ;
928  int nt = mg->get_nt(0) ;
929  int np = mg->get_np(0) ;
930 
931  if (etat != ETATZERO) {
932  assert (etat == ETATQCQ) ;
933  ylm() ;
934  for (int lz=nz1; lz<=nz2; lz++) {
935  if (c_cf->operator()(lz).get_etat() != ETATZERO)
936  for (int k=0; k<np+1; k++)
937  for (int j=0; j<nt; j++) {
938  int l_q, m_q, base_r ;
939  if (nullite_plm(j, nt, k, np, base) == 1) {
940  donne_lm(nz, lz, j, k, base, m_q, l_q,base_r) ;
941  if (l_q > nn)
942  for (int i=0; i<nr; i++)
943  c_cf->set(lz, k, j, i) = 0. ;
944  }
945  }
946  }
947  if (c != 0x0) {
948  delete c ;
949  c = 0x0 ;
950  }
951  }
952 
953 }
954 }
Lorene::Valeur::p_dsdp
Valeur * p_dsdp
Pointer on .
Definition: valeur.h:323
Lorene::Base_val
Bases of the spectral expansions.
Definition: base_val.h:322
Lorene::Mg3d::get_np
int get_np(int l) const
Returns the number of points in the azimuthal direction ( ) in domain no. l.
Definition: grilles.h:462
Lorene::Valeur::set_etat_cf_qcq
void set_etat_cf_qcq()
Sets the logical state to ETATQCQ (ordinary state) for values in the configuration space (Mtbl_cf c_c...
Definition: valeur.C:712
Lorene::Mtbl_cf::val_point
double val_point(int l, double x, double theta, double phi) const
Computes the value of the field represented by *this at an arbitrary point, by means of the spectral ...
Definition: mtbl_cf_val_point.C:115
Lorene::Valeur
Values and coefficients of a (real-value) function.
Definition: valeur.h:287
Lorene::Valeur::val_point
double val_point(int l, double x, double theta, double phi) const
Computes the value of the field represented by *this at an arbitrary point, by means of the spectral ...
Definition: valeur.C:882
Lorene::Valeur::set_etat_c_qcq
void set_etat_c_qcq()
Sets the logical state to ETATQCQ (ordinary state) for values in the configuration space (Mtbl c ).
Definition: valeur.C:701
Lorene::Valeur::p_d2sdt2
Valeur * p_d2sdt2
Pointer on .
Definition: valeur.h:317
Lorene::Base_val::set_base_t
void set_base_t(int base_t)
Sets the expansion basis for functions in all domains.
Definition: base_val.C:195
Lorene::Valeur::p_sx
Valeur * p_sx
Pointer on .
Definition: valeur.h:312
Lorene::Base_val::set_base_r
void set_base_r(int l, int base_r)
Sets the expansion basis for r ( ) functions in a given domain.
Definition: base_val.C:185
Lorene::Valeur::get_etat
int get_etat() const
Returns the logical state.
Definition: valeur.h:726
Lorene::Valeur::del_deriv
void del_deriv()
Logical destructor of the derivatives.
Definition: valeur.C:638
Lorene::Valeur::del_t
void del_t()
Logical destructor.
Definition: valeur.C:626
Lorene::Mg3d
Multi-domain grid.
Definition: grilles.h:273
Lorene::Valeur::c_cf
Mtbl_cf * c_cf
Coefficients of the spectral expansion of the function.
Definition: valeur.h:302
Lorene::Mtbl_cf::val_point_jk
double val_point_jk(int l, double x, int j, int k) const
Computes the value of the field represented by *this at an arbitrary point in , but collocation point...
Definition: mtbl_cf_val_point.C:258
Lorene::Valeur::p_dsdx
Valeur * p_dsdx
Pointer on .
Definition: valeur.h:310
Lorene::Mg3d::get_nt
int get_nt(int l) const
Returns the number of points in the co-latitude direction ( ) in domain no. l.
Definition: grilles.h:457
Lorene
Lorene prototypes.
Definition: app_hor.h:64
Lorene::Valeur::sauve
void sauve(FILE *) const
Save in a file.
Definition: valeur.C:475
Lorene::Mtbl::annule
void annule(int l_min, int l_max)
Sets the Mtbl to zero in some domains.
Definition: mtbl.C:329
Lorene::Mg3d::std_base_scal_odd
Base_val std_base_scal_odd() const
Returns the standard odd spectral bases for a scalar.
Definition: mg3d_std_base.C:131
Lorene::Mtbl
Multi-domain array.
Definition: mtbl.h:118
Lorene::Valeur::filtre_tp
void filtre_tp(int nn, int nz1, int nz2)
Sets the n lasts coefficients in to 0 in the domain nz1 to nz2 when expressed in spherical harmonics...
Definition: valeur.C:924
Lorene::Valeur::p_ssint
Valeur * p_ssint
Pointer on .
Definition: valeur.h:318
Lorene::Valeur::p_mult_sp
Valeur * p_mult_sp
Pointer on .
Definition: valeur.h:327
Lorene::Mtbl::get_mg
const Mg3d * get_mg() const
Gives the Mg3d on which the Mtbl is defined.
Definition: mtbl.h:274
Lorene::Valeur::Valeur
Valeur(const Mg3d &mgrid)
Constructor.
Definition: valeur.C:200
Lorene::Mtbl_cf::get_etat
int get_etat() const
Returns the logical state.
Definition: mtbl_cf.h:456
Lorene::Valeur::get_mg
const Mg3d * get_mg() const
Returns the Mg3d on which the this is defined.
Definition: valeur.h:729
Lorene::Valeur::annule
void annule(int l)
Sets the Valeur to zero in a given domain.
Definition: valeur.C:744
Lorene::Valeur::set_base_r
void set_base_r(int l, int base_r)
Sets the expansion basis for r ( ) functions in a given domain.
Definition: valeur.C:836
Lorene::Mtbl_cf::set
Tbl & set(int l)
Read/write of the Tbl containing the coefficients in a given domain.
Definition: mtbl_cf.h:294
Lorene::Valeur::set_base_p
void set_base_p(int base_p)
Sets the expansion basis for functions in all domains.
Definition: valeur.C:862
Lorene::Valeur::c
Mtbl * c
Values of the function at the points of the multi-grid
Definition: valeur.h:299
Lorene::Valeur::coef
void coef() const
Computes the coeffcients of *this.
Definition: valeur_coef.C:148
Lorene::Mtbl::sauve
void sauve(FILE *) const
Save in a file.
Definition: mtbl.C:209
Lorene::Valeur::annule_hard
void annule_hard()
Sets the Valeur to zero in a hard way.
Definition: valeur.C:723
Lorene::Valeur::affiche_seuil
void affiche_seuil(ostream &ostr, int type=0, int precision=4, double threshold=1.e-7) const
Prints only the values greater than a given threshold.
Definition: valeur.C:558
Lorene::Valeur::mg
const Mg3d * mg
Multi-grid Mgd3 on which this is defined.
Definition: valeur.h:292
Lorene::Mtbl_cf::annule_hard
void annule_hard()
Sets the Mtbl_cf to zero in a hard way.
Definition: mtbl_cf.C:312
Lorene::Valeur::p_scost
Valeur * p_scost
Pointer on
Definition: valeur.h:319
Lorene::Mtbl::affiche_seuil
void affiche_seuil(ostream &ostr, int precision=4, double threshold=1.e-7) const
Prints only the values greater than a given threshold.
Definition: mtbl.C:393
Lorene::Base_val::sauve
void sauve(FILE *) const
Save in a file.
Definition: base_val.C:228
Lorene::Valeur::std_base_scal
void std_base_scal()
Sets the bases for spectral expansions (member base ) to the standard ones for a scalar.
Definition: valeur.C:824
Lorene::Valeur::etat
int etat
Logical state (ETATNONDEF , ETATQCQ or ETATZERO ).
Definition: valeur.h:295
Lorene::Valeur::set_etat_nondef
void set_etat_nondef()
Sets the logical state to ETATNONDEF (undefined).
Definition: valeur.C:695
Lorene::Valeur::p_mult_x
Valeur * p_mult_x
Pointer on .
Definition: valeur.h:314
Lorene::Mtbl::get_etat
int get_etat() const
Gives the logical state.
Definition: mtbl.h:277
Lorene::Mg3d::sauve
void sauve(FILE *fd, bool save_base=false) const
Saves into a file.
Definition: mg3d.C:371
Lorene::Mg3d::get_nzone
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:448
Lorene::Base_val::set_base_nondef
void set_base_nondef()
Sets the spectral bases to NONDEF.
Definition: base_val.C:326
Lorene::Mtbl_cf::affiche_seuil
void affiche_seuil(ostream &ostr, int precision=4, double threshold=1.e-7) const
Prints only the values greater than a given threshold.
Definition: mtbl_cf.C:397
Lorene::Mg3d::get_nr
int get_nr(int l) const
Returns the number of points in the radial direction ( ) in domain no. l.
Definition: grilles.h:452
Lorene::Valeur::~Valeur
~Valeur()
Destructor.
Definition: valeur.C:297
Lorene::fread_be
int fread_be(int *aa, int size, int nb, FILE *fich)
Reads integer(s) from a binary file according to the big endian convention.
Definition: fread_be.C:69
Lorene::Valeur::set_base
void set_base(const Base_val &)
Sets the bases for spectral expansions (member base )
Definition: valeur.C:810
Lorene::Base_val::set_base_p
void set_base_p(int base_p)
Sets the expansion basis for functions in all domains.
Definition: base_val.C:204
Lorene::Valeur::p_mult_st
Valeur * p_mult_st
Pointer on .
Definition: valeur.h:321
Lorene::Valeur::p_sx2
Valeur * p_sx2
Pointer on .
Definition: valeur.h:313
Lorene::Valeur::p_d2sdx2
Valeur * p_d2sdx2
Pointer on .
Definition: valeur.h:311
Lorene::Valeur::p_stdsdp
Valeur * p_stdsdp
Pointer on .
Definition: valeur.h:324
Lorene::Mtbl_cf::get_mg
const Mg3d * get_mg() const
Returns the Mg3d on which the Mtbl_cf is defined.
Definition: mtbl_cf.h:453
Lorene::Valeur::set_base_t
void set_base_t(int base_t)
Sets the expansion basis for functions in all domains.
Definition: valeur.C:849
Lorene::Mtbl_cf::base
Base_val base
Bases of the spectral expansions.
Definition: mtbl_cf.h:200
Lorene::Mtbl_cf::annule
void annule(int l_min, int l_max)
Sets the Mtbl_cf to zero in some domains.
Definition: mtbl_cf.C:332
Lorene::Valeur::val_point_jk
double val_point_jk(int l, double x, int j, int k) const
Computes the value of the field represented by *this at an arbitrary point in , but collocation point...
Definition: valeur.C:900
Lorene::Valeur::p_d2sdp2
Valeur * p_d2sdp2
Pointer on .
Definition: valeur.h:325
Lorene::Valeur::std_base_scal_odd
void std_base_scal_odd()
Sets the bases for spectral expansions (member base ) to the standard odd ones for a scalar.
Definition: valeur.C:830
Lorene::Valeur::p_dsdt
Valeur * p_dsdt
Pointer on .
Definition: valeur.h:316
Lorene::Valeur::operator=
void operator=(const Valeur &a)
Assignement to another Valeur.
Definition: valeur.C:327
Lorene::Mtbl::annule_hard
void annule_hard()
Sets the Mtbl to zero in a hard way.
Definition: mtbl.C:312
Lorene::fwrite_be
int fwrite_be(const int *aa, int size, int nb, FILE *fich)
Writes integer(s) into a binary file according to the big endian convention.
Definition: fwrite_be.C:70
Lorene::Valeur::p_lapang
Valeur * p_lapang
Pointer on the angular Laplacian.
Definition: valeur.h:329
Lorene::Valeur::coef_i
void coef_i() const
Computes the physical value of *this.
Definition: valeur_coef_i.C:140
Lorene::Valeur::display_coef
void display_coef(double threshold=1.e-7, int precision=4, ostream &ostr=cout) const
Displays the spectral coefficients and the associated basis functions.
Definition: valeur.C:536
Lorene::Valeur::base
Base_val base
Bases on which the spectral expansion is performed.
Definition: valeur.h:305
Lorene::Valeur::set_etat_zero
void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition: valeur.C:689
Lorene::Valeur::nouveau
void nouveau()
Memory allocation.
Definition: valeur.C:617
Lorene::Valeur::p_mult_ct
Valeur * p_mult_ct
Pointer on .
Definition: valeur.h:320
Lorene::Valeur::set_der_0x0
void set_der_0x0()
Sets the pointers for derivatives to 0x0.
Definition: valeur.C:665
Lorene::Valeur::ylm
void ylm()
Computes the coefficients of *this.
Definition: valeur_ylm.C:138
Lorene::Mtbl_cf::display
void display(double threshold=1.e-7, int precision=4, ostream &ostr=cout) const
Prints the coefficients whose values are greater than a given threshold, as well as the corresponding...
Definition: mtbl_cf_display.C:56
Lorene::Mg3d::std_base_scal
Base_val std_base_scal() const
Returns the standard spectral bases for a scalar.
Definition: mg3d_std_base.C:101
Lorene::Mtbl_cf
Coefficients storage for the multi-domain spectral method.
Definition: mtbl_cf.h:186
Lorene::Valeur::p_mult_cp
Valeur * p_mult_cp
Pointer on .
Definition: valeur.h:326