Télécharger la plaquette Cast3M


fr en

Télécharger gacul.dgibi

Retour à la liste

Numérotation des lignes : 

  1. * fichier : gacul.dgibi
  2. ************************************************************************
  3. ************************************************************************
  4. GRAPH = FAUX ;
  5. ******************** CAS TEST : gacul.dgibi ************************
  6. *
  7. *
  8. * Test de fonctionnement de DARCYSAT en 2D avec effet de gravite
  9. * ==============================================================
  10. * Infiltration d'eau dans une colonne verticale de sable uniformément
  11. * désaturé.
  12. *
  13. * -  condition initiale : désaturation uniforme correspondant à une
  14. *    succion de 1,15 m ;
  15. * -  a l'instant initial, une extrémité de la colonne est noyée.
  16. *    La frontière est supposée rester a pression nulle ;
  17. *
  18. * ===================================================================
  19. *    Les options de modélisation declarées dans la table transmise à
  20. *    la procédure DARCYSAT sont les suivantes :
  21. *
  22. *
  23. * -  les effets gravitationnels sont pris en compte (indice GRAVITE de
  24. *    valeur VRAI);
  25. *
  26. * -  Une liste de temps de sauvegarde est fournie en valeur de l'indice
  27. *    TEMPS_SAUVES ;
  28. *
  29. * -  Le pas de temps est d'abord automatique (indice TEMPS_CALCULES absent)
  30. *    L'utilisateur fournit
  31. *    > le pas de temps initial (indice 'DT_INITIAL'),
  32. *    > le nombre d'itérations recherché par pas de temps (indice 'NITER')
  33. *    > le nombre de pas de temps (indice 'NPAS')
  34. *
  35. * -  Une sous-relaxation de 0,8 est prise en compte (indice SOUSRELAX)
  36. *
  37. * -  L'homogenéisation spatiale des propriétes physiques s'effectue
  38. *    par moyenne arithmétique des valeurs obtenues aux faces (indice
  39. *    HOMOGENEISATION de valeur DECENTRE)
  40. *
  41. *!!!!!!!!!!!!!!!!!!!!!!!!
  42. * -  La précision de convergence demandée est de 0,5 mm (indice RESIDU_MAX)
  43. * cas test tiré du rapport DMT 97/25 :
  44. * "Implémentation dans CASTEM 2000 d'un modèle de transfert hydrique
  45. *  en milieu poreux non saturé"
  46. *
  47. *
  48. ********************* CAS TEST : gacul.dgibi ************************
  49.  
  50. 'OPTION' 'ECHO' 1 ;
  51. 'SAUTER' 'PAGE';
  52. *
  53. 'TITRE' 'infiltration verticale dans le sable : gacul.dgibi';
  54. 'OPTION' 'DIME' 2 'ELEM' 'QUA4' ;
  55. 'OPTION' 'ISOV' 'LIGN' 'TRACER' 'PSC' ;
  56.  
  57. *--------------------------------------------------------------------
  58. * Création maillage
  59. *
  60. *- Discrétisation :
  61. ENX = 1      ;
  62. ENY = 200    ;
  63. 'DENS' (1./eny) ;
  64.  
  65. *- Création des points et des droites
  66. *
  67. A0 = -.5   1.0D0;  B0 = .5  1.0D0;
  68. A1 = -.5   0.5D0;  B1 = .5  0.5D0;
  69. A2 = -.5   0.25D0;  B2 = .5  0.25D0 ;
  70. A3 = -.5   0.0D0;  B3 = .5  0.0D0 ;
  71.  
  72. *- Création des droites
  73. *
  74. AB0 = 'DROIT' ENX A0 B0 ;
  75. AB1 = 'DROIT' ENX A1 B1 ;
  76. AB2 = 'DROIT' ENX A2 B2 ;
  77. AB3 = 'DROIT' ENX A3 B3 ;
  78.  
  79. *- Creation des surfaces
  80. *
  81. MASSIF0 = AB3 'REGLER' AB2 'REGLER' AB1 'REGLER' AB0 ;
  82. ENY  = 'NBEL' MASSIF0 ;
  83. ENTR = 'COULEUR' ('INVERSE' AB0) 'ROUGE' ;
  84. SORT = 'COULEUR' (AB3) 'ROUGE' ;
  85.  
  86. 'SI' GRAPH  ;
  87.   'TRACER' (MASSIF0 'ET' ENTR 'ET' SORT) ;
  88. 'FINSI' ;
  89.  
  90. *- Creation des maillages contenant tous les points
  91. *
  92. QFTOT  = 'CHANGER' MASSIF0 'QUAF' ;
  93. QFSORT = 'CHANGER' SORT    'QUAF' ;
  94. QFENTR = 'CHANGER' ENTR    'QUAF' ;
  95.  
  96.  
  97. 'ELIMINATION' 0.00001 (QFTOT 'ET' QFSORT 'ET' QFENTR) ;
  98.  
  99. *- Modèles
  100. *
  101. MODHYB  = 'MODELE' QFTOT  'DARCY' 'ISOTROPE'        ;
  102. MODENTR = 'MODELE' QFENTR 'DARCY' 'ISOTROPE'        ;
  103. MODSORT = 'MODELE' QFSORT 'DARCY' 'ISOTROPE'        ;
  104. CEENTR  = 'DOMA'  MODENTR 'CENTRE' ;
  105. CESORT  = 'DOMA'  MODSORT 'CENTRE' ;
  106. HYCEN   = 'DOMA'  MODHYB  'CENTRE' ;
  107. HYFAC   = 'DOMA'  MODHYB  'FACE';
  108.  
  109. *- ordonnées
  110. *
  111. NDIME = 'VALEUR' 'DIME' ;
  112. ZCC   = 'COOR' HYCEN NDIME ;
  113. ZFF   = 'COOR' HYFAC NDIME ;
  114.  
  115. *- Création ligne de suivi pour le post-traitement et le test
  116. * ligne des points centres (cas 1D)
  117. *
  118. 'REPETER' BCL (ENY - 1) ;
  119.     IP = &BCL ;
  120.     JP = IP + 1 ;
  121.     PI = 'POINT' HYCEN IP ;
  122.     PJ = 'POINT' HYCEN JP ;
  123.    'SI' (IP 'EGA' 1);
  124.       LCENC = ('MANU' 'SEG2' PI PJ) ;
  125.    'SINON' ;
  126.       LCENC = LCENC 'ET' ('MANU' 'SEG2' PI PJ) ;
  127.    'FINSI' ;
  128.  FIN BCL ;
  129.  
  130. *--------------------------------------------------------------------
  131. *- pression initiale (metre d'eau) dans le sable
  132. HN    = -1.15 ;
  133.  
  134. *--------------------------------------------------------------------
  135. *- Conditions aux limites
  136.  
  137. *- frontière en limite du domaine de calcul (milieu désaturé)
  138. ESORT = 'NOMC' (('REDU' ZFF CESORT) + HN) 'TH' 'NATURE' 'DISCRET';
  139.  
  140.  
  141. *- frontière mouillée
  142. EENTR = 'NOMC' ('REDU' ZFF CEENTR) 'TH' 'NATURE' 'DISCRET';
  143.  
  144. *- chargement des CLs
  145. LICALC = 'PROG' 0.D0  1.e20 ;
  146. LIST1  = 'PROG' 2 * 1.D0 ;
  147. VALI0  = 'CHAR' 'THIM' (ESORT et EENTR)  ('EVOL' 'MANU' LICALC LIST1) ;
  148.  
  149. *--------------------------------------------------------------------
  150. *- initialisation des inconnues
  151. *  (doit être compatible avec les conditions aux limites)
  152. *
  153.  
  154. *- trace de charge d'eau
  155. z      = 'REDU' (ZFF + HN) ('DIFF' CEENTR HYFAC) ;
  156.  
  157.  
  158. *- charge d'eau
  159. H0     = 'NOMC' (ZCC + HN) 'H' 'NATURE' 'DISCRET';
  160.  
  161. *- flux
  162. QFACE0 = 'MANU' 'CHPO' HYFAC 1 'FLUX' 0.D0 'NATURE' 'DISCRET' ;
  163.  
  164.  
  165. *                                         ---------------------------
  166. *                                         = Table DARCY_TRANSITOIRE =
  167. *                                         ---------------------------
  168. *- initialisation table
  169. SATUR                     = 'TABLE' ;
  170. SATUR. 'TEMPS'            = 'TABLE' ;
  171. SATUR. 'CHARGE'           = 'TABLE' ;
  172. SATUR. 'FLUX'             = 'TABLE' ;
  173.  
  174.  
  175. SATUR . 'ITMAX' = 20;
  176.  
  177. *- données géommétriques
  178. SATUR. 'SOUSTYPE'         = 'DARCY_TRANSATUR' ;
  179. SATUR. 'MODELE'           = MODHYB              ;
  180.  
  181. *- instant initial
  182. SATUR. 'TEMPS' . 0        = 0.                  ;
  183. SATUR. 'CHARGE' . 0       = H0                  ;
  184. SATUR. 'FLUX' . 0         = QFACE0              ;
  185.  
  186. *- conditions aux limites et chargements
  187. *SATUR. 'BLOCAGES_DARCY'   = BSORT 'ET' BENTR      ;
  188. *SATUR. 'CHARGEMENT'       = VALI0  ;
  189. SATUR . 'TRACE_IMPOSE' = VALI0 ;
  190.  
  191.  
  192. SATUR . 'LUMP'   = FAUX;
  193. SATUR . 'TYPDISCRETISATION'   = 'EFMH' ;
  194.  
  195. * GBM MODIFIE
  196.  
  197. TABRES = table METHINV;
  198. TABRES . 'TYPINV' = 2;
  199. TABRES . 'PRECOND' = 5;
  200.  
  201. SATUR . 'METHINV' = TABRES;
  202. SATUR . 'METHINV' . RESID = 1.D-20;
  203.  
  204.  
  205. *- données physiques
  206. * loi de perméabilité
  207. LoiP                      = 'TABLE' 'PUISSANCE';
  208. LoiP.'ALPHA'              = 7.3 ;
  209. LoiP.'PERMSAT'            = 5.85E-2 / 3600. ;
  210. SATUR.'LOI_PERMEABILITE'  = 'TABLE' LoiP ;
  211.  
  212. * loi de succion
  213. LoiS                      = 'TABLE' 'VAN_GENUCHTEN';
  214. LoiS. 'PORO'              = 0.3 ;
  215. LoiS. 'TERESIDU'          = 0.055 ;
  216. LoiS. 'NEXP'              = 2.0304  ;
  217. LoiS. 'MEXP'              = 0.5075 ;
  218. LoiS. 'BHETA'             = 0.029227 * 100. ;
  219. SATUR.'LOI_SATURATION'    = 'TABLE' LoiS ;
  220.  
  221. * gravité (rho_w g,  en mètres d'eau, vaut 1)
  222. SATUR. 'FORCE_GRAVITE'    = 1.;
  223.  
  224. *- données numériques
  225. SATUR. 'SOUS_RELAXATION'  = 1.0D0;
  226. SATUR. 'TEMPS_FINAL'      = 3000. ;
  227. SATUR. 'HOMOGENEISATION'  = 'CHAINE' 'DECENTRE' ;
  228. SATUR. 'NPAS'             = 1000 ;
  229. SATUR. 'RESIDU_MAX'       = 1.D-2 ;
  230. SATUR. 'NITER'            = 12 ;
  231. SATUR. 'DT_INITIAL'       = 5. ;
  232. SATUR. 'TEMPS_SAUVES'     = ('PROG' 1 'PAS' 1 6.)*500. ;
  233. *SATUR. 'MESSAGES'         = 0 ;
  234.  
  235.  
  236. 'SI' GRAPH ;
  237. *--------------------------------------------------------------------
  238. *--- Visualisation de la loi de capacité capillaire.
  239. * droite support des variables
  240.   dx  = 'DROIT' (0. -1.15 ) 1000  (0. 0.) ;
  241.   zc  = 'COOR' dx 2 ;
  242.   ev2 = 'EVOL' 'BLEU' 'CHPO' zc 'SCAL' dx ;
  243.   px  = 'EXTR' ev2 'ORDO' ;
  244. * calcul de la teneur en eau pour la pression zc
  245.   TAB = 'TABLE'  ;
  246.   s0 t0 cap = HT_PRO (SATUR.'LOI_SATURATION') ZC  ;
  247.   ev0 = 'EVOL' 'TURQ' 'CHPO' s0 'SCAL' dx ;
  248.   evt = 'EVOL' 'VERT' 'MANU' px (100. * ('EXTR' ev0 'ORDO')) ;
  249.   'DESSIN' evt 'TITX' 'Pc(m)' 'TITY' 'S(%)'
  250.                'TITRE' 'Loi capillaire S(Pc)' ;
  251. * calcul de la teneur en eau pour la pression zc - dp
  252.   dp  = 1.e-4 ;
  253.   s1 t1 cap = HT_PRO (SATUR.'LOI_SATURATION')
  254.                ('KOPS' zc '-'  dp) ;
  255. * représentation de la capacité
  256.   c1  = (t0 - t1) / dp;
  257.   ev1 = 'EVOL' 'ROUGE' 'CHPO' c1 'SCAL' dx ;
  258.   evc = 'EVOL' 'JAUNE' 'MANU' px ('EXTR' ev1 'ORDO') ;
  259.   'DESSIN' evc 'TITX' 'Pc(m)' 'TITY' 'Capa(1/m)'
  260.                'TITRE' 'Capacite capillaire' ;
  261. 'FINSI' ;
  262.  
  263. *                                           ===========================
  264. *                                           |         CALCUL          |
  265. *                                           ===========================
  266. DARCYSAT SATUR ;
  267.  
  268. * Post-traitement
  269. * ===============
  270. *-- Test de non régression
  271. * Le test est effectué en vérifiant la solution à 1 heure de temps physique
  272.  
  273. *- solution a 1 heure de temps physique
  274. lp = 'PROG' 35*-1.1500 -1.1501 -1.1496 -1.1516 -0.95203 -0.40387
  275.             -0.22041 -0.15351 -0.11533 -8.84629E-02 -6.74391E-02
  276.             -4.98338E-02 -3.43395E-02 -2.01365E-02 -6.63732E-03 ;
  277.  
  278. *- construction de la pression aux centres
  279. SATUR.'PRESSION' = 'TABLE'  ;
  280. 'REPETER' bb (dime SATUR.TEMPS) ;
  281.    id = &bb - 1 ;
  282.    HN = SATUR.'CHARGE'.id ;
  283.    SATUR.'PRESSION'.id = (nomc HN SCAL 'NATURE' 'DISCRET') - ZCC ;
  284. 'FIN' bb ;
  285.  
  286. *ev   = 'EVOL' 'CHPO' SATUR.'PRESSION'. 4 'SCAL' LCENC ;
  287. *lp2  = 'EXTR' ev 'ORDO' ;
  288. *err1 = 'MAXIMUM' ('ABS' (lp - lp2)) ;
  289. *mess 'régression : ' err1 ;
  290.  
  291. *-- Tracé (tous les temps)
  292.   LT        = 'LECT' 0  PAS 1 ((dime SATUR.TEMPS) - 1)  ;
  293.   liopt     = 'MOTS' 'MIMA' 'AXES';
  294.   TDES      = TRACHIS SATUR 'CHARGE' LCENC LT 'PREF' ' ' 'UNIT' 's' ;
  295. 'SI' GRAPH ;
  296.   DESTRA TDES liopt 'TITX' 'z (m)' 'TITY' 'Pw (m)' ;
  297. 'FINSI' ;
  298.  
  299. toto = 'PRIM' TDES . 6 . VALEUR;
  300. inta = 'EXTRAIRE' toto 'ORDO';
  301. inta = 'EXTRAIRE'  ('DIME' inta) inta;
  302.  
  303. err = 'ABS' (inta '+' 0.235 '/' 0.235);
  304. 'LISTE' err;
  305.  
  306. 'SI' (err > 2.D-2) ;
  307.    'ERREUR'  5;
  308. 'FINSI'  ;
  309.  
  310.  
  311. 'FIN';
  312.  
  313.  
  314.  
  315.  
  316.  
  317.  
  318.  
  319.  
  320.  
  321.  

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