Download nlin_lapnssphe_axi.dgibi

Back to the list

* fichier : nlin_lapnssphe_axi.dgibi
'OPTION' 'ECHO' 0 ;
*BEGINPROCEDUR gmass
************************************************************************
* NOM         : GMASS
* DESCRIPTION : Une matrice de masse
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v2, 14/03/2006, mise à jour NLIN évolué
* VERSION    : v1, 13/05/2004, version initiale
* HISTORIQUE : v1, 13/05/2004, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* Prière de PRENDRE LE TEMPS de compléter les commentaires
* en cas de modification de ce sous-programme afin de faciliter
* la maintenance !
************************************************************************
*
*
'DEBPROC' GMASS ;
'ARGUMENT'     _mt*'MAILLAGE' ;
'ARGUMENT'     _smt/'MAILLAGE' ;
'ARGUMENT'   gdisc*'MOT     ' ;
'ARGUMENT'    nomt*'MOT     ' ;
'ARGUMENT'   disct*'MOT     ' ;
'ARGUMENT'    nomq*'MOT     ' ;
'ARGUMENT'   discq*'MOT     ' ;
'ARGUMENT' coef/'FLOTTANT' ;
   'SI' ('NON' ('EXISTE' coef)) ;
      'ARGUMENT' coef2/'CHPOINT ' ;
      'SI' ('NON' ('EXISTE' coef2)) ;
         'ERREUR' 'Il faut donner un coef FLOTTANT ou CHPOINT' ;
      'SINON' ;
         coef = coef2 ;
         'ARGUMENT' discc*'MOT     ' ;
      'FINSI' ;                  
   'SINON' ;
      discc = 'CSTE' ;
   'FINSI' ;
'ARGUMENT' methgau/'MOT     ' ;
   'SI' ('NON' ('EXISTE' methgau)) ;
      methgau = 'GAU7' ;
   'FINSI' ;
'ARGUMENT' chpop/'CHPOINT' ;
'ARGUMENT' chpod/'CHPOINT' ;
*
   vdim = 'VALEUR' 'DIME' ;
   vmod = 'VALEUR' 'MODE' ;
   idim = 0 ;
   'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'PLANDEFO')) ;
      idim = 2 ;
      iaxi = FAUX ;
   'FINSI' ;
   'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'AXIS')) ;
      idim = 2 ;
      iaxi = VRAI ;
   'FINSI' ;
   'SI' ('ET' ('EGA' vdim 3) ('EGA' vmod 'TRID')) ;
      idim = 3 ;
      iaxi = FAUX ;
   'FINSI' ;
   'SI' ('EGA' vdim 1) ;
      idim = 1 ;
      iaxi = FAUX ;
   'FINSI' ;
   'SI' ('EGA' idim 0) ;
      'ERREUR' ('CHAINE' 'vmod=' vmod ' et vdim=' vdim ' non prevu') ;
   'FINSI' ;
   'SI' iaxi ;
      dprmt = '*' ('COORDONNEE' 1 _mt) ('*' PI 2.D0) ;
   'FINSI' ;      
   numop  = 1 ;
   numder = idim ;
   mmt = 'MOTS' nomt ;
   mmq = 'MOTS' nomq ;
   numvar = 1 ;
   numdat = 1 ;
   numcof = 1 ;
*
   A = ININLIN numop numvar numdat numcof numder ;   
   A . 'VAR' . 1 . 'NOMDDL' = mmt ;
   A . 'VAR' . 1 . 'DISC'   = disct ;
   'SI' ('EXISTE' chpop) ;
      A . 'VAR' . 1 . 'VALEUR'   = chpop ;
   'FINSI' ;   
   A . 'DAT' . 1 . 'NOMDDL' = 'MOTS' 'SCAL' ;
   A . 'DAT' . 1 . 'DISC'   = discc ;
   A . 'DAT' . 1 . 'VALEUR' = coef ;
   A . 'COF' . 1 . 'COMPOR' = 'IDEN' ;
   A . 'COF' . 1 . 'LDAT'   = 'LECT' 1 ;
*
   A . 1 . 1 . 0 = 'LECT' 1 ;
