Télécharger la plaquette Cast3M


fr en

Télécharger nlin_lapnssphe_3d.dgibi

Retour à la liste

Numérotation des lignes : 

  1. * fichier : nlin_lapnssphe_3d.dgibi
  2. ************************************************************************
  3. ************************************************************************
  4. 'OPTION' 'ECHO' 0 ;
  5. *BEGINPROCEDUR gmass
  6. ************************************************************************
  7. * NOM         : GMASS
  8. * DESCRIPTION : Une matrice de masse
  9. *
  10. *
  11. *
  12. * LANGAGE     : GIBIANE-CAST3M
  13. * AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  14. *               mél : gounand@semt2.smts.cea.fr
  15. **********************************************************************
  16. * VERSION    : v2, 14/03/2006, mise à jour NLIN évolué
  17. * VERSION    : v1, 13/05/2004, version initiale
  18. * HISTORIQUE : v1, 13/05/2004, création
  19. * HISTORIQUE :
  20. * HISTORIQUE :
  21. ************************************************************************
  22. * Prière de PRENDRE LE TEMPS de compléter les commentaires
  23. * en cas de modification de ce sous-programme afin de faciliter
  24. * la maintenance !
  25. ************************************************************************
  26. *
  27. *
  28. 'DEBPROC' GMASS ;
  29. 'ARGUMENT'     _mt*'MAILLAGE' ;
  30. 'ARGUMENT'     _smt/'MAILLAGE' ;
  31. 'ARGUMENT'   gdisc*'MOT     ' ;
  32. 'ARGUMENT'    nomt*'MOT     ' ;
  33. 'ARGUMENT'   disct*'MOT     ' ;
  34. 'ARGUMENT'    nomq*'MOT     ' ;
  35. 'ARGUMENT'   discq*'MOT     ' ;
  36. 'ARGUMENT' coef/'FLOTTANT' ;
  37.    'SI' ('NON' ('EXISTE' coef)) ;
  38.       'ARGUMENT' coef2/'CHPOINT ' ;
  39.       'SI' ('NON' ('EXISTE' coef2)) ;
  40.          'ERREUR' 'Il faut donner un coef FLOTTANT ou CHPOINT' ;
  41.       'SINON' ;
  42.          coef = coef2 ;
  43.          'ARGUMENT' discc*'MOT     ' ;
  44.       'FINSI' ;                  
  45.    'SINON' ;
  46.       discc = 'CSTE' ;
  47.    'FINSI' ;
  48. 'ARGUMENT' methgau/'MOT     ' ;
  49.    'SI' ('NON' ('EXISTE' methgau)) ;
  50.       methgau = 'GAU7' ;
  51.    'FINSI' ;
  52. 'ARGUMENT' chpop/'CHPOINT' ;
  53. 'ARGUMENT' chpod/'CHPOINT' ;
  54. *
  55.    vdim = 'VALEUR' 'DIME' ;
  56.    vmod = 'VALEUR' 'MODE' ;
  57.    idim = 0 ;
  58.    'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'PLANDEFO')) ;
  59.       idim = 2 ;
  60.       iaxi = FAUX ;
  61.    'FINSI' ;
  62.    'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'AXIS')) ;
  63.       idim = 2 ;
  64.       iaxi = VRAI ;
  65.    'FINSI' ;
  66.    'SI' ('ET' ('EGA' vdim 3) ('EGA' vmod 'TRID')) ;
  67.       idim = 3 ;
  68.       iaxi = FAUX ;
  69.    'FINSI' ;
  70.    'SI' ('EGA' vdim 1) ;
  71.       idim = 1 ;
  72.       iaxi = FAUX ;
  73.    'FINSI' ;
  74.    'SI' ('EGA' idim 0) ;
  75.       'ERREUR' ('CHAINE' 'vmod=' vmod ' et vdim=' vdim ' non prevu') ;
  76.    'FINSI' ;
  77.    'SI' iaxi ;
  78.       dprmt = '*' ('COORDONNEE' 1 _mt) ('*' PI 2.D0) ;
  79.    'FINSI' ;      
  80.    numop  = 1 ;
  81.    numder = idim ;
  82.    mmt = 'MOTS' nomt ;
  83.    mmq = 'MOTS' nomq ;
  84.    numvar = 1 ;
  85.    numdat = 1 ;
  86.    numcof = 1 ;
  87. *
  88.    A = ININLIN numop numvar numdat numcof numder ;   
  89.    A . 'VAR' . 1 . 'NOMDDL' = mmt ;
  90.    A . 'VAR' . 1 . 'DISC'   = disct ;
  91.    'SI' ('EXISTE' chpop) ;
  92.       A . 'VAR' . 1 . 'VALEUR'   = chpop ;
  93.    'FINSI' ;   
  94.    A . 'DAT' . 1 . 'NOMDDL' = 'MOTS' 'SCAL' ;
  95.    A . 'DAT' . 1 . 'DISC'   = discc ;
  96.    A . 'DAT' . 1 . 'VALEUR' = coef ;
  97.    A . 'COF' . 1 . 'COMPOR' = 'IDEN' ;
  98.    A . 'COF' . 1 . 'LDAT'   = 'LECT' 1 ;
  99. *
  100.    A . 1 . 1 . 0 = 'LECT' 1 ;
  101. *
  102.    'SI' iaxi ;
  103.       numdat = 1 ;
  104.       numcof = 1 ;
  105.    'SINON' ;      
  106.       numdat = 0 ;
  107.       numcof = 0 ;
  108.    'FINSI' ;      
  109.    B = ININLIN numop numvar numdat numcof numder ;   
  110.    B . 'VAR' . 1 . 'NOMDDL' = mmq ;
  111.    B . 'VAR' . 1 . 'DISC'   = discq ;
  112.    'SI' ('EXISTE' chpod) ;
  113.       B . 'VAR' . 1 . 'VALEUR'   = chpod ;
  114.    'FINSI' ;   
  115. *
  116.    'SI' iaxi ;
  117.       B . 'DAT' . 1 . 'NOMDDL' = 'MOTS' 'SCAL' ;
  118.       B . 'DAT' . 1 . 'DISC'   = gdisc ;
  119.       B . 'DAT' . 1 . 'VALEUR' = dprmt ;
  120.       B . 'COF' . 1 . 'COMPOR' = 'IDEN' ;
  121.       B . 'COF' . 1 . 'LDAT'   = 'LECT' 1 ;
  122.    'FINSI' ;
  123.    'SI' iaxi ;      
  124.       B . 1 . 1 . 0 = 'LECT' 1 ;
  125.    'SINON' ;      
  126.       B . 1 . 1 . 0 = 'LECT' ;
  127.    'FINSI' ;      
  128. *
  129.    'SI' ('EXISTE' _smt) ;
  130.        mgmass = 'NLIN' gdisc _mt _smt A B methgau ;
  131.    'SINON' ;       
  132.        mgmass = 'NLIN' gdisc _mt A B methgau ;
  133.    'FINSI' ;       
  134. *
  135. 'RESPRO' mgmass ;
  136. 'FINPROC' ;
  137. *
  138. * End of procedure file GMASS
  139. *
  140. *ENDPROCEDUR gmass
  141. *BEGINPROCEDUR glapn
  142. ************************************************************************
  143. * NOM         : GLAPN
  144. * DESCRIPTION : Un laplacien scalaire 
  145. *
  146. *
  147. *
  148. * LANGAGE     : GIBIANE-CAST3M
  149. * AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  150. *               mél : gounand@semt2.smts.cea.fr
  151. **********************************************************************
  152. * VERSION    : v2, 14/03/2006, mise à jour NLIN évolué
  153. * VERSION    : v1, 13/05/2004, version initiale
  154. * HISTORIQUE : v1, 13/05/2004, création
  155. * HISTORIQUE :
  156. * HISTORIQUE :
  157. ************************************************************************
  158. * Prière de PRENDRE LE TEMPS de compléter les commentaires
  159. * en cas de modification de ce sous-programme afin de faciliter
  160. * la maintenance !
  161. ************************************************************************
  162. *
  163. *
  164. 'DEBPROC' GLAPN ;
  165. 'ARGUMENT'     _mt*'MAILLAGE' ;
  166. 'ARGUMENT'   gdisc*'MOT     ' ;
  167. 'ARGUMENT'    nomt*'MOT     ' ;
  168. 'ARGUMENT'   disct*'MOT     ' ;
  169. 'ARGUMENT'    nomq*'MOT     ' ;
  170. 'ARGUMENT'   discq*'MOT     ' ;
  171. 'ARGUMENT' coef/'FLOTTANT' ;
  172.    'SI' ('NON' ('EXISTE' coef)) ;
  173.       'ARGUMENT' coef2/'CHPOINT ' ;
  174.       'SI' ('NON' ('EXISTE' coef2)) ;
  175.          'ERREUR' 'Il faut donner un coef FLOTTANT ou CHPOINT' ;
  176.       'SINON' ;
  177.          coef = coef2 ;
  178.          'ARGUMENT' discc*'MOT     ' ;
  179.       'FINSI' ;                  
  180.    'SINON' ;
  181.       discc = 'CSTE' ;
  182.    'FINSI' ;            
  183.    'ARGUMENT' methgau/'MOT     ' ;
  184.    'SI' ('NON' ('EXISTE' methgau)) ;
  185.       methgau = 'GAU7' ;
  186.    'FINSI' ;      
  187. 'ARGUMENT' chpop/'CHPOINT' ;
  188. 'ARGUMENT' chpod/'CHPOINT' ;
  189. *
  190.    vdim = 'VALEUR' 'DIME' ;
  191.    vmod = 'VALEUR' 'MODE' ;
  192.    idim = 0 ;
  193.    'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'PLANDEFO')) ;
  194.       idim = 2 ;
  195.       iaxi = FAUX ;
  196.    'FINSI' ;
  197.    'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'AXIS')) ;
  198.       idim = 2 ;
  199.       iaxi = VRAI ;
  200.    'FINSI' ;
  201.    'SI' ('ET' ('EGA' vdim 3) ('EGA' vmod 'TRID')) ;
  202.       idim = 3 ;
  203.       iaxi = FAUX ;
  204.    'FINSI' ;
  205.    'SI' ('EGA' vdim 1) ;
  206.       idim = 1 ;
  207.       iaxi = FAUX ;
  208.    'FINSI' ;
  209.    'SI' ('EGA' idim 0) ;
  210.       'ERREUR' ('CHAINE' 'vmod=' vmod ' et vdim=' vdim ' non prevu') ;
  211.    'FINSI' ;
  212.    'SI' iaxi ;
  213.       dprmt = '*' ('COORDONNEE' 1 _mt) ('*' PI 2.D0) ;
  214.    'FINSI' ;      
  215. * Test bête...
  216.    'SI' ('EGA' ('TYPE' coef) 'CHPOINT ') ;
  217.       mincoef = 'MINIMUM' coef ;
  218.    'SINON' ;
  219.       mincoef = coef ;
  220.    'FINSI' ;            
  221.    'SI' ('<' mincoef 0.D0) ;
  222.       'ERREUR' 'Le coef (une viscosité) doit etre positive' ;
  223.    'FINSI' ;      
  224. *
  225.    'SI' iaxi ;
  226.       lcoef = 'MOTS' 'NURR' 'NUZZ' ;
  227.    'SINON' ;
  228.       'SI' ('EGA' idim 2) ;
  229.          lcoef = 'MOTS' 'NUXX' 'NUYY' ;
  230.       'SINON' ;
  231.          lcoef = 'MOTS' 'NUXX' 'NUYY' 'NUZZ' ;
  232.       'FINSI' ;
  233.    'FINSI' ;      
  234. *
  235.    ltens = FAUX ;
  236.    'SI' ('EGA' ('TYPE' coef) 'CHPOINT') ;
  237.       ncomp = 'DIME' ('EXTRAIRE' coef 'COMP') ;
  238.       'SI' ('EGA' ncomp idim) ;
  239.          ltens = VRAI ;
  240.       'FINSI' ;
  241.    'FINSI' ;
  242. *
  243.    tcoef = 'TABLE' ;
  244.    'SI' ltens ;
  245.       'REPETER' iidim idim ;
  246.          tcoef . &iidim = 'EXCO' ('EXTRAIRE' lcoef &iidim) coef ;
  247.       'FIN' iidim ;
  248.    'SINON' ;               
  249.       'REPETER' iidim idim ;
  250.          tcoef . &iidim = coef ;
  251.       'FIN' iidim ;
  252.    'FINSI' ;
  253. *      
  254.    numop  = idim ;
  255.    numder = idim ;
  256.    mmt = 'MOTS' nomt ;
  257.    mmq = 'MOTS' nomq ;
  258.    numvar = 1 ;
  259.    numdat = idim ;
  260.    numcof = idim ;
  261. *
  262.    A = ININLIN numop numvar numdat numcof numder ;   
  263.    A . 'VAR' . 1 . 'NOMDDL' = mmt ;
  264.    A . 'VAR' . 1 . 'DISC'   = disct ;
  265.    'SI' ('EXISTE' chpop) ;
  266.       A . 'VAR' . 1 . 'VALEUR'   = chpop ;
  267.    'FINSI' ;   
  268.    'REPETER' iidim idim ;
  269.       A . 'DAT' . &iidim . 'NOMDDL' = 'MOTS' 'SCAL' ;
  270.       A . 'DAT' . &iidim . 'DISC'   = discc ;
  271.       A . 'DAT' . &iidim . 'VALEUR' = tcoef . &iidim ;
  272.       A . 'COF' . &iidim . 'COMPOR' = 'IDEN' ;
  273.       A . 'COF' . &iidim . 'LDAT'   = 'LECT' &iidim ;
  274.    'FIN' iidim ;
  275. *
  276.    'REPETER' iidim idim ;
  277.       A . &iidim . 1 . &iidim = 'LECT' &iidim ;
  278.    'FIN' iidim ;      
  279. *
  280.    'SI' iaxi ;
  281.       numdat = 1 ;
  282.       numcof = 1 ;
  283.    'SINON' ;      
  284.       numdat = 0 ;
  285.       numcof = 0 ;
  286.    'FINSI' ;      
  287. *
  288.    B = ININLIN numop numvar numdat numcof numder ;   
  289.    B . 'VAR' . 1 . 'NOMDDL' = mmq ;
  290.    B . 'VAR' . 1 . 'DISC'   = discq ;
  291.    'SI' ('EXISTE' chpod) ;
  292.       B . 'VAR' . 1 . 'VALEUR'   = chpod ;
  293.    'FINSI' ;   
  294.    'SI' iaxi ;
  295.       B . 'DAT' . 1 . 'NOMDDL' = 'MOTS' 'SCAL' ;
  296.       B . 'DAT' . 1 . 'DISC'   = gdisc ;
  297.       B . 'DAT' . 1 . 'VALEUR' = dprmt ;
  298.       B . 'COF' . 1 . 'COMPOR' = 'IDEN' ;
  299.       B . 'COF' . 1 . 'LDAT'   = 'LECT' 1 ;
  300.    'FINSI' ;
  301. *
  302.    'REPETER' iidim idim ;
  303.       'SI' iaxi ;      
  304.          B . &iidim . 1 . &iidim = 'LECT' 1 ;
  305.       'SINON' ;      
  306.          B . &iidim . 1 . &iidim = 'LECT' ;
  307.       'FINSI' ;      
  308.    'FIN' iidim ;      
  309. *
  310.    mglapn = 'NLIN' gdisc _mt A B methgau ;
  311. * Integration par parties
  312.    mglapn = '*' mglapn -1.D0 ;
  313. *
  314. 'RESPRO' mglapn ;
  315. 'FINPROC' ;
  316. *
  317. * End of procedure file GLAPN
  318. *
  319. *ENDPROCEDUR glapn
  320. *BEGINPROCEDUR formar
  321. ************************************************************************
  322. * NOM         : FORMAR
  323. * DESCRIPTION : formate un réel de facon courte
  324. *               pratique pour les noms de
  325. *               sauvegarde
  326. *               Exemples :
  327. * 'MESSAGE' ('CHAINE' (formar 2.9e5   1)) ;
  328. * 2.9E5
  329. * 'MESSAGE' ('CHAINE' (formar -2.9e5  1)) ;
  330. * -2.9E5
  331. * 'MESSAGE' ('CHAINE' (formar 2.9e-5  1)) ;
  332. * 2.9E-5
  333. * 'MESSAGE' ('CHAINE' (formar -2.9e-5 1)) ;
  334. * -2.9E-5
  335. * 'MESSAGE' ('CHAINE' (formar 2.9     1)) ;
  336. * 2.9
  337. * 'MESSAGE' ('CHAINE' (formar -2.9    1)) ;
  338. * -2.9
  339. * 'MESSAGE' ('CHAINE' (formar 0       1)) ;
  340. * 0
  341. * 'MESSAGE' ('CHAINE' (formar 0       1)) ;
  342. * 0
  343. * 'MESSAGE' ('CHAINE' (formar 2.9e5   0)) ;
  344. * 3E5
  345. * 'MESSAGE' ('CHAINE' (formar -2.9e5  0)) ;
  346. * -3E5
  347. * 'MESSAGE' ('CHAINE' (formar 2.9e-5  0)) ;
  348. * 3E-5
  349. * 'MESSAGE' ('CHAINE' (formar -2.9e-5 0)) ;
  350. * -3E-5
  351. * 'MESSAGE' ('CHAINE' (formar 2.9     0)) ;
  352. * 3
  353. * 'MESSAGE' ('CHAINE' (formar -2.9    0)) ;
  354. * -3
  355. * 'MESSAGE' ('CHAINE' (formar 0       0)) ;
  356. * 0
  357. * 'MESSAGE' ('CHAINE' (formar 0       0)) ;
  358. * 0
  359. *               
  360. *
  361. *
  362. * LANGAGE     : GIBIANE-CAST3M
  363. * AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  364. *               mél : gounand@semt2.smts.cea.fr
  365. **********************************************************************
  366. * VERSION    : v1, 18/02/2003, version initiale
  367. * HISTORIQUE : v1, 18/02/2003, création
  368. * HISTORIQUE :
  369. * HISTORIQUE :
  370. ************************************************************************
  371. * Prière de PRENDRE LE TEMPS de compléter les commentaires
  372. * en cas de modification de ce sous-programme afin de faciliter
  373. * la maintenance !
  374. ************************************************************************
  375. *
  376. *
  377. 'DEBPROC' FORMAR ;
  378. 'ARGUMENT' fl*'FLOTTANT' ;
  379. 'ARGUMENT' vir/'ENTIER  ' ;
  380. 'SI' ('NON' ('EXISTE' vir)) ;
  381.    vir = 1 ;
  382. 'SINON' ;   
  383.    'SI' ('<' vir 0) ;
  384.       'ERREUR' 'fournir un entier positif' ;
  385.    'FINSI' ;
  386. 'FINSI' ;         
  387. 'SI' ('<' ('ABS' fl) 10.D-100) ;
  388.    chfl = 'CHAINE' '0' ;
  389. 'SINON' ;
  390. *! sans le 1.D-10, ca ne fonctionne pas
  391. *! qd on entre pile poil une puissance de 10
  392.    lfl = LOG10 ('ABS' fl) ;
  393. *   lfl = '+' (LOG10 ('ABS' fl)) 1.D-10  ;
  394.    slfl = 'SIGNE' ('ENTIER' lfl) ;
  395.    'SI' ('EGA' slfl 1) ;
  396.       elfl = 'ENTIER' lfl ;
  397.    'SINON' ;
  398.       elfl = '-' ('ENTIER' lfl) 1 ;
  399.    'FINSI' ;
  400.    man = '/' fl ('**' 10.D0 elfl) ;
  401. *
  402. * Une verrue pour des histoires de précision...
  403. *   
  404.    'SI' ('EGA' man 10.D0 ('**' 10.D0 ('*' vir -1.D0))) ;
  405.       man  = '/' man 10.D0 ;
  406.       elfl = '+' elfl 1 ;
  407.    'FINSI' ;
  408. *   
  409.    sman = 'SIGNE' man ;
  410.    'SI' ('EGA' sman 1) ;   
  411.       fman = 'CHAINE' '(F' ('+' vir 2) '.0' vir ')' ;
  412.    'SINON' ;
  413.       fman = 'CHAINE' '(F' ('+' vir 3) '.0' vir ')' ;
  414.    'FINSI' ;
  415.    'SI' ('NEG' vir 0) ;   
  416.       'SI' ('NEG' elfl 0) ;
  417.          chfl = 'CHAINE' 'FORMAT' fman man 'E' elfl ;
  418.       'SINON' ;
  419.          chfl = 'CHAINE' 'FORMAT' fman man ;
  420.       'FINSI' ;
  421.    'SINON' ;
  422.       man2 = 'ENTIER' ('+' man ('*' 0.5D0 sman)) ;
  423.       'SI' ('NEG' elfl 0) ;
  424.          chfl = 'CHAINE' man2 'E' elfl ;
  425.       'SINON' ;
  426.          chfl = 'CHAINE' man2 ;
  427.       'FINSI' ;
  428.    'FINSI' ;            
  429. 'FINSI' ;   
  430. 'RESPRO' chfl ;
  431. *
  432. * End of procedure file FORMAR
  433. *
  434. 'FINPROC' ;
  435. *ENDPROCEDUR formar
  436. *BEGINPROCEDUR log10
  437. ************************************************************************
  438. * NOM         : LOG10
  439. * DESCRIPTION : Log_10
  440. *
  441. *
  442. *
  443. * LANGAGE     : GIBIANE-CAST3M
  444. * AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  445. *               mél : gounand@semt2.smts.cea.fr
  446. **********************************************************************
  447. * VERSION    : v1, 18/02/2003, version initiale
  448. * HISTORIQUE : v1, 18/02/2003, création
  449. * HISTORIQUE :
  450. * HISTORIQUE :
  451. ************************************************************************
  452. * Prière de PRENDRE LE TEMPS de compléter les commentaires
  453. * en cas de modification de ce sous-programme afin de faciliter
  454. * la maintenance !
  455. ************************************************************************
  456. *
  457. *
  458. 'DEBPROC' LOG10 ;
  459. 'REPETER' bouc ;
  460.    ok = FAUX ;
  461.    'ARGUMENT' fl/'FLOTTANT' ;
  462.    'ARGUMENT' lr/'LISTREEL' ;
  463.    'ARGUMENT' cp/'CHPOINT ' ;
  464.    'ARGUMENT' cm/'MCHAML  ' ;
  465.    'SI' ('EXISTE' fl) ;
  466.       ok = VRAI ;
  467.       'RESPRO' ('/' ('LOG' fl) ('LOG' 10.D0)) ;
  468.    'FINSI' ;
  469.    'SI' ('EXISTE' lr) ;
  470.       ok = VRAI ;
  471.       'RESPRO' ('/' ('LOG' lr) ('LOG' 10.D0)) ;
  472.    'FINSI' ;
  473.    'SI' ('EXISTE' cp) ;
  474.       ok = VRAI ;
  475.       'RESPRO' ('/' ('LOG' cp) ('LOG' 10.D0)) ;
  476.    'FINSI' ;
  477.    'SI' ('EXISTE' cm) ;
  478.       ok = VRAI ;
  479.       'RESPRO' ('/' ('LOG' cm) ('LOG' 10.D0)) ;
  480.    'FINSI' ;
  481.    'SI' ('NON' ok) ;
  482.       'QUITTER' bouc ;
  483.    'FINSI' ;
  484. 'FIN' bouc ;         
  485. *
  486. * End of procedure file LOG10
  487. *
  488. 'FINPROC' ;
  489. *ENDPROCEDUR log10
  490. 'OPTION' 'ECHO' 1 ;
  491. ************************************************************************
  492. * NOM         : NLIN_LAPNSSPHE_3D
  493. * DESCRIPTION : Vérification de NLIN : on résout un laplacien sur un
  494. *               bout de surface de sphère de rayon R en 3D. La
  495. *               solution exacte est une fonction du type :
  496. *               a * ln (tan (\theta / 2)) + b
  497. *               Il y a une singularité en theta = 0. Si on veut observer
  498. *               l'ordre correct pour theta petit, il faut bien raffiner.
  499. *
  500. *               c'est une primitive de a / sin (\theta).
  501. *               On vérifie les ordres de convergence sur cette solution
  502. *               et sur la surface calculée du bout de sphère.
  503. *
  504. * LANGAGE     : GIBIANE-CAST3M
  505. * AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  506. *               mél : gounand@semt2.smts.cea.fr
  507. **********************************************************************
  508. * VERSION    : v1, 23/10/2006, version initiale
  509. * HISTORIQUE : v1, 23/10/2006, création
  510. * HISTORIQUE :
  511. * HISTORIQUE :
  512. ************************************************************************
  513. * Prière de PRENDRE LE TEMPS de compléter les commentaires
  514. * en cas de modification de ce sous-programme afin de faciliter
  515. * la maintenance !
  516. ************************************************************************
  517. *
  518. *
  519. interact = FAUX ;
  520. graph    = FAUX ;
  521. *interact = VRAI ;
  522. *graph    = VRAI ;
  523. 'OPTION' 'DIME' 3 ;
  524. 'SI' ('NON' interact) ;
  525.   'OPTION' 'TRAC' 'PS' 'ECHO' 0 ;
  526. 'SINON' ;
  527.   'OPTION' 'TRAC' 'X' 'ECHO' 1 ;
  528. 'FINSI' ;  
  529. *
  530. * Paramètres
  531. *
  532. ldisc = 'MOTS' 'LINE' 'QUAI' ;
  533. ttamin  = 40. ;
  534. ttaminr = ttamin '*' ('/' PI 180.) ;
  535. tordsur = 'TABLE' ;
  536. tordsol = 'TABLE' ;
  537. ldens = 'PROG' 2.D-1 1.D-1 5.D-2 ;
  538. *ldens = 'PROG' 1.D-1 5.D-2 1.D-2 ;
  539. rayo  = '**' 2. 0.5D0 ;
  540. tabev = 'TABLE' ;
  541. tabt  = 'TABLE' ;
  542. itab = 0 ;
  543. ok = VRAI ;
  544. *
  545. * Boucle sur les discrétisations 
  546. *
  547. 'REPETER' iidisc ('DIME' ldisc) ;
  548.    idisc = &iidisc ;
  549.    disc  = 'EXTRAIRE' ldisc idisc ;
  550.    tordsur . disc = 'TABLE' ;
  551.    tordsol . disc = 'TABLE' ;
  552.    discg = disc ;
  553.    disct = disc ;
  554.    'SI' ('EGA' disc 'LINE') ;
  555.       'OPTION' 'ELEM' 'TRI3' ;
  556.    'SINON' ;      
  557.       'OPTION' 'ELEM' 'TRI6' ;
  558.    'FINSI' ;      
  559.    lerrsur = 'PROG' ;
  560.    lerrsol = 'PROG' ;
  561.    'REPETER' iidens ('DIME' ldens) ;
  562.       idens = &iidens ;
  563.       vdens = 'EXTRAIRE' ldens idens ;
  564.       'DENS' vdens ;
  565.       tit = 'CHAINE' 'Discretisation :' disc '  Densité = ' vdens ;
  566.       'SAUTER' 1 'LIGNE' ;
  567.       'MESSAGE' tit ;
  568. * Maillage
  569.       phi = 40. ;
  570.       phir = 40. '*' ('/' PI 180.) ;
  571.       p0 = 'POIN' 0. 0. 0. ;
  572.       p1 = 'POIN' 0. 0. ('*' ('COS' ttamin) rayo) ;
  573.       pA = POINTSPH rayo 0. 0. ;
  574.       pB = POINTSPH rayo phi 0. ;
  575.       pC = POINTSPH rayo phi ('-' 90. ttamin) ;
  576.       pD = POINTSPH rayo 0.  ('-' 90. ttamin) ;
  577.       lAB = 'CERCLE' pA p0 pB ;
  578.       lBC = 'CERCLE' pB p0 pC ;
  579.       lCD = 'CERCLE' pC p1 pD ;
  580.       lDA = 'CERCLE' pD p0 pA ;
  581.       cmt = lAB 'ET' lBC 'ET' lCD 'ET' lDA ;
  582. *      cmt = lCD ;
  583. *      xmt ymt zmt = 'COORDONNEE' cmt ;
  584. *      rmt = '**' (('**' xmt 2) '+' ('**' ymt 2) '+' ('**' zmt 2))
  585. *                 0.5D0 ;
  586. *      'LISTE' ('-' ('MAXIMUM' rmt) rayo) ;                 
  587. *      'LISTE' ('-' ('MINIMUM' rmt) rayo) ;
  588. *      'OPTION' 'DONN' 5 ;
  589.       mt = 'SURFACE' cmt 'SPHERIQUE' p0 ;
  590. *      'OPTION' 'DONN' 5 ;
  591. *      mt = 'PROJETER' mt 'CONI' p0 'SPHE' p0 pA ;
  592. *      lAB = 'PROJETER' lAB 'CONI' p0 'SPHE' p0 pA ;
  593. *      lBC = 'PROJETER' lBC 'CONI' p0 'SPHE' p0 pA ;
  594. *      lCD = 'PROJETER' lCD 'CONI' p0 'SPHE' p0 pA ;
  595. *      lDA = 'PROJETER' lDA 'CONI' p0 'SPHE' p0 pA ;
  596. *      'ELIMINATION' (mt 'ET' lAB 'ET' lBC 'ET' lCD 'ET' lDA)
  597. *                    ('/' vdens 100.) ;
  598.       'SI' graph ;
  599.          tit = 'CHAINE' 'Maillage' ;
  600.          'TRACER'  'CACH' mt 'TITR' tit ;
  601.       'FINSI' ;         
  602.       _mt = 'CHANGER' mt 'QUAF' ;
  603. * 1 : calcul de la surface du bout de sphère
  604.       ch1 = 'MANUEL' 'CHPO' mt 1 'T' 1. ;
  605.       gmas = GMASS _mt discg 'T' disct 'Q' disct 1. ;
  606.       surana = phir '*' rayo '*' rayo '*' ('COS' ttamin) ;
  607.       sur = xtmx ch1 gmas ;
  608.       'MESSAGE' ('CHAINE' 'Surface calculee =' sur) ;
  609.       'MESSAGE' ('CHAINE' 'Surface attendue =' surana) ;
  610.       errsur = 'ABS' ('-' sur surana) ;
  611.       'MESSAGE' ('CHAINE' 'Erreur sur la surface =' errsur) ;
  612.       lerrsur = 'ET' lerrsur ('PROG' errsur) ;
  613. * 2 : calcul d'une solution analytique
  614.       xmt ymt zmt = 'COORDONNEE' mt ;
  615.       rmt = '**' ('+' ('**' xmt 2) ('**' ymt 2)) 0.5D0  ;
  616.       tta = 'ATG' rmt zmt ;
  617.       ttamsur2 = '/' ttamin 2. ;
  618.       coef = 'LOG' ('/' ('SIN' ttamsur2) ('COS' ttamsur2)) ;
  619.       ttasur2 = '/' tta 2.D0 ;
  620.       solana = '/' ('LOG' ('/' ('SIN' ttasur2) ('COS' ttasur2))) coef ;
  621.       solana = 'NOMC' 'T' solana ;
  622.       mat = GLAPN _mt discg 'T' disct 'Q' disct 1. ;
  623.       mblo = 'BLOQUE' 'T' (lAB 'ET' lCD) ;
  624.       cblo = '+' ('MANUEL' 'CHPO'  lAB 1 'T' 0.)
  625.                  ('MANUEL' 'CHPO'  lCD 1 'T' 1.) ;
  626. *      mblo = 'BLOQUE' 'T' cmt ;
  627. *      cblo = 'REDU' solana cmt ;
  628.       fblo = 'DEPIMPOSE' mblo cblo ;   
  629. *
  630.       mtot = mat 'ET' mblo ;
  631.       ftot = fblo ;
  632.       solu = 'RESOUD' mtot ftot ;
  633.       sol  = 'EXCO' 'T' solu 'T' ;
  634.       matm = GMASS _mt discg 'T' disct 'Q' disct 1. ;
  635.       errsol = '**' (xtmx ('-' sol solana) matm) 0.5D0 ;
  636.       'MESSAGE' ('CHAINE' 'Erreur sur la sol. analytique =' errsol) ;
  637.       'SI' graph ;
  638.          tit = 'CHAINE' 'Sol Approchée' ;
  639.          'TRACER'  sol mt 'TITR' tit ;
  640.          tit = 'CHAINE' 'Sol Approchée moins Analytique' ;
  641.          'TRACER'  ('-' sol solana) mt 'TITR' tit ;
  642.       'FINSI' ;
  643. *      'OPTION' 'DONN' 5 ;
  644.       lerrsol = 'ET' lerrsol ('PROG' errsol) ;
  645.    'FIN' iidens ;
  646. * Calcul des ordres
  647.    lh lerrsur lerrsol = LOG10 ldens lerrsur lerrsol ;
  648.    evsur = 'EVOL' 'MANU' lh lerrsur ;
  649.    cpl dummy = @POMI evsur 1 'IDEM' ;
  650.    ord = cpl . 1 ;
  651.    tordsur . disc = ord ;
  652.    itab = '+' itab 1 ;
  653.    tabev . itab = evsur ;
  654.    tabt  . itab = 'CHAINE' 'sur;' disc ';ord=' (formar ord 1) ;
  655. *
  656.    evsol = 'EVOL' 'MANU' lh lerrsol ;
  657.    cpl dummy = @POMI evsol 1 'IDEM' ;
  658.    ord = cpl . 1 ;
  659.    tordsol . disc = ord ;
  660.    itab = '+' itab 1 ;
  661.    tabev . itab = evsol ;
  662.    tabt  . itab = 'CHAINE' 'sol;' disc ';ord=' (formar ord 1) ;
  663. 'FIN' iidisc ;
  664. * Tracés
  665. 'SI' graph ;
  666.    evt = @STBL tabev ;
  667.    tit = 'CHAINE' 'Ordre de convergence sur : surface, solution du '
  668.                   'laplacien'  ;
  669.    tix = 'CHAINE' 'Log10 h' ;   
  670.    tiy = 'CHAINE' 'Log10 err' ;
  671.    DESSEVOL evt tabt tit tix tiy ;
  672. 'FINSI' ;
  673. * Vérification des ordres
  674. valvoul = 1.5 ;
  675. valobt  = tordsur . 'LINE' ;
  676. test = valobt '>' valvoul ;
  677. ok = ok 'ET' test ;
  678. 'MESSAGE' ('CHAINE' 'Ordre de convergence sur la mesure de la surface'
  679.                     '(discretisation lineaire)=' valobt) ;
  680. 'SI' ('NON' test) ;
  681.    'MESSAGE' ('CHAINE' 'On aurait voulu avoir :' valvoul) ;
  682. 'FINSI' ;
  683. *
  684. valvoul = 1.5 ;
  685. valobt  = tordsol . 'LINE' ;
  686. test = valobt '>' valvoul ;
  687. ok = ok 'ET' test ;
  688. 'MESSAGE' ('CHAINE' 'Ordre de convergence sur la sol. du laplacien '
  689.                     '(discretisation lineaire)=' valobt) ;
  690. 'SI' ('NON' test) ;
  691.    'MESSAGE' ('CHAINE' 'On aurait voulu avoir :' valvoul) ;
  692. 'FINSI' ;
  693. *
  694. valvoul = 2.5 ;
  695. valobt  = tordsur . 'QUAI' ;
  696. test = valobt '>' valvoul ;
  697. ok = ok 'ET' test ;
  698. 'MESSAGE' ('CHAINE' 'Ordre de convergence sur la mesure de la surface'
  699.                     '(discretisation quadratique)=' valobt) ;
  700. 'SI' ('NON' test) ;
  701.    'MESSAGE' ('CHAINE' 'On aurait voulu avoir :' valvoul) ;
  702. 'FINSI' ;
  703. *
  704. valvoul = 2.5 ;
  705. valobt  = tordsol . 'QUAI' ;
  706. test = valobt '>' valvoul ;
  707. ok = ok 'ET' test ;
  708. 'MESSAGE' ('CHAINE' 'Ordre de convergence sur la sol. du laplacien '
  709.                     '(discretisation quadratique)=' valobt) ;
  710. 'SI' ('NON' test) ;
  711.    'MESSAGE' ('CHAINE' 'On aurait voulu avoir :' valvoul) ;
  712. 'FINSI' ;
  713. *
  714. 'SI' ('NON' ok) ;
  715.    'MESSAGE' ('CHAINE' 'Il y a eu des erreurs') ;
  716.    'ERREUR' 5 ;
  717. 'SINON' ;
  718.    'MESSAGE' ('CHAINE' 'Tout sest bien passe !') ;
  719. 'FINSI' ;
  720. *
  721. 'SI' interact ;
  722.    'OPTION' 'DONN' 5 ;
  723. 'FINSI' ;      
  724. *
  725. * End of dgibi file NLIN_LAPNSSPHE_3D
  726. *
  727. 'FIN' ;
  728.  
  729.  
  730.  
  731.  
  732.  

© Cast3M 2003 - Tous droits réservés.
Mentions légales