Download satnsathoriz.dgibi

Back to the list

* fichier : satnsathoriz.dgibi
'OPTI' 'TRAC' PSC ;
GRAPH = FAUX  ;


********************************************
DEBPROC HHT_pro H1*FLOTTANT ;
T1 = H1 * BET ** N + 1. ** m * TSR + TR ;
FINPROC T1 ;
********************************************

********************************************
DEBPROC KKR_pro H0*FLOTTANT ;
*K0 = H0 * ALF EXP * KS  ;
tt = HHT_PRO H0;
K0 = (tt '-' tr) '/' TSR ** ALF * KS  ;
FINPROC K0 ;
********************************************

********************************************
DEBPROC TH_pro T1*FLOTTANT ;
H1 = T1 - TR / TSR ** M1 - 1. ** N1 / BET ;
FINPROC H1 ;
********************************************



********************************************
DEBPROC DHDT_pro H1*FLOTTANT T1*FLOTTANT  ;
DHDT0 = BET*H1 ** MOINSN + 1. * H1 / (T1 - TR) * N1 * M1 ;
*MESS 'DHDT0' DHDT0;
DHDT0 = 1./BET *N1*M1*1./(TSR) *
        (((T1 - TR / TSR) ** M1  - 1.)
        ** (N1 - 1.)) *((T1 - TR / TSR) ** (M1 - 1.)) ;
*MESS 'DHDT0' DHDT0;
FINPROC DHDT0 ;
********************************************



********************************************
DEBPROC KINPHI int00*FLOTTANT ;
  phi0 = 0. ;
  mess 'int00' int00;
phit = (2* KS * PSI0) / (DT * int00) ;
  'MESS' 'phit' phit ;
*  phi0 = 1.15;
  phi0 = phit;
  int0 = int00 - (phi0/2.);  ;
  Int_p0 = prog int0 ;
  phi_p0 = prog phi0  ;
  repeter b1 (nb-1 );
    d0 = extr D2_prog &b1 ;
    phi0 =  (2. * d0 / int0) + phi0 ;
    phi_p0 = phi_p0 et (prog phi0) ;
    int0 = int0 - phi0 ;
    Int_p0 = Int_p0 et (Prog int0) ;
  fin b1 ;
    d0 = extr D2_prog NB;
    phi0 =  2. * d0 / iNT0 + phi0 ;
    phi_p0 = phi_p0 et (prog phi0) ;
*    list phi_p0 ;
*    list int_p0 ;
     PNM1 = EXTR phi_p0   NB ;
     MESS 'PNM' PNM1 ;
     int = EXTR int_p0   NB ;
     MESS 'int' int ;
FINPROC phi_p0 Int_p0 ;
********************************************



********************************************
DEBPROC KPROG T00*FLOTTANT TN0*FLOTTANT NB0*ENTIER MO1*MOT ;
  DT = T00 - TN0 / NB0 ;

SI (EGA MO1 DEMI ) ;
MESS 'DEMI = ' MO1   ;
  T = T00 - (dt*0.5) ;
  h_pr0 = prog ;
  K_pr0 = prog ;
  D_pr0 = prog ;
  C_pr0 = prog ;
  T_pr0 = prog ;
    H = TH_pro T ;
    K = KKR_pro H ;
    DHDT = DHDT_pro H  T ;
    D = K*DHDT ;
    C = 1./DHDT ;
  h_pr0 = prog H;
  K_pr0 = prog K;
  D_pr0 = prog D;
  C_pr0 = prog C;
  T_pr0 = prog T;
  T = T - DT ;

SINON ;
MESS 'TOTA = ' MO1   ;
    T = T00 ;
    H = TH_pro T ;
    K = KKR_pro H ;
    DHDT = DHDT_pro H  T ;
    D = K*DHDT ;
    C = 1./DHDT ;
  h_pr0 = prog H;
  K_pr0 = prog K;
  D_pr0 = prog D;
  C_pr0 = prog C;
  T_pr0 = prog T;
  T = T - DT ;
