Télécharger stokes_lagaug.dgibi

Retour à la liste

Numérotation des lignes :

  1. * fichier stokes_lagaug.dgibi
  2. ************************************************************************
  3. ************************************************************************
  4. *
  5. 'OPTI' 'ECHO' 0 ;
  6. *
  7. 'SAUTER' 2 'LIGNE' ;
  8. 'MESSAGE' ' Execution de stokes_lagaug.dgibi' ;
  9. 'SAUTER' 2 'LIGNE' ;
  10. *
  11. graph = faux ;
  12. interact = faux ;
  13. ************************************************************************
  14. * NOM : STOKES_LAGAUG
  15. * DESCRIPTION : Cas-test équation de Stokes incompressible
  16. * Méthode directe et de Lagrangien augmenté pour les
  17. * contraintes de divergence.
  18. * Cavité carrée entrainee
  19. *
  20. *
  21. * LANGAGE : GIBIANE-CAST3M
  22. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  23. * mél : gounand@semt2.smts.cea.fr
  24. **********************************************************************
  25. * VERSION : v1, 10/06/2011, version initiale
  26. * HISTORIQUE : v1, 10/06/2011, création
  27. * HISTORIQUE :
  28. * HISTORIQUE :
  29. ************************************************************************
  30. 'OPTION' 'DIME' 2 'MODE' 'PLAN' 'ISOV' 'SURF' ;
  31. ************************************************************************
  32. *
  33. *
  34. * PROCEDURES
  35. *
  36. *
  37. ************************************************************************
  38. *
  39. *
  40. *
  41. 'DEBPROC' ERRREL ;
  42. 'ARGUMENT' val*'FLOTTANT' ;
  43. 'ARGUMENT' valref*'FLOTTANT' ;
  44. *
  45. 'SI' ('<' ('ABS' valref) 1.D-10) ;
  46. echref = 1.D0 ;
  47. 'SINON' ;
  48. echref = valref ;
  49. 'FINSI' ;
  50. *
  51. errabs = 'ABS' ('/' ('-' val valref) echref);
  52. *
  53. 'RESPRO' errabs ;
  54. *
  55. * End of procedure file ERRREL
  56. *
  57. 'FINPROC' ;
  58. *ENDPROCEDUR errrel
  59. *BEGINPROCEDUR quafme
  60. ************************************************************************
  61. * NOM : QUAFME
  62. * DESCRIPTION :
  63. *
  64. *
  65. *
  66. * LANGAGE : GIBIANE-CAST3M
  67. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  68. * mél : gounand@semt2.smts.cea.fr
  69. **********************************************************************
  70. * VERSION : v1, 01/12/2004, version initiale
  71. * HISTORIQUE : v1, 01/12/2004, création
  72. * HISTORIQUE :
  73. * HISTORIQUE :
  74. ************************************************************************
  75. * Prière de PRENDRE LE TEMPS de compléter les commentaires
  76. * en cas de modification de ce sous-programme afin de faciliter
  77. * la maintenance !
  78. ************************************************************************
  79. *
  80. *
  81. 'DEBPROC' QUAFME ;
  82. 'ARGUMENT' tmail/'TABLE' ;
  83. 'SI' ('NON' ('EXISTE' tmail)) ;
  84. tmail = 'TABLE' ;
  85. 'FINSI' ;
  86. 'REPETER' bcl ;
  87. 'ARGUMENT' mquad/'MAILLAGE' ;
  88. 'SI' ('EXISTE' mquad) ;
  89. 'REPETER' bbcl 1 ;
  90. 'ARGUMENT' nom*'MOT' ;
  91. tmail . nom = 'TABLE' ;
  92. dimm = DEADUTIL 'DIMM' mquad ;
  93. 'SI' ('EGA' dimm 0) ;
  94. mquaf = mquad ;
  95. tmail . nom = mquad ;
  96. 'QUITTER' bbcl ;
  97. 'FINSI' ;
  98. typm = DEADUTIL 'TYPM' mquad ;
  99. 'SI' ('EGA' typm 'QUAF') ;
  100. mquaf = mquad ;
  101. tmail . nom . 'QUAF' = mquaf ;
  102. 'SI' ('EGA' dimm 2) ;
  103. tmail . nom . 'QUAI' = mquaf ;
  104. tmail . nom . 'QUAD' = mquaf ;
  105. 'FINSI' ;
  106. 'QUITTER' bbcl ;
  107. 'FINSI' ;
  108. 'SI' ('EGA' typm 'QUAI') ;
  109. mquaf = 'CHANGER' mquad 'QUAF' ;
  110. tmail . nom . 'QUAF' = mquaf ;
  111. tmail . nom . 'QUAI' = mquad ;
  112. tmail . nom . 'LINE' = 'CHANGER' 'LINEAIRE' mquad ;
  113. 'QUITTER' bbcl ;
  114. 'FINSI' ;
  115. 'SI' ('EGA' typm 'LINE') ;
  116. mquai = 'CHANGER' 'QUADRATIQUE' mquad ;
  117. mquaf = 'CHANGER' mquai 'QUAF' ;
  118. tmail . nom . 'QUAF' = mquaf ;
  119. tmail . nom . 'QUAI' = mquai ;
  120. tmail . nom . 'LINE' = mquad ;
  121. 'QUITTER' bbcl ;
  122. 'FINSI' ;
  123. 'FIN' bbcl ;
  124. * ielim = '+' ielim 1 ;
  125. * telim . ielim = mquaf ;
  126. 'SINON' ;
  127. 'QUITTER' bcl ;
  128. 'FINSI' ;
  129. 'FIN' bcl ;
  130. *
  131. *
  132. telim = 'TABLE' 'ESCLAVE' ;
  133. ielim = 0 ;
  134. tidx = 'INDEX' tmail ;
  135. 'REPETER' iidx ('DIME' tidx) ;
  136. idx = tidx . &iidx ;
  137. val = tmail . idx ;
  138. tval = 'TYPE' val ;
  139. 'SI' ('EGA' tval 'MAILLAGE') ;
  140. ielim = '+' ielim 1 ;
  141. telim . ielim = val ;
  142. 'FINSI' ;
  143. 'SI' ('EGA' tval 'TABLE') ;
  144. 'SI' ('EXISTE' val 'QUAF') ;
  145. ielim = '+' ielim 1 ;
  146. telim . ielim = val . 'QUAF' ;
  147. 'FINSI' ;
  148. 'FINSI' ;
  149. 'FIN' iidx ;
  150. *
  151. *
  152. 'ARGUMENT' tol*'FLOTTANT' ;
  153. 'ELIMINATION' ('ET' telim) tol ;
  154. *
  155. 'RESPRO' tmail ;
  156. 'FINPROC' ;
  157. *
  158. * End of procedure file QUAFME
  159. *
  160. *ENDPROCEDUR quafme
  161. *BEGINPROCEDUR modenlin
  162. ************************************************************************
  163. * NOM : MODENLIN
  164. * DESCRIPTION :
  165. *
  166. *
  167. *
  168. * LANGAGE : GIBIANE-CAST3M
  169. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  170. * mél : gounand@semt2.smts.cea.fr
  171. **********************************************************************
  172. * VERSION : v1, ??/??/2007, version initiale
  173. * HISTORIQUE : v1, ??/??/2007, création
  174. * HISTORIQUE :
  175. * HISTORIQUE :
  176. ************************************************************************
  177. * Prière de PRENDRE LE TEMPS de compléter les commentaires
  178. * en cas de modification de ce sous-programme afin de faciliter
  179. * la maintenance !
  180. ************************************************************************
  181. *
  182. *
  183. 'DEBPROC' MODENLIN ;
  184. *
  185. vdim = 'VALEUR' 'DIME' ;
  186. vmod = 'VALEUR' 'MODE' ;
  187. *
  188. 'SI' ('ET' ('NEG' vdim 2) ('NEG' vdim 3)) ;
  189. 'ERREUR' ('CHAINE' 'vdim=' vdim) ;
  190. 'FINSI' ;
  191. *
  192. * Noms par défaut
  193. *
  194. 'SI' ('EGA' vdim 2) ;
  195. 'SI' ('NEG' vmod 'AXIS') ;
  196. nomvec = 'MOTS' 'UX' 'UY' ;
  197. nomfor = 'MOTS' 'FX' 'FY' ;
  198. 'SINON' ;
  199. nomvec = 'MOTS' 'UR' 'UZ' ;
  200. nomfor = 'MOTS' 'FR' 'FZ' ;
  201. 'FINSI' ;
  202. 'FINSI' ;
  203. 'SI' ('EGA' vdim 3) ;
  204. nomvec = 'MOTS' 'UX' 'UY' 'UZ' ;
  205. nomfor = 'MOTS' 'FX' 'FY' 'FZ' ;
  206. 'FINSI' ;
  207. *
  208. nomsca = 'MOTS' 'SCAL' ;
  209. *nomflu = 'MOTS' 'SCAL' ;
  210. *
  211. 'ARGUMENT' tmode/'TABLE' ;
  212. 'SI' ('NON' ('EXISTE' tmode)) ;
  213. tmode = 'TABLE' ;
  214. *
  215. * Initialisation des inconnues par défaut
  216. *
  217. lmdisc = 'MOTS' 'CSTE' 'LINE' 'QUAF' ;
  218. 'REPETER' iidisc ('DIME' lmdisc) ;
  219. mdisc = 'EXTRAIRE' lmdisc &iidisc ;
  220. tmode . mdisc = 'TABLE' ;
  221. tmode . mdisc . 'DISC' = mdisc ;
  222. tmode . mdisc . 'NOMPRI' = 'TABLE' ;
  223. tmode . mdisc . 'NOMPRI' . 1 = nomsca ;
  224. tmode . mdisc . 'NOMDUA' = 'TABLE' ;
  225. tmode . mdisc . 'NOMDUA' . 1 = nomsca ;
  226. mdiscv = 'CHAINE' mdisc 'V' ;
  227. tmode . mdiscv = 'TABLE' ;
  228. tmode . mdiscv . 'DISC' = mdisc ;
  229. tmode . mdiscv . 'NOMPRI' = 'TABLE' ;
  230. 'REPETER' idim vdim ;
  231. TMODE . mdiscv. 'NOMPRI' . &idim =
  232. 'MOTS' ('EXTRAIRE' nomvec &idim) ;
  233. 'FIN' idim ;
  234. tmode . mdiscv . 'NOMDUA' = 'TABLE' ;
  235. 'REPETER' idim vdim ;
  236. TMODE . mdiscv. 'NOMDUA' . &idim =
  237. 'MOTS' ('EXTRAIRE' nomfor &idim) ;
  238. 'FIN' idim ;
  239. 'FIN' iidisc ;
  240. 'FINSI' ;
  241. *
  242. * Lecture des mots clés et des inconnues
  243. *
  244. lmotcle = 'MOTS' 'GEOM' 'INCO' ;
  245. ldiscdd = 'MOTS' 'LINM' 'CUBI' ;
  246. ltypinc = 'MOTS' 'SCAL' 'VECT' ;
  247. 'REPETER' imotcle ;
  248. 'ARGUMENT' motcle/'MOT' ;
  249. 'SI' ('NON' ('EXISTE' motcle)) ; 'QUITTER' imotcle ; 'FINSI' ;
  250. 'SI' ('NON' ('EXISTE' lmotcle motcle)) ;
  251. cherr = 'CHAINE' 'Keyword ' motcle ' unknown.' ; 'ERREUR' cherr ;
  252. 'FINSI' ;
  253. 'SI' ('EGA' motcle 'GEOM') ;
  254. 'ARGUMENT' discg*'MOT' ;
  255. TMODE . 'GEOM' = 'TABLE' ;
  256. TMODE . 'GEOM' . 'DISC' = discg ;
  257. 'FINSI' ;
  258. 'SI' ('EGA' motcle 'INCO') ;
  259. 'ARGUMENT' nominc*'MOT' ;
  260. 'SI' ('EGA' nominc 'GEOM') ;
  261. 'ERREUR' 'GEOM nest pas un nom dinconnu acceptable' ;
  262. 'FINSI' ;
  263. TMODE . nominc = 'TABLE' ;
  264. 'ARGUMENT' disinc*'MOT' ;
  265. TMODE . nominc . 'DISC' = disinc ;
  266. ldd = 'EXISTE' ldiscdd disinc ;
  267. 'ARGUMENT' typinc*'MOT' ;
  268. 'SI' ('NON' ('EXISTE' ltypinc typinc)) ;
  269. cherr = 'CHAINE' 'Type=' typinc ' unknown.' ;
  270. 'ERREUR' cherr ;
  271. 'FINSI' ;
  272. 'SI' ('EGA' typinc 'SCAL') ;
  273. nbinc = 1 ;
  274. 'SINON' ;
  275. nbinc = vdim ;
  276. 'FINSI' ;
  277. TMODE . nominc . 'NOMPRI' = 'TABLE' ;
  278. 'REPETER' ibinc nbinc ;
  279. 'SI' ldd ;
  280. 'ARGUMENT' nomcom*'LISTMOTS' ;
  281. 'SINON' ;
  282. 'ARGUMENT' nomcom*'MOT' ;
  283. nomcom = 'MOTS' nomcom ;
  284. 'FINSI' ;
  285. TMODE . nominc . 'NOMPRI' . &ibinc = nomcom ;
  286. 'FIN' ibinc ;
  287. TMODE . nominc . 'NOMDUA' = 'TABLE' ;
  288. 'REPETER' ibinc nbinc ;
  289. 'SI' ldd ;
  290. 'ARGUMENT' nomcom*'LISTMOTS' ;
  291. 'SINON' ;
  292. 'ARGUMENT' nomcom*'MOT' ;
  293. nomcom = 'MOTS' nomcom ;
  294. 'FINSI' ;
  295. TMODE . nominc . 'NOMDUA' . &ibinc = nomcom ;
  296. 'FIN' ibinc ;
  297. 'FINSI' ;
  298. 'FIN' imotcle ;
  299. *
  300. 'RESPRO' tmode ;
  301. *
  302. * End of procedure file MODENLIN
  303. *
  304. 'FINPROC' ;
  305. *ENDPROCEDUR modenlin
  306. *BEGINPROCEDUR grig2
  307. ************************************************************************
  308. * NOM : GRIG2
  309. * DESCRIPTION :
  310. *
  311. *
  312. *
  313. * LANGAGE : GIBIANE-CAST3M
  314. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  315. * mél : gounand@semt2.smts.cea.fr
  316. **********************************************************************
  317. * VERSION : v1, ??/??/2007, version initiale
  318. * HISTORIQUE : v1, ??/??/2007, création
  319. * HISTORIQUE :
  320. * HISTORIQUE :
  321. ************************************************************************
  322. * Prière de PRENDRE LE TEMPS de compléter les commentaires
  323. * en cas de modification de ce sous-programme afin de faciliter
  324. * la maintenance !
  325. ************************************************************************
  326. *
  327. *
  328. 'DEBPROC' GRIG2 ;
  329. 'ARGUMENT' _mt*'MAILLAGE' ;
  330. 'ARGUMENT' tdisc*'TABLE' ;
  331. *
  332. * Lectures
  333. *
  334. debug = FAUX ;
  335. lmotcle = 'MOTS' 'NPRI' 'FPRI' 'CPRI' 'NDUA' 'FDUA' 'CDUA'
  336. 'NCOF' 'FCOF' 'CCOF' 'LAPN' 'GMBT' ;
  337. * Il faut initialiser valt et valq, sinon on peut capturer ceux de
  338. * la procédure appelante
  339. valt = 'valt' ; valq = 'valq' ;
  340. llapn = 0 ;
  341. 'REPETER' imotcle ;
  342. 'ARGUMENT' motcle/'MOT' ;
  343. 'SI' ('NON' ('EXISTE' motcle)) ; 'QUITTER' imotcle ; 'FINSI' ;
  344. 'SI' ('NON' ('EXISTE' lmotcle motcle)) ;
  345. cherr = 'CHAINE' 'Keyword ' motcle ' unknown.' ; 'ERREUR' cherr ;
  346. 'FINSI' ;
  347. 'SI' ('EGA' motcle 'NPRI') ; 'ARGUMENT' nomt*'MOT' ; 'FINSI' ;
  348. 'SI' ('EGA' motcle 'NDUA') ; 'ARGUMENT' nomq*'MOT' ; 'FINSI' ;
  349. 'SI' ('EGA' motcle 'NCOF') ; 'ARGUMENT' nomo*'MOT' ; 'FINSI' ;
  350. tst1 = 'EGA' motcle 'FPRI' ; tst2 = 'EGA' motcle 'FDUA' ;
  351. tst = tst1 'OU' tst2 ;
  352. 'SI' tst ;
  353. 'SI' tst1 ; tt = TDISC . nomt ; 'FINSI' ;
  354. 'SI' tst2 ; tt = TDISC . nomq ; 'FINSI' ;
  355. isvec = ('>' ('DIME' (tt. 'NOMPRI')) 1) ;
  356. 'SI' isvec ; 'ARGUMENT' valv*'LISTREEL' ; 'SINON' ;
  357. 'ARGUMENT' valv*'FLOTTANT' ;
  358. 'FINSI' ;
  359. 'SI' tst1 ; valt = valv ; 'FINSI' ;
  360. 'SI' tst2 ; valq = valv ; 'FINSI' ;
  361. 'FINSI' ;
  362. 'SI' ('EGA' motcle 'FCOF') ; 'ARGUMENT' valo*'FLOTTANT' ; 'FINSI' ;
  363. 'SI' ('EGA' motcle 'CPRI') ; 'ARGUMENT' valt*'CHPOINT' ; 'FINSI' ;
  364. 'SI' ('EGA' motcle 'CDUA') ; 'ARGUMENT' valq*'CHPOINT' ; 'FINSI' ;
  365. 'SI' ('EGA' motcle 'CCOF') ; 'ARGUMENT' valo*'CHPOINT' ; 'FINSI' ;
  366. 'SI' ('EGA' motcle 'LAPN') ; llapn = 1 ; 'FINSI' ;
  367. 'SI' ('EGA' motcle 'GMBT') ; llapn = 2 ; 'FINSI' ;
  368. 'FIN' imotcle ;
  369. *
  370. * Tests
  371. *
  372. discg = TDISC . 'GEOM' . 'DISC' ;
  373. 'SI' ('EXISTE' tdisc 'methgau') ;
  374. methgau = tdisc . 'methgau' . 'rigi' ;
  375. 'SINON' ;
  376. methgau = 'GAU7' ;
  377. 'FINSI' ;
  378. tnomt = TDISC . nomt ;
  379. lvalt = 'NEG' ('TYPE' valt) 'MOT' ;
  380. tnomq = TDISC . nomq ;
  381. lvalq = 'NEG' ('TYPE' valq) 'MOT' ;
  382. * Scalaire ou vecteur
  383. ninct = 'DIME' (tnomt . 'NOMPRI') ;
  384. nincq = 'DIME' (tnomq . 'NOMPRI') ;
  385. 'SI' ('NEG' ninct nincq) ;
  386. cherr = 'CHAINE'
  387. 'les primales et duales nont pas le meme nombre de composantes' ;
  388. 'ERREUR' cherr ;
  389. 'FINSI' ;
  390. 'SI' ('NEG' ninct ('VALEUR' 'DIME')) ;
  391. cherr = 'CHAINE'
  392. 'les inconnues doivent etre vectorielles' ;
  393. 'ERREUR' cherr ;
  394. 'FINSI' ;
  395. *
  396. ninc = ninct ;
  397. *
  398. lcof = 'EXISTE' TDISC nomo ;
  399. 'SI' lcof ; ncof = 1 ; tcof = TDISC . nomo ;
  400. 'SINON' ; ncof = 0 ;
  401. 'FINSI' ;
  402. *
  403. 'SI' debug ;
  404. 'SI' lcof ; 'MESSAGE' 'Un coef a ete detecte' ;
  405. 'SINON' ; 'MESSAGE' 'pas de coef detecte' ;
  406. 'FINSI' ;
  407. 'FINSI' ;
  408. *
  409. vdim = 'VALEUR' 'DIME' ;
  410. vmod = 'VALEUR' 'MODE' ;
  411. idim = 0 ;
  412. 'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'PLANDEFO')) ;
  413. idim = 2 ;
  414. iaxi = FAUX ;
  415. 'FINSI' ;
  416. 'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'AXIS')) ;
  417. idim = 2 ;
  418. iaxi = VRAI ;
  419. 'FINSI' ;
  420. 'SI' ('ET' ('EGA' vdim 3) ('EGA' vmod 'TRID')) ;
  421. idim = 3 ;
  422. iaxi = FAUX ;
  423. 'FINSI' ;
  424. 'SI' ('EGA' vdim 1) ;
  425. idim = 1 ;
  426. iaxi = FAUX ;
  427. 'FINSI' ;
  428. * 'MESSAGE' ('CHAINE' 'iaxi=' iaxi );
  429. 'SI' ('EGA' idim 0) ;
  430. 'ERREUR' ('CHAINE' 'vmod=' vmod ' et vdim=' vdim ' non prevu') ;
  431. 'FINSI' ;
  432. 'SI' iaxi ;
  433. rmt = 'COORDONNEE' 1 _mt ;
  434. deupi = '*' PI 2.D0 ;
  435. 'FINSI' ;
  436. *
  437. * Optimisation possible : construire la matrice par blocs
  438. * qd valt et valq ne sont pas donnés
  439. *
  440. *
  441. *Bug ? numop = ('**' ninc 2) '+' 1 ;
  442. numop = '**' ninc 2 ;
  443. 'SI' iaxi ; numop = '+' numop 1 ; 'FINSI' ;
  444. numder = idim ;
  445. numvar = ninc ;
  446. ncof = '+' ncof 1 ;
  447. *delete 'SI' iaxi ; ncof = '+' ncof 1 ; 'FINSI' ;
  448. numdat = ncof ;
  449. numcof = ncof ;
  450. *
  451. A = ININLIN numop numvar numdat numcof numder ;
  452. *
  453. lvt = 'EGA' ('TYPE' valt) 'LISTREEL' ;
  454. 'REPETER' iiinc ninc ;
  455. iinc = &iiinc ;
  456. A . 'VAR' . iinc . 'NOMDDL' = tnomt . 'NOMPRI' . iinc ;
  457. A . 'VAR' . iinc . 'DISC' = tnomt . 'DISC' ;
  458. 'SI' lvalt ;
  459. 'SI' lvt ;
  460. A . 'VAR' . iinc . 'VALEUR' = 'EXTRAIRE' valt iinc ;
  461. 'SINON' ;
  462. A . 'VAR' . iinc . 'VALEUR' = valt ;
  463. 'FINSI' ;
  464. 'FINSI' ;
  465. 'FIN' iiinc ;
  466. *
  467. icof = 0 ;
  468. ll = 'LECT' ;
  469. icof = '+' icof 1 ;
  470. A . 'DAT' . icof . 'NOMDDL' = 'MOTS' 'SCAL' ;
  471. A . 'DAT' . icof . 'DISC' = 'CSTE' ;
  472. A . 'DAT' . icof . 'VALEUR' = 2.D0 ;
  473. A . 'COF' . icof . 'COMPOR' = 'IDEN' ;
  474. A . 'COF' . icof . 'LDAT' = 'LECT' icof ;
  475. ll = 'LECT' ;
  476. ll2 = 'LECT' icof ;
  477. 'SI' lcof ;
  478. icof = '+' icof 1 ;
  479. A . 'DAT' . icof . 'NOMDDL' = tcof . 'NOMPRI' . 1 ;
  480. A . 'DAT' . icof . 'DISC' = tcof . 'DISC' ;
  481. A . 'DAT' . icof . 'VALEUR' = valo ;
  482. A . 'COF' . icof . 'COMPOR' = 'IDEN' ;
  483. A . 'COF' . icof . 'LDAT' = 'LECT' icof ;
  484. ll = 'ET' ll ('LECT' icof) ;
  485. ll2 = 'ET' ll2 ('LECT' icof) ;
  486. 'FINSI' ;
  487. *
  488. iop = 0 ;
  489. 'REPETER' iidim idim ;
  490. 'REPETER' jidim idim ;
  491. iop = '+' iop 1 ;
  492. 'SI' ('EGA' &iidim &jidim) ;
  493. 'SI' ('EGA' llapn 0) ;
  494. A . iop . &iidim . &jidim = ll2 ;
  495. 'SINON' ;
  496. A . iop . &iidim . &jidim = ll ;
  497. 'FINSI' ;
  498. 'SINON' ;
  499. 'SI' ('NEG' llapn 2) ;
  500. A . iop . &iidim . &jidim = ll ;
  501. 'FINSI' ;
  502. 'SI' ('NEG' llapn 1) ;
  503. A . iop . &jidim . &iidim = ll ;
  504. 'FINSI' ;
  505. 'FINSI' ;
  506. 'FIN' jidim ;
  507. 'FIN' iidim ;
  508. 'SI' iaxi ;
  509. iop = '+' iop 1 ;
  510. 'SI' ('EGA' llapn 0) ;
  511. A . iop . 1 . 0 = ll2 ;
  512. 'SINON' ;
  513. A . iop . 1 . 0 = ll ;
  514. 'FINSI' ;
  515. 'FINSI' ;
  516. *
  517. 'SI' iaxi ;
  518. numdat = 2 ;
  519. numcof = 2 ;
  520. 'SINON' ;
  521. numdat = 0 ;
  522. numcof = 0 ;
  523. 'FINSI' ;
  524. *
  525. B = ININLIN numop numvar numdat numcof numder ;
  526. *
  527. icof = 0 ;
  528. *
  529. 'SI' iaxi ;
  530. icof = '+' icof 1 ;
  531. B . 'DAT' . icof . 'NOMDDL' = 'MOTS' 'SCAL' ;
  532. B . 'DAT' . icof . 'DISC' = discg ;
  533. B . 'DAT' . icof . 'VALEUR' = rmt ;
  534. B . 'COF' . icof . 'COMPOR' = 'IDEN' ;
  535. B . 'COF' . icof . 'LDAT' = 'LECT' icof ;
  536. icof = '+' icof 1 ;
  537. B . 'DAT' . icof . 'NOMDDL' = 'MOTS' 'SCAL' ;
  538. B . 'DAT' . icof . 'DISC' = 'CSTE' ;
  539. B . 'DAT' . icof . 'VALEUR' = deupi ;
  540. B . 'COF' . icof . 'COMPOR' = 'IDEN' ;
  541. B . 'COF' . icof . 'LDAT' = 'LECT' icof ;
  542. ll = 'LECT' ('-' icof 1) icof ;
  543. mic = '*' ('-' icof 1) -1 ;
  544. llb = 'LECT' mic icof ;
  545. 'FINSI' ;
  546. *
  547. lvq = 'EGA' ('TYPE' valq) 'LISTREEL' ;
  548. 'REPETER' iiinc ninc ;
  549. iinc = &iiinc ;
  550. B . 'VAR' . iinc . 'NOMDDL' = tnomq . 'NOMDUA' . iinc ;
  551. B . 'VAR' . iinc . 'DISC' = tnomq . 'DISC' ;
  552. 'SI' lvalq ;
  553. 'SI' lvq ;
  554. B . 'VAR' . iinc . 'VALEUR' = 'EXTRAIRE' valq iinc ;
  555. 'SINON' ;
  556. B . 'VAR' . iinc . 'VALEUR' = valq ;
  557. 'FINSI' ;
  558. 'FINSI' ;
  559. 'FIN' iiinc ;
  560. *
  561. iop = 0 ;
  562. 'REPETER' iidim idim ;
  563. 'REPETER' jidim idim ;
  564. iop = '+' iop 1 ;
  565. B . iop . &iidim . &jidim = ll ;
  566. 'FIN' jidim ;
  567. 'FIN' iidim ;
  568. 'SI' iaxi ;
  569. iop = '+' iop 1 ;
  570. B . iop . 1 . 0 = llb ;
  571. * B . iop . 1 . 0 = ll ;
  572. 'FINSI' ;
  573. *
  574. mgrig = NLINP discg _mt A B methgau ;
  575. * Integration par parties
  576. * mgrig = '*' mgrig -1.D0 ;
  577. *
  578. 'RESPRO' mgrig ;
  579. 'FINPROC' ;
  580. *
  581. * End of procedure file GRIG2
  582. *
  583. *ENDPROCEDUR grig2
  584. *BEGINPROCEDUR gdiv2
  585. ************************************************************************
  586. * NOM : GDIV2
  587. * DESCRIPTION : Une matrice de masse
  588. *
  589. *
  590. *
  591. * LANGAGE : GIBIANE-CAST3M
  592. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  593. * mél : gounand@semt2.smts.cea.fr
  594. **********************************************************************
  595. * VERSION : v2, 14/03/2006, mise à jour NLIN évolué
  596. * VERSION : v1, 13/05/2004, version initiale
  597. * HISTORIQUE : v1, 13/05/2004, création
  598. * HISTORIQUE :
  599. * HISTORIQUE :
  600. ************************************************************************
  601. * Prière de PRENDRE LE TEMPS de compléter les commentaires
  602. * en cas de modification de ce sous-programme afin de faciliter
  603. * la maintenance !
  604. ************************************************************************
  605. *
  606. *
  607. 'DEBPROC' GDIV2 ;
  608. 'ARGUMENT' _mt*'MAILLAGE' ;
  609. 'ARGUMENT' _smt/'MAILLAGE' ;
  610. 'ARGUMENT' tdisc*'TABLE' ;
  611. *
  612. * Lectures
  613. *
  614. debug = FAUX ;
  615. lmotcle = 'MOTS' 'NPRI' 'FPRI' 'CPRI' 'NDUA' 'FDUA' 'CDUA'
  616. 'NCOF' 'FCOF' 'CCOF' 'GBBT' 'GMBT' ;
  617. * Il faut initialiser valt et valq, sinon on peut capturer ceux de
  618. * la procédure appelante
  619. valt = 'valt' ; valq = 'valq' ;
  620. lbbt = 0 ;
  621. *
  622. 'REPETER' imotcle ;
  623. 'ARGUMENT' motcle/'MOT' ;
  624. 'SI' ('NON' ('EXISTE' motcle)) ; 'QUITTER' imotcle ; 'FINSI' ;
  625. 'SI' ('NON' ('EXISTE' lmotcle motcle)) ;
  626. cherr = 'CHAINE' 'Keyword ' motcle ' unknown.' ; 'ERREUR' cherr ;
  627. 'FINSI' ;
  628. 'SI' ('EGA' motcle 'NPRI') ; 'ARGUMENT' nomt*'MOT' ; 'FINSI' ;
  629. 'SI' ('EGA' motcle 'NDUA') ; 'ARGUMENT' nomq*'MOT' ; 'FINSI' ;
  630. 'SI' ('EGA' motcle 'NCOF') ; 'ARGUMENT' nomo*'MOT' ; 'FINSI' ;
  631. 'SI' ('EGA' motcle 'FPRI') ; 'ARGUMENT' valt*'LISTREEL' ; 'FINSI' ;
  632. 'SI' ('EGA' motcle 'FDUA') ; 'ARGUMENT' valq*'FLOTTANT' ; 'FINSI' ;
  633. 'SI' ('EGA' motcle 'FCOF') ; 'ARGUMENT' valo*'FLOTTANT' ; 'FINSI' ;
  634. 'SI' ('EGA' motcle 'CPRI') ; 'ARGUMENT' valt*'CHPOINT' ; 'FINSI' ;
  635. 'SI' ('EGA' motcle 'CDUA') ; 'ARGUMENT' valq*'CHPOINT' ; 'FINSI' ;
  636. 'SI' ('EGA' motcle 'CCOF') ; 'ARGUMENT' valo*'CHPOINT' ; 'FINSI' ;
  637. 'SI' ('EGA' motcle 'GBBT') ; lbbt = 1 ; 'FINSI' ;
  638. 'SI' ('EGA' motcle 'GMBT') ; lbbt = 2 ; 'FINSI' ;
  639. 'FIN' imotcle ;
  640. *
  641. * Tests
  642. *
  643. discg = TDISC . 'GEOM' . 'DISC' ;
  644. 'SI' ('EXISTE' tdisc 'methgau') ;
  645. methgau = tdisc . 'methgau' . 'amor' ;
  646. 'SINON' ;
  647. methgau = 'GAU7' ;
  648. 'FINSI' ;
  649. tnomt = TDISC . nomt ;
  650. lvalt = 'NEG' ('TYPE' valt) 'MOT' ;
  651. tnomq = TDISC . nomq ;
  652. lvalq = 'NEG' ('TYPE' valq) 'MOT' ;
  653. *
  654. lcof = 'EXISTE' TDISC nomo ;
  655. 'SI' lcof ; ncof = 1 ; tcof = TDISC . nomo ;
  656. 'SINON' ; ncof = 0 ;
  657. 'FINSI' ;
  658. *
  659. 'SI' debug ;
  660. 'SI' lcof ; 'MESSAGE' 'Un coef a ete detecte' ;
  661. 'SINON' ; 'MESSAGE' 'pas de coef detecte' ;
  662. 'FINSI' ;
  663. 'FINSI' ;
  664. *
  665. vdim = 'VALEUR' 'DIME' ;
  666. vmod = 'VALEUR' 'MODE' ;
  667. idim = 0 ;
  668. 'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'PLANDEFO')) ;
  669. idim = 2 ;
  670. iaxi = FAUX ;
  671. 'FINSI' ;
  672. 'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'AXIS')) ;
  673. idim = 2 ;
  674. iaxi = VRAI ;
  675. 'FINSI' ;
  676. 'SI' ('ET' ('EGA' vdim 3) ('EGA' vmod 'TRID')) ;
  677. idim = 3 ;
  678. iaxi = FAUX ;
  679. 'FINSI' ;
  680. 'SI' ('EGA' vdim 1) ;
  681. idim = 1 ;
  682. iaxi = FAUX ;
  683. 'FINSI' ;
  684. 'SI' ('EGA' idim 0) ;
  685. 'ERREUR' ('CHAINE' 'vmod=' vmod ' et vdim=' vdim ' non prevu') ;
  686. 'FINSI' ;
  687. 'SI' iaxi ;
  688. dp = ('*' PI 2.D0) ;
  689. rmt = 'COORDONNEE' 1 _mt ;
  690. ncof = ncof '+' 2 ;
  691. 'FINSI' ;
  692. * Scalaire ou vecteur
  693. ninct = 'DIME' (tnomt . 'NOMPRI') ;
  694. nincq = 'DIME' (tnomq . 'NOMPRI') ;
  695. 'SI' ('NEG' ninct idim) ;
  696. cherr = 'CHAINE'
  697. 'la primale doit etre un vecteur' ;
  698. 'ERREUR' cherr ;
  699. 'FINSI' ;
  700. 'SI' ('NEG' nincq 1) ;
  701. cherr = 'CHAINE'
  702. 'la duale doit etre un scalaire' ;
  703. 'ERREUR' cherr ;
  704. 'FINSI' ;
  705. *
  706. numop = 1 ; numder = idim ; numvar = ninct ;
  707. numdat = ncof ; numcof = ncof ;
  708. A = ININLIN numop numvar numdat numcof numder ;
  709. *
  710. lvt = 'EGA' ('TYPE' valt) 'LISTREEL' ;
  711. 'REPETER' iiinct ninct ;
  712. iinct = &iiinct ;
  713. A . 'VAR' . iinct . 'NOMDDL' = tnomt . 'NOMPRI' . iinct ;
  714. A . 'VAR' . iinct . 'DISC' = tnomt . 'DISC' ;
  715. 'SI' lvalt ;
  716. 'SI' lvt ;
  717. A . 'VAR' . iinct . 'VALEUR' = 'EXTRAIRE' valt iinct ;
  718. 'SINON' ;
  719. A . 'VAR' . iinct . 'VALEUR' = valt ;
  720. 'FINSI' ;
  721. 'FINSI' ;
  722. 'FIN' iiinct ;
  723. *
  724. icof = 0 ;
  725. 'SI' lcof ;
  726. icof = '+' icof 1 ;
  727. A . 'DAT' . icof . 'NOMDDL' = tcof . 'NOMPRI' . 1 ;
  728. A . 'DAT' . icof . 'DISC' = tcof . 'DISC' ;
  729. A . 'DAT' . icof . 'VALEUR' = valo ;
  730. A . 'COF' . icof . 'COMPOR' = 'IDEN' ;
  731. A . 'COF' . icof . 'LDAT' = 'LECT' icof ;
  732. ll = 'LECT' 1 ;
  733. 'SINON' ;
  734. ll = 'LECT' ;
  735. 'FINSI' ;
  736. *
  737. 'SI' iaxi ;
  738. icof = '+' icof 1 ;
  739. A . 'DAT' . icof . 'NOMDDL' = 'MOTS' 'SCAL' ;
  740. A . 'DAT' . icof . 'DISC' = 'CSTE' ;
  741. A . 'DAT' . icof . 'VALEUR' = dp ;
  742. A . 'COF' . icof . 'COMPOR' = 'IDEN' ;
  743. A . 'COF' . icof . 'LDAT' = 'LECT' icof ;
  744. icof = '+' icof 1 ;
  745. A . 'DAT' . icof . 'NOMDDL' = 'MOTS' 'SCAL' ;
  746. A . 'DAT' . icof . 'DISC' = discg ;
  747. A . 'DAT' . icof . 'VALEUR' = rmt ;
  748. A . 'COF' . icof . 'COMPOR' = 'IDEN' ;
  749. A . 'COF' . icof . 'LDAT' = 'LECT' icof ;
  750. lldpr = ll 'ET' ('LECT' ('-' icof 1) icof) ;
  751. lldp = ll 'ET' ('LECT' ('-' icof 1)) ;
  752. 'FINSI' ;
  753. *
  754. 'SI' iaxi ;
  755. 'REPETER' iidim idim ;
  756. A . 1 . &iidim . &iidim = lldpr ;
  757. 'FIN' iidim ;
  758. A . 1 . 1 . 0 = lldp ;
  759. 'SINON' ;
  760. 'REPETER' iidim idim ;
  761. A . 1 . &iidim . &iidim = ll ;
  762. 'FIN' iidim ;
  763. 'FINSI' ;
  764. *
  765. numvar = 1 ;
  766. numdat = 0 ;
  767. numcof = 0 ;
  768. *
  769. B = ININLIN numop numvar numdat numcof numder ;
  770. B . 'VAR' . 1 . 'NOMDDL' = tnomq . 'NOMDUA' . 1 ;
  771. B . 'VAR' . 1 . 'DISC' = tnomq . 'DISC' ;
  772. 'SI' lvalq ;
  773. B . 'VAR' . 1 . 'VALEUR' = valq ;
  774. 'FINSI' ;
  775. B . 1 . 1 . 0 = 'LECT' ;
  776. *
  777. 'SI' ('OU' ('EGA' lbbt 0) ('EGA' lbbt 1)) ;
  778. 'SI' ('EXISTE' _smt) ;
  779. mgdiv2 = 'NLIN' discg _mt _smt A B methgau ;
  780. 'SINON' ;
  781. mgdiv2 = NLINP discg _mt A B methgau ;
  782. 'FINSI' ;
  783. 'FINSI' ;
  784. 'SI' ('OU' ('EGA' lbbt 1) ('EGA' lbbt 2)) ;
  785. B . 'VAR' . 1 . 'NOMDDL' = tnomq . 'NOMPRI' . 1 ;
  786. 'REPETER' iiinct ninct ;
  787. iinct = &iiinct ;
  788. A . 'VAR' . iinct . 'NOMDDL' = tnomt . 'NOMDUA' . iinct ;
  789. 'FIN' iiinct ;
  790. 'SI' ('EXISTE' _smt) ;
  791. mgdiv3 = 'NLIN' discg _mt _smt B A methgau ;
  792. 'SINON' ;
  793. mgdiv3 = NLINP discg _mt B A methgau ;
  794. 'FINSI' ;
  795. 'FINSI' ;
  796. 'SI' ('EGA' lbbt 0) ;
  797. mgdiv = mgdiv2 ;
  798. 'FINSI' ;
  799. 'SI' ('EGA' lbbt 1) ;
  800. mgdiv = mgdiv2 'ET' mgdiv3 ;
  801. 'FINSI' ;
  802. 'SI' ('EGA' lbbt 2) ;
  803. mgdiv = mgdiv3 ;
  804. 'FINSI' ;
  805. 'RESPRO' mgdiv ;
  806. 'FINPROC' ;
  807. *
  808. * End of procedure file GDIV2
  809. *
  810. *ENDPROCEDUR gdiv2
  811. *BEGINPROCEDUR gnor2
  812. ************************************************************************
  813. * NOM : GNOR2
  814. * DESCRIPTION : Calcule le champ de normales à une surface.
  815. * Peut servir à calculer une pression, un potentiel
  816. * lié à la gravité, un volume contenu dans une surface.
  817. * Attention à l'orientation de la surface !
  818. *
  819. * Computes a field of normal to a surface.
  820. * Also useful to compute a pressure field,
  821. * a gravity potential field, a volume enclosed
  822. * by a surface.
  823. * WARNING : The orientation of the surface matters !
  824. *
  825. *
  826. *
  827. *
  828. * LANGAGE : GIBIANE-CAST3M
  829. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  830. * mél : gounand@semt2.smts.cea.fr
  831. **********************************************************************
  832. * VERSION : v1, 22/04/2011
  833. * HISTORIQUE : v1, 22/04/2011, création
  834. * HISTORIQUE :
  835. * HISTORIQUE :
  836. ************************************************************************
  837. *
  838. *
  839. 'DEBPROC' GNOR2 ;
  840. 'ARGUMENT' _mt*'MAILLAGE' ;
  841. 'ARGUMENT' tdisc*'TABLE' ;
  842. *
  843. * Lectures
  844. *
  845. dim = 'VALEUR' 'DIME' ;
  846. debug = FAUX ;
  847. lmotcle = 'MOTS' 'NPRI' 'FPRI' 'CPRI' 'NDUA' 'FDUA' 'CDUA'
  848. 'NCOF' 'FCOF' 'CCOF' ;
  849. * Il faut initialiser valt et valq, sinon on peut capturer ceux de
  850. * la procédure appelante
  851. valt = 'valt' ; valq = 'valq' ;
  852. 'REPETER' imotcle ;
  853. 'ARGUMENT' motcle/'MOT' ;
  854. 'SI' ('NON' ('EXISTE' motcle)) ; 'QUITTER' imotcle ; 'FINSI' ;
  855. 'SI' ('NON' ('EXISTE' lmotcle motcle)) ;
  856. cherr = 'CHAINE' 'Keyword ' motcle ' unknown.' ; 'ERREUR' cherr ;
  857. 'FINSI' ;
  858. 'SI' ('EGA' motcle 'NPRI') ; 'ARGUMENT' nomt*'MOT' ; 'FINSI' ;
  859. 'SI' ('EGA' motcle 'NDUA') ; 'ARGUMENT' nomq*'MOT' ; 'FINSI' ;
  860. 'SI' ('EGA' motcle 'NCOF') ; 'ARGUMENT' nomo*'MOT' ; 'FINSI' ;
  861. tst1 = 'EGA' motcle 'FPRI' ; tst2 = 'EGA' motcle 'FDUA' ;
  862. tst3 = 'EGA' motcle 'FCOF' ;
  863. tst = tst1 'OU' tst2 'OU' tst3 ;
  864. 'SI' tst ;
  865. 'SI' tst1 ; tt = TDISC . nomt ; 'FINSI' ;
  866. 'SI' tst2 ; tt = TDISC . nomq ; 'FINSI' ;
  867. 'SI' tst3 ; tt = TDISC . nomo ; 'FINSI' ;
  868. isvec = ('>' ('DIME' (tt. 'NOMPRI')) 1) ;
  869. 'SI' isvec ; 'ARGUMENT' valv*'LISTREEL' ; 'SINON' ;
  870. 'ARGUMENT' valv*'FLOTTANT' ;
  871. 'FINSI' ;
  872. 'SI' tst1 ; valt = valv ; 'FINSI' ;
  873. 'SI' tst2 ; valq = valv ; 'FINSI' ;
  874. 'SI' tst3 ; valo = valv ; 'FINSI' ;
  875. 'FINSI' ;
  876. 'SI' ('EGA' motcle 'CPRI') ; 'ARGUMENT' valt*'CHPOINT' ; 'FINSI' ;
  877. 'SI' ('EGA' motcle 'CDUA') ; 'ARGUMENT' valq*'CHPOINT' ; 'FINSI' ;
  878. 'SI' ('EGA' motcle 'CCOF') ; 'ARGUMENT' valo*'CHPOINT' ; 'FINSI' ;
  879. 'FIN' imotcle ;
  880. *
  881. * Tests
  882. *
  883. discg = TDISC . 'GEOM' . 'DISC' ;
  884. 'SI' ('EXISTE' tdisc 'methgau') ;
  885. methgau = tdisc . 'methgau' . 'mass' ;
  886. 'SINON' ;
  887. methgau = 'GAU7' ;
  888. 'FINSI' ;
  889. tnomt = TDISC . nomt ;
  890. lvalt = 'NEG' ('TYPE' valt) 'MOT' ;
  891. tnomq = TDISC . nomq ;
  892. lvalq = 'NEG' ('TYPE' valq) 'MOT' ;
  893. * Scalaire ou vecteur
  894. ninct = 'DIME' (tnomt . 'NOMPRI') ;
  895. nincq = 'DIME' (tnomq . 'NOMDUA') ;
  896. 'SI' ('ET' ('NEG' ninct 1) ('NEG' ninct dim)) ;
  897. cherr = 'CHAINE'
  898. 'la primale doit etre un scalaire ou un vecteur' ;
  899. 'ERREUR' cherr ;
  900. 'FINSI' ;
  901. 'SI' ('NEG' nincq dim) ;
  902. cherr = 'CHAINE'
  903. 'la duale doit etre un vecteur' ;
  904. 'ERREUR' cherr ;
  905. 'FINSI' ;
  906. *ninc = ninct ;
  907. *
  908. lcof = 'EXISTE' TDISC nomo ;
  909. 'SI' lcof ; tcof = TDISC . nomo ;
  910. ncof = 'DIME' (tcof . 'NOMPRI') ;
  911. 'SINON' ; ncof = 0 ;
  912. 'FINSI' ;
  913. *
  914. 'SI' debug ;
  915. 'SI' lcof ; 'MESSAGE' 'Un coef a ete detecte' ;
  916. 'SINON' ; 'MESSAGE' 'pas de coef detecte' ;
  917. 'FINSI' ;
  918. 'FINSI' ;
  919. *
  920. vdim = 'VALEUR' 'DIME' ;
  921. vmod = 'VALEUR' 'MODE' ;
  922. idim = 0 ;
  923. 'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'PLANDEFO')) ;
  924. idim = 2 ;
  925. iaxi = FAUX ;
  926. 'FINSI' ;
  927. 'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'AXIS')) ;
  928. idim = 2 ;
  929. iaxi = VRAI ;
  930. 'FINSI' ;
  931. 'SI' ('ET' ('EGA' vdim 3) ('EGA' vmod 'TRID')) ;
  932. idim = 3 ;
  933. iaxi = FAUX ;
  934. 'FINSI' ;
  935. 'SI' ('EGA' vdim 1) ;
  936. idim = 1 ;
  937. iaxi = FAUX ;
  938. 'FINSI' ;
  939. * 'MESSAGE' ('CHAINE' 'iaxi=' iaxi );
  940. 'SI' ('EGA' idim 0) ;
  941. 'ERREUR' ('CHAINE' 'vmod=' vmod ' et vdim=' vdim ' non prevu') ;
  942. 'FINSI' ;
  943. 'SI' iaxi ;
  944. dprmt = '*' ('COORDONNEE' 1 _mt) ('*' PI 2.D0) ;
  945. 'FINSI' ;
  946. *
  947. * Optimisation possible : construire la matrice par blocs
  948. * qd valt et valq ne sont pas donnés
  949. *
  950. numop = nincq ; numder = idim ; numvar = ninct ;
  951. numdat = ncof ; numcof = ncof ;
  952. A = ININLIN numop numvar numdat numcof numder ;
  953. 'SI' lcof ;
  954. lvo = 'EGA' ('TYPE' valo) 'LISTREEL' ;
  955. 'REPETER' iicof ncof ;
  956. icof = &iicof ;
  957. A . 'DAT' . icof . 'NOMDDL' = tcof . 'NOMPRI' . icof ;
  958. A . 'DAT' . icof . 'DISC' = tcof . 'DISC' ;
  959. 'SI' lvo ;
  960. A . 'DAT' . icof . 'VALEUR' = 'EXTRAIRE' valo icof ;
  961. 'SINON' ;
  962. A . 'DAT' . icof . 'VALEUR' = valo ;
  963. 'FINSI' ;
  964. A . 'COF' . icof . 'COMPOR' = 'IDEN' ;
  965. A . 'COF' . icof . 'LDAT' = 'LECT' icof ;
  966. 'FIN' iicof ;
  967. 'SINON' ;
  968. ll = 'LECT' ;
  969. 'FINSI' ;
  970. lvt = 'EGA' ('TYPE' valt) 'LISTREEL' ;
  971. 'REPETER' iiincq nincq ;
  972. iincq = &iiincq ;
  973. iinct = 'MINIMUM' ('LECT' iincq ninct) ;
  974. A . 'VAR' . iinct . 'NOMDDL' = tnomt . 'NOMPRI' . iinct ;
  975. A . 'VAR' . iinct . 'DISC' = tnomt . 'DISC' ;
  976. 'SI' lvalt ;
  977. 'SI' lvt ;
  978. A . 'VAR' . iinct . 'VALEUR' = 'EXTRAIRE' valt iinct ;
  979. 'SINON' ;
  980. A . 'VAR' . iinct . 'VALEUR' = valt ;
  981. 'FINSI' ;
  982. 'FINSI' ;
  983. 'SI' lcof ;
  984. icof = 'MINIMUM' ('LECT' iincq ncof) ;
  985. A . iincq . iinct . 0 = 'LECT' icof ;
  986. 'SINON' ;
  987. A . iincq . iinct . 0 = ll ;
  988. 'FINSI' ;
  989. 'FIN' iiincq ;
  990. *
  991. 'SI' iaxi ;
  992. numdat = 1 ;
  993. numcof = dim '+' 1 ;
  994. 'SINON' ;
  995. numdat = 0 ;
  996. numcof = dim ;
  997. 'FINSI' ;
  998. numvar = nincq ;
  999. B = ININLIN numop numvar numdat numcof numder ;
  1000. icof = 0 ;
  1001. 'REPETER' iiidim idim ;
  1002. icof = '+' icof 1 ;
  1003. B . 'COF' . icof . 'COMPOR' = 'CHAINE' 'VNOR' &iiidim ;
  1004. B . 'COF' . icof . 'LDAT' = 'LECT' ;
  1005. 'FIN' iiidim ;
  1006. *
  1007. 'SI' iaxi ;
  1008. icof = '+' icof 1 ;
  1009. B . 'DAT' . 1 . 'NOMDDL' = 'MOTS' 'SCAL' ;
  1010. B . 'DAT' . 1 . 'DISC' = discg ;
  1011. B . 'DAT' . 1 . 'VALEUR' = dprmt ;
  1012. B . 'COF' . icof . 'COMPOR' = 'IDEN' ;
  1013. B . 'COF' . icof . 'LDAT' = 'LECT' 1 ;
  1014. ll = 'LECT' icof ;
  1015. 'SINON' ;
  1016. ll = 'LECT' ;
  1017. 'FINSI' ;
  1018. lvq = 'EGA' ('TYPE' valq) 'LISTREEL' ;
  1019. 'REPETER' iiincq nincq ;
  1020. iincq = &iiincq ;
  1021. B . 'VAR' . iincq . 'NOMDDL' = tnomq . 'NOMDUA' . iincq ;
  1022. B . 'VAR' . iincq . 'DISC' = tnomq . 'DISC' ;
  1023. 'SI' lvalq ;
  1024. 'SI' lvq ;
  1025. B . 'VAR' . iincq . 'VALEUR' = 'EXTRAIRE' valq iincq ;
  1026. 'SINON' ;
  1027. B . 'VAR' . iincq . 'VALEUR' = valq ;
  1028. 'FINSI' ;
  1029. 'FINSI' ;
  1030. B . iincq . iincq . 0 = ('LECT' iincq) 'ET' ll ;
  1031. 'FIN' iiincq ;
  1032. *
  1033. mgnor2 = NLIN discg _mt A B methgau ;
  1034. *
  1035. 'RESPRO' mgnor2 ;
  1036. 'FINPROC' ;
  1037. *
  1038. * End of procedure file GNOR2
  1039. *
  1040. *ENDPROCEDUR gnor2
  1041. *BEGINPROCEDUR gmass2
  1042. ************************************************************************
  1043. * NOM : GMASS2
  1044. * DESCRIPTION : Une matrice de masse
  1045. *
  1046. *
  1047. *
  1048. * LANGAGE : GIBIANE-CAST3M
  1049. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  1050. * mél : gounand@semt2.smts.cea.fr
  1051. **********************************************************************
  1052. * VERSION : v2, 14/03/2006, mise à jour NLIN évolué
  1053. * VERSION : v1, 13/05/2004, version initiale
  1054. * HISTORIQUE : v1, 13/05/2004, création
  1055. * HISTORIQUE :
  1056. * HISTORIQUE :
  1057. ************************************************************************
  1058. * Prière de PRENDRE LE TEMPS de compléter les commentaires
  1059. * en cas de modification de ce sous-programme afin de faciliter
  1060. * la maintenance !
  1061. ************************************************************************
  1062. *
  1063. *
  1064. 'DEBPROC' GMASS2 ;
  1065. 'ARGUMENT' _mt*'MAILLAGE' ;
  1066. 'ARGUMENT' _smt/'MAILLAGE' ;
  1067. 'ARGUMENT' tdisc*'TABLE' ;
  1068. *
  1069. * Lectures
  1070. *
  1071. debug = FAUX ;
  1072. lmotcle = 'MOTS' 'NPRI' 'FPRI' 'CPRI' 'NDUA' 'FDUA' 'CDUA'
  1073. 'NCOF' 'FCOF' 'CCOF' ;
  1074. * Il faut initialiser valt et valq, sinon on peut capturer ceux de
  1075. * la procédure appelante
  1076. valt = 'valt' ; valq = 'valq' ;
  1077. 'REPETER' imotcle ;
  1078. 'ARGUMENT' motcle/'MOT' ;
  1079. 'SI' ('NON' ('EXISTE' motcle)) ; 'QUITTER' imotcle ; 'FINSI' ;
  1080. 'SI' ('NON' ('EXISTE' lmotcle motcle)) ;
  1081. cherr = 'CHAINE' 'Keyword ' motcle ' unknown.' ; 'ERREUR' cherr ;
  1082. 'FINSI' ;
  1083. 'SI' ('EGA' motcle 'NPRI') ; 'ARGUMENT' nomt*'MOT' ; 'FINSI' ;
  1084. 'SI' ('EGA' motcle 'NDUA') ; 'ARGUMENT' nomq*'MOT' ; 'FINSI' ;
  1085. 'SI' ('EGA' motcle 'NCOF') ; 'ARGUMENT' nomo*'MOT' ; 'FINSI' ;
  1086. tst1 = 'EGA' motcle 'FPRI' ; tst2 = 'EGA' motcle 'FDUA' ;
  1087. tst = tst1 'OU' tst2 ;
  1088. 'SI' tst ;
  1089. 'SI' tst1 ; tt = TDISC . nomt ; 'FINSI' ;
  1090. 'SI' tst2 ; tt = TDISC . nomq ; 'FINSI' ;
  1091. isvec = ('>' ('DIME' (tt. 'NOMPRI')) 1) ;
  1092. 'SI' isvec ; 'ARGUMENT' valv*'LISTREEL' ; 'SINON' ;
  1093. 'ARGUMENT' valv*'FLOTTANT' ;
  1094. 'FINSI' ;
  1095. 'SI' tst1 ; valt = valv ; 'FINSI' ;
  1096. 'SI' tst2 ; valq = valv ; 'FINSI' ;
  1097. 'FINSI' ;
  1098. 'SI' ('EGA' motcle 'FCOF') ; 'ARGUMENT' valo*'FLOTTANT' ; 'FINSI' ;
  1099. 'SI' ('EGA' motcle 'CPRI') ; 'ARGUMENT' valt*'CHPOINT' ; 'FINSI' ;
  1100. 'SI' ('EGA' motcle 'CDUA') ; 'ARGUMENT' valq*'CHPOINT' ; 'FINSI' ;
  1101. 'SI' ('EGA' motcle 'CCOF') ; 'ARGUMENT' valo*'CHPOINT' ; 'FINSI' ;
  1102. 'FIN' imotcle ;
  1103. *
  1104. * Tests
  1105. *
  1106. discg = TDISC . 'GEOM' . 'DISC' ;
  1107. 'SI' ('EXISTE' tdisc 'methgau') ;
  1108. methgau = tdisc . 'methgau' . 'mass' ;
  1109. 'SINON' ;
  1110. methgau = 'GAU7' ;
  1111. 'FINSI' ;
  1112. tnomt = TDISC . nomt ;
  1113. lvalt = 'NEG' ('TYPE' valt) 'MOT' ;
  1114. tnomq = TDISC . nomq ;
  1115. lvalq = 'NEG' ('TYPE' valq) 'MOT' ;
  1116. * Scalaire ou vecteur
  1117. ninct = 'DIME' (tnomt . 'NOMPRI') ;
  1118. nincq = 'DIME' (tnomq . 'NOMPRI') ;
  1119. 'SI' ('NEG' ninct nincq) ;
  1120. cherr = 'CHAINE'
  1121. 'les primales et duales nont pas le meme nombre de composantes' ;
  1122. 'ERREUR' cherr ;
  1123. 'FINSI' ;
  1124. ninc = ninct ;
  1125. *
  1126. lcof = 'EXISTE' TDISC nomo ;
  1127. 'SI' lcof ; ncof = 1 ; tcof = TDISC . nomo ;
  1128. 'SINON' ; ncof = 0 ;
  1129. 'FINSI' ;
  1130. *
  1131. 'SI' debug ;
  1132. 'SI' lcof ; 'MESSAGE' 'Un coef a ete detecte' ;
  1133. 'SINON' ; 'MESSAGE' 'pas de coef detecte' ;
  1134. 'FINSI' ;
  1135. 'FINSI' ;
  1136. *
  1137. vdim = 'VALEUR' 'DIME' ;
  1138. vmod = 'VALEUR' 'MODE' ;
  1139. idim = 0 ;
  1140. 'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'PLANDEFO')) ;
  1141. idim = 2 ;
  1142. iaxi = FAUX ;
  1143. 'FINSI' ;
  1144. 'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'AXIS')) ;
  1145. idim = 2 ;
  1146. iaxi = VRAI ;
  1147. 'FINSI' ;
  1148. 'SI' ('ET' ('EGA' vdim 3) ('EGA' vmod 'TRID')) ;
  1149. idim = 3 ;
  1150. iaxi = FAUX ;
  1151. 'FINSI' ;
  1152. 'SI' ('EGA' vdim 1) ;
  1153. idim = 1 ;
  1154. iaxi = FAUX ;
  1155. 'FINSI' ;
  1156. * 'MESSAGE' ('CHAINE' 'iaxi=' iaxi );
  1157. 'SI' ('EGA' idim 0) ;
  1158. 'ERREUR' ('CHAINE' 'vmod=' vmod ' et vdim=' vdim ' non prevu') ;
  1159. 'FINSI' ;
  1160. 'SI' iaxi ;
  1161. dprmt = '*' ('COORDONNEE' 1 _mt) ('*' PI 2.D0) ;
  1162. 'FINSI' ;
  1163. *
  1164. * Optimisation possible : construire la matrice par blocs
  1165. * qd valt et valq ne sont pas donnés
  1166. *
  1167. numop = ninc ; numder = idim ; numvar = ninc ;
  1168. numdat = ncof ; numcof = ncof ;
  1169. A = ININLIN numop numvar numdat numcof numder ;
  1170. 'SI' lcof ;
  1171. A . 'DAT' . 1 . 'NOMDDL' = tcof . 'NOMPRI' . 1 ;
  1172. A . 'DAT' . 1 . 'DISC' = tcof . 'DISC' ;
  1173. A . 'DAT' . 1 . 'VALEUR' = valo ;
  1174. A . 'COF' . 1 . 'COMPOR' = 'IDEN' ;
  1175. A . 'COF' . 1 . 'LDAT' = 'LECT' 1 ;
  1176. ll = 'LECT' 1 ;
  1177. 'SINON' ;
  1178. ll = 'LECT' ;
  1179. 'FINSI' ;
  1180. lvt = 'EGA' ('TYPE' valt) 'LISTREEL' ;
  1181. 'REPETER' iiinc ninc ;
  1182. iinc = &iiinc ;
  1183. A . 'VAR' . iinc . 'NOMDDL' = tnomt . 'NOMPRI' . iinc ;
  1184. A . 'VAR' . iinc . 'DISC' = tnomt . 'DISC' ;
  1185. 'SI' lvalt ;
  1186. 'SI' lvt ;
  1187. A . 'VAR' . iinc . 'VALEUR' = 'EXTRAIRE' valt iinc ;
  1188. 'SINON' ;
  1189. A . 'VAR' . iinc . 'VALEUR' = valt ;
  1190. 'FINSI' ;
  1191. 'FINSI' ;
  1192. A . iinc . iinc . 0 = ll ;
  1193. 'FIN' iiinc ;
  1194. *
  1195. 'SI' iaxi ;
  1196. numdat = 1 ;
  1197. numcof = 1 ;
  1198. 'SINON' ;
  1199. numdat = 0 ;
  1200. numcof = 0 ;
  1201. 'FINSI' ;
  1202. B = ININLIN numop numvar numdat numcof numder ;
  1203. 'SI' iaxi ;
  1204. B . 'DAT' . 1 . 'NOMDDL' = 'MOTS' 'SCAL' ;
  1205. B . 'DAT' . 1 . 'DISC' = discg ;
  1206. B . 'DAT' . 1 . 'VALEUR' = dprmt ;
  1207. B . 'COF' . 1 . 'COMPOR' = 'IDEN' ;
  1208. B . 'COF' . 1 . 'LDAT' = 'LECT' 1 ;
  1209. ll = 'LECT' 1 ;
  1210. 'SINON' ;
  1211. ll = 'LECT' ;
  1212. 'FINSI' ;
  1213. lvq = 'EGA' ('TYPE' valq) 'LISTREEL' ;
  1214. 'REPETER' iiinc ninc ;
  1215. iinc = &iiinc ;
  1216. B . 'VAR' . iinc . 'NOMDDL' = tnomq . 'NOMDUA' . iinc ;
  1217. B . 'VAR' . iinc . 'DISC' = tnomq . 'DISC' ;
  1218. 'SI' lvalq ;
  1219. 'SI' lvq ;
  1220. B . 'VAR' . iinc . 'VALEUR' = 'EXTRAIRE' valq iinc ;
  1221. 'SINON' ;
  1222. B . 'VAR' . iinc . 'VALEUR' = valq ;
  1223. 'FINSI' ;
  1224. 'FINSI' ;
  1225. B . iinc . iinc . 0 = ll ;
  1226. 'FIN' iiinc ;
  1227. *
  1228. 'SI' ('EXISTE' _smt) ;
  1229. mgmass2 = 'NLIN' discg _mt _smt A B methgau ;
  1230. 'SINON' ;
  1231. mgmass2 = NLINP discg _mt A B methgau ;
  1232. 'FINSI' ;
  1233. *
  1234. 'RESPRO' mgmass2 ;
  1235. 'FINPROC' ;
  1236. *
  1237. * End of procedure file GMASS2
  1238. *
  1239. *ENDPROCEDUR gmass2
  1240. *
  1241. *
  1242. *************************************************************************
  1243. *
  1244. *
  1245. * MAIN : 1) MESH
  1246. * 2) COMPUTATIONAL LOOP
  1247. * 3) TESTs
  1248. *
  1249. ************************************************************************
  1250. *
  1251. *
  1252. * Maillage et modèle
  1253. *
  1254. *nx = 1 ; ny =1 ;
  1255. *nx = 1 ; ny =2 ;
  1256. nx = 10 ; ny= 10 ;
  1257. pA = 0. 0. ; pB = 1. 0. ; pC = 1. 1. ; pD = 0. 1. ;
  1258. ang = 32.1 ;
  1259. pA pB pC pD = 'TOURNER' pA pB pC pD ang pA ;
  1260. discu = 'QUAF' ; discp = 'LINE' ;
  1261. 'OPTI' 'ELEM' 'QUA4' ;
  1262. tol = 1.D-6 ;
  1263. bas = 'DROIT' nx pA pB ;
  1264. dro = 'DROIT' ny pB pC ;
  1265. hau = 'DROIT' nx pC pD ;
  1266. gau = 'DROIT' ny pD pA ;
  1267. mt = 'DALLER' bas dro hau gau ;
  1268. *
  1269. tmail = QUAFME mt 'mt' tol ;
  1270. tmail = QUAFME tmail bas 'bas' dro 'dro' hau 'hau' gau 'gau' tol ;
  1271. _mt = tmail . 'mt' . 'QUAF' ;
  1272. mt = tmail . 'mt' . discu ;
  1273. mtp = tmail . 'mt' . discp ;
  1274. cmt = 'CONTOUR' mt ;
  1275. bas = tmail . 'bas' . discu ;
  1276. dro = tmail . 'dro' . discu ;
  1277. hau = tmail . 'hau' . discu ;
  1278. gau = tmail . 'gau' . discu ;
  1279. _bas = tmail . 'bas' . 'QUAF' ;
  1280. _dro = tmail . 'dro' . 'QUAF' ;
  1281. _hau = tmail . 'hau' . 'QUAF' ;
  1282. _gau = tmail . 'gau' . 'QUAF' ;
  1283. *
  1284. tmode = MODENLIN 'GEOM' 'LINE'
  1285. 'INCO' 'UN' discu 'VECT' 'UX' 'UY' 'FX' 'FY'
  1286. 'INCO' 'PN' discp 'SCAL' 'LXP' 'LXP'
  1287. ;
  1288. nompp = 'EXTRAIRE' (tmode . 'PN' .'NOMPRI' . 1) ('LECT' 1) ;
  1289. nompd = 'EXTRAIRE' (tmode . 'PN' .'NOMDUA' . 1) ('LECT' 1) ;
  1290. nomvp = 'ET' (tmode . 'UN' .'NOMPRI' . 1)
  1291. (tmode . 'UN' .'NOMPRI' . 2) ;
  1292. nomvd = 'ET' (tmode . 'UN' .'NOMDUA' . 1)
  1293. (tmode . 'UN' .'NOMDUA' . 2) ;
  1294. discg = tmode . 'GEOM' . 'DISC' ;
  1295. *
  1296. * Matrice de rigidité
  1297. *
  1298. *rig = GRIG2 _mt tmode 'NPRI' 'UN' 'NDUA' 'UN' 'LAPN' ;
  1299. rig = GRIG2 _mt tmode 'NPRI' 'UN' 'NDUA' 'UN' ;
  1300. *
  1301. * Matrice de la contrainte de divergence nulle et sa transposée
  1302. *
  1303. cnt1 = GDIV2 _mt tmode 'NPRI' 'UN' 'NDUA' 'PN' ;
  1304. cnt2 = GDIV2 _mt tmode 'NPRI' 'UN' 'NDUA' 'PN' 'GMBT' ;
  1305. mcnt = cnt1 'ET' cnt2 ;
  1306. *
  1307. * Conditions aux limites
  1308. *
  1309. * On relache la contrainte sur la vitesse normale en 1 point (pmil) pour
  1310. * fixer le mode de pression constante
  1311. pmil = * ('PLUS' pC pD) 0.5 ;
  1312. pmil = 'POIN' hau 'PROC' pmil ;
  1313. haur = 'DIFF' ('CHANGER' 'POI1' hau)
  1314. (pmil 'ET' pC 'ET' pD) ;
  1315. gaur = 'DIFF' ('CHANGER' 'POI1' gau)
  1316. ('ET' pA pD) ;
  1317. dror = 'DIFF' ('CHANGER' 'POI1' dro)
  1318. ('ET' pB pC) ;
  1319. xbas ybas = 'COORDONNEE' bas ;
  1320. sbas = '**' ('+' ('**' xbas 2.) ('**' ybas 2)) 0.5 ;
  1321. sbas = '**' ('SIN' ('*' sbas 180.)) 2. ;
  1322. ssbas = 'SIN' ('*' 90. sbas) ;
  1323. vAB = 'MOIN' pB pA ;
  1324. xab = 'COORDONNEE' 1 vAB ;
  1325. yab = 'COORDONNEE' 2 vAB ;
  1326. *
  1327. ssbas = '+' ('NOMC' 'UX' ('*' ssbas xab))
  1328. ('NOMC' 'UY' ('*' ssbas yab)) ;
  1329. *
  1330. * Construction de la normale
  1331. *
  1332. vnorh = GNOR2 _hau tmode 'NPRI' discg 'FPRI' 1. 'NDUA' 'UN' ;
  1333. vnorg = GNOR2 _gau tmode 'NPRI' discg 'FPRI' 1. 'NDUA' 'UN' ;
  1334. vnord = GNOR2 _dro tmode 'NPRI' discg 'FPRI' 1. 'NDUA' 'UN' ;
  1335. vnorh = 'NOMC' nomvd nomvp vnorh ;
  1336. vnorg = 'NOMC' nomvd nomvp vnorg ;
  1337. vnord = 'NOMC' nomvd nomvp vnord ;
  1338. *
  1339. mblh = 'BLOQUE' 'DEPL' 'DIRE' vnorh haur ;
  1340. mblg = 'BLOQUE' 'DEPL' 'DIRE' vnorg gaur ;
  1341. mbld = 'BLOQUE' 'DEPL' 'DIRE' vnord dror ;
  1342. mblb = 'BLOQUE' 'DEPL' bas ;
  1343. mblp = 'BLOQUE' 'DEPL' (pC 'ET' pD) ;
  1344. bltot = mblh 'ET' mblg 'ET' mbld 'ET' mblb 'ET' mblp ;
  1345. fbltot = 'DEPIMPOSE' mblb ssbas ;
  1346. *
  1347. * Matrice et second membre totaux
  1348. * Résolution
  1349. *
  1350. mtot = rig 'ET' mcnt 'ET' bltot ;
  1351. *mtot = rig 'ET' bltot ;
  1352. ftot = fbltot ;
  1353. *'LISTE' ('EXTRAIRE' mtot 'COMP') ;
  1354. *'LISTE' ('EXTRAIRE' mtot 'COMP' 'DUAL') ;
  1355. *'LISTE' ('EXTRAIRE' ftot 'COMP') ;
  1356. * opti impi 1 ;
  1357. lpass = VRAI ;
  1358. sol = 'KRES' mtot ftot 'LDEPE' FAUX 'IMPR' 2 ;
  1359. sol2 = 'KRES' mtot ftot 'LDEPE' VRAI 'IMPR' 2 ;
  1360. vit = 'EXCO' ('MOTS' 'UX' 'UY') sol ;
  1361. 'SI' graph ;
  1362. vvit = 'VECT' vit 'DEPL' 'JAUN' ;
  1363. 'TRACER' vvit cmt 'TITR' 'Vitesses' ;
  1364. 'FINSI' ;
  1365. *
  1366. * Comme l'élimination des dépendances est un préconditionnement,
  1367. * on peut vérifier qu'il marche
  1368. *
  1369. sol3 = 'KRES' mtot ftot 'LDEPE' FAUX 'IMPR' 2 ;
  1370. sol4 = 'KRES' mtot ftot 'LDEPE' VRAI 'IMPR' 2 ;
  1371. errv = 1.D-14 ;
  1372. err1 = 'MAXIMUM' ('-' sol sol3 'ABS') ;
  1373. err2 = 'MAXIMUM' ('-' sol2 sol4 'ABS') ;
  1374. tst1 = '<' err1 errv ; tst2 = '<' err2 errv ;
  1375. tst = tst1 'ET' tst2 ;
  1376. 'MESSAGE' ('CHAINE' 'Test 1 sur la reproductibilité de lelimination') ;
  1377. 'MESSAGE' ' err1 = ' err1 ;
  1378. 'MESSAGE' ' err2 = ' err2 ;
  1379. 'SI' tst ;
  1380. 'MESSAGE' 'Test 1 OK' ;
  1381. 'SINON' ;
  1382. 'MESSAGE' '!!! Test 1 not passed ' ;
  1383. 'FINSI' ;
  1384. lpass = lpass 'ET' tst ;
  1385. *
  1386. * Vérification de CONDENSE et EVAPORE
  1387. *
  1388. mtotb = '*' mtot 1.D0 ;
  1389. mtotc ftotc ftot1 = 'KOPS' 'CONDENSE' mtot ftot ;
  1390. solc = 'KRES' mtotc ftotc 'LDEPE' FAUX 'IMPR' 2 ;
  1391. * Vérification de la résolution du système réduit (cf. fiche ano 8186)
  1392. err0 = '/' ('MAXIMUM' ('-' ftotc ('*' mtotc solc)) 'ABS')
  1393. ('MAXIMUM' ftotc 'ABS') ;
  1394. sol5 = 'KOPS' 'EVAPORE' solc mtot ftot ftot1 ;
  1395. mtotc ftotc ftot1 = 'KOPS' 'CONDENSE' mtotb ftot ;
  1396. solc = 'KRES' mtotc ftotc 'LDEPE' FAUX 'IMPR' 2 ;
  1397. sol6 = 'KOPS' 'EVAPORE' solc mtot ftot ftot1 ;
  1398. errv = 1.D-14 ;
  1399. err1 = 'MAXIMUM' ('-' sol2 sol5 'ABS') ;
  1400. err2 = 'MAXIMUM' ('-' sol2 sol6 'ABS') ;
  1401. tst0 = '<' err0 errv ;
  1402. tst1 = '<' err1 errv ; tst2 = '<' err2 errv ;
  1403. tst = tst0 'ET' tst1 'ET' tst2 ;
  1404. 'MESSAGE' ('CHAINE' 'Test 2 sur CONDENSE et EVAPORE') ;
  1405. 'MESSAGE' ' err0 = ' err0 ;
  1406. 'MESSAGE' ' err1 = ' err1 ;
  1407. 'MESSAGE' ' err2 = ' err2 ;
  1408. 'SI' tst ;
  1409. 'MESSAGE' 'Test 2 OK' ;
  1410. 'SINON' ;
  1411. 'MESSAGE' '!!! Test 2 not passed ' ;
  1412. 'FINSI' ;
  1413. lpass = lpass 'ET' tst ;
  1414. *
  1415. * On vérifie que l'élimination des conditions aux limites sur
  1416. * la vitesse s'est bien faite.
  1417. * Les résultats avec et sans élimination des conditions aux limites ne
  1418. * sont pas si proches car la matrice Vitesse-Pression est mal
  1419. * conditionnée (l'erreur est plus élevée en pression)
  1420. *
  1421. errv = 1.D-13 ; errp = 1.D-11 ; errlx = 4.D-13 ;
  1422. err1 = 'MAXIMUM' ('-' sol sol2) 'ABS' 'AVEC' nomvp ;
  1423. err2 = 'MAXIMUM' ('-' sol sol2) 'ABS' 'AVEC' nompp ;
  1424. err3 = 'MAXIMUM' ('-' sol sol2) 'ABS' 'SANS' ('ET' nomvp nompp) ;
  1425. tst1 = '<' err1 errv ;
  1426. tst2 = '<' err2 errp ;
  1427. tst3 = '<' err3 errlx ;
  1428. tst = tst1 'ET' tst2 'ET' tst3 ;
  1429. 'MESSAGE' ('CHAINE' 'Test 3 sur la solution avec OU sans elimination'
  1430. ' des conditions aux limites') ;
  1431. 'MESSAGE' ' err UX = ' err1 ;
  1432. 'MESSAGE' ' err P = ' err2 ;
  1433. 'MESSAGE' ' err LX = ' err3 ;
  1434. 'SI' tst ;
  1435. 'MESSAGE' 'Test 3 OK' ;
  1436. 'SINON' ;
  1437. 'MESSAGE' '!!! Test 3 not passed ' ;
  1438. 'FINSI' ;
  1439. lpass = lpass 'ET' tst ;
  1440. *
  1441. *
  1442. * On essaie le raffinement itératif pour rapprocher
  1443. * les solutions avec et sans élimination
  1444. * On choisit nit=1 car sinon on dégrade à nouveau la précision
  1445. *
  1446. nit = 1 ;
  1447. sol7 = sol ;
  1448. 'REPETER' it nit ;
  1449. smb = '-' ftot ('*' mtot sol7) ;
  1450. dsol = 'KRES' mtot smb 'LDEPE' FAUX ;
  1451. sol7 = '+' sol7 dsol ;
  1452. 'FIN' it ;
  1453. sol8 = sol2 ;
  1454. 'REPETER' it nit ;
  1455. smb = '-' ftot ('*' mtot sol8) ;
  1456. dsol = 'KRES' mtot smb 'LDEPE' VRAI ;
  1457. sol8 = '+' sol8 dsol ;
  1458. 'FIN' it ;
  1459. errv = 2.D-14 ; errp = 2.D-12 ; errlx = 8.D-14 ;
  1460. err1 = 'MAXIMUM' ('-' sol7 sol8) 'ABS' 'AVEC' nomvp ;
  1461. err2 = 'MAXIMUM' ('-' sol7 sol8) 'ABS' 'AVEC' nompp ;
  1462. err3 = 'MAXIMUM' ('-' sol7 sol8) 'ABS' 'SANS' ('ET' nomvp nompp) ;
  1463. tst1 = '<' err1 errv ;
  1464. tst2 = '<' err2 errp ;
  1465. tst3 = '<' err3 errlx ;
  1466. tst = tst1 'ET' tst2 'ET' tst3 ;
  1467. 'MESSAGE' ('CHAINE' 'Test 4 sur la solution avec OU sans elimination'
  1468. ' 'ET' raffinement iteratif') ;
  1469. 'MESSAGE' ' err UX = ' err1 ;
  1470. 'MESSAGE' ' err P = ' err2 ;
  1471. 'MESSAGE' ' err LX = ' err3 ;
  1472. 'SI' tst ;
  1473. 'MESSAGE' 'Test 4 OK' ;
  1474. 'SINON' ;
  1475. 'MESSAGE' '!!! Test 4 not passed ' ;
  1476. 'FINSI' ;
  1477. lpass = lpass 'ET' tst ;
  1478. * On choisit ensuite sol8 comme solution de référence
  1479. solref = sol8 ;
  1480. *
  1481. * Méthode de lagrangien augmenté itérée
  1482. *
  1483. * Eliminons d'abord le plus de contraintes possible
  1484. mtotc ftotc ftot1 = 'KOPS' 'CONDENSE' mtot ftot ;
  1485. * On vérifie que les 'FLX' ont disparus dans la matrice condensée et de
  1486. * la force (en fait elle existe toujours dans la force)
  1487. * MAJ: suite à la fiche ano 8186, on peut rétablir le test
  1488. lp1 = 'EXTRAIRE' mtotc 'COMP' 'DUAL' ;
  1489. lp2 = 'EXTRAIRE' ftotc 'COMP' ;
  1490. tst1 = 'NON' ('EXISTE' lp1 'FLX') ;
  1491. tst2 = 'NON' ('EXISTE' lp2 'FLX') ;
  1492. tst = tst1 'ET' tst2 ;
  1493. *tst = tst1 ;
  1494. 'MESSAGE' ('CHAINE' 'Test 5 sur condensation des LX et FLX') ;
  1495. 'MESSAGE' ('CHAINE' 'lp1=') ; 'LISTE' lp1 ;
  1496. 'MESSAGE' ('CHAINE' 'lp2=') ; 'LISTE' lp2 ;
  1497. 'MESSAGE' ('CHAINE' ' tst1 = ' tst1) ;
  1498. 'MESSAGE' ('CHAINE' ' tst2 = ' tst2) ;
  1499. 'SI' tst ;
  1500. 'MESSAGE' 'Test 5 OK' ;
  1501. 'SINON' ;
  1502. 'MESSAGE' '!!! Test 5 not passed ' ;
  1503. 'FINSI' ;
  1504. lpass = lpass 'ET' tst ;
  1505. * Ensuite extraction rigidité et contraintes
  1506. krig = 'EXTRAIRE' mtotc nomvp nomvd ;
  1507. kcnt = 'EXTRAIRE' mtotc nomvp ('MOTS' 'LXP') ;
  1508. kcnt2 = 'EXTRAIRE' mtotc ('MOTS' 'LXP') nomvd ;
  1509. * Matrice masse (+ diagonalisé pour la pression)
  1510. kmasp = GMASS2 _mt tmode 'NPRI' 'PN' 'NDUA' 'PN' ;
  1511. kmaspd = GMASS2 _mt tmode 'NPRI' 'PN' 'FPRI' 1.D0 'NDUA' 'PN' ;
  1512. mas = 'MAXIMUM' ('RESULT' kmaspd) ;
  1513. dkmasp = 'EXTRAIRE' kmasp 'DIAG' ;
  1514. masd = 'MAXIMUM' ('RESULT' dkmasp) ;
  1515. 'MESSAGE' ('CHAINE' 'mas=' mas) ;
  1516. 'MESSAGE' ('CHAINE' 'masd=' masd) ;
  1517. * Cela ne paraît pas forcément utile
  1518. *dkmass = '*' dkmasp ('/' mas masd) ;
  1519. dkmass = dkmasp ;
  1520. dkmassi = 'INVERSE' dkmass ;
  1521. *
  1522. * Matrice de pénalisation
  1523. *
  1524. mlagd = 'KOPS' 'CMCT' kcnt2 dkmassi kcnt2 ;
  1525. * Scaling
  1526. drig = 'EXTRAIRE' krig 'DIAG' ;
  1527. dlagd = 'EXTRAIRE' mlagd 'DIAG' ;
  1528. dsca = '/' drig dlagd nomvp nomvp nomvp ;
  1529. * Petite vérif de la division nouvellement créée
  1530. drrig = '*' dsca dlagd nomvp nomvp nomvp ;
  1531. errv = 1.D-14 ;
  1532. err1 = 'MAXIMUM' ('-' drrig drig 'ABS') ;
  1533. tst1 = '<' err1 errv ;
  1534. tst = tst1 ;
  1535. 'MESSAGE' ('CHAINE' 'Test 6 sur /') ;
  1536. 'MESSAGE' ' err1 = ' err1 ;
  1537. 'SI' tst ;
  1538. 'MESSAGE' 'Test 6 OK' ;
  1539. 'SINON' ;
  1540. 'MESSAGE' '!!! Test 6 not passed ' ;
  1541. 'FINSI' ;
  1542. lpass = lpass 'ET' tst ;
  1543. *'OPTI' 'ECHO' 1 ;
  1544. *'LISTE' ('MAXIMUM' dsca) ;
  1545. *'LISTE' ('MINIMUM' dsca) ;
  1546. *'OPTION' 'DONN' 5 ;
  1547. * Matrice de pénalisation scalée : on scale kcnt
  1548. **mdsca = 'MANUEL' 'RIGI' ('**' dsca 0.5) ;
  1549. **kcnts = 'KOPS' 'CMCT' kcnt mdsca ;
  1550. **kcnt2s = 'KOPS' 'TRANSPOS' kcnts ;
  1551. *mlagds = '*' mlagd ('MINIMUM' dsca) ;
  1552. *dlagds = 'EXTRAIRE' mlagds 'DIAG' ;
  1553. *dscas = '/' drig dlagds nomvp nomvp nomvp ;
  1554. *'LISTE' ('MAXIMUM' dscas) ;
  1555. *'LISTE' ('MINIMUM' dscas) ;
  1556. *
  1557. lalphad = 'PROG' 1. 10. 100. ;
  1558. lerrv = 'PROG' 9.D-2 9.D-4 9.D-6 ;
  1559. lerrp = 'PROG' 6.D-1 6.D-2 6.D-4 ;
  1560. lerrlx = 'PROG' 5.D-2 5.D-4 5.D-6 ;
  1561. ncas = 'DIME' lalphad ;
  1562. 'REPETER' iicas ncas ;
  1563. icas = &iicas ;
  1564. alphadd = 'EXTRAIRE' lalphad icas ;
  1565. alphad = '*' alphadd ('MINIMUM' dsca) ;
  1566. mlagds = '*' mlagd alphad ;
  1567. * Directe
  1568. solc = 'KRES' mtotc ftotc 'LDEPE' FAUX ;
  1569. *
  1570. nit = 2 ;
  1571. sol9c = '*' solc 0. ;
  1572. ftot9c = ftotc ;
  1573. mtot9c = mtotc ;
  1574. maug = krig 'ET' mlagds ;
  1575. 'REPETER' it nit ;
  1576. smb = '-' ftot9c ('*' mtot9c sol9c) ;
  1577. smbf = 'EXCO' nomvd smb ;
  1578. smbg = 'EXCO' nompd smb ;
  1579. smbf2 = '*' ('*' ('*' dkmassi smbg nompp nompd nompp) kcnt2)
  1580. alphad ;
  1581. smbaug = '+' smbf smbf2 ;
  1582. dsolu tt = 'KRES' maug smbaug 'LDEPE' FAUX
  1583. 'LTIME' VRAI 'IMPR' 1 'TYPINV' 1 ;
  1584. * 'LISTE' tt ;
  1585. dsolp = '*' ('*' dkmassi ('-' ('*' kcnt dsolu) smbg)
  1586. nompp nompd nompp)
  1587. ('*' alphad +1.) ;
  1588. dsol = dsolu '+' dsolp ;
  1589. sol9c = sol9c '+' dsol ;
  1590. 'FIN' it ;
  1591. errv = 'EXTRAIRE' lerrv icas ;
  1592. errp = 'EXTRAIRE' lerrp icas ;
  1593. errlx = 'EXTRAIRE' lerrlx icas ;
  1594. sol9 = 'KOPS' 'EVAPORE' sol9c mtot ftot ftot1 ;
  1595. 'SI' graph ;
  1596. vvit2 = 'VECT' sol9 'DEPL' 'ROUG' ;
  1597. tit = 'CHAINE' 'Vitesses ' 'jaun=direct ' 'roug=lagrangien alpha='
  1598. alphadd ;
  1599. 'TRACER' ('ET' vvit vvit2) cmt 'TITR' tit ;
  1600. 'FINSI' ;
  1601. err1 = 'MAXIMUM' ('-' sol8 sol9) 'ABS' 'AVEC' nomvp ;
  1602. err2 = 'MAXIMUM' ('-' sol8 sol9) 'ABS' 'AVEC' nompp ;
  1603. err3 = 'MAXIMUM' ('-' sol8 sol9) 'ABS' 'SANS' ('ET' nomvp nompp) ;
  1604. tst1 = '<' err1 errv ;
  1605. tst2 = '<' err2 errp ;
  1606. tst3 = '<' err3 errlx ;
  1607. tst = tst1 'ET' tst2 'ET' tst3 ;
  1608. 'MESSAGE' ('CHAINE' 'Test 7 Lagrangien augmenté ' nit ' iterations '
  1609. 'alpha=' alphadd) ;
  1610. 'MESSAGE' ('CHAINE' ' err UX = ' err1) ;
  1611. 'MESSAGE' ('CHAINE' ' err P = ' err2) ;
  1612. 'MESSAGE' ('CHAINE' ' err LX = ' err3) ;
  1613. 'SI' tst ;
  1614. 'MESSAGE' 'Test 7 OK' ;
  1615. 'SINON' ;
  1616. 'MESSAGE' '!!! Test 7 not passed ' ;
  1617. 'FINSI' ;
  1618. lpass = lpass 'ET' tst ;
  1619. 'FIN' iicas ;
  1620. *
  1621. 'SAUTER' 2 'LIGNE' ;
  1622. 'SI' lpass ;
  1623. 'MESSAGE' 'Tout sest bien passe' ;
  1624. 'SINON' ;
  1625. 'MESSAGE' 'Il y a eu des erreurs' ;
  1626. 'FINSI' ;
  1627. 'SAUTER' 2 'LIGNE' ;
  1628. 'SI' ('NON' lpass) ;
  1629. 'ERREUR' 5 ;
  1630. 'FINSI' ;
  1631. ************************************************************************
  1632. *
  1633. * END OF TEST PART
  1634. *
  1635. ************************************************************************
  1636. 'FINSI' ;
  1637. 'SI' interact ;
  1638. 'OPTION' 'DONN' 5 'ECHO' 1 ;
  1639. 'FINSI' ;
  1640. *
  1641. 'FIN' ;
  1642.  
  1643.  
  1644.  
  1645.  
  1646.  
  1647.  
  1648.  
  1649.  
  1650.  
  1651.  
  1652.  
  1653.  
  1654.  
  1655.  
  1656.  
  1657.  

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