Télécharger transsat.dgibi

Retour à la liste

Numérotation des lignes :

  1. * fichier : transsat.dgibi
  2. GRAPH = vrai ;
  3. ******************** CAS TEST : transsat.dgibi ************************
  4. *
  5. *
  6. * Test de fonctionnement de DARCYSAT sur un probleme multizone.
  7. * ===================================================================
  8. * Infiltration d'eau dans une barrière ouvragée dans son site d'acceuil.
  9. *
  10. * - condition initiale :
  11. * > BO : désaturation uniforme correspondant a une succion de 7827,42 m ;
  12. * > SITE saturé sous 500m d'eau
  13. * - les paramètres des lois de perméabilité diffèrent entre le SITE et la BO
  14. *
  15. * ===================================================================
  16. * Les options de modélisation déclarées dans la table transmise à
  17. * la procedure DARCYSAT sont les suivantes :
  18. *
  19. * - les effets gravitationnels ne sont pas pris en compte : pression =
  20. * charge (indice GRAVITE absent ; valeur par defaut : FAUX);
  21. *
  22. * - Le pas de temps est automatique (indice TEMPS_CALCULES absent)
  23. * sont fournis :
  24. * > le pas de temps initial (indice 'DT_INITIAL'),
  25. * > le nombre d'itérations recherchées par pas de temps (indice 'NITER')
  26. * > le nombre de pas de temps maximum (indice 'NPAS')
  27. * > le nombre de CFL initial (indice 'CFL')
  28. *
  29. * - La méthode de pénalisation est désactivée (indice XI absent) ;
  30. *
  31. * - Homogénéisation : les propriétés physiques sont calculées à partir
  32. * des charges aux éléments (indice HOMOGENEISATION de valeur CENTRE)
  33. *
  34. *
  35. ************************************************************************
  36. 'OPTION' 'ECHO' 1 ;
  37. 'SAUTER' 'PAGE';
  38. *
  39.  
  40. TITRE 'Saturation d une barriere ouvragee : transsat.dgibi' ;
  41. 'OPTION' 'DIME' 2 'ELEM' 'QUA4';
  42. 'OPTION' 'ISOV' 'LIGN' 'TRACER' 'PSC' ;
  43.  
  44. *--------------------------------------------------------------------
  45. * Création maillage
  46. *
  47. *- Discrétisation :
  48. ENX = 1 ;
  49. ENY = 100 ;
  50. 'DENS' (3./eny) ;
  51.  
  52. *- Création des points et des droites
  53. *
  54. A0 = -.5 6.D0 ; B0 = .5 6.D0 ;
  55. A1 = -.5 4.D0 ; B1 = .5 4.D0 ;
  56. A2 = -.5 1.D0 ; B2 = .5 1.D0 ;
  57. A3 = -.5 0.D0 ; B3 = .5 0.D0 ;
  58.  
  59. *- Création des droites
  60. *
  61. AB0 = 'DROIT' ENX A0 B0 ;
  62. AB1 = 'DROIT' ENX A1 B1 ;
  63. AB2 = 'DROIT' ENX A2 B2 ;
  64. AB3 = 'DROIT' ENX A3 B3 ;
  65.  
  66. *- Création des surfaces
  67. *
  68. BO = 'REGLER' AB3 AB2 ;
  69. BO = 'COULEUR' BO 'VERT' ;
  70. BO = 'INVERSE' BO ;
  71. SITE = AB2 'REGLER' DINI (3./eny) DFIN 0.2 AB1 'REGLER' 10 AB0 ;
  72. SITE = 'COULEUR' SITE 'BLEU' ;
  73. SITE = 'INVERSE' SITE ;
  74. MASSIF0 = BO 'ET' SITE;
  75.  
  76. * entrée
  77. ENTR = 'COULEUR' ('INVERSE' AB0) 'ROUGE' ;
  78.  
  79. 'SI' GRAPH ;
  80. 'TRACER' (MASSIF0 'ET' ENTR) ;
  81. 'FINSI' ;
  82.  
  83. *- Création des maillages contenant tous les points
  84. *************************************************************************
  85. * IMPORTANT : ne pas concaténer de maillages QUAF *
  86. *************************************************************************
  87. *
  88. QMASSIF = 'CHANGER' MASSIF0 'QUAF' ;
  89. QSITE = 'CHANGER' SITE 'QUAF' ;
  90. QBO = 'CHANGER' BO 'QUAF' ;
  91. QFENTR = 'CHANGER' ENTR 'QUAF' ;
  92.  
  93. 'ELIMINATION' (QMASSIF 'ET' QFENTR 'ET' QBO 'ET' QSITE) 0.00001 ;
  94.  
  95. *- Modèles
  96. *************************************************************************
  97. * IMPORTANT : Il faut concaténer les modèles et ne surtout pas créer un *
  98. * modèle total à une seule zone *
  99. *************************************************************************
  100. MBO = 'MODELE' QBO 'DARCY' 'ISOTROPE' ;
  101. MSITE = 'MODELE' QSITE 'DARCY' 'ISOTROPE' ;
  102. MENTR = 'MODELE' QFENTR 'DARCY' 'ISOTROPE' ;
  103. MODHYB = MSITE 'ET' MBO;
  104.  
  105. CEENTR = 'DOMA' MENTR 'CENTRE' ;
  106. HYCEN = 'DOMA' MODHYB 'CENTRE' ;
  107. HYFAC = 'DOMA' MODHYB 'FACE';
  108.  
  109. *- Création ligne de suivi pour le post-traitement
  110. * ligne des points centres (cas 1D)
  111. *
  112. * on ordonnance de bas en haut les points centres
  113. orig = 0. 0. ;
  114. lcen1 = 'ORDO' (orig 'ET' hycen) ;
  115. lcen2 = 'DIFF' lcen1 ('MANU' 'POI1' orig) ;
  116. * et on crée la ligne de seg2
  117. lcenc = 'QUELCONQUE' lcen2 'SEG2' ;
  118.  
  119. *- ordonnées
  120. *
  121. NDIME = 'VALEUR' 'DIME' ;
  122. ZCC = 'COOR' HYCEN NDIME ;
  123. ZFF = 'COOR' HYFAC NDIME ;
  124.  
  125. *--------------------------------------------------------------------
  126. *
  127. * Paramètres
  128. *
  129. unjour = 24 * 3600. ;
  130. unmois = 30 * unjour ;
  131. unan = unjour * 365.25 ;
  132. Tini = 0. ;
  133. Tfin = 1000 * unan ;
  134.  
  135. *--- profondeur du site
  136. Profsite = 500. ;
  137. T = 50 + 273.15 ;
  138. hc_bo = -7827.42 ;
  139. RHOg = 9810. ;
  140. poro_bo = 0.3 ;
  141.  
  142. *--- coefficient d'emmagasinement
  143. * C = poro * rhoeau * g *
  144. * (compressibilite_eau - compressibilite_solide + alpha/poro)
  145. * avec alpha = 3 (1 - 2 Nupoisson) / Eyoung
  146. * pour FOCA : nu = 0.25 ; E = 500 MPa
  147. * poro = 0.3 ; RHOg = 9810 ; compres_eau = 5E-10 Pa-1 ; compres_solide = 0
  148. emma_sit = 3.83E-6 ;
  149. nu_bo = 0.4 ;
  150. Young_bo = 150.E6 * 2 * (1 + nu_bo) ;
  151. alph_bo = 3 * (1 - (2 * nu_bo)) / Young_bo ;
  152. emma_bo = 1.D0 *poro_bo * RHOg * (5E-10 - 0 + (alph_bo/poro_bo)) ;
  153. emma0 = ('MANU' 'CHPO' ('DOMA' msite 'CENTRE') 1
  154. 'SCAL' emma_sit 'NATURE' 'DISCRET')
  155. 'ET' ('MANU' 'CHPO' ('DOMA' mbo 'CENTRE') 1
  156. 'SCAL' emma_bo 'NATURE' 'DISCRET') ;
  157. emma0 = 'KCHT' modhyb 'SCAL' 'CENTRE' 'COMP' 'SCAL' emma0 ;
  158.  
  159.  
  160. *--------------------------------------------------------------------
  161. * Conditions aux limites
  162. * 1) flux nul sur la frontiere du colis => rien a faire
  163. * 2) flux nuls en latéral => rien a faire
  164. * 3.1) pression imposee en limite de site Pliq = h avec h =500m
  165. BENTR = 'BLOQUE' CEENTR 'TH' ;
  166. EENTR = 'MANUEL' 'CHPO' CEENTR 1 'TH' Profsite ;
  167.  
  168. * chargement des CLs
  169. LICALC = 'PROG' 0.D0 1.e11 ;
  170. LIST1 = 'PROG' 2 * 1.D0 ;
  171. VALI0 = 'CHAR' 'THIM' EENTR ('EVOL' 'MANU' LICALC LIST1) ;
  172.  
  173. *--------------------------------------------------------------------
  174. *- initialisation des inconnues
  175. * (doit être compatible avec les conditions aux limites)
  176. *
  177. Thi_site = 'MANU' 'CHPO' ('DOMA' Msite 'FACE') 1
  178. 'TH' Profsite 'NATURE' 'DISCRET' ;
  179. Thi_bo = 'MANU' 'CHPO' ('DOMA' Mbo 'FACE') 1
  180. 'TH' hc_bo 'NATURE' 'DISCRET' ;
  181. Th_ini = 'KCHT' modHYB 'SCAL' 'FACE' 'COMP' 'TH' Thi_site Thi_bo;
  182. *
  183. Hi_site = 'MANU' 'CHPO' ('DOMA' Msite 'CENTRE') 1
  184. 'H' Profsite 'NATURE' 'DISCRET' ;
  185. Hi_bo = 'MANU' 'CHPO' ('DOMA' Mbo 'CENTRE') 1
  186. 'H' hc_bo 'NATURE' 'DISCRET' ;
  187. H_ini = Hi_site 'ET' Hi_bo ;
  188. Flux_ini = 'MANU' 'CHPO' ('DOMA' modhyb 'FACE') 1
  189. 'FLUX' 0. 'NATURE' 'DISCRET' ;
  190.  
  191. * ---------------------------
  192. * = Table DARCY_TRANSITOIRE =
  193. * ---------------------------
  194. *- initialisation table
  195. SATUR = 'TABLE' ;
  196. SATUR. 'TEMPS' = 'TABLE' ;
  197. SATUR. 'TRACE_CHARGE' = 'TABLE' ;
  198. SATUR. 'CHARGE' = 'TABLE' ;
  199. SATUR. 'FLUX' = 'TABLE' ;
  200.  
  201. *- données géométriques
  202. SATUR. 'SOUSTYPE' = 'DARCY_TRANSATUR' ;
  203. SATUR. 'MODELE' = MODHYB ;
  204.  
  205. *- conditions initiales
  206. SATUR. 'TEMPS' . 0 = Tini ;
  207. SATUR. 'TRACE_CHARGE' . 0 = Th_ini ;
  208. SATUR. 'CHARGE' . 0 = H_ini ;
  209. SATUR. 'FLUX' . 0 = Flux_ini ;
  210.  
  211. *- conditions aux limites et chargements
  212. SATUR. 'TRACE_IMPOSE' = VALI0 ;
  213.  
  214.  
  215. *- définition de la loi de perméabilité
  216. SATUR.'LOI_PERMEABILITE' = 'TABLE' 'MULTIZONE' ;
  217.  
  218. * pour le site
  219. LoiP = 'TABLE' 'PUISSANCE';
  220. LoiP . 'ALPHA' = RHOg / 1.D7 ;
  221. LoiP . 'PERMSAT' = 1.D-12 ;
  222. LoiP . 'MODELE' = MSITE ;
  223. SATUR. 'LOI_PERMEABILITE'. 'SITE' = LoiP ;
  224.  
  225. * pour la BO
  226. LoiP = 'TABLE' 'PUISSANCE';
  227. LoiP .'ALPHA' = RHOg / (7.D6);
  228. LoiP .'PERMSAT' = 6.D-14 ;
  229. LoiP .'MODELE' = MBO ;
  230. SATUR.'LOI_PERMEABILITE'. 'BO' = LoiP ;
  231.  
  232. *- définition de la teneuren equ¼/span>
  233. SATUR.'LOI_SATURATION' = 'TABLE' 'MULTIZONE' ;
  234.  
  235. * pour le site
  236. LoiS = 'TABLE' 'VAN_GENUCHTEN';
  237. LoiS . 'PORO' = 0.05 ;
  238. LoiS . 'TERESIDU' = 0. ;
  239. LoiS . 'NEXP' = 1.7 ;
  240. LoiS . 'MEXP' = 1 - (1 / 1.7) ;
  241. LoiS . 'BHETA' = RHOg / 10.D6 ;
  242. LoiS . 'MODELE' = MSITE ;
  243. SATUR.'LOI_SATURATION' . 'SITE' = LoiS ;
  244.  
  245. 'SI' GRAPH ;
  246. *- Vérification du choix du dp minimum pour le calcul de la capacité.
  247. * droite support des variables
  248. dx = 'DROIT' (0. -7827.4 ) 1000 (0. 0.) ;
  249. zc = 'COOR' dx 2 ;
  250. ev2 = 'EVOL' 'BLEU' 'CHPO' zc 'SCAL' dx ;
  251. px = 'EXTR' ev2 'ORDO' ;
  252. * calcul de la teneur en eau pour la pression zc
  253. s0 t0 dum = HT_PRO (SATUR.'LOI_SATURATION'.'SITE') ZC ;
  254. ev0 = 'EVOL' 'TURQ' 'CHPO' s0 'SCAL' dx ;
  255. evt = 'EVOL' 'VERT' 'MANU' px (100. * ('EXTR' ev0 'ORDO')) ;
  256. 'DESSIN' evt 'TITX' 'Pc(m)' 'TITY' 'S(%)'
  257. 'TITRE' 'Loi capillaire S(Pc)' ;
  258. * calcul de la teneur en eau pour la pression zc - dp
  259. dp = 1. ;
  260. s1 t1 dum = HT_PRO (SATUR.'LOI_SATURATION'.'SITE') ('KOPS' zc '-' dp);
  261. * représentation de la capacité
  262. c1 = (t0 - t1) / dp;
  263. ev1 = 'EVOL' 'ROUGE' 'CHPO' c1 'SCAL' dx ;
  264. evc = 'EVOL' 'TURQ' 'MANU' px ('EXTR' ev1 'ORDO') ;
  265. 'DESSIN' evc 'TITX' 'Pc(m)' 'TITY' 'Capa(1/m)'
  266. 'TITRE' 'Capacite capillaire SITE' ;
  267. FINSI ;
  268.  
  269. * pour la BO
  270. LoiS = TABLE 'VAN_GENUCHTEN';
  271. LoiS . 'PORO' = poro_bo ;
  272. LoiS . 'TERESIDU' = 0. ;
  273. LoiS . 'NEXP' = 1. / (1. - 0.27) ;
  274. LoiS . 'MEXP' = 0.27 ;
  275. LoiS . 'BHETA' = RHOg / 17.E6 ;
  276. LoiS . 'MODELE' = MBO ;
  277. SATUR.'LOI_SATURATION' . 'BO' = LoiS ;
  278.  
  279. 'SI' GRAPH ;
  280. *- Vérification du choix du dp minimum pour le calcul de la capacité.
  281. *--- calcul de la teneur en eau pour p = zc
  282. * calcul de la teneur en eau pour la pression zc
  283. s0 t0 dum = HT_PRO (SATUR.'LOI_SATURATION'.'BO') ZC ;
  284. ev0 = 'EVOL' 'TURQ' 'CHPO' s0 'SCAL' dx ;
  285. evt = 'EVOL' 'VERT' 'MANU' px (100. * ('EXTR' ev0 'ORDO')) ;
  286. 'DESSIN' evt 'TITX' 'Pc(m)' 'TITY' 'S(%)'
  287. 'TITRE' 'Loi capillaire S(Pc)' ;
  288. * calcul de la teneur en eau pour la pression zc - dp
  289. dp = 1. ;
  290. s1 t1 dum = HT_PRO (SATUR.'LOI_SATURATION'.'BO') ('KOPS' zc '-' dp);
  291. * représentation de la capacité
  292. c1 = (t0 - t1) / dp;
  293. ev1 = 'EVOL' 'ROUGE' 'CHPO' c1 'SCAL' dx ;
  294. evc = 'EVOL' 'TURQ' 'MANU' px ('EXTR' ev1 'ORDO') ;
  295. 'DESSIN' evc 'TITX' 'Pc(m)' 'TITY' 'Capa(1/m)'
  296. 'TITRE' 'Capacite capillaire B.O.' ;
  297. FINSI ;
  298.  
  299. SATUR. 'COEF_EMMAGASINEMENT' = emma0 ;
  300.  
  301. * données numériques
  302. SATUR . 'ITMAX' = 40;
  303. SATUR. 'TEMPS_FINAL' = 1.e10 ;
  304. SATUR. 'HOMOGENEISATION' = 'CENTRE' ;
  305. SATUR. 'NPAS' = 100 ;
  306. SATUR . 'RESIDU_MAX' = 1.D-4;
  307. SATUR. 'DT_INITIAL' = 12 * 3600. ;
  308. SATUR. 'MESSAGES' = 1 ;
  309.  
  310. SATUR . 'LUMP' = FAUX ;
  311. SATUR . 'TYPDISCRETISATION'= 'EFMH' ;
  312. SATUR . 'DECENTR' = FAUX;
  313. TABRES = table METHINV;
  314. TABRES . 'TYPINV' = 3;
  315. TABRES . 'PRECOND' = 3;
  316.  
  317. SATUR . 'METHINV' = TABRES;
  318. SATUR . 'METHINV' . RESID = 1.D-20;
  319.  
  320. SATUR. 'TEMPS_SAUVES' = ('PROG' 1 'PAS' 1 41.)*36000. ;
  321. SATUR. 'TEMPS_CALCULES' = ('PROG' 1 'PAS' 1 41.)*36000. ;
  322.  
  323.  
  324. * ===============
  325. * Calcul
  326. * ===============
  327. DARCYSAT SATUR ;
  328.  
  329. * ===============
  330. * Post-traitement
  331. * ===============
  332. * Le test est effectué en vérifiant la solution
  333.  
  334. *- solution après 40 pas de temps
  335. lp = 'PROG' -7827.4 -7827.4 -7827.4 -7827.4 -7827.4 -7827.4 -7827.4
  336. -7827.4 -7827.4 -7827.4 -7827.3 -7827.2 -7826.9 -7826.3
  337. -7825.2 -7823.1 -7819.2 -7812.2 -7800.1 -7779.9 -7747.0
  338. -7695.1 -7615.8 -7498.0 -7328.3 -7090.4 -6765.7 -6334.1
  339. -5773.9 -5064.7 -4189.3 -3141.8 -1944.8 -1272.1 -1188.6
  340. -1101.0 -1009.7 -915.24 -818.60 -720.85 ;
  341. lp = lp 'ET' ('PROG'
  342. -623.36 -527.60 -435.04 -346.88 -263.88 -186.13 -113.58
  343. -47.363 14.082 73.448 130.79 185.41 236.60 283.74
  344. 326.27 363.79 396.07 423.06 444.92 462.03 474.90
  345. 484.16 490.51 494.64 497.16 498.60 499.36 499.72
  346. 499.88 499.95 499.98 499.99 500.00 500.00 500.00
  347. 500.00);
  348.  
  349. id = 40 ;
  350. ev = 'EVOL' 'CHPO' SATUR.'CHARGE'.id 'H' LCENC ;
  351. lp2 = 'EXTR' ev 'ORDO' ;
  352. err1 = 'MAXIMUM' ('ABS' (lp - lp2));
  353. mess 'régression : ' err1 ;
  354.  
  355. 'LISTE' lp2 ;
  356.  
  357.  
  358. 'SI' (err1 > 5.D-1);
  359. mess 'erreur anormale après ' (@ARR id 0) ' pas : ' err1 ;
  360. 'ERREUR' 5;
  361. 'SINON' ;
  362. 'ERREUR' 0;
  363. 'FINSI' ;
  364.  
  365. * ===============
  366. * Post-traitement
  367. * ===============
  368.  
  369. 'SI' GRAPH ;
  370. *--- liste des indices des temps postraites
  371. NN = (('DIME' SATUR.TEMPS) - 1) / 3 ;
  372. LT = 'LECT' 0 'PAS' 3 (3 * NN) ;
  373. TDES = TRACHIS SATUR 'SATURATION' LCENC LT
  374. ('MOTS' 'SCAL') ('MOTS' ' ') 'PREF' ' ' ;
  375. DESTRA TDES 'MIMA' 'AXES' 'TITX' 'z (m)' 'TITY' 'S (s.d.)' ;
  376. TDES = TRACHIS SATUR 'CHARGE' LCENC LT
  377. ('MOTS' 'H') ('MOTS' ' ') 'PREF' ' ' ;
  378. DESTRA TDES 'MIMA' 'AXES' 'TITX' 'z (m)' 'TITY' 'Pw (m)' ;
  379. 'FINSI' ;
  380.  
  381. 'FIN';
  382.  
  383.  
  384.  
  385.  
  386.  
  387.  
  388.  
  389.  
  390.  

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