LORENE
param_elliptic.C
1/*
2 * Copyright (c) 2004 Philippe Grandclement
3 *
4 * This file is part of LORENE.
5 *
6 * LORENE is free software; you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License version 2
8 * as published by the Free Software Foundation.
9 *
10 * LORENE is distributed in the hope that it will be useful,
11 * but WITHOUT ANY WARRANTY; without even the implied warranty of
12 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13 * GNU General Public License for more details.
14 *
15 * You should have received a copy of the GNU General Public License
16 * along with LORENE; if not, write to the Free Software
17 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
18 *
19 */
20
21char param_elliptic_C[] = "$Header: /cvsroot/Lorene/C++/Source/Param_elliptic/param_elliptic.C,v 1.21 2014/10/13 08:53:37 j_novak Exp $" ;
22
23/*
24 * $Id: param_elliptic.C,v 1.21 2014/10/13 08:53:37 j_novak Exp $
25 * $Log: param_elliptic.C,v $
26 * Revision 1.21 2014/10/13 08:53:37 j_novak
27 * Lorene classes and functions now belong to the namespace Lorene.
28 *
29 * Revision 1.20 2014/10/06 15:13:15 j_novak
30 * Modified #include directives to use c++ syntax.
31 *
32 * Revision 1.19 2007/05/06 10:48:14 p_grandclement
33 * Modification of a few operators for the vorton project
34 *
35 * Revision 1.18 2007/04/24 09:04:13 p_grandclement
36 * Addition of an operator for the vortons
37 *
38 * Revision 1.17 2005/05/12 09:49:44 j_novak
39 * Temptative treatment of the case where the source is null in the CED (putting
40 * dzpuis to 4). May be a bad idea...
41 *
42 * Revision 1.16 2005/02/15 15:43:17 j_novak
43 * First version of the block inversion for the vector Poisson equation (method 6).
44 *
45 * Revision 1.15 2004/12/23 16:30:16 j_novak
46 * New files and class for the solution of the rr component of the tensor Poisson
47 * equation.
48 *
49 * Revision 1.14 2004/08/24 09:14:49 p_grandclement
50 * Addition of some new operators, like Poisson in 2d... It now requieres the
51 * GSL library to work.
52 *
53 * Also, the way a variable change is stored by a Param_elliptic is changed and
54 * no longer uses Change_var but rather 2 Scalars. The codes using that feature
55 * will requiere some modification. (It should concern only the ones about monopoles)
56 *
57 * Revision 1.13 2004/06/22 13:46:52 j_novak
58 * Deplacement du conte++ dans set_pois_vect_r
59 *
60 * Revision 1.12 2004/06/22 08:49:59 p_grandclement
61 * Addition of everything needed for using the logarithmic mapping
62 *
63 * Revision 1.11 2004/06/14 15:07:12 j_novak
64 * New methods for the construction of the elliptic operator appearing in
65 * the vector Poisson equation (acting on eta).
66 *
67 * Revision 1.10 2004/05/17 15:50:54 j_novak
68 * Removed unused nr variables
69 *
70 * Revision 1.9 2004/05/17 15:42:22 j_novak
71 * The method 1 of vector Poisson eq. solves directly for F^r.
72 * Some bugs were corrected in the operator poisson_vect.
73 *
74 * Revision 1.8 2004/05/14 08:51:02 p_grandclement
75 * *** empty log message ***
76 *
77 * Revision 1.7 2004/05/10 15:28:22 j_novak
78 * First version of functions for the solution of the r-component of the
79 * vector Poisson equation.
80 *
81 * Revision 1.6 2004/03/05 09:18:49 p_grandclement
82 * Addition of operator sec_order_r2
83 *
84 * Revision 1.5 2004/01/15 09:15:39 p_grandclement
85 * Modification and addition of the Helmholtz operators
86 *
87 * Revision 1.4 2004/01/07 14:36:38 p_grandclement
88 * Modif mineure in Param_elliptic.set_variable
89 *
90 * Revision 1.3 2003/12/11 16:11:38 e_gourgoulhon
91 * Changed #include <iostream.h> to #include "headcpp.h".
92 *
93 * Revision 1.2 2003/12/11 15:57:27 p_grandclement
94 * include stdlib.h encore ...
95 *
96 * Revision 1.1 2003/12/11 14:48:51 p_grandclement
97 * Addition of ALL (and that is a lot !) the files needed for the general elliptic solver ... UNDER DEVELOPEMENT...
98 *
99 *
100 * $Header: /cvsroot/Lorene/C++/Source/Param_elliptic/param_elliptic.C,v 1.21 2014/10/13 08:53:37 j_novak Exp $
101 *
102 */
103
104#include "headcpp.h"
105
106#include <cmath>
107#include <cstdlib>
108
109#include "map.h"
110#include "ope_elementary.h"
111#include "param_elliptic.h"
112#include "base_val.h"
113#include "scalar.h"
114#include "proto.h"
115
116
117namespace Lorene {
119
120 switch (type_map) {
121 case MAP_AFF:
122 return *mp_af ;
123 break ;
124 case MAP_LOG:
125 return *mp_log ;
126 break ;
127 default:
128 cout << "Unknown mapping in Param_elliptic" << endl ;
129 abort() ;
130 return *mp_af ;
131 }
132}
133
134double Param_elliptic::get_alpha(int l) const {
135
136 switch (type_map) {
137 case MAP_AFF:
138 return mp_af->get_alpha()[l] ;
139 break ;
140 case MAP_LOG:
141 return mp_log->get_alpha(l) ;
142 break ;
143 default:
144 cout << "Unknown mapping in Param_elliptic" << endl ;
145 abort() ;
146 return 1 ;
147 }
148}
149
150double Param_elliptic::get_beta(int l) const {
151
152 switch (type_map) {
153 case MAP_AFF:
154 return mp_af->get_beta()[l] ;
155 break ;
156 case MAP_LOG:
157 return mp_log->get_beta(l) ;
158 break ;
159 default:
160 cout << "Unknown mapping in Param_elliptic" << endl ;
161 abort() ;
162 return 1 ;
163 }
164}
165
166int Param_elliptic::get_type(int l) const {
167
168 switch (type_map) {
169 case MAP_AFF:
170 return AFFINE ;
171 break ;
172 case MAP_LOG:
173 return mp_log->get_type(l) ;
174 break ;
175 default:
176 cout << "Unknown mapping in Param_elliptic" << endl ;
177 abort() ;
178 return 1 ;
179 }
180}
181
182// Construction (By default it is set to Poisson with appropriate dzpuis...)
183Param_elliptic::Param_elliptic(const Scalar& so) : var_F(so.get_mp()), var_G(so.get_mp()),
184 done_F (so.get_mp().get_mg()->get_nzone(),
185 so.get_mp().get_mg()->get_np(0) + 1,
186 so.get_mp().get_mg()->get_nt(0)),
187 done_G (so.get_mp().get_mg()->get_nzone()),
188 val_F_plus (so.get_mp().get_mg()->get_nzone(),
189 so.get_mp().get_mg()->get_np(0) + 1,
190 so.get_mp().get_mg()->get_nt(0)),
191 val_F_minus (so.get_mp().get_mg()->get_nzone(),
192 so.get_mp().get_mg()->get_np(0) + 1,
193 so.get_mp().get_mg()->get_nt(0)),
194 val_dF_plus (so.get_mp().get_mg()->get_nzone(),
195 so.get_mp().get_mg()->get_np(0) + 1,
196 so.get_mp().get_mg()->get_nt(0)),
197 val_dF_minus (so.get_mp().get_mg()->get_nzone(),
198 so.get_mp().get_mg()->get_np(0) + 1,
199 so.get_mp().get_mg()->get_nt(0)),
200 val_G_plus (so.get_mp().get_mg()->get_nzone()),
201 val_G_minus (so.get_mp().get_mg()->get_nzone()),
202 val_dG_plus (so.get_mp().get_mg()->get_nzone()),
203 val_dG_minus (so.get_mp().get_mg()->get_nzone())
204
205
206{
207
208 // On passe en Ylm
209 Scalar auxi(so) ;
210 auxi.set_spectral_va().coef() ;
211 auxi.set_spectral_va().ylm() ;
212
213 Base_val base (auxi.get_spectral_va().base) ;
214 int dzpuis = (auxi.dz_nonzero() ? auxi.get_dzpuis() : 4) ; //## to be modified??
215
216 // Right now, only applicable with affine mapping
217 const Map_af* map_affine = dynamic_cast <const Map_af*> (&so.get_mp()) ;
218 const Map_log* map_log = dynamic_cast <const Map_log*> (&so.get_mp()) ;
219
220
221 if ((map_affine == 0x0) && (map_log == 0x0)) {
222 cout << "Param_elliptic not yet defined on this type of mapping" << endl ;
223 abort() ;
224 }
225 else {
226
227 if (map_affine != 0x0) {
228 type_map = MAP_AFF ;
229 mp_af = map_affine ;
230 mp_log = 0x0 ;
231 }
232 if (map_log != 0x0) {
233 type_map = MAP_LOG ;
234 mp_af = 0x0 ;
235 mp_log = map_log ;
236 }
237 int nz = get_mp().get_mg()->get_nzone() ;
238 int nbr = 0 ;
239 for (int l=0 ; l<nz ; l++)
240 nbr += (get_mp().get_mg()->get_np(l)+1)*
241 get_mp().get_mg()->get_nt(l) ;
242
243 operateurs = new Ope_elementary* [nbr] ;
244
245 for (int l=0 ; l<nbr ; l++)
246 operateurs[l] = 0x0 ;
247
248 int nr ;
249 int base_r, m_quant, l_quant ;
250
251 int conte = 0 ;
252 for (int l=0 ; l<nz ; l++) {
253
254 nr = get_mp().get_mg()->get_nr(l) ;
255
256 for (int k=0 ; k<get_mp().get_mg()->get_np(l)+1 ; k++)
257 for (int j=0 ; j<get_mp().get_mg()->get_nt(l) ; j++) {
258
259 so.get_spectral_va().base.give_quant_numbers
260 (l, k, j, m_quant, l_quant, base_r) ;
261
262 switch (type_map) {
263 case MAP_AFF:
264 operateurs[conte] = new
265 Ope_poisson(nr, base_r, get_alpha(l),
266 get_beta(l), l_quant, dzpuis) ;
267 break ;
268 case MAP_LOG:
269 if (mp_log->get_type(l) == AFFINE)
270 operateurs[conte] = new
271 Ope_poisson(nr, base_r, get_alpha(l),
272 get_beta(l), l_quant, dzpuis) ;
273 else
274 operateurs[conte] = new
275 Ope_sec_order (nr, base_r, get_alpha(l), get_beta(l),
276 1., 2. , l_quant) ;
277 break ;
278 default :
279 cout << "Unknown mapping in Param_elliptic::Param_elliptic"
280 << endl ;
281
282 }
283 conte ++ ;
284 }
285 }
286 }
287
288 // STD VARIABLE CHANGE :
290 var_F.set_spectral_va().set_base (so.get_spectral_va().get_base()) ;
291
293 var_G = 1 ;
295
298
303
308}
309
311
312 int nbr = 0 ;
313 for (int l=0 ; l<get_mp().get_mg()->get_nzone() ; l++)
314 nbr += (get_mp().get_mg()->get_np(l)+1)*get_mp().get_mg()->get_nt(l) ;
315
316 for (int i=0 ; i<nbr ; i++)
317 if (operateurs[i] != 0x0)
318 delete operateurs[i] ;
319
320 delete [] operateurs ;
321}
322
324
325 int np, nt ;
326
327 int conte = 0 ;
328 for (int l=0 ; l<get_mp().get_mg()->get_nzone() ; l++) {
329
330 np = get_mp().get_mg()->get_np(l) ;
331 nt = get_mp().get_mg()->get_nt(l) ;
332
333 for (int k=0 ; k<np+1 ; k++)
334 for (int j=0 ; j<nt ; j++) {
335 if ((operateurs[conte] != 0x0) && (l==zone))
337 conte ++ ;
338 }
339 }
340}
341
342
344
345 source.set_spectral_va().coef() ;
346 source.set_spectral_va().ylm() ;
347
348 int nz = get_mp().get_mg()->get_nzone() ;
349 int nr ;
350 int conte = 0 ;
351 int m_quant, base_r_1d, l_quant ;
352
353 switch (type_map) {
354
355 case MAP_AFF:
356
357 for (int l=0 ; l<nz ; l++) {
358
359 nr = get_mp().get_mg()->get_nr(l) ;
360
361 for (int k=0 ; k<get_mp().get_mg()->get_np(l)+1 ; k++)
362 for (int j=0 ; j<get_mp().get_mg()->get_nt(l) ; j++) {
363 if (l==zone) {
364 if (operateurs[conte] != 0x0) {
365 delete operateurs[conte] ;
366 source.get_spectral_va().base.give_quant_numbers
367 (l, k, j, m_quant, l_quant, base_r_1d) ;
368 operateurs[conte] = new Ope_helmholtz_minus (nr, base_r_1d, l_quant, get_alpha(l),
369 get_beta(l) , masse) ;
370 }
371 }
372 conte ++ ;
373 }
374 }
375 break ;
376
377 case MAP_LOG :
378 if (mp_log->get_type(zone) != AFFINE) {
379 cout << "Operator not define with LOG mapping..." << endl ;
380 abort() ;
381 }
382 else {
383 for (int l=0 ; l<nz ; l++) {
384
385 nr = get_mp().get_mg()->get_nr(l) ;
386
387 for (int k=0 ; k<get_mp().get_mg()->get_np(l)+1 ; k++)
388 for (int j=0 ; j<get_mp().get_mg()->get_nt(l) ; j++) {
389 if (l==zone) {
390 if (operateurs[conte] != 0x0) {
391 delete operateurs[conte] ;
392 source.get_spectral_va().base.give_quant_numbers
393 (l, k, j, m_quant, l_quant, base_r_1d) ;
394 operateurs[conte] = new Ope_helmholtz_minus (nr, base_r_1d, l_quant,
395 get_alpha(l), get_beta(l), masse) ;
396 }
397 }
398 conte ++ ;
399 }
400 }
401 }
402 break ;
403
404 default :
405 cout << "Unkown mapping in set_helmhotz_minus" << endl ;
406 abort() ;
407 break ;
408 }
409}
410
412
413 source.set_spectral_va().coef() ;
414 source.set_spectral_va().ylm() ;
415
416 int nz = get_mp().get_mg()->get_nzone() ;
417 int nr ;
418 int conte = 0 ;
419 int m_quant, base_r_1d, l_quant ;
420
421 switch (type_map) {
422
423 case MAP_AFF:
424
425 for (int l=0 ; l<nz ; l++) {
426
427 nr = get_mp().get_mg()->get_nr(l) ;
428
429 for (int k=0 ; k<get_mp().get_mg()->get_np(l)+1 ; k++)
430 for (int j=0 ; j<get_mp().get_mg()->get_nt(l) ; j++) {
431 if (l==zone) {
432 if (operateurs[conte] != 0x0) {
434 // PROVISOIRE, DANS LE NOYAU SEUL LE CAS SPHERIQUE EST IMPLEMENTE
435 if ((old_base != R_CHEBI)) {
436 delete operateurs[conte] ;
437 source.get_spectral_va().base.give_quant_numbers
438 (l, k, j, m_quant, l_quant, base_r_1d) ;
439 operateurs[conte] = new Ope_helmholtz_plus (nr, base_r_1d, l_quant, get_alpha(l),
440 get_beta(l) , masse) ;
441 }
442 }
443 }
444 conte ++ ;
445 }
446 }
447 break ;
448
449 case MAP_LOG :
450 if (mp_log->get_type(zone) != AFFINE) {
451 cout << "Operator not define with LOG mapping..." << endl ;
452 abort() ;
453 }
454 else {
455 for (int l=0 ; l<nz ; l++) {
456
457 nr = get_mp().get_mg()->get_nr(l) ;
458
459 for (int k=0 ; k<get_mp().get_mg()->get_np(l)+1 ; k++)
460 for (int j=0 ; j<get_mp().get_mg()->get_nt(l) ; j++) {
461 if (l==zone) {
462 if (operateurs[conte] != 0x0) {
464 // PROVISOIRE, DANS LE NOYAU SEUL LE CAS SPHERIQUE EST IMPLEMENTE
465 if ((old_base != R_CHEBI)) {
466 delete operateurs[conte] ;
467 source.get_spectral_va().base.give_quant_numbers
468 (l, k, j, m_quant, l_quant, base_r_1d) ;
469 operateurs[conte] = new Ope_helmholtz_plus (nr, base_r_1d, l_quant,
470 get_alpha(l), get_beta(l), masse) ;
471 }
472 }
473 }
474 conte ++ ;
475 }
476 }
477 }
478 break ;
479
480 default :
481 cout << "Unkown mapping in set_helmhotz_minus" << endl ;
482 abort() ;
483 break ;
484 }
485}
486
488
489 if (type_map != MAP_AFF) {
490 cout << "set_poisson_vect_r only defined for an affine mapping..." << endl ;
491 abort() ;
492 }
493 else {
494
495 int nz = get_mp().get_mg()->get_nzone() ;
496
497 int nr ;
498
499 int conte = 0 ;
500 for (int l=0 ; l<nz ; l++) {
501
502 nr = get_mp().get_mg()->get_nr(l) ;
503
504 for (int k=0 ; k<get_mp().get_mg()->get_np(l)+1 ; k++)
505 for (int j=0 ; j<get_mp().get_mg()->get_nt(l) ; j++) {
506 if ((operateurs[conte] != 0x0) && (l==zone)) {
509 dynamic_cast<Ope_poisson*>(operateurs[conte]) ;
510 assert (op_pois !=0x0) ;
511 int lq_old = op_pois->get_lquant() ;
512 int dzp = op_pois->get_dzpuis() ;
513
514 delete operateurs[conte] ;
515 if (type_map == MAP_AFF) {
516 if ((!only_l_zero)||(lq_old == 0)) {
517 operateurs[conte] = new Ope_pois_vect_r(nr, old_base,get_alpha(l),
518 get_beta(l), lq_old, dzp) ;
519 }
520 else {
521 operateurs[conte] = 0x0 ;
522 }
523 }
524 else
525 operateurs[conte] = 0x0 ;
526 }
527 conte ++ ;
528 }
529 }
530 }
531}
532
534
535 int nz = get_mp().get_mg()->get_nzone() ;
536
537 int conte = 0 ;
538 for (int l=0 ; l<nz ; l++) {
539
540 bool ced = (get_mp().get_mg()->get_type_r(l) == UNSURR ) ;
541 for (int k=0 ; k<get_mp().get_mg()->get_np(l)+1 ; k++)
542 for (int j=0 ; j<get_mp().get_mg()->get_nt(l) ; j++) {
543 if ((operateurs[conte] != 0x0) && (l==zone)) {
545 dynamic_cast<Ope_poisson*>(operateurs[conte]) ;
546 assert (op_pois !=0x0) ;
547 int lq_old = op_pois->get_lquant() ;
548 if (lq_old == 0) {
549 delete operateurs[conte] ;
550 operateurs[conte] = 0x0 ;
551 }
552 else
553 ced ? op_pois->inc_l_quant() : op_pois->dec_l_quant() ;
554 }
555 conte ++ ;
556 }
557 }
558}
559
561
562 if (type_map != MAP_AFF) {
563 cout << "set_poisson_tens_rr only defined for an affine mapping..."
564 << endl ;
565 abort() ;
566 }
567 else {
568
569 int nz = get_mp().get_mg()->get_nzone() ;
570
571 int nr ;
572
573 int conte = 0 ;
574 for (int l=0 ; l<nz ; l++) {
575
576 nr = get_mp().get_mg()->get_nr(l) ;
577
578 for (int k=0 ; k<get_mp().get_mg()->get_np(l)+1 ; k++)
579 for (int j=0 ; j<get_mp().get_mg()->get_nt(l) ; j++) {
580 if ((operateurs[conte] != 0x0) && (l==zone)) {
583 dynamic_cast<Ope_poisson*>(operateurs[conte]) ;
584 assert (op_pois !=0x0) ;
585 int lq_old = op_pois->get_lquant() ;
586 int dzp = op_pois->get_dzpuis() ;
587
588 delete operateurs[conte] ;
589 if (lq_old >= 2) {
591 (nr, old_base,get_alpha(l), get_beta(l), lq_old, dzp) ;
592 }
593 else
594 operateurs[conte] = 0x0 ;
595 }
596 conte ++ ;
597 }
598 }
599 }
600}
601
602void Param_elliptic::set_sec_order_r2 (int zone, double a, double b, double c){
603
604
605 int nz = get_mp().get_mg()->get_nzone() ;
606 int nr ;
607 int conte = 0 ;
608
609 switch (type_map) {
610
611 case MAP_AFF:
612
613 for (int l=0 ; l<nz ; l++) {
614
615 nr = get_mp().get_mg()->get_nr(l) ;
616
617 for (int k=0 ; k<get_mp().get_mg()->get_np(l)+1 ; k++)
618 for (int j=0 ; j<get_mp().get_mg()->get_nt(l) ; j++) {
619 if ((operateurs[conte] != 0x0) && (l==zone)) {
621 // PROVISOIRE, DANS LE NOYAU SEUL LE CAS SPHERIQUE EST IMPLEMENTE
622 if ((old_base != R_CHEBI)) {
623 delete operateurs[conte] ;
625 get_alpha(l),
626 get_beta(l), a, b, c) ;
627 }
628 }
629 conte ++ ;
630 }
631 }
632 break ;
633
634 case MAP_LOG :
635 if (mp_log->get_type(zone) != AFFINE) {
636 cout << "Operator not define with LOG mapping..." << endl ;
637 abort() ;
638 }
639 else {
640 for (int l=0 ; l<nz ; l++) {
641
642 nr = get_mp().get_mg()->get_nr(l) ;
643
644 for (int k=0 ; k<get_mp().get_mg()->get_np(l)+1 ; k++)
645 for (int j=0 ; j<get_mp().get_mg()->get_nt(l) ; j++) {
646 if ((operateurs[conte] != 0x0) && (l==zone)) {
648 // PROVISOIRE, DANS LE NOYAU SEUL LE CAS SPHERIQUE EST IMPLEMENTE
649 if ((old_base != R_CHEBI)) {
650 delete operateurs[conte] ;
652 get_alpha(l),
653 get_beta(l), a, b, c) ;
654 }
655 }
656 conte ++ ;
657 }
658 }
659 }
660 break ;
661
662 default :
663 cout << "Unkown mapping in set_sec_order_r2" << endl ;
664 abort() ;
665 break ;
666 }
667}
668
669void Param_elliptic::set_sec_order (int zone, double a, double b, double c){
670
671 if ((type_map == MAP_AFF) || (mp_log->get_type(zone) == AFFINE)) {
672 cout << "set_sec_order only defined for a log mapping" << endl ;
673 abort() ;
674 }
675 else {
676
677 int nz = get_mp().get_mg()->get_nzone() ;
678
679 int nr ;
680
681 int conte = 0 ;
682 for (int l=0 ; l<nz ; l++) {
683
684 nr = get_mp().get_mg()->get_nr(l) ;
685
686 for (int k=0 ; k<get_mp().get_mg()->get_np(l)+1 ; k++)
687 for (int j=0 ; j<get_mp().get_mg()->get_nt(l) ; j++) {
688 if ((operateurs[conte] != 0x0) && (l==zone)) {
689
691 // PROVISOIRE, DANS LE NOYAU SEUL LE CAS SPHERIQUE EST IMPLEMENTE
692 if (old_base != R_CHEBI) {
693 delete operateurs[conte] ;
695 get_alpha(l),
696 get_beta(l), a, b, c) ;
697 }
698 }
699 conte ++ ;
700 }
701 }
702 }
703}
704
706
707 source.set_spectral_va().coef() ;
708 source.set_spectral_va().ylm() ;
709
710 int nz = get_mp().get_mg()->get_nzone() ;
711 int nr ;
712 int conte = 0 ;
713 int m_quant, base_r_1d, l_quant ;
714 int dzpuis = source.get_dzpuis() ;
715
716 switch (type_map) {
717
718 case MAP_AFF:
719
720 for (int l=0 ; l<nz ; l++) {
721 int dz = (l==nz-1) ? dzpuis : 0 ;
722 nr = get_mp().get_mg()->get_nr(l) ;
723
724 for (int k=0 ; k<get_mp().get_mg()->get_np(l)+1 ; k++)
725 for (int j=0 ; j<get_mp().get_mg()->get_nt(l) ; j++) {
726 if (l==zone) {
727 if (operateurs[conte] != 0x0) {
728 delete operateurs[conte] ;
729 source.get_spectral_va().base.give_quant_numbers
730 (l, k, j, m_quant, l_quant, base_r_1d) ;
731 operateurs[conte] = new Ope_vorton (nr, base_r_1d, get_alpha(l),
732 get_beta(l), l_quant, dz) ;
733 }
734 }
735 conte ++ ;
736 }
737 }
738 break ;
739
740 case MAP_LOG :
741 if (mp_log->get_type(zone) != AFFINE) {
742 cout << "Operator not define with LOG mapping..." << endl ;
743 abort() ;
744 }
745 else {
746 for (int l=0 ; l<nz ; l++) {
747 int dz = (l==nz-1) ? dzpuis : 0 ;
748 nr = get_mp().get_mg()->get_nr(l) ;
749
750 for (int k=0 ; k<get_mp().get_mg()->get_np(l)+1 ; k++)
751 for (int j=0 ; j<get_mp().get_mg()->get_nt(l) ; j++) {
752 if (l==zone) {
753 if (operateurs[conte] != 0x0) {
754 delete operateurs[conte] ;
755 source.get_spectral_va().base.give_quant_numbers
756 (l, k, j, m_quant, l_quant, base_r_1d) ;
758 get_alpha(l), get_beta(l), l_quant, dz) ;
759 }
760 }
761 conte ++ ;
762 }
763 }
764 }
765 break ;
766
767 default :
768 cout << "Unkown mapping in set_ope_vorton" << endl ;
769 abort() ;
770 break ;
771 }
772}
773
775
776 assert (so.get_etat() != ETATNONDEF) ;
777 assert (so.get_mp() == get_mp()) ;
778
779 var_F = so ;
781}
782
784
785 assert (so.get_etat() != ETATNONDEF) ;
786 assert (so.get_mp() == get_mp()) ;
787
788 var_G = so ;
790}
791}
Bases of the spectral expansions.
Definition base_val.h:322
Time evolution with partial storage (*** under development ***).
Definition evolution.h:371
void annule_hard()
Sets the Itbl to zero in a hard way.
Definition itbl.C:272
Affine radial mapping.
Definition map.h:2027
Logarithmic radial mapping.
Definition map.h:3583
Base class for pure radial mappings.
Definition map.h:1536
Basic class for elementary elliptic operators.
int get_base_r() const
Returns base_r}.
virtual void inc_l_quant()=0
Increases the quatum number l by one unit.
Class for the Helmholtz operator ( ).
Class for the Helmholtz operator (m > 0).
Class for the operator of the rr component of the divergence-free tensor Poisson equation.
Class for the operator of the r component of the vector Poisson equation.
Class for the operator of the Poisson equation (i.e.
Class for operator of the type .
Class for operator of the type .
Class for the operator appearing for the vortons.
Scalar var_G
Multiplicative variable change that must be sphericaly symetric !
Itbl done_F
Stores what has been computed for F.
Tbl val_dG_minus
Values of the derivative of G at the inner boundaries of the various domains.
void set_sec_order_r2(int zone, double a, double b, double c)
Set the operator to in one domain (only in the shells).
void inc_l_quant(int zone)
Increases the quantum number l in the domain zone .
void set_helmholtz_minus(int zone, double mas, Scalar &so)
Set the operator to in one domain (not in the nucleus).
Tbl val_dF_minus
Values of the derivative of F at the inner boundaries of the various domains.
void set_poisson_vect_eta(int zone)
Sets the operator to be a regular elliptic operator to solve for the component of the vector Poisson...
const Map_log * mp_log
The mapping if log type.
Param_elliptic(const Scalar &)
Standard constructor from a Scalar.
const Map_radial & get_mp() const
Returns the mapping.
Tbl val_F_minus
Values of F at the inner boundaries of the various domains.
Tbl val_F_plus
Values of F at the outer boundaries of the various domains.
void set_sec_order(int zone, double a, double b, double c)
Set the operator to in one domain (only in the shells).
Tbl val_G_minus
Values of G at the inner boundaries of the various domains.
Tbl val_dG_plus
Values of the derivative of G at the outer boundaries of the various domains.
Ope_elementary ** operateurs
Array on the elementary operators.
void set_poisson_vect_r(int zone, bool only_l_zero=false)
Sets the operator to in all domains, for ; and to in all domains otherwise.
const Map_af * mp_af
The mapping, if affine.
void set_poisson_tens_rr(int zone)
Sets the operator to in all domains, for only.
void set_variable_F(const Scalar &)
Changes the variable function F.
void set_helmholtz_plus(int zone, double mas, Scalar &so)
Set the operator to in one domain (only in the shells).
void set_variable_G(const Scalar &)
Changes the variable function G.
int type_map
Type of mapping either MAP_AFF or MAP_LOG.
Tbl val_dF_plus
Values of the derivative of F at the outer boundaries of the various domains.
Tbl val_G_plus
Values of G at the outer boundaries of the various domains.
void set_ope_vorton(int zone, Scalar &so)
Set the operator to in one domain (not implemented in the nucleus).
Itbl done_G
Stores what has been computed for G.
Scalar var_F
Additive variable change function.
Tensor field of valence 0 (or component of a tensorial field).
Definition scalar.h:387
virtual void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition scalar.C:353
virtual void std_spectral_base()
Sets the spectral bases of the Valeur va to the standard ones for a scalar field.
Definition scalar.C:784
Valeur & set_spectral_va()
Returns va (read/write version)
Definition scalar.h:604
void annule_hard()
Sets the Scalar to zero in a hard way.
Definition scalar.C:380
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition tbl.C:361
void set_base(const Base_val &)
Sets the bases for spectral expansions (member base )
Definition valeur.C:810
#define R_CHEBI
base de Cheb. impaire (rare) seulement
Lorene prototypes.
Definition app_hor.h:64