*
   'SI' iaxi ;
      numdat = 1 ;
      numcof = 1 ;
   'SINON' ;      
      numdat = 0 ;
      numcof = 0 ;
   'FINSI' ;      
   B = ININLIN numop numvar numdat numcof numder ;   
   B . 'VAR' . 1 . 'NOMDDL' = mmq ;
   B . 'VAR' . 1 . 'DISC'   = discq ;
   'SI' ('EXISTE' chpod) ;
      B . 'VAR' . 1 . 'VALEUR'   = chpod ;
   'FINSI' ;   
*
   'SI' iaxi ;
      B . 'DAT' . 1 . 'NOMDDL' = 'MOTS' 'SCAL' ;
      B . 'DAT' . 1 . 'DISC'   = gdisc ;
      B . 'DAT' . 1 . 'VALEUR' = dprmt ;
      B . 'COF' . 1 . 'COMPOR' = 'IDEN' ;
      B . 'COF' . 1 . 'LDAT'   = 'LECT' 1 ;
   'FINSI' ;
   'SI' iaxi ;      
      B . 1 . 1 . 0 = 'LECT' 1 ;
   'SINON' ;      
      B . 1 . 1 . 0 = 'LECT' ;
   'FINSI' ;      
*
   'SI' ('EXISTE' _smt) ;
       mgmass = 'NLIN' gdisc _mt _smt A B methgau ;
   'SINON' ;       
       mgmass = 'NLIN' gdisc _mt A B methgau ;
   'FINSI' ;       
*
'RESPRO' mgmass ;
'FINPROC' ;
*
* End of procedure file GMASS
*
*ENDPROCEDUR gmass
*BEGINPROCEDUR glapn
************************************************************************
* NOM         : GLAPN
* DESCRIPTION : Un laplacien scalaire 
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v2, 14/03/2006, mise à jour NLIN évolué
* VERSION    : v1, 13/05/2004, version initiale
* HISTORIQUE : v1, 13/05/2004, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* Prière de PRENDRE LE TEMPS de compléter les commentaires
* en cas de modification de ce sous-programme afin de faciliter
* la maintenance !
************************************************************************
*
*
'DEBPROC' GLAPN ;
'ARGUMENT'     _mt*'MAILLAGE' ;
'ARGUMENT'   gdisc*'MOT     ' ;
'ARGUMENT'    nomt*'MOT     ' ;
'ARGUMENT'   disct*'MOT     ' ;
'ARGUMENT'    nomq*'MOT     ' ;
'ARGUMENT'   discq*'MOT     ' ;
'ARGUMENT' coef/'FLOTTANT' ;
   'SI' ('NON' ('EXISTE' coef)) ;
      'ARGUMENT' coef2/'CHPOINT ' ;
      'SI' ('NON' ('EXISTE' coef2)) ;
         'ERREUR' 'Il faut donner un coef FLOTTANT ou CHPOINT' ;
      'SINON' ;
         coef = coef2 ;
         'ARGUMENT' discc*'MOT     ' ;
      'FINSI' ;                  
   'SINON' ;
      discc = 'CSTE' ;
   'FINSI' ;            
   'ARGUMENT' methgau/'MOT     ' ;
   'SI' ('NON' ('EXISTE' methgau)) ;
      methgau = 'GAU7' ;
   'FINSI' ;      
'ARGUMENT' chpop/'CHPOINT' ;
'ARGUMENT' chpod/'CHPOINT' ;
*
   vdim = 'VALEUR' 'DIME' ;
   vmod = 'VALEUR' 'MODE' ;
   idim = 0 ;
   'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'PLANDEFO')) ;
      idim = 2 ;
      iaxi = FAUX ;
   'FINSI' ;
   'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'AXIS')) ;
      idim = 2 ;
      iaxi = VRAI ;
   'FINSI' ;
   'SI' ('ET' ('EGA' vdim 3) ('EGA' vmod 'TRID')) ;
      idim = 3 ;
      iaxi = FAUX ;
   'FINSI' ;
   'SI' ('EGA' vdim 1) ;
      idim = 1 ;
      iaxi = FAUX ;
   'FINSI' ;
   'SI' ('EGA' idim 0) ;
      'ERREUR' ('CHAINE' 'vmod=' vmod ' et vdim=' vdim ' non prevu') ;
   'FINSI' ;
   'SI' iaxi ;
      dprmt = '*' ('COORDONNEE' 1 _mt) ('*' PI 2.D0) ;
   'FINSI' ;      