FINSI;

  repe b0 nb0 ;
    H = TH_pro T ;
    K = KKR_pro H ;
    DHDT = DHDT_pro H  T ;
    D = K*DHDT ;
    C = 1./DHDT ;
    h_pr0 = h_pr0 ET (prog H) ;
    K_pr0 = K_pr0 ET (prog K) ;
    D_pr0 = D_pr0 ET (prog D) ;
    C_pr0 = C_pr0 ET (prog C) ;
    T_pr0 = T_pr0 ET (prog T);
    T = T - DT ;
  FIN B0 ;
FINPROC h_pr0 K_pr0 D_pr0 C_pr0 T_pr0 ;
********************************************




********************************************
DEBPROC A  x*FLOTTANT  ;
si (x < 2.9 ) ;
X2 = (-1.)*X*X ;
A0 = X2 EXP * X * (1. - (ERF X) ** (-1.)) * RA2PI + (2.*X2) ;
sinon ;
X2 = X ** (-2.) * -0.5 ;
A0 = 108830.*X2+8162.*X2+706.*X2+74.*X2+10.*X2+2.*X2+1.;
*A0 = 8162.*X2+706.*X2+74.*X2+10.*X2+2.*X2+1.;
*A0 = 706.*X2+74.*X2+10.*X2+2.*X2+1.;
*A0 = 74.*X2+10.*X2+2.*X2+1.;
FINSI ;
FINPROC A0 ;
********************************************



********************************************
DEBPROC INTN  phi*LISTREEL ;
PNM1 = EXTR phi   NB ;
DD0 = EXTR D2_PROG NB ;
XX = PNM1 * 0.5D0 * (DD0 ** -0.5D0) ;
IC0 = PNM1*0.5D0 + (2.D0*DD0/PNM1*(A XX)) ;
FINPROC IC0 ;
********************************************



'DEBPROC' CALPHILI ;

h_prog K_prog D_prog C_prog T_prog = KPROG T0 TN NB 'TOTA';

dinf = extr D_prog 1 ;
mess 'dinf' dinf;
dinf = extr D_prog 2 ;
mess 'dinf2' dinf;
T1 = extr T_prog 1 ;
mess 'T1' T1;
TN1 = extr T_prog (DIME T_prog) ;
T01 = extr T_prog 1 ;
db0 = 2. / ((T01 - TN1) ** 2.) ;

Y = T_prog - (prog (DIME T_prog)*TN1) * K_prog ;
KfH_evol = evol manu H_prog Y ;
*-------------------------- Dmoy par K
Db = (-1.) * ((somm KfH_evol) extr 1) * db0 ;
*-------------------------- calcul de I 1/2    (j = 1)
*????????????????????????????????????
Int = Db/pi**0.5*2.*nb ; 'MESS' 'Int' Int ;
*Int = 1.e-3/(ts - tn)*nb  ;list Int ;
*Int = int * 2. ;
*opti donn 5 ;

*-------------------------- calcul de D (I + 1/2)
h2_prog K2_prog D2_prog C2_prog T2_prog = KPROG T0 TN nb DEMI ;


*-------------------------- calcul de PHI (I + 1/2)
phi_prog Int_prog = KINPHI Int ;


I2_evol = evol manu T_prog phi_prog ;
I2mi = I2_evol somm  * (-1.) / dt ; LIST I2mi ;

*I2_evol = evol manu phi_prog T_prog  ;
*DESS I2_evol ;
*opti donn 5 ;

*-------------------------- calcul de ICHA (N - 1/2)
ICHA = INTN phi_prog ;

*-------------------------- DELTA (J=1)
DELTA = (EXTR Int_prog NB) - ICHA ;

*-------------------------- calcul de I 1/2    (j = 2)
int1 = int ;
Int = int1 - (DELTA * 0.5d0) ;
deltI = Int - Int1 ;

