LORENE
map_star_fait.C
1/*
2 * This file is part of LORENE.
3 *
4 * LORENE is free software; you can redistribute it and/or modify
5 * it under the terms of the GNU General Public License as published by
6 * the Free Software Foundation; either version 2 of the License, or
7 * (at your option) any later version.
8 *
9 * LORENE is distributed in the hope that it will be useful,
10 * but WITHOUT ANY WARRANTY; without even the implied warranty of
11 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12 * GNU General Public License for more details.
13 *
14 * You should have received a copy of the GNU General Public License
15 * along with LORENE; if not, write to the Free Software
16 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 *
18 */
19
20
21// Includes
22#include <cassert>
23#include <cstdlib>
24#include <cmath>
25
26#include "mtbl.h"
27#include "map.h"
28#include "proto.h"
29#include "valeur.h"
30
31
32
33
34
35
36
37
38 //----------------//
39 // Coord. radiale //
40 //----------------//
41
42namespace Lorene {
43Mtbl* map_star_fait_r(const Map* cvi) {
44
45 // recup du changement de variable
46 const Map_star* cv = static_cast<const Map_star*>(cvi) ;
47 const Mg3d* mg = cv->get_mg() ;
48 int nz = mg->get_nzone() ;
49
50 // Le resultat
51 Mtbl* mti = new Mtbl(mg) ;
52 mti->set_etat_qcq() ;
53
54 // Pour le confort
55 const Valeur& alpha = cv->get_alpha() ;
56 const Valeur& beta = cv->get_beta() ;
57
58 int i, j, k ;
59 for (int l=0 ; l<nz ; l++) {
60 int ir = mg->get_nr(l);
61 int it = mg->get_nt(l) ;
62 int ip = mg->get_np(l) ;
63 const Grille3d* g = mg->get_grille3d(l) ;
64 Tbl* tb = (mti->t)[l] ;
65 tb->set_etat_qcq() ;
66 double* p_r = tb->t ;
67
68 switch(mg->get_type_r(l)) {
69 case RARE: case FIN:
70 for (k=0 ; k<ip ; k++) {
71 for (j=0 ; j<it ; j++) {
72 for (i=0 ; i<ir ; i++) {
73 *p_r = alpha(l,k,j,0) * (g->x)[i] + beta(l,k,j,0) ;
74 p_r++ ;
75 } // Fin de boucle sur r
76 } // Fin de boucle sur theta
77 } // Fin de boucle sur phi
78 break ;
79
80 // case FIN:{
81 // cout << "map_star_fait_r: Shells not implemented yet..." << endl;
82 // abort() ;
83 // break ;
84 // }
85
86 case UNSURR:
87 for (k=0 ; k<ip ; k++) {
88 for (j=0 ; j<it ; j++) {
89 for (i=0 ; i<ir ; i++) {
90 cout << "map_star_fait_r: Warning ! No compactified zone allowed !" << endl;
91 abort() ;
92 } // Fin de boucle sur r
93 } // Fin de boucle sur theta
94 } // Fin de boucle sur phi
95 break ;
96
97 default:
98 cout << "map_star_fait_r: unknown type_r !\n" ;
99 abort () ;
100 exit(-1) ;
101
102 } // Fin du switch
103 } // Fin de boucle sur zone
104
105 // Termine
106 return mti ;
107}
108
109 //--------------//
110 // Coord. Theta //
111 //--------------//
112
113Mtbl* map_star_fait_tet(const Map* cvi) {
114
115 // recup du changement de variable
116 const Map_star* cv = static_cast<const Map_star*>(cvi) ;
117 const Mg3d* mg = cv->get_mg() ;
118 int nz = mg->get_nzone() ;
119
120 // Le resultat
121 Mtbl* mti = new Mtbl(mg) ;
122 mti->set_etat_qcq() ;
123
124 int i, j, k ;
125 for (int l=0 ; l<nz ; l++) {
126 int ir = mg->get_nr(l);
127 int it = mg->get_nt(l);
128 int ip = mg->get_np(l);
129 const Grille3d* g = mg->get_grille3d(l) ;
130 Tbl* tb = (mti->t)[l] ;
131 tb->set_etat_qcq() ;
132 double* p_r = tb->t ;
133 for (k=0 ; k<ip ; k++) {
134 for (j=0 ; j<it ; j++) {
135 for (i=0 ; i<ir ; i++) {
136 *p_r = (g->tet)[j] ;
137 p_r++ ;
138 } // Fin de boucle sur r
139 } // Fin de boucle sur theta
140 } // Fin de boucle sur phi
141 } // Fin de boucle sur zone
142
143 // Termine
144 return mti ;
145}
146
147 //------------//
148 // Coord. Phi //
149 //------------//
150
151Mtbl* map_star_fait_phi(const Map* cvi) {
152
153 // recup du changement de variable
154 const Map_star* cv = static_cast<const Map_star*>(cvi) ;
155 const Mg3d* mg = cv->get_mg() ;
156 int nz = mg->get_nzone() ;
157
158 // Le resultat
159 Mtbl* mti = new Mtbl(mg) ;
160 mti->set_etat_qcq() ;
161
162 int i, j, k ;
163 for (int l=0 ; l<nz ; l++) {
164 int ir = mg->get_nr(l);
165 int it = mg->get_nt(l);
166 int ip = mg->get_np(l);
167 const Grille3d* g = mg->get_grille3d(l) ;
168 Tbl* tb = (mti->t)[l] ;
169 tb->set_etat_qcq() ;
170 double* p_r = tb->t ;
171 for (k=0 ; k<ip ; k++) {
172 for (j=0 ; j<it ; j++) {
173 for (i=0 ; i<ir ; i++) {
174 *p_r = (g->phi)[k] ;
175 p_r++ ;
176 } // Fin de boucle sur r
177 } // Fin de boucle sur theta
178 } // Fin de boucle sur phi
179 } // Fin de boucle sur zone
180
181 // Termine
182 return mti ;
183}
184
185 //----------//
186 // Coord. X //
187 //----------//
188
189Mtbl* map_star_fait_x(const Map* cvi) {
190
191 // recup de la grille
192 const Mg3d* mg = cvi->get_mg() ;
193
194 // Le resultat
195 Mtbl* mti = new Mtbl(mg) ;
196
197 *mti = (cvi->r) * (cvi->sint) * (cvi->cosp) ;
198
199 // Termine
200 return mti ;
201}
202
203 //----------//
204 // Coord. Y //
205 //----------//
206
207Mtbl* map_star_fait_y(const Map* cvi) {
208
209 // recup de la grille
210 const Mg3d* mg = cvi->get_mg() ;
211
212 // Le resultat
213 Mtbl* mti = new Mtbl(mg) ;
214
215 *mti = (cvi->r) * (cvi->sint) * (cvi->sinp) ;
216
217 // Termine
218 return mti ;
219}
220
221 //----------//
222 // Coord. Z //
223 //----------//
224
225Mtbl* map_star_fait_z(const Map* cvi) {
226
227 // recup de la grille
228 const Mg3d* mg = cvi->get_mg() ;
229
230 // Le resultat
231 Mtbl* mti = new Mtbl(mg) ;
232
233 *mti = (cvi->r) * (cvi->cost) ;
234
235 // Termine
236 return mti ;
237}
238
239 //--------------------//
240 // Coord. X "absolue" //
241 //--------------------//
242
243Mtbl* map_star_fait_xa(const Map* cvi) {
244
245 // recup de la grille
246 const Mg3d* mg = cvi->get_mg() ;
247
248 // Le resultat
249 Mtbl* mti = new Mtbl(mg) ;
250
251 double r_phi = cvi->get_rot_phi() ;
252 double t_x = cvi->get_ori_x() ;
253
254 if ( fabs(r_phi) < 1.e-14 ) {
255 *mti = (cvi->x) + t_x ;
256 }
257 else if ( fabs(r_phi - M_PI) < 1.e-14 ) {
258 *mti = - (cvi->x) + t_x ;
259 }
260 else {
261 Mtbl phi = cvi->phi + r_phi ;
262 *mti = (cvi->r) * (cvi->sint) * cos(phi) + t_x ;
263 }
264
265 // Termine
266 return mti ;
267}
268
269 //--------------------//
270 // Coord. Y "absolue" //
271 //--------------------//
272
273Mtbl* map_star_fait_ya(const Map* cvi) {
274
275 // recup de la grille
276 const Mg3d* mg = cvi->get_mg() ;
277
278 // Le resultat
279 Mtbl* mti = new Mtbl(mg) ;
280
281 double r_phi = cvi->get_rot_phi() ;
282 double t_y = cvi->get_ori_y() ;
283
284 if ( fabs(r_phi) < 1.e-14 ) {
285 *mti = (cvi->y) + t_y ;
286 }
287 else if ( fabs(r_phi - M_PI) < 1.e-14 ) {
288 *mti = - (cvi->y) + t_y ;
289 }
290 else {
291 Mtbl phi = cvi->phi + r_phi ;
292 *mti = (cvi->r) * (cvi->sint) * sin(phi) + t_y ;
293 }
294
295 // Termine
296 return mti ;
297}
298
299 //--------------------//
300 // Coord. Z "absolue" //
301 //--------------------//
302
303Mtbl* map_star_fait_za(const Map* cvi) {
304
305 // recup de la grille
306 const Mg3d* mg = cvi->get_mg() ;
307
308 double t_z = cvi->get_ori_z() ;
309
310 // Le resultat
311 Mtbl* mti = new Mtbl(mg) ;
312
313 *mti = (cvi->r) * (cvi->cost) + t_z ;
314
315 // Termine
316 return mti ;
317}
318
319 //---------------//
320 // Trigonometrie //
321 //---------------//
322
323Mtbl* map_star_fait_sint(const Map* cvi) {
324
325 // recup du changement de variable
326 const Map_star* cv = static_cast<const Map_star*>(cvi) ;
327 const Mg3d* mg = cv->get_mg() ;
328 int nz = mg->get_nzone() ;
329
330 // Le resultat
331 Mtbl* mti = new Mtbl(mg) ;
332 mti->set_etat_qcq() ;
333
334 int i, j, k ;
335 for (int l=0 ; l<nz ; l++) {
336 int ir = mg->get_nr(l);
337 int it = mg->get_nt(l);
338 int ip = mg->get_np(l);
339 const Grille3d* g = mg->get_grille3d(l) ;
340 Tbl* tb = (mti->t)[l] ;
341 tb->set_etat_qcq() ;
342 double* p_r = tb->t ;
343 for (k=0 ; k<ip ; k++) {
344 for (j=0 ; j<it ; j++) {
345 for (i=0 ; i<ir ; i++) {
346 *p_r = sin(g->tet[j]) ;
347 p_r++ ;
348 } // Fin de boucle sur r
349 } // Fin de boucle sur theta
350 } // Fin de boucle sur phi
351 } // Fin de boucle sur zone
352
353 // Termine
354 return mti ;
355}
356
357Mtbl* map_star_fait_cost(const Map* cvi) {
358
359 // recup du changement de variable
360 const Map_star* cv = static_cast<const Map_star*>(cvi) ;
361 const Mg3d* mg = cv->get_mg() ;
362 int nz = mg->get_nzone() ;
363
364 // Le resultat
365 Mtbl* mti = new Mtbl(mg) ;
366 mti->set_etat_qcq() ;
367
368 int i, j, k ;
369 for (int l=0 ; l<nz ; l++) {
370 int ir = mg->get_nr(l);
371 int it = mg->get_nt(l);
372 int ip = mg->get_np(l);
373 const Grille3d* g = mg->get_grille3d(l) ;
374 Tbl* tb = (mti->t)[l] ;
375 tb->set_etat_qcq() ;
376 double* p_r = tb->t ;
377 for (k=0 ; k<ip ; k++) {
378 for (j=0 ; j<it ; j++) {
379 for (i=0 ; i<ir ; i++) {
380 *p_r = cos(g->tet[j]) ;
381 p_r++ ;
382 } // Fin de boucle sur r
383 } // Fin de boucle sur theta
384 } // Fin de boucle sur phi
385 } // Fin de boucle sur zone
386
387 // Termine
388 return mti ;
389}
390
391Mtbl* map_star_fait_sinp(const Map* cvi) {
392
393 // recup du changement de variable
394 const Map_star* cv = static_cast<const Map_star*>(cvi) ;
395 const Mg3d* mg = cv->get_mg() ;
396 int nz = mg->get_nzone() ;
397
398 // Le resultat
399 Mtbl* mti = new Mtbl(mg) ;
400 mti->set_etat_qcq() ;
401
402 int i, j, k ;
403 for (int l=0 ; l<nz ; l++) {
404 int ir = mg->get_nr(l);
405 int it = mg->get_nt(l);
406 int ip = mg->get_np(l);
407 const Grille3d* g = mg->get_grille3d(l) ;
408 Tbl* tb = (mti->t)[l] ;
409 tb->set_etat_qcq() ;
410 double* p_r = tb->t ;
411 for (k=0 ; k<ip ; k++) {
412 for (j=0 ; j<it ; j++) {
413 for (i=0 ; i<ir ; i++) {
414 *p_r = sin(g->phi[k]) ;
415 p_r++ ;
416 } // Fin de boucle sur r
417 } // Fin de boucle sur theta
418 } // Fin de boucle sur phi
419 } // Fin de boucle sur zone
420
421 // Termine
422 return mti ;
423}
424
425Mtbl* map_star_fait_cosp(const Map* cvi) {
426
427 // recup du changement de variable
428 const Map_star* cv = static_cast<const Map_star*>(cvi) ;
429 const Mg3d* mg = cv->get_mg() ;
430 int nz = mg->get_nzone() ;
431
432 // Le resultat
433 Mtbl* mti = new Mtbl(mg) ;
434 mti->set_etat_qcq() ;
435
436 int i, j, k ;
437 for (int l=0 ; l<nz ; l++) {
438 int ir = mg->get_nr(l);
439 int it = mg->get_nt(l);
440 int ip = mg->get_np(l);
441 const Grille3d* g = mg->get_grille3d(l) ;
442 Tbl* tb = (mti->t)[l] ;
443 tb->set_etat_qcq() ;
444 double* p_r = tb->t ;
445 for (k=0 ; k<ip ; k++) {
446 for (j=0 ; j<it ; j++) {
447 for (i=0 ; i<ir ; i++) {
448 *p_r = cos(g->phi[k]) ;
449 p_r++ ;
450 } // Fin de boucle sur r
451 } // Fin de boucle sur theta
452 } // Fin de boucle sur phi
453 } // Fin de boucle sur zone
454
455 // Termine
456 return mti ;
457}
458
459/*
460 ************************************************************************
461 * x/R dans le noyau, 1/R dans les coquilles, (x-1)/U dans la ZEC
462 ************************************************************************
463 */
464
465Mtbl* map_star_fait_xsr(const Map* cvi) {
466
467 // recup du changement de variable
468 const Map_star* cv = static_cast<const Map_star*>(cvi) ;
469 const Mg3d* mg = cv->get_mg() ;
470 int nz = mg->get_nzone() ;
471
472 // Le resultat
473 Mtbl* mti = new Mtbl(mg) ;
474 mti->set_etat_qcq() ;
475
476 // Pour le confort
477 const Valeur& alpha = cv->get_alpha() ;
478 const Valeur& beta = cv->get_beta() ;
479
480 int i, j, k ;
481 for (int l=0 ; l<nz ; l++) {
482 int ir = mg->get_nr(l);
483 int it = mg->get_nt(l) ;
484 int ip = mg->get_np(l) ;
485 const Grille3d* g = mg->get_grille3d(l) ;
486 Tbl* tb = (mti->t)[l] ;
487 tb->set_etat_qcq() ;
488 double* p_r = tb->t ;
489
490 switch(mg->get_type_r(l)) {
491
492 case RARE:
493 assert(beta(l).get_etat()==ETATZERO) ;
494 for (k=0 ; k<ip ; k++) {
495 for (j=0 ; j<it ; j++) {
496 for (i=0 ; i<ir ; i++) {
497 *p_r = 1. / alpha(l,k,j,0) ;
498 p_r++ ;
499 } // Fin de boucle sur r
500 } // Fin de boucle sur theta
501 } // Fin de boucle sur phi
502 break ;
503
504 case FIN:
505 for (k=0 ; k<ip ; k++) {
506 for (j=0 ; j<it ; j++) {
507 if (ir == 1) { //Some hack for angular grid case...
508 *p_r = 1. / beta(l,k,j,0) ;
509 p_r++ ;
510 }
511 else
512 for (i=0 ; i<ir ; i++) {
513 *p_r = 1. / ( alpha(l,k,j,0) * (g->x)[i] + beta(l,k,j,0) ) ;
514 p_r++ ;
515 } // Fin de boucle sur r
516 } // Fin de boucle sur theta
517 } // Fin de boucle sur phi
518 break ;
519
520 case UNSURR:
521 cout << "map_star_fait_xsr: Compactified domain not allowed !" << endl;
522 abort() ;
523 break ;
524
525 default:
526 cout << "map_star_fait_xsr: unknown type_r !" << endl ;
527 abort() ;
528
529 } // Fin du switch
530 } // Fin de boucle sur zone
531
532 // Termine
533 return mti ;
534
535}
536
537/*
538 ************************************************************************
539 * 1/(dR/dx) ( -1/(dU/dx) ds la ZEC )
540 ************************************************************************
541 */
542
543Mtbl* map_star_fait_dxdr(const Map* cvi) {
544
545 // recup du changement de variable
546 const Map_star* cv = static_cast<const Map_star*>(cvi) ;
547 const Mg3d* mg = cv->get_mg() ;
548 int nz = mg->get_nzone() ;
549
550 // Le resultat
551 Mtbl* mti = new Mtbl(mg) ;
552 mti->set_etat_qcq() ;
553
554 // Pour le confort
555 const Valeur& alpha = cv->get_alpha() ;
556
557 int i, j, k ;
558 for (int l=0 ; l<nz ; l++) {
559 int ir = mg->get_nr(l);
560 int it = mg->get_nt(l) ;
561 int ip = mg->get_np(l) ;
562 Tbl* tb = (mti->t)[l] ;
563 tb->set_etat_qcq() ;
564 double* p_r = tb->t ;
565
566 switch(mg->get_type_r(l)) {
567
568 case RARE: case FIN:
569 for (k=0 ; k<ip ; k++) {
570 for (j=0 ; j<it ; j++) {
571 for (i=0 ; i<ir ; i++) {
572 *p_r = 1. / alpha(l,k,j,0) ;
573 p_r++ ;
574 } // Fin de boucle sur r
575 } // Fin de boucle sur theta
576 } // Fin de boucle sur phi
577 break ;
578
579 // case FIN:{
580 // cout << "map_star_fait_dxdr: Shells not implemented yet..." << endl;
581 // abort() ;
582 // break ;
583 // }
584
585 case UNSURR:
586 cout << "map_star_fait_dxdr: Compactified domain not allowed !" << endl;
587 abort() ;
588 break ;
589
590
591 default:
592 cout << "map_star_fait_dxdr: unknown type_r !" << endl ;
593 abort() ;
594 } // Fin du switch
595 } // Fin de boucle sur zone
596
597 // Termine
598 return mti ;
599}
600
601/*
602 ************************************************************************
603 * dR/dtheta
604 ************************************************************************
605 */
606
607Mtbl* map_star_fait_drdt(const Map* cvi) {
608
609 // recup du changement de variable
610 const Map_star* cv = static_cast<const Map_star*>(cvi) ;
611 const Mg3d* mg = cv->get_mg() ;
612 int nz = mg->get_nzone() ;
613
614 // Le resultat
615 Mtbl* mti = new Mtbl(mg) ;
616 mti->set_etat_qcq() ;
617
618 // Pour le confort
619 const Valeur& alpha = cv->get_alpha() ;
620 const Valeur& dalphadt = alpha.dsdt() ;
621 const Valeur& beta = cv->get_beta() ;
622 const Valeur& dbetadt = beta.dsdt() ;
623
624
625 for (int l=0 ; l<nz ; l++) {
626
627 int nr = mg->get_nr(l);
628 int nt = mg->get_nt(l) ;
629 int np = mg->get_np(l) ;
630 const Grille3d* g = mg->get_grille3d(l) ;
631 Tbl* tb = (mti->t)[l] ;
632 tb->set_etat_qcq() ;
633 double* p_r = tb->t ;
634
635 switch(mg->get_type_r(l)) {
636
637
638 case RARE: case FIN: {
639 for (int k=0 ; k<np ; k++) {
640 for (int j=0 ; j<nt ; j++) {
641 for (int i=0 ; i<nr ; i++) {
642 *p_r = dalphadt(l,k,j,0) * (g->x)[i] + dbetadt(l,k,j,0) ;
643 p_r++ ;
644 } // Fin de boucle sur r
645 } // Fin de boucle sur theta
646 } // Fin de boucle sur phi
647 break ;
648 }
649 // case FIN:{
650 // cout << "map_star_fait_drdt: Shells not implemented yet..." << endl;
651 // abort() ;
652 // break ;
653 // }
654
655 case UNSURR: {
656 cout << "map_star_fait_drdt: Compactified domain not allowed !" << endl;
657 abort() ;
658 break ;
659 }
660
661 default: {
662 cout << "map_star_fait_drdt: unknown type_r !" << endl ;
663 abort() ;
664 }
665 } // Fin du switch
666 } // Fin de boucle sur zone
667
668 // Termine
669
670 return mti ;
671}
672
673/*
674 ************************************************************************
675 * 1/sin(theta) dR/dphi
676 ************************************************************************
677 */
678
679Mtbl* map_star_fait_stdrdp(const Map* cvi) {
680
681 // recup du changement de variable
682 const Map_star* cv = static_cast<const Map_star*>(cvi) ;
683 const Mg3d* mg = cv->get_mg() ;
684 int nz = mg->get_nzone() ;
685
686 // Le resultat
687 Mtbl* mti = new Mtbl(mg) ;
688 mti->set_etat_qcq() ;
689
690 // Pour le confort
691 const Valeur& alpha = cv->get_alpha() ;
692 const Valeur& stalphadp = alpha.stdsdp() ;
693 const Valeur& beta = cv->get_beta() ;
694 const Valeur& stbetadp = beta.stdsdp() ;
695
696
697 for (int l=0 ; l<nz ; l++) {
698
699 int nr = mg->get_nr(l);
700 int nt = mg->get_nt(l) ;
701 int np = mg->get_np(l) ;
702 const Grille3d* g = mg->get_grille3d(l) ;
703 Tbl* tb = (mti->t)[l] ;
704 tb->set_etat_qcq() ;
705 double* p_r = tb->t ;
706
707 switch(mg->get_type_r(l)) {
708
709
710 case RARE: case FIN: {
711 for (int k=0 ; k<np ; k++) {
712 for (int j=0 ; j<nt ; j++) {
713 for (int i=0 ; i<nr ; i++) {
714 *p_r = stalphadp(l,k,j,0) * (g->x)[i] + stbetadp(l,k,j,0) ;
715 p_r++ ;
716 } // Fin de boucle sur r
717 } // Fin de boucle sur theta
718 } // Fin de boucle sur phi
719 break ;
720 }
721
722 // case FIN:{
723 // cout << "map_star_fait_stdrdp: Shells not implemented yet..." << endl;
724 // abort() ;
725 // break ;
726 // }
727
728 case UNSURR: {
729 cout << "map_star_fait_stdrdp: Compactified domain not allowed !" << endl;
730 abort() ;
731 break ;
732 }
733
734 default: {
735 cout << "map_star_fait_stdrdp: unknown type_r !" << endl ;
736 abort() ;
737 }
738 } // Fin du switch
739 } // Fin de boucle sur zone
740
741 // Termine
742
743 return mti ;
744}
745
746/*
747 ************************************************************************
748 * 1/R dR/dtheta
749 ************************************************************************
750 */
751
752Mtbl* map_star_fait_srdrdt(const Map* cvi) {
753
754 // recup du changement de variable
755 const Map_star* cv = static_cast<const Map_star*>(cvi) ;
756 const Mg3d* mg = cv->get_mg() ;
757 int nz = mg->get_nzone() ;
758
759 // Le resultat
760 Mtbl* mti = new Mtbl(mg) ;
761 mti->set_etat_qcq() ;
762
763 // Pour le confort
764 const Valeur& alpha = cv->get_alpha() ;
765 const Valeur& dalphadt = alpha.dsdt() ;
766 const Valeur& beta = cv->get_beta() ;
767 const Valeur& dbetadt = beta.dsdt() ;
768
769
770 for (int l=0 ; l<nz ; l++) {
771
772 int nr = mg->get_nr(l);
773 int nt = mg->get_nt(l) ;
774 int np = mg->get_np(l) ;
775 const Grille3d* g = mg->get_grille3d(l) ;
776 Tbl* tb = (mti->t)[l] ;
777 tb->set_etat_qcq() ;
778 double* p_r = tb->t ;
779
780 switch(mg->get_type_r(l)) {
781
782
783 case RARE : {
784 for (int k=0 ; k<np ; k++) {
785 for (int j=0 ; j<nt ; j++) {
786 for (int i=0 ; i<nr ; i++) {
787 *p_r = dalphadt(l,k,j,0) / alpha(l,k,j,0) ;
788 p_r++ ;
789 } // Fin de boucle sur r
790 } // Fin de boucle sur theta
791 } // Fin de boucle sur phi
792 break ;
793 }
794
795 case FIN:{
796 for (int k=0 ; k<np ; k++) {
797 for (int j=0 ; j<nt ; j++) {
798 for (int i=0 ; i<nr ; i++) {
799 *p_r = (dalphadt(l,k,j,0)*(g->x)[i] + dbetadt(l,k,j,0)) / (alpha(l,k,j,0)*(g->x)[i] + beta(l,k,j,0)) ;
800 p_r++ ;
801 } // Fin de boucle sur r
802 } // Fin de boucle sur theta
803 } // Fin de boucle sur phi
804 break ;
805 }
806
807 case UNSURR: {
808 cout << "map_star_fait_srdrdt: Compactified domain not allowed !" << endl;
809 abort() ;
810 break ;
811 }
812
813 default: {
814 cout << "map_star_fait_srdrdt: unknown type_r !" << endl ;
815 abort() ;
816 }
817 } // Fin du switch
818 } // Fin de boucle sur zone
819
820 // Termine
821
822 return mti ;
823}
824
825/*
826 ************************************************************************
827 * 1/(R sin(theta)) dR/dphi
828 ************************************************************************
829 */
830
831Mtbl* map_star_fait_srstdrdp(const Map* cvi) {
832
833 // recup du changement de variable
834 const Map_star* cv = static_cast<const Map_star*>(cvi) ;
835 const Mg3d* mg = cv->get_mg() ;
836 int nz = mg->get_nzone() ;
837
838 // Le resultat
839 Mtbl* mti = new Mtbl(mg) ;
840 mti->set_etat_qcq() ;
841
842 // Pour le confort
843 const Valeur& alpha = cv->get_alpha() ;
844 const Valeur& stalphadp = alpha.stdsdp() ;
845 const Valeur& beta = cv->get_beta() ;
846 const Valeur& stbetadp = beta.stdsdp() ;
847
848
849 for (int l=0 ; l<nz ; l++) {
850
851 int nr = mg->get_nr(l);
852 int nt = mg->get_nt(l) ;
853 int np = mg->get_np(l) ;
854 const Grille3d* g = mg->get_grille3d(l) ;
855 Tbl* tb = (mti->t)[l] ;
856 tb->set_etat_qcq() ;
857 double* p_r = tb->t ;
858
859 switch(mg->get_type_r(l)) {
860
861
862 case RARE : {
863 for (int k=0 ; k<np ; k++) {
864 for (int j=0 ; j<nt ; j++) {
865 for (int i=0 ; i<nr ; i++) {
866 *p_r = stalphadp(l,k,j,0) / alpha(l,k,j,0) ;
867 p_r++ ;
868 } // Fin de boucle sur r
869 } // Fin de boucle sur theta
870 } // Fin de boucle sur phi
871 break ;
872 }
873
874 case FIN:{
875 for (int k=0 ; k<np ; k++) {
876 for (int j=0 ; j<nt ; j++) {
877 for (int i=0 ; i<nr ; i++) {
878 *p_r = (stalphadp(l,k,j,0)*(g->x)[i] + stbetadp(l,k,j,0)) / (alpha(l,k,j,0)*(g->x)[i] + beta(l,k,j,0)) ;
879 p_r++ ;
880 } // Fin de boucle sur r
881 } // Fin de boucle sur theta
882 } // Fin de boucle sur phi
883 break ;
884 }
885
886 case UNSURR: {
887 cout << "map_star_fait_srstdrdp: Compactified domain not allowed !" << endl;
888 abort() ;
889 break ;
890 }
891
892 default: {
893 cout << "map_star_fait_srstdrdp: unknown type_r !" << endl ;
894 abort() ;
895 }
896 } // Fin du switch
897 } // Fin de boucle sur zone
898
899 // Termine
900
901 return mti ;
902}
903
904/*
905 ************************************************************************
906 * 1/R^2 dR/dtheta
907 ************************************************************************
908 */
909
910Mtbl* map_star_fait_sr2drdt(const Map* cvi) {
911
912 // recup du changement de variable
913 const Map_star* cv = static_cast<const Map_star*>(cvi) ;
914 const Mg3d* mg = cv->get_mg() ;
915 int nz = mg->get_nzone() ;
916
917 // Le resultat
918 Mtbl* mti = new Mtbl(mg) ;
919 mti->set_etat_qcq() ;
920
921 // Pour le confort
922 const Valeur& alpha = cv->get_alpha() ;
923 const Valeur& dalphadt = alpha.dsdt() ;
924
925
926 for (int l=0 ; l<nz ; l++) {
927
928 int nr = mg->get_nr(l);
929 int nt = mg->get_nt(l) ;
930 int np = mg->get_np(l) ;
931 const Grille3d* g = mg->get_grille3d(l) ;
932 Tbl* tb = (mti->t)[l] ;
933 tb->set_etat_qcq() ;
934 double* p_r = tb->t ;
935
936 switch(mg->get_type_r(l)) {
937
938
939 case RARE : {
940 for (int k=0 ; k<np ; k++) {
941 for (int j=0 ; j<nt ; j++) {
942 for (int i=0 ; i<nr ; i++) {
943 double ww = alpha(l,k,j,0) ;
944 *p_r = dalphadt(l,k,j,0) / (ww*ww * (g->x)[i]) ;
945 p_r++ ;
946 } // Fin de boucle sur r
947 } // Fin de boucle sur theta
948 } // Fin de boucle sur phi
949 break ;
950 }
951
952 case FIN:{
953 cout << "map_star_fait_sr2drdt: Shells not implemented yet..." << endl;
954 abort() ;
955 break ;
956 }
957
958 case UNSURR: {
959 cout << "map_star_fait_sr2drdt: Compactified domain not allowed !" << endl;
960 abort() ;
961 break ;
962 }
963
964 default: {
965 cout << "map_star_fait_sr2drdt: unknown type_r !" << endl ;
966 abort() ;
967 }
968 } // Fin du switch
969 } // Fin de boucle sur zone
970
971 // Termine
972
973 return mti ;
974}
975
976/*
977 ************************************************************************
978 * 1/(R^2 sin(theta)) dR/dphi
979 ************************************************************************
980 */
981
982Mtbl* map_star_fait_sr2stdrdp(const Map* cvi) {
983
984 // recup du changement de variable
985 const Map_star* cv = static_cast<const Map_star*>(cvi) ;
986 const Mg3d* mg = cv->get_mg() ;
987 int nz = mg->get_nzone() ;
988
989 // Le resultat
990 Mtbl* mti = new Mtbl(mg) ;
991 mti->set_etat_qcq() ;
992
993 // Pour le confort
994 const Valeur& alpha = cv->get_alpha() ;
995 const Valeur& stalphadp = alpha.stdsdp() ;
996
997
998 for (int l=0 ; l<nz ; l++) {
999
1000 int nr = mg->get_nr(l);
1001 int nt = mg->get_nt(l) ;
1002 int np = mg->get_np(l) ;
1003 const Grille3d* g = mg->get_grille3d(l) ;
1004 Tbl* tb = (mti->t)[l] ;
1005 tb->set_etat_qcq() ;
1006 double* p_r = tb->t ;
1007
1008 switch(mg->get_type_r(l)) {
1009
1010
1011 case RARE : {
1012 for (int k=0 ; k<np ; k++) {
1013 for (int j=0 ; j<nt ; j++) {
1014 for (int i=0 ; i<nr ; i++) {
1015 double ww = alpha(l,k,j,0) ;
1016 *p_r = stalphadp(l,k,j,0) / (ww*ww * (g->x)[i]) ;
1017 p_r++ ;
1018 } // Fin de boucle sur r
1019 } // Fin de boucle sur theta
1020 } // Fin de boucle sur phi
1021 break ;
1022 }
1023
1024 case FIN:{
1025 cout << "map_star_fait_sr2stdrdp: Shells not implemented yet..." << endl;
1026 abort() ;
1027 break ;
1028 }
1029
1030 case UNSURR: {
1031 cout << "map_star_fait_sr2stdrdp: Compactified domain not allowed !" << endl;
1032 abort() ;
1033 break ;
1034 }
1035
1036 default: {
1037 cout << "map_star_fait_sr2stdrdp: unknown type_r !" << endl ;
1038 abort() ;
1039 }
1040 } // Fin du switch
1041 } // Fin de boucle sur zone
1042
1043 // Termine
1044
1045 return mti ;
1046}
1047
1048/*
1049 ************************************************************************
1050 * d^2R/dx^2
1051 ************************************************************************
1052 */
1053
1054Mtbl* map_star_fait_d2rdx2(const Map* cvi) {
1055
1056 // recup de la grille
1057 const Mg3d* mg = cvi->get_mg() ;
1058
1059 // Le resultat est nul :
1060 Mtbl* mti = new Mtbl(mg) ;
1061 mti->set_etat_zero() ;
1062
1063 return mti ;
1064}
1065
1066/*
1067 *****************************************************************************
1068 * 1/R^2 ( 1/sin(th) d/dth( sin(th) dR/dth ) + 1/sin(th)^2 d^2R/dphi^2 )
1069 *****************************************************************************
1070 */
1071
1072Mtbl* map_star_fait_lapr_tp(const Map* cvi) {
1073
1074 // recup du changement de variable
1075 const Map_star* cv = static_cast<const Map_star*>(cvi) ;
1076 const Mg3d* mg = cv->get_mg() ;
1077 int nz = mg->get_nzone() ;
1078
1079 // Le resultat
1080 Mtbl* mti = new Mtbl(mg) ;
1081 mti->set_etat_qcq() ;
1082
1083 // Pour le confort
1084 const Valeur& alpha = cv->get_alpha() ;
1085
1086 Valeur alpha_tmp = alpha ;
1087 alpha_tmp.ylm() ;
1088 const Valeur& lapalpha = alpha_tmp.lapang() ;
1089
1090
1091 for (int l=0 ; l<nz ; l++) {
1092
1093 const Grille3d* g = mg->get_grille3d(l) ;
1094
1095 int nr = mg->get_nr(l);
1096 int nt = mg->get_nt(l) ;
1097 int np = mg->get_np(l) ;
1098
1099 Tbl* tb = (mti->t)[l] ;
1100 tb->set_etat_qcq() ;
1101 double* p_r = tb->t ;
1102
1103 switch(mg->get_type_r(l)) {
1104
1105 case RARE: {
1106
1107 for (int k=0 ; k<np ; k++) {
1108 for (int j=0 ; j<nt ; j++) {
1109 for (int i=0 ; i<nr ; i++) {
1110 *p_r = (g->x)[i]*lapalpha(l,k,j,0) ;
1111 p_r++ ;
1112 } // Fin de boucle sur r
1113 } // Fin de boucle sur theta
1114 } // Fin de boucle sur phi
1115 break ;
1116 }
1117
1118 case FIN:{
1119 cout << "map_star_fait_lapr_tp: Shells not implemented yet..." << endl;
1120 abort() ;
1121 break ;
1122 }
1123
1124 case UNSURR: {
1125 cout << "map_star_fait_lapr_tp: Compactified domain not allowed !" << endl;
1126 abort() ;
1127 break ;
1128 }
1129
1130 default: {
1131 cout << "map_star_fait_lapr_tp: unknown type_r !" << endl ;
1132 abort() ;
1133 }
1134 } // Fin du switch
1135 } // Fin de boucle sur zone
1136
1137 // Termine
1138
1139 return mti ;
1140}
1141
1142/*
1143 ************************************************************************
1144 * d^2R/dthdx
1145 ************************************************************************
1146 */
1147
1148Mtbl* map_star_fait_d2rdtdx(const Map* cvi) {
1149
1150 // recup du changement de variable
1151 const Map_star* cv = static_cast<const Map_star*>(cvi) ;
1152 const Mg3d* mg = cv->get_mg() ;
1153 int nz = mg->get_nzone() ;
1154
1155 // Le resultat
1156 Mtbl* mti = new Mtbl(mg) ;
1157 mti->set_etat_qcq() ;
1158
1159 // Pour le confort
1160 const Valeur& alpha = cv->get_alpha() ;
1161 const Valeur& d2alphadxdt = alpha.dsdt().dsdx() ;
1162
1163
1164 for (int l=0 ; l<nz ; l++) {
1165
1166 int nr = mg->get_nr(l);
1167 int nt = mg->get_nt(l) ;
1168 int np = mg->get_np(l) ;
1169
1170 Tbl* tb = (mti->t)[l] ;
1171 tb->set_etat_qcq() ;
1172 double* p_r = tb->t ;
1173
1174 switch(mg->get_type_r(l)) {
1175
1176
1177 case RARE : {
1178 for (int k=0 ; k<np ; k++) {
1179 for (int j=0 ; j<nt ; j++) {
1180 for (int i=0 ; i<nr ; i++) {
1181 *p_r = d2alphadxdt(l,k,j,0)/alpha(l,k,j,0) ;
1182 p_r++ ;
1183 } // Fin de boucle sur r
1184 } // Fin de boucle sur theta
1185 } // Fin de boucle sur phi
1186 break ;
1187 }
1188 case FIN:{
1189 cout << "map_star_fait_d2rdtdx: Shells not implemented yet..." << endl;
1190 abort() ;
1191 break ;
1192 }
1193
1194 case UNSURR: {
1195 cout << "map_star_fait_d2rdtdx: Compactified domain not allowed !" << endl;
1196 abort() ;
1197 break ;
1198 }
1199
1200 default: {
1201 cout << "map_star_fait_d2rdtdx: unknown type_r !" << endl ;
1202 abort() ;
1203 }
1204 } // Fin du switch
1205 } // Fin de boucle sur zone
1206
1207 // Termine
1208
1209 return mti ;
1210}
1211
1212/*
1213 ************************************************************************
1214 * 1/sin(th) d^2R/dphidx
1215 ************************************************************************
1216 */
1217
1218Mtbl* map_star_fait_sstd2rdpdx(const Map* cvi) {
1219
1220 // recup du changement de variable
1221 const Map_star* cv = static_cast<const Map_star*>(cvi) ;
1222 const Mg3d* mg = cv->get_mg() ;
1223 int nz = mg->get_nzone() ;
1224
1225 // Le resultat
1226 Mtbl* mti = new Mtbl(mg) ;
1227 mti->set_etat_qcq() ;
1228
1229 // Pour le confort
1230 const Valeur& alpha = cv->get_alpha() ;
1231 const Valeur& stalphadp = alpha.stdsdp() ;
1232
1233
1234 for (int l=0 ; l<nz ; l++) {
1235
1236 int nr = mg->get_nr(l);
1237 int nt = mg->get_nt(l) ;
1238 int np = mg->get_np(l) ;
1239 const Grille3d* g = mg->get_grille3d(l) ;
1240 Tbl* tb = (mti->t)[l] ;
1241 tb->set_etat_qcq() ;
1242 double* p_r = tb->t ;
1243
1244 switch(mg->get_type_r(l)) {
1245
1246
1247 case RARE : {
1248 for (int k=0 ; k<np ; k++) {
1249 for (int j=0 ; j<nt ; j++) {
1250 for (int i=0 ; i<nr ; i++) {
1251 *p_r = stalphadp(l,k,j,0) ;
1252 p_r++ ;
1253 } // Fin de boucle sur r
1254 } // Fin de boucle sur theta
1255 } // Fin de boucle sur phi
1256 break ;
1257 }
1258
1259 case FIN:{
1260 cout << "map_star_fait_sstd2rdpdx: Shells not implemented yet..." << endl;
1261 abort() ;
1262 break ;
1263 }
1264
1265 case UNSURR: {
1266 cout << "map_star_fait_sstd2rdpdx: Compactified domain not allowed !" << endl;
1267 abort() ;
1268 break ;
1269 }
1270
1271 default: {
1272 cout << "map_star_fait_sstd2rdpdx: unknown type_r !" << endl ;
1273 abort() ;
1274 }
1275 } // Fin du switch
1276 } // Fin de boucle sur zone
1277
1278 // Termine
1279
1280
1281 return mti ;
1282}
1283
1284/*
1285 ************************************************************************
1286 * 1/R^2 d^2R/dtheta^2
1287 ************************************************************************
1288 */
1289
1290Mtbl* map_star_fait_sr2d2rdt2(const Map* cvi) {
1291
1292 // recup du changement de variable
1293 const Map_star* cv = static_cast<const Map_star*>(cvi) ;
1294 const Mg3d* mg = cv->get_mg() ;
1295 int nz = mg->get_nzone() ;
1296
1297 // Le resultat
1298 Mtbl* mti = new Mtbl(mg) ;
1299 mti->set_etat_qcq() ;
1300
1301 // Pour le confort
1302 const Valeur& alpha = cv->get_alpha() ;
1303 const Valeur& d2alphadt2 = alpha.d2sdt2() ;
1304
1305
1306 for (int l=0 ; l<nz ; l++) {
1307
1308 int nr = mg->get_nr(l);
1309 int nt = mg->get_nt(l) ;
1310 int np = mg->get_np(l) ;
1311 const Grille3d* g = mg->get_grille3d(l) ;
1312 Tbl* tb = (mti->t)[l] ;
1313 tb->set_etat_qcq() ;
1314 double* p_r = tb->t ;
1315
1316 switch(mg->get_type_r(l)) {
1317
1318
1319 case RARE : {
1320 for (int k=0 ; k<np ; k++) {
1321 for (int j=0 ; j<nt ; j++) {
1322 for (int i=0 ; i<nr ; i++) {
1323 double ww = alpha(l,k,j,0) ;
1324 *p_r = d2alphadt2(l,k,j,0) / (ww*ww * (g->x)[i]) ;
1325 p_r++ ;
1326 } // Fin de boucle sur r
1327 } // Fin de boucle sur theta
1328 } // Fin de boucle sur phi
1329 break ;
1330 }
1331 case FIN:{
1332 cout << "map_star_fait_drdt: Shells not implemented yet..." << endl;
1333 abort() ;
1334 break ;
1335 }
1336
1337 case UNSURR: {
1338 cout << "map_star_fait_drdt: Compactified domain not allowed !" << endl;
1339 abort() ;
1340 break ;
1341 }
1342
1343 default: {
1344 cout << "map_star_fait_drdt: unknown type_r !" << endl ;
1345 abort() ;
1346 }
1347 } // Fin du switch
1348 } // Fin de boucle sur zone
1349
1350 // Termine
1351
1352 return mti ;
1353}
1354
1355}
Lorene::Grille3d
3D grid class in one domain.
Definition grilles.h:200
Lorene::Grille3d::x
double * x
Array of values of at the nr collocation points.
Definition grilles.h:215
Lorene::Grille3d::tet
double * tet
Array of values of at the nt collocation points.
Definition grilles.h:217
Lorene::Grille3d::get_nr
int get_nr() const
Returns nr.
Definition grilles.h:236
Lorene::Map_star
Definition map.h:3917
Lorene::Map_star::get_alpha
const Valeur & get_alpha() const
Returns the reference on the Tbl alpha.
Definition map_star.C:237
Lorene::Map_star::Map_star
Map_star(const Mg3d &mgrille, const double *r_limits)
Standard Constructor.
Definition map_star.C:52
Lorene::Map_star::map_star_fait_r
friend Mtbl * map_star_fait_r(const Map *)
< Not implemented
Definition map_star_fait.C:43
Lorene::Map_star::alpha
Valeur alpha
Array (size: mg->nzone*Nt*Np ) of the values of in each domain.
Definition map.h:3924
Lorene::Mg3d
Multi-domain grid.
Definition grilles.h:279
Lorene::Mg3d::get_grille3d
const Grille3d * get_grille3d(int l) const
Returns a pointer on the 3D mono-grid for domain no. l.
Definition grilles.h:517
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:479
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:474
Lorene::Mg3d::get_nzone
int get_nzone() const
Returns the number of domains.
Definition grilles.h:465
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:469
Lorene::Mg3d::get_type_r
int get_type_r(int l) const
Returns the type of sampling in the radial direction in domain no.
Definition grilles.h:491
Lorene::Mtbl
Multi-domain array.
Definition mtbl.h:118
Lorene::Mtbl::t
Tbl ** t
Array (size nzone ) of pointers on the Tbl 's.
Definition mtbl.h:132
Lorene::Mtbl::set_etat_qcq
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition mtbl.C:302
Lorene::Tbl
Basic array class.
Definition tbl.h:161
Lorene::Tbl::set_etat_qcq
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition tbl.C:364
Lorene::Tbl::t
double * t
The array of double.
Definition tbl.h:173
Lorene::Valeur
Values and coefficients of a (real-value) function.
Definition valeur.h:297
Lorene::Valeur::d2sdt2
const Valeur & d2sdt2() const
Returns of *this.
Definition valeur_d2sdt2.C:110
Lorene::Valeur::stdsdp
const Valeur & stdsdp() const
Returns of *this.
Definition valeur_stdsdp.C:63
Lorene::Valeur::ylm
void ylm()
Computes the coefficients of *this.
Definition valeur_ylm.C:141
Lorene::Valeur::dsdt
const Valeur & dsdt() const
Returns of *this.
Definition valeur_dsdt.C:115
Lorene::Valeur::dsdx
const Valeur & dsdx() const
Returns of *this.
Definition valeur_dsdx.C:114
Lorene::Valeur::lapang
const Valeur & lapang() const
Returns the angular Laplacian of *this.
Definition valeur_lapang.C:75
Lorene::sin
Cmp sin(const Cmp &)
Sine.
Definition cmp_math.C:72
Lorene::cos
Cmp cos(const Cmp &)
Cosine.
Definition cmp_math.C:97
Lorene
Lorene prototypes.
Definition app_hor.h:67
Lorene::Map
Map(const Mg3d &)
Constructor from a multi-domain 3D grid.
Definition map.C:142
Lorene::phi
Coord phi
coordinate centered on the grid
Definition map.h:732