* Test bête...
   'SI' ('EGA' ('TYPE' coef) 'CHPOINT ') ;
      mincoef = 'MINIMUM' coef ;
   'SINON' ;
      mincoef = coef ;
   'FINSI' ;            
   'SI' ('<' mincoef 0.D0) ;
      'ERREUR' 'Le coef (une viscosité) doit etre positive' ;
   'FINSI' ;      
*
   'SI' iaxi ;
      lcoef = 'MOTS' 'NURR' 'NUZZ' ;
   'SINON' ;
      'SI' ('EGA' idim 2) ;
         lcoef = 'MOTS' 'NUXX' 'NUYY' ;
      'SINON' ;
         lcoef = 'MOTS' 'NUXX' 'NUYY' 'NUZZ' ;
      'FINSI' ;
   'FINSI' ;      
*
   ltens = FAUX ;
   'SI' ('EGA' ('TYPE' coef) 'CHPOINT') ;
      ncomp = 'DIME' ('EXTRAIRE' coef 'COMP') ;
      'SI' ('EGA' ncomp idim) ;
         ltens = VRAI ;
      'FINSI' ;
   'FINSI' ;
*
   tcoef = 'TABLE' ;
   'SI' ltens ;
      'REPETER' iidim idim ;
         tcoef . &iidim = 'EXCO' ('EXTRAIRE' lcoef &iidim) coef ;
      'FIN' iidim ;
   'SINON' ;               
      'REPETER' iidim idim ;
         tcoef . &iidim = coef ;
      'FIN' iidim ;
   'FINSI' ;
*      
   numop  = idim ;
   numder = idim ;
   mmt = 'MOTS' nomt ;
   mmq = 'MOTS' nomq ;
   numvar = 1 ;
   numdat = idim ;
   numcof = idim ;
*
   A = ININLIN numop numvar numdat numcof numder ;   
   A . 'VAR' . 1 . 'NOMDDL' = mmt ;
   A . 'VAR' . 1 . 'DISC'   = disct ;
   'SI' ('EXISTE' chpop) ;
      A . 'VAR' . 1 . 'VALEUR'   = chpop ;
   'FINSI' ;   
   'REPETER' iidim idim ;
      A . 'DAT' . &iidim . 'NOMDDL' = 'MOTS' 'SCAL' ;
      A . 'DAT' . &iidim . 'DISC'   = discc ;
      A . 'DAT' . &iidim . 'VALEUR' = tcoef . &iidim ;
      A . 'COF' . &iidim . 'COMPOR' = 'IDEN' ;
      A . 'COF' . &iidim . 'LDAT'   = 'LECT' &iidim ;
   'FIN' iidim ;
*
   'REPETER' iidim idim ;
      A . &iidim . 1 . &iidim = 'LECT' &iidim ;
   'FIN' iidim ;      
*
   'SI' iaxi ;
      numdat = 1 ;
      numcof = 1 ;
   'SINON' ;      
      numdat = 0 ;
      numcof = 0 ;
   'FINSI' ;      
*
   B = ININLIN numop numvar numdat numcof numder ;   
   B . 'VAR' . 1 . 'NOMDDL' = mmq ;
   B . 'VAR' . 1 . 'DISC'   = discq ;
   'SI' ('EXISTE' chpod) ;
      B . 'VAR' . 1 . 'VALEUR'   = chpod ;
   'FINSI' ;   
   'SI' iaxi ;
      B . 'DAT' . 1 . 'NOMDDL' = 'MOTS' 'SCAL' ;
      B . 'DAT' . 1 . 'DISC'   = gdisc ;
      B . 'DAT' . 1 . 'VALEUR' = dprmt ;
      B . 'COF' . 1 . 'COMPOR' = 'IDEN' ;
      B . 'COF' . 1 . 'LDAT'   = 'LECT' 1 ;
   'FINSI' ;