repeter b2 20 ;
*-------------------------- calcul de PHI (I + 1/2)
  phi_prog Int_prog = KINPHI Int ;

  'MESSAGE' 'valeur de xphi';
  'LISTE' ('EXTRAIRE' 1 phi_prog);

  I2_evol = evol manu T_prog phi_prog ;
  I2mi = I2_evol somm * (-1.) / dt ;
*-------------------------- calcul de ICHA (N - 1/2)
  ICHA = INTN phi_prog ;
mess
'iter :' &b2 'icha' ICHA 'DELTA' DELTA 'Itot' (maxi I2mi) 'int' int;

*-------------------------- DELTA (J=1)
delta1 = DELTA ;
DELTA = (EXTR Int_prog NB) - ICHA ;

*-------------------------- calcul de I 1/2
int1 = int ;
*int = int1 - (DELTA * delta1 / (DELTA + delta1)) ;
*int = int1 - (DELTA * deltI / (DELTA + delta1)) ;
Int = int1 - (DELTA * 0.5d0 ) ;
deltI = Int - Int1 ;
si ( (DELTA abs) < (ICHA*1.e-4))
quitter  b2 ;
finsi;
menage ;
fin b2 ;


'FINPROC' phi_prog h_prog T_prog ;










'MESS' 'Calcul de la solution de Philip en milieu saturé non saturé';
xf = 0.0;
'MESS' 'xf (frontiére saturé / non saturé) ' xf;
PSI0 = 1.D0;

'MESS' 'Quelle est la pression d injection dans le milieu saturé?' ;
mess 'PSI0'  PSI0;

ALF   = 0.08 * 100. ;
KS    = 5.85E-2 / 3600. ;
TS    = 0.3 ;
TR    = 0.055 ;
TSR   = TS - TR ;
N     = 2.0304  ;   N1 = 1./N ;   MOINSN = (-1.)*N ;
M     = -0.5075 ;   M1 = 1./M ;
BET   = -0.029227 * 100. ;

nb  = 500 ;
H00 = -1.E-20 ;
*H00 = 0;
Hinf = -1.15D0 ;

RA2PI = PI ** (-0.5) * 2. ;

phi0 = 0. ;

HN = HINF ;

TN = HHT_pro HN ; LIST TN ;

HN = TH_pro TN ; LIST HN ;

KN = KKR_pro HN ; LIST KN ;

DHDTN = DHDT_pro HN  TN ; LIST DHDTN ;

DN = KN*DHDTN ;

CN = 1./DHDTN ;

T0 = HHT_pro H00 ;

LIST (HN - (TH_pro (TN - 1.E-6)) / 1.E-6) ;

DT = T0 - TN / NB ;




LTEMPS   = 'PROG' 1000. 'PAS' 1000. (1000. * 3.) ;
ntemps   = 'DIME' LTEMPS ;
lastemps = 0. ;
TETSAT   = TS ;
TETNSAT  = HHT_pro HINF ;
OLDXDTET = (TETSAT - TETNSAT) * xf ;

*
*
* CALCUL DE PHI, H et TETA
*
  phi_prog h_prog t_prog = CALPHILI ;



   xfini = xf;
'REPETER' BTEMPS ntemps ;
   ttemps = 'EXTR' LTEMPS &BTEMPS ;
   DELTEMP = ttemps - lastemps ;
* profile de la charge dans la zone saturée
   LXSAT = 'PROG' 0. ;
   LHSAT = 'PROG' PSI0 ;
   LTETSAT = 'PROG' TETSAT ;


   LXNSAT   = phi_prog * ((TTEMPS)** 0.5) ;

   ID       = 'PROG' ('DIME' LXNSAT)*1. ;
   phif = EXTR phi_prog 1;
   mess 'phif' phif;
   xfr = xfini + (phif*(ttemps**0.5));
   mess 'xfr= 'xfr;
   LXNSAT   = LXNSAT + (ID*xfini) ;
   LHNSAT   = h_prog ;
   LTETNSAT = t_prog ;
