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
32char 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
199namespace Lorene {
200Valeur::Valeur(const Mg3d& mgi) : mg(&mgi), base(mgi.get_nzone()) {
201
202 // C'est nouveau
203 nouveau() ;
204}
205
206Valeur::Valeur(const Mg3d* mgi) : mg(mgi), base(mgi->get_nzone()) {
207
208 // C'est nouveau
209 nouveau() ;
210}
211
212// Copie
213Valeur::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:
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// --------------------------------------------------
263Valeur::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// -----------
307void 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
318
319 // Les bases
321
322 *c = x ;
323}
324
325// From Valeur
326// -----------
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// ---------
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// ------------
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
475void 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// -----------
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// ------------
536void 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
558void 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
744void Valeur::annule(int l) {
745
746 annule(l, l) ;
747}
748
749
750
751// Sets the Valeur to zero in several domains
752// ------------------------------------------
753
754void 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
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
836void 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
850
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
863
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
882double 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
900double 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
924void 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}
Bases of the spectral expansions.
Definition base_val.h:322
void sauve(FILE *) const
Save in a file.
Definition base_val.C:228
void set_base_t(int base_t)
Sets the expansion basis for functions in all domains.
Definition base_val.C:195
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
void set_base_nondef()
Sets the spectral bases to NONDEF.
Definition base_val.C:326
void set_base_p(int base_p)
Sets the expansion basis for functions in all domains.
Definition base_val.C:204
Time evolution with partial storage (*** under development ***).
Definition evolution.h:371
Multi-domain grid.
Definition grilles.h:273
Coefficients storage for the multi-domain spectral method.
Definition mtbl_cf.h:186
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...
Base_val base
Bases of the spectral expansions.
Definition mtbl_cf.h:200
int get_etat() const
Returns the logical state.
Definition mtbl_cf.h:456
void annule(int l_min, int l_max)
Sets the Mtbl_cf to zero in some domains.
Definition mtbl_cf.C:332
Tbl & set(int l)
Read/write of the Tbl containing the coefficients in a given domain.
Definition mtbl_cf.h:294
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
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...
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 ...
void annule_hard()
Sets the Mtbl_cf to zero in a hard way.
Definition mtbl_cf.C:312
Multi-domain array.
Definition mtbl.h:118
void sauve(FILE *) const
Save in a file.
Definition mtbl.C:209
void annule(int l_min, int l_max)
Sets the Mtbl to zero in some domains.
Definition mtbl.C:329
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
void annule_hard()
Sets the Mtbl to zero in a hard way.
Definition mtbl.C:312
Values and coefficients of a (real-value) function.
Definition valeur.h:287
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
void del_t()
Logical destructor.
Definition valeur.C:626
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
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
~Valeur()
Destructor.
Definition valeur.C:297
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
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
Valeur * p_dsdp
Pointer on .
Definition valeur.h:323
Valeur * p_mult_sp
Pointer on .
Definition valeur.h:327
Valeur * p_mult_x
Pointer on .
Definition valeur.h:314
Valeur * p_scost
Pointer on
Definition valeur.h:319
Valeur * p_sx2
Pointer on .
Definition valeur.h:313
Valeur * p_d2sdx2
Pointer on .
Definition valeur.h:311
void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition valeur.C:689
void set_base(const Base_val &)
Sets the bases for spectral expansions (member base )
Definition valeur.C:810
Valeur * p_dsdx
Pointer on .
Definition valeur.h:310
Valeur(const Mg3d &mgrid)
Constructor.
Definition valeur.C:200
void nouveau()
Memory allocation.
Definition valeur.C:617
void ylm()
Computes the coefficients of *this.
Definition valeur_ylm.C:138
Valeur * p_d2sdt2
Pointer on .
Definition valeur.h:317
Mtbl * c
Values of the function at the points of the multi-grid
Definition valeur.h:299
Valeur * p_mult_cp
Pointer on .
Definition valeur.h:326
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
void annule(int l)
Sets the Valeur to zero in a given domain.
Definition valeur.C:744
const Mg3d * mg
Multi-grid Mgd3 on which this is defined.
Definition valeur.h:292
void operator=(const Valeur &a)
Assignement to another Valeur.
Definition valeur.C:327
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
Valeur * p_sx
Pointer on .
Definition valeur.h:312
void coef_i() const
Computes the physical value of *this.
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
Mtbl_cf * c_cf
Coefficients of the spectral expansion of the function.
Definition valeur.h:302
void coef() const
Computes the coeffcients of *this.
Valeur * p_ssint
Pointer on .
Definition valeur.h:318
void set_der_0x0()
Sets the pointers for derivatives to 0x0.
Definition valeur.C:665
Valeur * p_stdsdp
Pointer on .
Definition valeur.h:324
void set_base_r(int l, int base_r)
Sets the expansion basis for r ( ) functions in a given domain.
Definition valeur.C:836
Valeur * p_dsdt
Pointer on .
Definition valeur.h:316
Valeur * p_lapang
Pointer on the angular Laplacian.
Definition valeur.h:329
Base_val base
Bases on which the spectral expansion is performed.
Definition valeur.h:305
void std_base_scal()
Sets the bases for spectral expansions (member base ) to the standard ones for a scalar.
Definition valeur.C:824
void annule_hard()
Sets the Valeur to zero in a hard way.
Definition valeur.C:723
Valeur * p_mult_ct
Pointer on .
Definition valeur.h:320
int etat
Logical state (ETATNONDEF , ETATQCQ or ETATZERO ).
Definition valeur.h:295
void set_base_t(int base_t)
Sets the expansion basis for functions in all domains.
Definition valeur.C:849
void set_base_p(int base_p)
Sets the expansion basis for functions in all domains.
Definition valeur.C:862
void set_etat_nondef()
Sets the logical state to ETATNONDEF (undefined).
Definition valeur.C:695
Valeur * p_d2sdp2
Pointer on .
Definition valeur.h:325
Valeur * p_mult_st
Pointer on .
Definition valeur.h:321
void del_deriv()
Logical destructor of the derivatives.
Definition valeur.C:638
void sauve(FILE *) const
Save in a file.
Definition valeur.C:475
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
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 prototypes.
Definition app_hor.h:64