*
   'REPETER' iidim idim ;
      'SI' iaxi ;      
         B . &iidim . 1 . &iidim = 'LECT' 1 ;
      'SINON' ;      
         B . &iidim . 1 . &iidim = 'LECT' ;
      'FINSI' ;      
   'FIN' iidim ;      
*
   mglapn = 'NLIN' gdisc _mt A B methgau ;
* Integration par parties
   mglapn = '*' mglapn -1.D0 ;
*
'RESPRO' mglapn ;
'FINPROC' ;
*
* End of procedure file GLAPN
*
*ENDPROCEDUR glapn
*BEGINPROCEDUR formar
************************************************************************
* NOM         : FORMAR
* DESCRIPTION : formate un réel de facon courte
*               pratique pour les noms de
*               sauvegarde
*               Exemples :
* 'MESSAGE' ('CHAINE' (formar 2.9e5   1)) ;
* 2.9E5
* 'MESSAGE' ('CHAINE' (formar -2.9e5  1)) ;
* -2.9E5
* 'MESSAGE' ('CHAINE' (formar 2.9e-5  1)) ;
* 2.9E-5
* 'MESSAGE' ('CHAINE' (formar -2.9e-5 1)) ;
* -2.9E-5
* 'MESSAGE' ('CHAINE' (formar 2.9     1)) ;
* 2.9
* 'MESSAGE' ('CHAINE' (formar -2.9    1)) ;
* -2.9
* 'MESSAGE' ('CHAINE' (formar 0       1)) ;
* 0
* 'MESSAGE' ('CHAINE' (formar 0       1)) ;
* 0
* 'MESSAGE' ('CHAINE' (formar 2.9e5   0)) ;
* 3E5
* 'MESSAGE' ('CHAINE' (formar -2.9e5  0)) ;
* -3E5
* 'MESSAGE' ('CHAINE' (formar 2.9e-5  0)) ;
* 3E-5
* 'MESSAGE' ('CHAINE' (formar -2.9e-5 0)) ;
* -3E-5
* 'MESSAGE' ('CHAINE' (formar 2.9     0)) ;
* 3
* 'MESSAGE' ('CHAINE' (formar -2.9    0)) ;
* -3
* 'MESSAGE' ('CHAINE' (formar 0       0)) ;
* 0
* 'MESSAGE' ('CHAINE' (formar 0       0)) ;
* 0
*               
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, 18/02/2003, version initiale
* HISTORIQUE : v1, 18/02/2003, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* Prière de PRENDRE LE TEMPS de compléter les commentaires
* en cas de modification de ce sous-programme afin de faciliter
* la maintenance !
************************************************************************
*
*
'DEBPROC' FORMAR ;
'ARGUMENT' fl*'FLOTTANT' ;
'ARGUMENT' vir/'ENTIER  ' ;
'SI' ('NON' ('EXISTE' vir)) ;
   vir = 1 ;
'SINON' ;   
   'SI' ('<' vir 0) ;
      'ERREUR' 'fournir un entier positif' ;
   'FINSI' ;
'FINSI' ;         
'SI' ('<' ('ABS' fl) 10.D-100) ;
   chfl = 'CHAINE' '0' ;
'SINON' ;
*! sans le 1.D-10, ca ne fonctionne pas
*! qd on entre pile poil une puissance de 10
   lfl = LOG10 ('ABS' fl) ;
*   lfl = '+' (LOG10 ('ABS' fl)) 1.D-10  ;
   slfl = 'SIGNE' ('ENTIER' lfl) ;
   'SI' ('EGA' slfl 1) ;
      elfl = 'ENTIER' lfl ;
   'SINON' ;
      elfl = '-' ('ENTIER' lfl) 1 ;
   'FINSI' ;
   man = '/' fl ('**' 10.D0 elfl) ;
*
* Une verrue pour des histoires de précision...
*   
   'SI' ('EGA' man 10.D0 ('**' 10.D0 ('*' vir -1.D0))) ;
      man  = '/' man 10.D0 ;
      elfl = '+' elfl 1 ;
   'FINSI' ;