*
   LX = LXSAT ET LXNSAT ;
   LTET = LTETSAT ET LTETNSAT ;
   LH = LHSAT ET LHNSAT ;

   'SI' (&BTEMPS 'EGA' 1) ;
*      (LH = 'ORDO' LH CROISSANT);
      uun = 'PROG' ('DIME' LX) * 1.D0;
      courb = ('EVOL' 'MANU'(uun '-'  LX) LH);
   'SINON' ;
*       (LH = 'ORDO' LH CROISSANT) ;
      courb = courb 'ET' ('EVOL' 'MANU' (uun '-' LX) LH);
   'FINSI' ;
*
'FIN' BTEMPS ;


TAB1 = 'TABLE' ;
TAB1 . 'TITRE' = 'TABLE' ;
TAB1 . 'TITRE' . 1 = 'exacte 1000' ;
TAB1 . 'TITRE' . 2 = 'exacte 2000' ;
TAB1 . 'TITRE' . 3 = 'exacte 3000' ;
TAB1 . 1 = 'TIRR ';
TAB1 . 2 = 'TIRR ';
TAB1 . 3 = 'TIRR ';


'SI' (GRAPH);
dess (courb) LEGE TAB1;
'FINSI' ;








******************** CAS TEST : cacul.dgibi ************************
*
*
* Test de fonctionnement de DARCYSAT en 2D sans effet de gravité.
* ===================================================================
* Infiltration d'eau dans une colonne horizontale de sable uniformément
* désaturé.
*
* -  condition initiale : désaturation uniforme correspondant à une
*    succion de 1,15 m ;
* -  a l'instant initial, une extrémité de la colonne est noyée.
*    La frontière est supposée rester a pression nulle ;
*
* ===================================================================
*    Les options de modélisation declarées dans la table transmise à
*    la procédure DARCYSAT sont les suivantes :
*
*
* -  les effets gravitationnels ne sont pas pris en compte : pression =
*    charge (indice GRAVITE absent ; valeur par défaut : FAUX);
*
* -  Une liste de temps de sauvegarde est fournie en valeur de l'indice
*    TEMPS_SAUVES ;
*
* -  Le pas de temps est d'abord automatique (indice TEMPS_CALCULES absent)
*    L'utilisateur fournit
*    > le pas de temps initial (indice 'DT_INITIAL'),
*    > le nombre d'itérations recherché par pas de temps (indice 'NITER')
*    > le nombre de pas de temps (indice 'NPAS')
*
* -  Le calcul est ensuite refais avec une liste de temps de calcul
*    donnée à l'indice TEMPS_CALCULES ;
*
* -  L'homogenéisation spatiale des propriétes physiques s'effectue
*    par moyenne arithmétique des valeurs obtenues aux faces (indice
*    HOMOGENEISATION de valeur DECENTRE)
*
*!!!!!!!!!!!!!!!!!!!!!!!!
* -  La précision de convergence demandée est de 0,5 mm (indice RESIDU_MAX)
* cas test tiré du rapport DMT 97/25 :
* "Implémentation dans CASTEM 2000 d'un modèle de transfert hydrique
*  en milieu poreux non saturé"
*
*
********************* CAS TEST : cacul.dgibi ************************

'OPTION' 'ECHO' 1 ;
'SAUTER' 'PAGE';
*

'TITRE' 'infiltration horizontale dans le sable : cacul.dgibi'   ;
'OPTION' 'DIME' 2 'ELEM' 'QUA4' ;
*'OPTION' 'ISOV' 'LIGN' ;
*'TRACER' 'PSC' ;

