Télécharger test_kres_lapn.dgibi
* fichier : test_kres_lapn.dgibi ************************************************************************ * NOM : TEST_KRES_LAPN * DESCRIPTION : Test d'options de KRES sur un laplacien * - Scaling ou pas * - solveurs itératifs : CG; BiCGStab, GMRES * - préconditionneurs ilu0, ilut, ilut + pivoting * * LANGAGE : GIBIANE-CAST3M * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF) * mél : gounand@semt2.smts.cea.fr ********************************************************************** * VERSION : v1, 14/02/2005, version initiale * HISTORIQUE : v1, 14/02/2005, 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 ; complet = FAUX ; graph = FAUX ; ************************************************************************ * NOM : MODULO * DESCRIPTION : Calcule un entier modulo un autre... * * 'DEBPROC' MODULO ; 'ARGUMENT' i*'ENTIER' j*'ENTIER' ; 'SI' ('EGA' j 0) ; 'MESSAGE' 'Impossible de faire modulo 0' ; 'ERREUR' 5 ; 'SINON' ; k=i '/' j ; mod=i '-' ( k '*'j ) ; 'RESPRO' mod ; 'FINSI' ; * * End of procedure file MODULO * 'FINPROC' ; ************************************************************************ * NOM : LOG10 * DESCRIPTION : Log_10 * 'DEBPROC' LOG10 ; 'ARGUMENT' fl/'FLOTTANT' ; 'ARGUMENT' lr/'LISTREEL' ; 'ARGUMENT' cp/'CHPOINT ' ; 'ARGUMENT' cm/'MCHAML ' ; 'SI' ('EXISTE' fl) ; 'RESPRO' ('/' ('LOG' fl) ('LOG' 10.D0)) ; 'FINSI' ; 'SI' ('EXISTE' lr) ; 'RESPRO' ('/' ('LOG' lr) ('LOG' 10.D0)) ; 'FINSI' ; 'SI' ('EXISTE' cp) ; 'RESPRO' ('/' ('LOG' cp) ('LOG' 10.D0)) ; 'FINSI' ; 'SI' ('EXISTE' cm) ; 'RESPRO' ('/' ('LOG' cm) ('LOG' 10.D0)) ; 'FINSI' ; * * End of procedure file LOG10 * 'FINPROC' ; ************************************************************************ * NOM : FORMAR * DESCRIPTION : formate un réel de facon courte * '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' ; * * Procédure CMAIL * 'DEBPROC' CMAIL ; 'ARGUMENT' imail*'ENTIER' ; 'ARGUMENT' iquad*'ENTIER' ; 'SI' ('EGA' imail 1) ; 'SI' complet ; nmail = 75 ; 'SINON' ; nmail = 10 ; 'FINSI' ; 'SI' ('EGA' iquad 1) ; nmail = '*' nmail 2 ; 'FINSI' ; pA = 0. 0. ; pB = 1. 0. ; pC = 1. 1. ; pD = 0. 1. ; l1 = 'DROIT' nmail pA pB ; l2 = 'DROIT' nmail pB pC ; l3 = 'DROIT' nmail pC pD ; l4 = 'DROIT' nmail pD pA ; mt = 'DALLER' l1 l2 l3 l4 ; cmt = l1 'ET' l2 'ET' l3 'ET' l4 ; _mt = 'CHANGER' mt 'QUAF' ; _cmt = 'CHANGER' cmt 'QUAF' ; 'FINSI' ; 'SI' ('EGA' imail 2) ; dini = 0.2D0 ; dfin = 0.004D0 ; 'SI' complet ; dini = '/' dini 4.D0 ; dfin = '/' dfin 4.D0 ; 'FINSI' ; 'SI' ('EGA' iquad 1) ; dini = '/' dini 2.D0 ; dfin = '/' dfin 2.D0 ; 'FINSI' ; pA = 0. 0. ; pB = 1. 0. ; pC = 1. 1. ; pD = 0. 1. ; l1 = 'DROIT' pA pB 'DINI' dini 'DFIN' dini ; l2 = 'DROIT' pB pC 'DINI' dini 'DFIN' dfin ; l3 = 'DROIT' pC pD 'DINI' dfin 'DFIN' dini ; l4 = 'DROIT' pD pA 'DINI' dini 'DFIN' dini ; cmt = l1 'ET' l2 'ET' l3 'ET' l4 ; mt = 'SURFACE' cmt ; _mt = 'CHANGER' mt 'QUAF' ; _cmt = 'CHANGER' cmt 'QUAF' ; 'FINSI' ; 'SI' ('EGA' imail 3) ; 'SI' complet ; nmail = 15 ; 'SINON' ; nmail = 5 ; 'FINSI' ; 'SI' ('EGA' iquad 1) ; nmail = '*' nmail 2 ; 'FINSI' ; pA = 0. 0. 0. ; pB = 1. 0. 0. ; pC = 1. 1. 0. ; pD = 0. 1. 0. ; l1 = 'DROIT' nmail pA pB ; l2 = 'DROIT' nmail pB pC ; l3 = 'DROIT' nmail pC pD ; l4 = 'DROIT' nmail pD pA ; smt = 'DALLER' l1 l2 l3 l4 'PLAN' ; cmt = 'ENVELOPPE' mt ; _mt = 'CHANGER' mt 'QUAF' ; _cmt = 'CHANGER' cmt 'QUAF' ; 'FINSI' ; 'ELIMINATION' ('ET' _mt _cmt) 1.D-5 ; 'SI' graph ; mtot = 'ET' _mt ('COULEUR' _cmt 'ROUG') ; 'TRACER' 'CACH' mtot 'TITR' titch ; 'FINSI' ; tabmod = 'TABLE' ; tabmod . 'MT' = 'TABLE' ; tabmod . 'MT' . 'QUAF' = _mt ; tabmod . 'MT' . 'LINE' = mt ; tabmod . 'CMT' = 'TABLE' ; tabmod . 'CMT' . 'QUAF' = _cmt ; tabmod . 'CMT' . 'LINE' = cmt ; 'RESPRO' tabmod ; 'FINPROC' ; * 'DEBPROC' RESOU ; 'ARGUMENT' mat/'RIGIDITE' ; 'SI' ('NON' ('EXISTE' mat) ); 'ARGUMENT' mat*'MATRIK' ; 'SINON' ; mat = 'KOPS' 'CHANINCO' lispri lispri lisdua lispri 'FINSI' ; 'ARGUMENT' clim*'CHPOINT' ; 'ARGUMENT' pcmlag*'ENTIER' ; 'ARGUMENT' scaling*'ENTIER' ; 'ARGUMENT' typinv*'ENTIER' ; 'ARGUMENT' typrec*'ENTIER' ; * rvc = 'TABLE' ; rvm = rv . 'METHINV' ; rvm . 'TYPINV' = typinv '+' 1 ; rvm . 'IMPINV' = 0 ; rvm . 'PCMLAG' = 'CHAINE' 'APR' ('+' pcmlag 1) ; rvm . 'SCALING' = '-' scaling 1 ; 'SI' ('EGA' typrec 1) ; rvm . 'PRECOND' = 3 ; 'SINON' ; 'SI' ('ET' ('>EG' typrec 2) ('<EG' typrec 8)) ; rvm . 'ILUTLFIL' = 2. ; 'SINON' ; rvm . 'ILUTLFIL' = 3. ; 'FINSI' ; typrec = '-' typrec 1 ; typrec = MODULO ('-' typrec 1) 7 ; 'SI' ('EGA' typrec 0) ; rvm . 'PRECOND' = 5 ; 'SINON' ; typrec = '-' typrec 1 ; typrec2 = modulo typrec 2 ; rvm . 'PRECOND' = '+' typrec2 7 ; typrec3 = '/' typrec 3 ; 'SI' ('EGA' typrec3 0) ; rvm . 'ILUTPPIV' = 0. ; 'FINSI' ; 'SI' ('EGA' typrec3 1) ; rvm . 'ILUTPPIV' = 1.D-3 ; 'FINSI' ; 'SI' ('EGA' typrec3 2) ; rvm . 'ILUTPPIV' = 1.D-1 ; 'FINSI' ; 'FINSI' ; 'FINSI' ; * 'TEMPS' 'ZERO' ; tcpu = TABTPS.'TEMPS_CPU'.'INITIAL'; nit = 'EXTRAIRE' (rvm . 'NMATVEC') nit ; * 'RESPRO' solu nit tcpu ; 'FINPROC' ; * * Procédure CALCUL * 'DEBPROC' CALCUL ; 'ARGUMENT' imail*'ENTIER' ; 'ARGUMENT' iquad*'ENTIER' ; 'ARGUMENT' pcmlag*'ENTIER' ; 'ARGUMENT' scaling*'ENTIER' ; 'ARGUMENT' typinv*'ENTIER' ; 'ARGUMENT' typrec*'ENTIER' ; tabmod = CMAIL imail iquad ; 'SI' ('EGA' iquad 2) ; disc = 'QUAF' ; 'SINON' ; disc = 'LINE' ; 'FINSI' ; * * Calcul * * Modèle * Solution exacte mt = tabmod . 'MT' . disc ; cmt = tabmod . 'CMT' . disc ; xmt = 'COORDONNEE' 1 mt ; ymt = 'COORDONNEE' 2 mt ; 'SI' dim3 ; zmt = 'COORDONNEE' 3 mt ; 'FINSI' ; solex = '+' ('*' xmt ('**' 2.D0 0.5D0)) ('*' ymt PI) ; 'SI' dim3 ; solex = '+' solex ('*' zmt ('**' PI 0.5D0)) ; 'FINSI' ; * Conditions aux limites _mt = tabmod . 'MT' . 'QUAF' ; * Laplacien mlap = GLAPN _mt 'LINE' 'T' disc 'Q' disc 1.D0 ; * Masse mmas = GMASS _mt 'LINE' 'T' disc 'Q' disc 1.D0 ; * Resolution solu nit tcpu = RESOU mlap clim pcmlag scaling typinv typrec ; 'RESPRO' nit tcpu resit ; 'FINPROC' ; * * Procédure CALCUL2 (opérateurs JPM) * 'DEBPROC' CALCUL2 ; 'ARGUMENT' imail*'ENTIER' ; 'ARGUMENT' iquad*'ENTIER' ; 'ARGUMENT' pcmlag*'ENTIER' ; 'ARGUMENT' scaling*'ENTIER' ; 'ARGUMENT' typinv*'ENTIER' ; 'ARGUMENT' typrec*'ENTIER' ; tabmod = CMAIL imail iquad ; 'SI' ('EGA' iquad 2) ; disc = 'QUAF' ; 'SINON' ; disc = 'LINE' ; 'FINSI' ; * * Calcul * * Modèle * Solution exacte mt = tabmod . 'MT' . disc ; cmt = tabmod . 'CMT' . disc ; xmt = 'COORDONNEE' 1 mt ; ymt = 'COORDONNEE' 2 mt ; 'SI' dim3 ; zmt = 'COORDONNEE' 3 mt ; 'FINSI' ; solex = '+' ('*' xmt ('**' 2.D0 0.5D0)) ('*' ymt PI) ; 'SI' dim3 ; solex = '+' solex ('*' zmt ('**' PI 0.5D0)) ; 'FINSI' ; * Conditions aux limites _mt = tabmod . 'MT' . 'QUAF' ; $mt = 'MODELISER' _mt 'NAVIER_STOKES' disc ; rv = 'EQEX' 'OPTI' 'EF' 'IMPL' 'OPTI' 'EF' 'IMPL' ; rv . 'INCO' = 'TABLE' 'INCO' ; * Resolution solu nit tcpu = RESOU mlap clim pcmlag scaling typinv typrec ; *resit = '**' (XTMX ('-' solex solu) mmas) 0.5D0 ; dif = '-' solex solu ; 'RESPRO' nit tcpu resit ; 'FINPROC' ; * * Test scaling * 'DEBPROC' TESTSCAL ; tabres = 'TABLE' ; pcmlag = 1 ; typinv = 2 ; typrec = 1 ; 'REPETER' iscal 3 ; tabres . &iscal = 'TABLE' ; scaling = &iscal ; 'REPETER' imail 3 ; tabres . &iscal . &imail = 'TABLE' ; 'REPETER' iquad 2 ; tabres . &iscal . &imail . &iquad = 'TABLE' ; nit tcpu resit = CALCUL2 &imail &iquad pcmlag scaling typinv typrec ; tabres . &iscal . &imail . &iquad . 'NIT' = nit ; tabres . &iscal . &imail . &iquad . 'TCPU' = tcpu ; tabres . &iscal . &imail . &iquad . 'RESIT' = resit ; 'FIN' iquad ; 'FIN' imail ; 'FIN' iscal ; * optecho = 'VALEUR' 'ECHO' ; 'OPTION' 'ECHO' 0 ; cc1 = 'CHAINE' 'Paramètres fixés : ' (rcvt . 'TYPINV' . typinv) ' ' (rcvt . 'TYPREC' . typrec) ; 'MESSAGE' cc1 ; 'REPETER' imail 3 ; 'REPETER' iquad 2 ; 'SAUTER' 'LIGNE' ; cc2 = 'CHAINE' (rcvt . 'MAILLAGE' . &imail) ' ' 'MESSAGE' cc2 ; cc3 = 'CHAINE' 'Nb iter.'/20 'Temps CPU'/40 'Residu reel'/60 ; 'MESSAGE' cc3 ; 'REPETER' iscal 3 ; nit = tabres . &iscal . &imail . &iquad . 'NIT' ; tcpu = tabres . &iscal . &imail . &iquad . 'TCPU' ; resit = tabres . &iscal . &imail . &iquad . 'RESIT' ; cc4 = 'CHAINE' 'SCALING : ' (rcvt . 'SCALING' . &iscal) nit*30 (formar tcpu 2)*50 (formar resit 4)*70 ; 'MESSAGE' cc4 ; 'FIN' iscal ; 'FIN' iquad ; 'FIN' imail ; 'OPTION' 'ECHO' optecho ; 'FINPROC' ; * * Test solveur * 'DEBPROC' TESTSOLV ; tabres = 'TABLE' ; pcmlag = 1 ; scaling = 1 ; typrec = 1 ; 'REPETER' itypinv 5 ; tabres . &itypinv = 'TABLE' ; typinv = &itypinv ; 'REPETER' imail 3 ; tabres . &itypinv . &imail = 'TABLE' ; 'REPETER' iquad 2 ; tabres . &itypinv . &imail . &iquad = 'TABLE' ; nit tcpu resit = CALCUL2 &imail &iquad pcmlag scaling typinv typrec ; tabres . &itypinv . &imail . &iquad . 'NIT' = nit ; tabres . &itypinv . &imail . &iquad . 'TCPU' = tcpu ; tabres . &itypinv . &imail . &iquad . 'RESIT' = resit ; 'FIN' iquad ; 'FIN' imail ; 'FIN' itypinv ; * optecho = 'VALEUR' 'ECHO' ; 'OPTION' 'ECHO' 0 ; cc1 = 'CHAINE' 'Paramètres fixés : ' (rcvt . 'SCALING' . scaling) ' ' (rcvt . 'TYPREC' . typrec) ; 'MESSAGE' cc1 ; 'REPETER' imail 3 ; 'REPETER' iquad 2 ; 'SAUTER' 'LIGNE' ; cc2 = 'CHAINE' (rcvt . 'MAILLAGE' . &imail) ' ' 'MESSAGE' cc2 ; cc3 = 'CHAINE' 'Nb iter.'/20 'Temps CPU'/40 'Residu reel'/60 ; 'MESSAGE' cc3 ; 'REPETER' itypinv 5 ; nit = tabres . &itypinv . &imail . &iquad . 'NIT' ; tcpu = tabres . &itypinv . &imail . &iquad . 'TCPU' ; resit = tabres . &itypinv . &imail . &iquad . 'RESIT' ; cc4 = 'CHAINE' 'Solveur : ' (rcvt . 'TYPINV' . &itypinv) nit*30 (formar tcpu 2)*50 (formar resit 4)*70 ; 'MESSAGE' cc4 ; 'FIN' itypinv ; 'FIN' iquad ; 'FIN' imail ; 'OPTION' 'ECHO' optecho ; 'FINPROC' ; * * Test preconditionneur * 'DEBPROC' TESTPREC ; tabres = 'TABLE' ; pcmlag = 1 ; scaling = 1 ; typinv = 2 ; 'REPETER' ityprec 15 ; tabres . &ityprec = 'TABLE' ; typrec = &ityprec ; 'REPETER' imail 3 ; tabres . &ityprec . &imail = 'TABLE' ; 'REPETER' iquad 2 ; * 'MESSAGE' ('CHAINE' 'typrec=' typrec) ; tabres . &ityprec . &imail . &iquad = 'TABLE' ; nit tcpu resit = CALCUL2 &imail &iquad pcmlag scaling typinv typrec ; * * Correction pour le gradient conjugué * tabres . &ityprec . &imail . &iquad . 'NIT' = nit ; tabres . &ityprec . &imail . &iquad . 'TCPU' = tcpu ; tabres . &ityprec . &imail . &iquad . 'RESIT' = resit ; 'FIN' iquad ; 'FIN' imail ; 'FIN' ityprec ; * optecho = 'VALEUR' 'ECHO' ; 'OPTION' 'ECHO' 0 ; cc1 = 'CHAINE' 'Paramètres fixés : ' (rcvt . 'SCALING' . scaling) ' ' (rcvt . 'TYPINV' . typinv) ; 'MESSAGE' cc1 ; 'REPETER' imail 3 ; 'REPETER' iquad 2 ; 'SAUTER' 'LIGNE' ; cc2 = 'CHAINE' (rcvt . 'MAILLAGE' . &imail) ' ' 'MESSAGE' cc2 ; cc3 = 'CHAINE' 'Nb iter.'/20 'Temps CPU'/40 'Residu reel'/60 ; 'MESSAGE' cc3 ; 'REPETER' ityprec 15 ; nit = tabres . &ityprec . &imail . &iquad . 'NIT' ; tcpu = tabres . &ityprec . &imail . &iquad . 'TCPU' ; resit = tabres . &ityprec . &imail . &iquad . 'RESIT' ; cc4 = 'CHAINE' 'Précond. : ' (rcvt . 'TYPREC' . &ityprec) nit*30 (formar tcpu 2)*50 (formar resit 4)*70 ; 'MESSAGE' cc4 ; 'FIN' ityprec ; 'FIN' iquad ; 'FIN' imail ; 'OPTION' 'ECHO' optecho ; 'FINPROC' ; * * * Fin des procédures * Début du programme * * * * Table de conversion texte * rcvt = 'TABLE' ; rcvt . 'MAILLAGE' = 'TABLE' ; rcvt . 'MAILLAGE' . 1 = 'CHAINE' 'Carre reg ' ; rcvt . 'MAILLAGE' . 2 = 'CHAINE' 'Tri. irreg' ; rcvt . 'MAILLAGE' . 3 = 'CHAINE' 'Cube reg' ; rcvt . 'SCALING' = 'TABLE' ; rcvt . 'SCALING' . 1 = 'CHAINE' 'No scal.' ; rcvt . 'SCALING' . 2 = 'CHAINE' 'L2 Scaling ' ; rcvt . 'SCALING' . 3 = 'CHAINE' 'L1 Scaling ' ; rcvt . 'TYPINV' = 'TABLE' ; rcvt . 'TYPINV' . 2 = 'CHAINE' 'BiCGStab ' ; rcvt . 'TYPINV' . 3 = 'CHAINE' 'BiCGStab(4)' ; rcvt . 'TYPINV' . 4 = 'CHAINE' 'GMRES(50) ' ; rcvt . 'TYPINV' . 5 = 'CHAINE' 'CGS ' ; rcvt . 'TYPREC' = 'TABLE' ; rcvt . 'TYPREC' . 1 = 'CHAINE' 'ILU(0) ' ; rcvt . 'TYPREC' . 2 = 'CHAINE' 'ILUT(2) ' ; rcvt . 'TYPREC' . 3 = 'CHAINE' 'ILUTP(2,0) ' ; rcvt . 'TYPREC' . 4 = 'CHAINE' 'ILUTP+0(2,0) ' ; rcvt . 'TYPREC' . 5 = 'CHAINE' 'ILUTP(2,-3) ' ; rcvt . 'TYPREC' . 6 = 'CHAINE' 'ILUTP+0(2,-3)' ; rcvt . 'TYPREC' . 7 = 'CHAINE' 'ILUTP(2,-1) ' ; rcvt . 'TYPREC' . 8 = 'CHAINE' 'ILUTP+0(2,-1)' ; rcvt . 'TYPREC' . 9 = 'CHAINE' 'ILUT(3) ' ; rcvt . 'TYPREC' . 10 = 'CHAINE' 'ILUTP(3,0) ' ; rcvt . 'TYPREC' . 11 = 'CHAINE' 'ILUTP+0(3,0) ' ; rcvt . 'TYPREC' . 12 = 'CHAINE' 'ILUTP(3,-3) ' ; rcvt . 'TYPREC' . 13 = 'CHAINE' 'ILUTP+0(3,-3)' ; rcvt . 'TYPREC' . 14 = 'CHAINE' 'ILUTP(3,-1) ' ; rcvt . 'TYPREC' . 15 = 'CHAINE' 'ILUTP+0(3,-1)' ; * TESTSCAL ; 'SAUTER' 10 'LIGNE' ; TESTSOLV ; 'SAUTER' 10 'LIGNE' ; TESTPREC ; 'SI' interact ; 'OPTION' 'DONN' 5 ; 'FINSI' ; * * End of dgibi file TEST_KRES * 'FIN' ;
© Cast3M 2003 - Tous droits réservés.
Mentions légales