*   
   sman = 'SIGNE' man ;
   'SI' ('EGA' sman 1) ;   
      fman = 'CHAINE' '(F' ('+' vir 2) '.0' vir ')' ;
   'SINON' ;
      fman = 'CHAINE' '(F' ('+' vir 3) '.0' vir ')' ;
   'FINSI' ;
   'SI' ('NEG' vir 0) ;   
      'SI' ('NEG' elfl 0) ;
         chfl = 'CHAINE' 'FORMAT' fman man 'E' elfl ;
      'SINON' ;
         chfl = 'CHAINE' 'FORMAT' fman man ;
      'FINSI' ;
   'SINON' ;
      man2 = 'ENTIER' ('+' man ('*' 0.5D0 sman)) ;
      'SI' ('NEG' elfl 0) ;
         chfl = 'CHAINE' man2 'E' elfl ;
      'SINON' ;
         chfl = 'CHAINE' man2 ;
      'FINSI' ;
   'FINSI' ;            
'FINSI' ;   
'RESPRO' chfl ;
*
* End of procedure file FORMAR
*
'FINPROC' ;
*ENDPROCEDUR formar
*BEGINPROCEDUR log10
************************************************************************
* NOM         : LOG10
* DESCRIPTION : Log_10
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, 18/02/2003, version initiale
* HISTORIQUE : v1, 18/02/2003, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* Prière de PRENDRE LE TEMPS de compléter les commentaires
* en cas de modification de ce sous-programme afin de faciliter
* la maintenance !
************************************************************************
*
*
'DEBPROC' LOG10 ;
'REPETER' bouc ;
   ok = FAUX ;
   'ARGUMENT' fl/'FLOTTANT' ;
   'ARGUMENT' lr/'LISTREEL' ;
   'ARGUMENT' cp/'CHPOINT ' ;
   'ARGUMENT' cm/'MCHAML  ' ;
   'SI' ('EXISTE' fl) ;
      ok = VRAI ;
      'RESPRO' ('/' ('LOG' fl) ('LOG' 10.D0)) ;
   'FINSI' ;
   'SI' ('EXISTE' lr) ;
      ok = VRAI ;
      'RESPRO' ('/' ('LOG' lr) ('LOG' 10.D0)) ;
   'FINSI' ;
   'SI' ('EXISTE' cp) ;
      ok = VRAI ;
      'RESPRO' ('/' ('LOG' cp) ('LOG' 10.D0)) ;
   'FINSI' ;
   'SI' ('EXISTE' cm) ;
      ok = VRAI ;
      'RESPRO' ('/' ('LOG' cm) ('LOG' 10.D0)) ;
   'FINSI' ;
   'SI' ('NON' ok) ;
      'QUITTER' bouc ;
   'FINSI' ;
'FIN' bouc ;         
*
* End of procedure file LOG10
*
'FINPROC' ;
*ENDPROCEDUR log10
'OPTION' 'ECHO' 1 ;
************************************************************************
* NOM         : NLIN_LAPNSSPHE_AXI
* DESCRIPTION : Vérification de NLIN : on résout un laplacien sur un
*               bout de surface de sphère de rayon R en 2D axi. La
*               solution exacte est une fonction du type :
*               a * ln (tan (\theta / 2)) + b
*               Il y a une singularité en theta = 0. Si on veut observer
*               l'ordre correct pour theta petit, il faut bien raffiner.
*
*               c'est une primitive de a / sin (\theta).
*               On vérifie les ordres de convergence sur cette solution
*               et sur la surface calculée du bout de sphère.
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, 23/10/2006, version initiale
* HISTORIQUE : v1, 23/10/2006, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* Prière de PRENDRE LE TEMPS de compléter les commentaires
* en cas de modification de ce sous-programme afin de faciliter
* la maintenance !
************************************************************************
*
*
interact = FAUX ;
graph    = FAUX ;
*interact = VRAI ;
*graph    = VRAI ;
'OPTION' 'DIME' 2 'MODE' 'AXIS' ;
'SI' ('NON' interact) ;
  'OPTION' 'TRAC' 'PS' 'ECHO' 0 ;
'SINON' ;
  'OPTION' 'TRAC' 'X' 'ECHO' 1 ;