************************************************************************
*  Fonction qui calcule la perméabilité en fonction de la
*  saturation réduite.
*
*  Loi puissance :
*  Perméabilité : K = Ks S^B
*
*  Paramètres physiques à définir dans la table LOI :
*  --------------------------------------------------
*    ALPHA   : coef. B (s.d.)
*    PERMSAT : coef. Ks, perméabilité à saturation (m/s)
*
*  Entrée :
*  --------
*    LOI     : table de données décrite ci-dessus
*    SAT     : saturation réduite
*    TEST    : mot facultatif 'NOTEST' qui permet de shunter les tests.
*
*  Sortie :
*  --------
*    K1      : perméabilité totale en eau (m/s)
*
************************************************************************


*--------------------------------------------------------------------
* Création maillage
*
*- Discrétisation :
ENX = 1      ;
ENY = 50    ;
'DENS' (1./eny) ;

*- Création des points et des droites
*
A0 = -0.5   1.0D0;   B0 = 0.5   1.0D0;
'DENS' (1./eny) ;
A1 = -0.5   0.5D0;   B1 = 0.5   0.5D0;
A2 = -0.5   0.25D0 ;  B2 = 0.5   0.25D0 ;
A3 = -0.5   0.0D0;   B3 = 0.5   0.0D0 ;

*- Création des droites
*
AB0 = 'DROIT' ENX A0 B0 ;
AB1 = 'DROIT' ENX A1 B1 ;
AB2 = 'DROIT' ENX A2 B2 ;
AB3 = 'DROIT' ENX A3 B3 ;

*- Creation des surfaces
*
MASSIF0 = AB3 'REGLER' AB2 'REGLER' AB1 'REGLER' AB0 ;
ENY  = 'NBEL' MASSIF0 ;
ENTR = 'COULEUR' ('INVERSE' AB0) 'ROUGE' ;
SORT = 'COULEUR' (AB3) 'ROUGE' ;

'SI' GRAPH  ;
  'TRACER' (MASSIF0 'ET' ENTR 'ET' SORT) ;
'FINSI' ;

*- Creation des maillages contenant tous les points
*
QFTOT  = 'CHANGER' MASSIF0 'QUAF' ;
QFSORT = 'CHANGER' SORT    'QUAF' ;
QFENTR = 'CHANGER' ENTR    'QUAF' ;


'ELIMINATION' 0.00001 (QFTOT 'ET' QFSORT 'ET' QFENTR) ;

*- Modèles
*
MODHYB  = 'MODELE' QFTOT  'DARCY' 'ISOTROPE'        ;
MODENTR = 'MODELE' QFENTR 'DARCY' 'ISOTROPE'        ;
MODSORT = 'MODELE' QFSORT 'DARCY' 'ISOTROPE'        ;
CEENTR  = 'DOMA'  MODENTR 'CENTRE' ;
CESORT  = 'DOMA'  MODSORT 'CENTRE' ;
HYCEN   = 'DOMA'  MODHYB  'CENTRE' ;
HYFAC   = 'DOMA'  MODHYB  'FACE';

*- Création ligne de suivi pour le post-traitement et le test
* ligne des points centres (cas 1D)
*
'REPETER' BCL (ENY - 1) ;
    IP = &BCL ;
    JP = IP + 1 ;
    PI = 'POINT' HYCEN IP ;
    PJ = 'POINT' HYCEN JP ;
   'SI' (IP 'EGA' 1);
      LCENC = ('MANU' 'SEG2' PI PJ) ;
   'SINON' ;
      LCENC = LCENC 'ET' ('MANU' 'SEG2' PI PJ) ;
   'FINSI' ;
 FIN BCL ;

*--------------------------------------------------------------------
*- pression initiale (metre d'eau) dans le sable
HN    = -1.15 ;

*--------------------------------------------------------------------
*- Conditions aux limites


* GBM modifie supprime bloque
*- frontière en limite du domaine de calcul (milieu désaturé)
ESORT = 'MANU' 'CHPO' CESORT 1 'TH' HN 'NATURE' 'DISCRET';

*- frontière mouillée
EENTR =  'MANU' 'CHPO' CEENTR 1 'TH' 1. 'NATURE' 'DISCRET';


*- chargement des CLs
LICALC = 'PROG' 0.D0  1.e20 ;
LIST1  = 'PROG' 2 * 1.D0 ;
VALI0  = 'CHAR' 'THIM' (ESORT et EENTR)  ('EVOL' 'MANU' LICALC LIST1) ;

*--------------------------------------------------------------------
*- initialisation des inconnues
*  (doit être compatible avec les conditions aux limites)
*

*- trace de charge d'eau
*TH0  = 'MANU' 'CHPO' ('DIFF' CEENTR HYFAC) 1 'TH' HN 'NATURE' 'DISCRET';
*TH0  = TH0 + SIMPE ;

*- charge d'eau
H0   = 'MANU' 'CHPO' HYCEN 1 'H' HN 'NATURE' 'DIFFUS';
H0 = 'MASQUE' ('COORDONNEE' 2 ('DOMA' modhyb centre))
             SUPERIEUR 0.98D0;
H0 = 1 '-' H0;
H0 = (HN * H0) + ((1.D0 - H0) * 1.D0);
H0 = 'NOMC' 'H' H0 NATURE DISCRET;
'TRACER' ('KCHA' MODHYB H0 'CHAM') modhyb;


*- flux
QFACE0 = 'MANU' 'CHPO' HYFAC 1 'FLUX' 0.D0 'NATURE' 'DISCRET' ;


*                                         ---------------------------
*                                         = Table DARCY_TRANSITOIRE =
*                                         ---------------------------
*- initialisation table
SATUR                     = 'TABLE' ;
SATUR. 'TEMPS'            = 'TABLE' ;
SATUR. 'CHARGE'           = 'TABLE' ;
SATUR. 'FLUX'             = 'TABLE' ;

SATUR . 'ITMAX' = 45;


*- données géommétriques
SATUR. 'SOUSTYPE'         = 'DARCY_TRANSATUR' ;
SATUR. 'MODELE'           = MODHYB              ;

*- instant initial
SATUR. 'TEMPS' . 0        = 0.                  ;
SATUR. 'CHARGE' . 0       = 'COPIER' H0                  ;
SATUR. 'FLUX' . 0         = 'COPIER' QFACE0              ;

*********************************************
* GBM
*- conditions aux limites et chargements
*SATUR. 'BLOCAGES_DARCY'   = BSORT 'ET' BENTR      ;
*SATUR. 'CHARGEMENT'       = VALI0  ;
SATUR . 'TRACE_IMPOSE' = VALI0 ;

SATUR . 'LUMP'   = FAUX;
SATUR . 'TYPDISCRETISATION'   = 'VF' ;



* GBM MODIFIE

TABRES = table METHINV;
TABRES . 'TYPINV' = 1;
TABRES . 'PRECOND' = 5;

SATUR . 'METHINV' = TABRES;
SATUR . 'METHINV' . RESID = 1.D-20;

*******************************************



*******************************************



*- données physiques
* loi de perméabilité
LoiP                      = 'TABLE' 'PUISSANCE';
LoiP. 'ALPHA'             = 8. ;
LoiP. 'PERMSAT'           = 5.85E-2 / 3600. ;
SATUR.'LOI_PERMEABILITE'  = 'TABLE' LoiP ;

* loi de succion
LoiS                      = 'TABLE' 'VAN_GENUCHTEN';
LoiS. 'PORO'              = 0.3 ;
LoiS. 'TERESIDU'          = 0.055 ;
LoiS. 'NEXP'              = 2.0304  ;
LoiS. 'MEXP'              = 0.5075 ;
LoiS. 'BHETA'             = 0.029227 * 100. ;
SATUR.'LOI_SATURATION'    = 'TABLE' LoiS ;

SATUR . 'COEF_EMMAGASINEMENT' = 0.D-6;
*SATUR . 'XI' = 1.D-5;

*- données numériques
SATUR. 'TEMPS_FINAL'     = 2000.D0 ;
SATUR. 'HOMOGENEISATION' = 'CHAINE' 'DECENTRE' ;


 SATUR.'SOUS_RELAXATION' = 1.;
SATUR. 'NPAS'            = 10000 ;
SATUR. 'RESIDU_MAX'      = 1.D-4 ;
SATUR. 'NITER'           = 10 ;
SATUR. 'DT_INITIAL'      = 0.1D0 ;
SATUR. 'TEMPS_SAUVES'    = ('PROG' 1 'PAS' 1 2.)*1000. ;
*SATUR. 'TEMPS_CALCULES'    = ('PROG' 1 'PAS' 1 3.)*1000. ;
SATUR. 'MESSAGES'        = 2 ;



*--------------------------------------------------------------------
*-  Vérification du choix du dp minimum pour le calcul de la capacité.
* droite support des variables
dx  = 'DROIT' (0. -1.15 ) 1000  (0. 0.) ;
zc  = 'COOR' dx 2 ;
ev2 = 'EVOL' 'BLEU' 'CHPO' zc 'SCAL' dx ;
px  = 'EXTR' ev2 'ORDO' ;
* calcul de la teneur en eau pour la pression zc
s0 t0 cap = HT_PRO (SATUR.'LOI_SATURATION') ZC  ;
ev0 = 'EVOL' 'TURQ' 'CHPO' s0 'SCAL' dx ;
evt = 'EVOL' 'VERT' 'MANU' px (100. * ('EXTR' ev0 'ORDO')) ;

'SI' GRAPH ;
  'DESSIN' evt 'TITX' 'Pc(m)' 'TITY' 'S(%)'
               'TITRE' 'Loi capillaire S(Pc)' ;
'FINSI'  ;

* calcul de la teneur en eau pour la pression zc - dp
dp  = 1.e-4 ;
s1 t1 cap = HT_PRO (SATUR.'LOI_SATURATION')
              ('KOPS' zc '-'  dp) ;
* représentation de la capacité
c1  = (t0 - t1) / dp;
ev1 = 'EVOL' 'ROUGE' 'CHPO' c1 'SCAL' dx ;
evc = 'EVOL' 'JAUNE' 'MANU' px ('EXTR' ev1 'ORDO') ;

'SI' GRAPH ;
  'DESSIN' evc 'TITX' 'Pc(m)' 'TITY' 'Capa(1/m)'
               'TITRE' 'Capacite capillaire' ;
'FINSI' ;

*                                           ===========================
*                                           |       1er CALCUL        |
*                                           ===========================
*- fonctionnement avec temps automatiques



DARCYSAT SATUR ;


* Post-traitement
* ===============

  LT        = 'LECT' 0  PAS 1 ((dime SATUR.TEMPS) - 1)  ;
  liopt     = 'MOTS' 'MIMA' 'AXES';
  TDES      = TRACHIS SATUR 'CHARGE' LCENC LT 'PREF' ' ' 'UNIT' 's' ;




* le fichier lu est produit par resanalyphilhorsatnsat.dgibi


  LT        = 'LECT' 1 'PAS' 1 (('DIME' SATUR.'TEMPS') - 1)  ;
  liopt     = 'MOTS' 'MIMA' 'AXES';
  TDES      = TRACHIS SATUR 'CHARGE' LCENC LT 'PREF' ' ' 'UNIT' 's' ;
  TDES . 3 = 'TABLE' ;
  TDES . 3 . VALEUR = 'EXTRAIRE'  courb 'COUR' 1;
  TDES . 3 . LEGEND2 = TAB1 . TITRE . 1;
  TDES . 3 . LEGEND1 = ' ' ;
  TDES . 4 = 'TABLE' ;
  TDES . 4 . VALEUR = 'EXTRAIRE'  courb 'COUR' 2;
  TDES . 4 . LEGEND2 = TAB1 . TITRE . 2;
  TDES . 4 . LEGEND1 = ' ' ;
*  TDES . 5 = 'TABLE' ;
*  TDES . 5 . VALEUR = 'EXTRAIRE'  courb 'COUR' 3;
*  TDES . 5 . LEGEND2 = TAB1 . TITRE . 3;
*  TDES . 5 . LEGEND1 = ' ' ;

'SI' GRAPH ;
  DESTRA TDES liopt 'TITX' 'z (m)' 'TITY' 'Pw (m)' ;
'FINSI' ;


toto = 'PRIM' TDES . 2 . VALEUR;
inta1 = 'EXTRAIRE' toto 'ORDO';
inta1 = 'EXTRAIRE'  ('DIME' inta1) inta1;

'LISTE' inta1;

err1 = 'ABS' (inta1);
'LISTE' err1;

'SI' ((err1 '-' 2.95D-2) > 1.D-3) ;
   'ERREUR'  5;
'FINSI'  ;


*opti donn 5;

*                                           ===========================
*                                           |       2nd  CALCUL       |
*                                           ===========================
*- fonctionnement avec une liste de temps calculés
SATUR. 'TEMPS'            = 'TABLE' ;
SATUR. 'TRACE_CHARGE'     = 'TABLE' ;
SATUR. 'CHARGE'           = 'TABLE' ;
SATUR. 'FLUX'             = 'TABLE' ;
SATUR. 'TEMPS' . 0        = 0.                  ;
SATUR. 'TRACE_CHARGE'  = TABLE                 ;
SATUR. 'CHARGE' . 0       = H0                  ;
SATUR. 'FLUX' . 0         = QFACE0              ;
SATUR . 'TYPDISCRETISATION'   = 'EFMH' ;




DARCYSAT SATUR ;

* Post-traitement
* ===============
*-- Test de non régression
* Le test est effectué en vérifiant la solution

*- solution après 10 pas de temps


*-- Tracé (tous les temps)

  LT        = 'LECT' 0 'PAS' 1 (('DIME' SATUR.'TEMPS') - 1)  ;
  liopt     = 'MOTS' 'MIMA' 'AXES';
  TDES      = TRACHIS SATUR 'CHARGE' LCENC LT 'PREF' ' ' 'UNIT' 's' ;
    LT        = 'LECT' 1 'PAS' 1 (('DIME' SATUR.'TEMPS') - 1)  ;
  liopt     = 'MOTS' 'MIMA' 'AXES';
  TDES      = TRACHIS SATUR LCENC 'CHARGE' LT 'PREF' ' ' 'UNIT' 's' ;
  TDES . 3 = 'TABLE' ;
  TDES . 3 . VALEUR = 'EXTRAIRE'  courb 'COUR' 1;
  TDES . 3 . LEGEND2 = TAB1 . TITRE . 1;
  TDES . 3 . LEGEND1 = ' ' ;
  TDES . 4 = 'TABLE' ;
  TDES . 4 . VALEUR = 'EXTRAIRE'  courb 'COUR' 2;
  TDES . 4 . LEGEND2 = TAB1 . TITRE . 2;
  TDES . 4 . LEGEND1 = ' ' ;
'SI' GRAPH ;
   DESTRA TDES liopt 'TITX' 'z (m)' 'TITY' 'Pw (m)' ;
'FINSI' ;



toto = 'PRIM' TDES . 2 . VALEUR;
inta2 = 'EXTRAIRE' toto 'ORDO';
inta2 = 'EXTRAIRE'  ('DIME' inta2) inta2;

'MESSAGE' inta1 inta2;

err2 = 'ABS' (inta2);
'LISTE' err2;

'SI' ((err2 '-' 1.22D-2) > 1.D-3) ;
   'ERREUR'  5;
'FINSI'  ;

'MESSAGE' err1 err2;

'FIN';







 

© Cast3M 2003 - All rights reserved.
Disclaimer