'FINSI' ;  
*
* Paramètres
*
ldisc = 'MOTS' 'LINE' 'QUAI' ;
ttamin  = 5. ;
ttaminr = ttamin '*' ('/' PI 180.) ;
tordsur = 'TABLE' ;
tordsol = 'TABLE' ;
ldens = 'PROG' 2.D-1 1.D-1 5.D-2 ;
rayo  = '**' 2. 0.5D0 ;
tabev = 'TABLE' ;
tabt  = 'TABLE' ;
itab = 0 ;
ok = VRAI ;
*
* Boucle sur les discrétisations 
*
'REPETER' iidisc ('DIME' ldisc) ;
   idisc = &iidisc ;
   disc  = 'EXTRAIRE' ldisc idisc ;
   tordsur . disc = 'TABLE' ;
   tordsol . disc = 'TABLE' ;
   discg = disc ;
   disct = disc ;
   'SI' ('EGA' disc 'LINE') ;
      'OPTION' 'ELEM' 'QUA4' ;
   'SINON' ;      
      'OPTION' 'ELEM' 'QUA8' ;
   'FINSI' ;      
   lerrsur = 'PROG' ;
   lerrsol = 'PROG' ;
   'REPETER' iidens ('DIME' ldens) ;
      idens = &iidens ;
      vdens = 'EXTRAIRE' ldens idens ;
      'DENS' vdens ;
      tit = 'CHAINE' 'Discretisation :' disc '  Densité = ' vdens ;
      'SAUTER' 1 'LIGNE' ;
      'MESSAGE' tit ;
* Maillage
      p0 = 'POIN' 0. 0. ;
      pA = POINTCYL rayo 0. ;
      pB = POINTCYL rayo ('-' 90. ttamin) ;
      mt = 'CERCLE' pB p0 pA ;
      _mt = 'CHANGER' mt 'QUAF' ;
* 1 : calcul de la surface du bout de sphère
      ch1 = 'MANUEL' 'CHPO' mt 1 'T' 1. ;
      gmas = GMASS _mt discg 'T' disct 'Q' disct 1. ;
      surana = 2. '*' PI '*' rayo '*' rayo '*' ('COS' ttamin) ;
      sur = xtmx ch1 gmas ;
      'MESSAGE' ('CHAINE' 'Surface calculee =' sur) ;
      'MESSAGE' ('CHAINE' 'Surface attendue =' surana) ;
      errsur = 'ABS' ('-' sur surana) ;
      'MESSAGE' ('CHAINE' 'Erreur sur la surface =' errsur) ;
      lerrsur = 'ET' lerrsur ('PROG' errsur) ;
* 2 : calcul d'une solution analytique
      rmt zmt = 'COORDONNEE' mt ;
      tta = 'ATG' rmt zmt ;
      ttamsur2 = '/' ttamin 2. ;
      coef = 'LOG' ('/' ('SIN' ttamsur2) ('COS' ttamsur2)) ;
      ttasur2 = '/' tta 2.D0 ;
      solana = '/' ('LOG' ('/' ('SIN' ttasur2) ('COS' ttasur2))) coef ;
      solana = 'NOMC' 'T' solana ;
      mat = GLAPN _mt discg 'T' disct 'Q' disct 1. ;
      mblo = 'BLOQUE' 'T' (pA 'ET' pB) ;
      cblo = '+' ('MANUEL' 'CHPO'  pA 1 'T' 0.)
                 ('MANUEL' 'CHPO'  pB 1 'T' 1.) ;
      fblo = 'DEPIMPOSE' mblo cblo ;   
*
      mtot = mat 'ET' mblo ;
      ftot = fblo ;
      solu = 'RESOUD' mtot ftot ;
      sol  = 'EXCO' 'T' solu 'T' ;
      matm = GMASS _mt discg 'T' disct 'Q' disct 1. ;
      errsol = '**' (xtmx ('-' sol solana) matm) 0.5D0 ;
      'MESSAGE' ('CHAINE' 'Erreur sur la sol. analytique =' errsol) ;
      'SI' graph ;
         tt = 'TABLE' ;
         evana = 'EVOL' 'CHPO' solana 'T' mt ;
         tt . 1 = 'CHAINE' 'Exact' ;
         ev    = 'EVOL' 'CHPO' sol 'T' mt ;
         tt . 2 = 'CHAINE' 'Approchée' ;
         evt = evana 'ET' ev ;
         tit = 'CHAINE' 'Solution du laplacien sur le bout de '
                        'surface de sphere' ;
         tix = 'CHAINE' 'theta' ;
         tiy = 'CHAINE' 'T' ;                        
         DESSEVOL evt tt tit tix tiy ;
      'FINSI' ;                
      lerrsol = 'ET' lerrsol ('PROG' errsol) ;
   'FIN' iidens ;
* Calcul des ordres
   lh lerrsur lerrsol = LOG10 ldens lerrsur lerrsol ;
   evsur = 'EVOL' 'MANU' lh lerrsur ;
   cpl dummy = @POMI evsur 1 'IDEM' ;
   ord = cpl . 1 ;
   tordsur . disc = ord ;
   itab = '+' itab 1 ;
   tabev . itab = evsur ;
   tabt  . itab = 'CHAINE' 'sur;' disc ';ord=' (formar ord 1) ;
*
   evsol = 'EVOL' 'MANU' lh lerrsol ;
   cpl dummy = @POMI evsol 1 'IDEM' ;
   ord = cpl . 1 ;
   tordsol . disc = ord ;
   itab = '+' itab 1 ;
   tabev . itab = evsol ;
   tabt  . itab = 'CHAINE' 'sol;' disc ';ord=' (formar ord 1) ;
'FIN' iidisc ;
* Tracés
'SI' graph ;
   evt = @STBL tabev ;
   tit = 'CHAINE' 'Ordre de convergence sur : surface, solution du '
                  'laplacien'  ;
   tix = 'CHAINE' 'Log10 h' ;   
   tiy = 'CHAINE' 'Log10 err' ;
   DESSEVOL evt tabt tit tix tiy ;
'FINSI' ;
* Vérification des ordres
valvoul = 1.5 ;
valobt  = tordsur . 'LINE' ;
test = valobt '>' valvoul ;
ok = ok 'ET' test ;
'MESSAGE' ('CHAINE' 'Ordre de convergence sur la mesure de la surface'
                    '(discretisation lineaire)=' valobt) ;
'SI' ('NON' test) ;
   'MESSAGE' ('CHAINE' 'On aurait voulu avoir :' valvoul) ;
'FINSI' ;
*
valvoul = 1.5 ;
valobt  = tordsol . 'LINE' ;
test = valobt '>' valvoul ;
ok = ok 'ET' test ;
'MESSAGE' ('CHAINE' 'Ordre de convergence sur la sol. du laplacien '
                    '(discretisation lineaire)=' valobt) ;
'SI' ('NON' test) ;
   'MESSAGE' ('CHAINE' 'On aurait voulu avoir :' valvoul) ;
'FINSI' ;
*
valvoul = 2.5 ;
valobt  = tordsur . 'QUAI' ;
test = valobt '>' valvoul ;
ok = ok 'ET' test ;
'MESSAGE' ('CHAINE' 'Ordre de convergence sur la mesure de la surface'
                    '(discretisation quadratique)=' valobt) ;
'SI' ('NON' test) ;
   'MESSAGE' ('CHAINE' 'On aurait voulu avoir :' valvoul) ;
'FINSI' ;
*
valvoul = 2.5 ;
valobt  = tordsol . 'QUAI' ;
test = valobt '>' valvoul ;
ok = ok 'ET' test ;
'MESSAGE' ('CHAINE' 'Ordre de convergence sur la sol. du laplacien '
                    '(discretisation quadratique)=' valobt) ;
'SI' ('NON' test) ;
   'MESSAGE' ('CHAINE' 'On aurait voulu avoir :' valvoul) ;
'FINSI' ;
*
'SI' ('NON' ok) ;
   'MESSAGE' ('CHAINE' 'Il y a eu des erreurs') ;
   'ERREUR' 5 ;
'SINON' ;
   'MESSAGE' ('CHAINE' 'Tout sest bien passe !') ;
'FINSI' ;
*
'SI' interact ;
   'OPTION' 'DONN' 5 ;
'FINSI' ;      
*
* End of dgibi file NLIN_LAPNSSPHE_AXI
*
'FIN' ;
 

© Cast3M 2003 - All rights reserved.
Disclaimer