Télécharger drop.dgibi

Retour à la liste

Numérotation des lignes :

  1. * fichier : drop.dgibi
  2. ************************************************************************
  3. ************************************************************************
  4. *
  5. 'OPTI' 'ECHO' 0 ;
  6. *
  7. 'SAUTER' 2 'LIGNE' ;
  8. 'MESSAGE' ' Execution de drop.dgibi' ;
  9. 'SAUTER' 2 'LIGNE' ;
  10. *
  11. graph = faux ;
  12. complet = faux ;
  13. interact = faux ;
  14. lmatrik = faux ;
  15. ************************************************************************
  16. * NOM : DROP
  17. * DESCRIPTION : Une goutte plane ou axi soumise à la gravité et à
  18. * la tension de surface.
  19. * Contraintes : on fixe le Delta P ou le Volume
  20. * on fixe la position des points ou l'angle
  21. * (cf. table tclim)
  22. *
  23. * A plane or axisymmetric drop subject to surface tension
  24. * and gravity and to the following constraints or forces
  25. * 1) Constant volume or constant pressure difference
  26. * between the interior and exterior
  27. * 2) On the triple line: given contact angle or given
  28. * position
  29. * (see tclim TABLE)
  30. *
  31. * Reference solution :
  32. * -> No gravity : the drop is spherical => analytical solution
  33. * -> upward gravity :
  34. *@Article{sumesh,
  35. * author = {P.T. Sumesh and Rama Govindrajan},
  36. * title = {The possible equilibrium shapes of static pendant drops},
  37. * journal = {The Journal of Chemical Physics},
  38. * year = {2010},
  39. * key = {144707},
  40. * volume = {133},
  41. * pages = {1--8},
  42. *}
  43. *
  44. *
  45. * LANGAGE : GIBIANE-CAST3M
  46. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  47. * mél : stephane.gounand@cea.fr
  48. **********************************************************************
  49. * VERSION : v1, 22/04/2011, version initiale
  50. * HISTORIQUE : v1, 22/04/2011, création
  51. * HISTORIQUE :
  52. * HISTORIQUE :
  53. ************************************************************************
  54. *
  55. 'OPTION' 'DIME' 2 'ELEM' 'QUA4' 'MODE' 'PLAN' ;
  56. ************************************************************************
  57. *
  58. *
  59. * PROCEDURES
  60. *
  61. *
  62. ************************************************************************
  63. *BEGINPROCEDUR affvar
  64. ************************************************************************
  65. * NOM : AFFVAR
  66. * DESCRIPTION : Affiche des variables
  67. *
  68. *
  69. *
  70. * LANGAGE : GIBIANE-CAST3M
  71. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  72. * mél : gounand@semt2.smts.cea.fr
  73. **********************************************************************
  74. *
  75. *
  76. 'DEBPROC' AFFVAR ;
  77. 'REPETER' bcl ;
  78. 'ARGUMENT' x/'FLOTTANT' ;
  79. 'SI' ('EXISTE' x) ;
  80. 'ARGUMENT' lx*'MOT' ;
  81. 'MESSAGE' ('CHAINE' lx '=' x) ;
  82. 'SINON' ;
  83. 'QUITTER' bcl ;
  84. 'FINSI' ;
  85. 'FIN' bcl ;
  86. 'FINPROC' ;
  87. *
  88. * End of procedure file AFFVAR
  89. *
  90. *ENDPROCEDUR affvar
  91. *BEGINPROCEDUR append
  92. ************************************************************************
  93. * NOM : APPEND
  94. * DESCRIPTION : Rajoute :
  95. * - un entier à un listentier
  96. * - un réel à un listreel
  97. * - un objet (liste, evolution, matrice ou chpoint)
  98. * à un indice de table ('MOT' ou 'ENTIER')
  99. * * si l'indice n'existe pas
  100. * * 'ET' si l'indice existe
  101. *
  102. *
  103. *
  104. * LANGAGE : GIBIANE-CAST3M
  105. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  106. * mél : gounand@semt2.smts.cea.fr
  107. **********************************************************************
  108. * VERSION : v1, 10/09/2004, version initiale
  109. * HISTORIQUE : v1, 10/09/2004, création
  110. * HISTORIQUE :
  111. * HISTORIQUE :
  112. ************************************************************************
  113. * Prière de PRENDRE LE TEMPS de compléter les commentaires
  114. * en cas de modification de ce sous-programme afin de faciliter
  115. * la maintenance !
  116. ************************************************************************
  117. *
  118. *
  119. 'DEBPROC' APPEND ;
  120. 'ARGUMENT' tab/'TABLE' ;
  121. 'SI' ('EXISTE' tab) ;
  122. 'ARGUMENT' itab/'MOT' ;
  123. 'SI' ('NON' ('EXISTE' itab)) ;
  124. 'ARGUMENT' itab*'ENTIER' ;
  125. 'FINSI' ;
  126. lobj = FAUX ;
  127. 'SI' ('NON' lobj) ;
  128. 'ARGUMENT' lr/'LISTREEL' ;
  129. 'SI' ('EXISTE' lr) ;
  130. obj = lr ; lobj = VRAI ;
  131. 'FINSI' ;
  132. 'FINSI' ;
  133. 'SI' ('NON' lobj) ;
  134. 'ARGUMENT' le/'LISTENTI' ;
  135. 'SI' ('EXISTE' le) ;
  136. obj = le ; lobj = VRAI ;
  137. 'FINSI' ;
  138. 'FINSI' ;
  139. 'SI' ('NON' lobj) ;
  140. 'ARGUMENT' lev/'EVOLUTION' ;
  141. 'SI' ('EXISTE' lev) ;
  142. obj = lev ; lobj = VRAI ;
  143. 'FINSI' ;
  144. 'FINSI' ;
  145. 'SI' ('NON' lobj) ;
  146. 'ARGUMENT' lm/'MAILLAGE' ;
  147. 'SI' ('EXISTE' lm) ;
  148. obj = lm ; lobj = VRAI ;
  149. 'FINSI' ;
  150. 'FINSI' ;
  151. 'SI' ('NON' lobj) ;
  152. 'ARGUMENT' chpo/'CHPOINT' ;
  153. 'SI' ('EXISTE' chpo) ;
  154. obj = chpo ; lobj = VRAI ;
  155. 'FINSI' ;
  156. 'FINSI' ;
  157. 'SI' ('NON' lobj) ;
  158. 'ARGUMENT' rig/'RIGIDITE' ;
  159. 'SI' ('EXISTE' rig) ;
  160. obj = rig ; lobj = VRAI ;
  161. 'FINSI' ;
  162. 'FINSI' ;
  163. 'SI' ('NON' lobj) ;
  164. 'ARGUMENT' matk/'MATRIK' ;
  165. 'SI' ('EXISTE' matk) ;
  166. obj = matk ; lobj = VRAI ;
  167. 'FINSI' ;
  168. 'FINSI' ;
  169. 'SI' ('NON' lobj) ;
  170. cherr = 'CHAINE'
  171. 'Il faut fournir un objet liste, evolution, matrice ou chpoint.'
  172. ;
  173. 'ERREUR' cherr ;
  174. 'FINSI' ;
  175. 'SI' ('EXISTE' tab itab) ;
  176. 'SI' ('EGA' ('TYPE' obj) 'CHPOINT') ;
  177. tab . itab = '+' (tab . itab) obj ;
  178. 'SINON' ;
  179. tab . itab = 'ET' (tab . itab) obj ;
  180. 'FINSI' ;
  181. 'SINON' ;
  182. tab . itab = obj ;
  183. 'FINSI' ;
  184. 'RESPRO' tab ;
  185. 'FINSI' ;
  186. 'ARGUMENT' lenti/'LISTENTI' ;
  187. 'ARGUMENT' lreel/'LISTREEL' ;
  188. 'SI' ('EXISTE' lenti) ;
  189. 'ARGUMENT' enti*'ENTIER' ;
  190. lenti = 'ET' lenti ('LECT' enti) ;
  191. 'RESPRO' lenti ;
  192. 'FINSI' ;
  193. 'SI' ('EXISTE' lreel) ;
  194. 'ARGUMENT' reel*'FLOTTANT' ;
  195. lreel = 'ET' lreel ('PROG' reel) ;
  196. 'RESPRO' lreel ;
  197. 'FINSI' ;
  198. *
  199. * End of procedure file APPEND
  200. *
  201. 'FINPROC' ;
  202. *ENDPROCEDUR append
  203. *BEGINPROCEDUR errrel
  204. ************************************************************************
  205. * NOM : ERRREL
  206. * DESCRIPTION : Calcul d'une erreur relative
  207. *
  208. *
  209. *
  210. * LANGAGE : GIBIANE-CAST3M
  211. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  212. * mél : gounand@semt2.smts.cea.fr
  213. **********************************************************************
  214. * VERSION : v1, 23/04/2003, version initiale
  215. * HISTORIQUE : v1, 23/04/2003, création
  216. * HISTORIQUE :
  217. * HISTORIQUE :
  218. ************************************************************************
  219. * Prière de PRENDRE LE TEMPS de compléter les commentaires
  220. * en cas de modification de ce sous-programme afin de faciliter
  221. * la maintenance !
  222. ************************************************************************
  223. *
  224. *
  225. 'DEBPROC' ERRREL ;
  226. 'ARGUMENT' val*'FLOTTANT' ;
  227. 'ARGUMENT' valref*'FLOTTANT' ;
  228. *
  229. 'SI' ('<' ('ABS' valref) 1.D-10) ;
  230. echref = 1.D0 ;
  231. 'SINON' ;
  232. echref = valref ;
  233. 'FINSI' ;
  234. *
  235. errabs = 'ABS' ('/' ('-' val valref) echref);
  236. *
  237. 'RESPRO' errabs ;
  238. *
  239. * End of procedure file ERRREL
  240. *
  241. 'FINPROC' ;
  242. *ENDPROCEDUR errrel
  243. *BEGINPROCEDUR exmomod
  244. ************************************************************************
  245. * NOM : EXMOMOD
  246. * DESCRIPTION : Extraction d'un mot d'un listmots
  247. *
  248. *
  249. *
  250. * LANGAGE : GIBIANE-CAST3M
  251. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  252. * mél : gounand@semt2.smts.cea.fr
  253. **********************************************************************
  254. * VERSION : v1, 23/06/2003, version initiale
  255. * HISTORIQUE : v1, 23/06/2003, création
  256. * HISTORIQUE :
  257. * HISTORIQUE :
  258. ************************************************************************
  259. * Prière de PRENDRE LE TEMPS de compléter les commentaires
  260. * en cas de modification de ce sous-programme afin de faciliter
  261. * la maintenance !
  262. ************************************************************************
  263. *
  264. *
  265. 'DEBPROC' EXMOMOD ;
  266. 'ARGUMENT' lm*'LISTMOTS' i*'ENTIER' ;
  267. j = 'DIME' lm ;
  268. k = '+' (MODULO ('-' i 1) j) 1 ;
  269. lemot = 'EXTRAIRE' lm k ;
  270. * Usage de l'opérateur text pour éviter que lemot
  271. * ne soit interprété comme un opérateur
  272. 'RESPRO' 'TEXTE' lemot ;
  273. *
  274. * End of procedure file EXMOMOD
  275. *
  276. 'FINPROC' ;
  277. *ENDPROCEDUR exmomod
  278. *BEGINPROCEDUR formar
  279. ************************************************************************
  280. * NOM : FORMAR
  281. * DESCRIPTION : formate un réel de facon courte
  282. * pratique pour les noms de
  283. * sauvegarde
  284. * Exemples :
  285. * 'MESSAGE' ('CHAINE' (formar 2.9e5 1)) ;
  286. * 2.9E5
  287. * 'MESSAGE' ('CHAINE' (formar -2.9e5 1)) ;
  288. * -2.9E5
  289. * 'MESSAGE' ('CHAINE' (formar 2.9e-5 1)) ;
  290. * 2.9E-5
  291. * 'MESSAGE' ('CHAINE' (formar -2.9e-5 1)) ;
  292. * -2.9E-5
  293. * 'MESSAGE' ('CHAINE' (formar 2.9 1)) ;
  294. * 2.9
  295. * 'MESSAGE' ('CHAINE' (formar -2.9 1)) ;
  296. * -2.9
  297. * 'MESSAGE' ('CHAINE' (formar 0 1)) ;
  298. * 0
  299. * 'MESSAGE' ('CHAINE' (formar 0 1)) ;
  300. * 0
  301. * 'MESSAGE' ('CHAINE' (formar 2.9e5 0)) ;
  302. * 3E5
  303. * 'MESSAGE' ('CHAINE' (formar -2.9e5 0)) ;
  304. * -3E5
  305. * 'MESSAGE' ('CHAINE' (formar 2.9e-5 0)) ;
  306. * 3E-5
  307. * 'MESSAGE' ('CHAINE' (formar -2.9e-5 0)) ;
  308. * -3E-5
  309. * 'MESSAGE' ('CHAINE' (formar 2.9 0)) ;
  310. * 3
  311. * 'MESSAGE' ('CHAINE' (formar -2.9 0)) ;
  312. * -3
  313. * 'MESSAGE' ('CHAINE' (formar 0 0)) ;
  314. * 0
  315. * 'MESSAGE' ('CHAINE' (formar 0 0)) ;
  316. * 0
  317. *
  318. *
  319. *
  320. * LANGAGE : GIBIANE-CAST3M
  321. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  322. * mél : gounand@semt2.smts.cea.fr
  323. **********************************************************************
  324. * VERSION : v1, 18/02/2003, version initiale
  325. * HISTORIQUE : v1, 18/02/2003, création
  326. * HISTORIQUE :
  327. * HISTORIQUE :
  328. ************************************************************************
  329. * Prière de PRENDRE LE TEMPS de compléter les commentaires
  330. * en cas de modification de ce sous-programme afin de faciliter
  331. * la maintenance !
  332. ************************************************************************
  333. *
  334. *
  335. 'DEBPROC' FORMAR ;
  336. 'ARGUMENT' fl*'FLOTTANT' ;
  337. 'ARGUMENT' vir/'ENTIER ' ;
  338. 'SI' ('NON' ('EXISTE' vir)) ;
  339. vir = 1 ;
  340. 'SINON' ;
  341. 'SI' ('<' vir 0) ;
  342. 'ERREUR' 'fournir un entier positif' ;
  343. 'FINSI' ;
  344. 'FINSI' ;
  345. 'SI' ('<' ('ABS' fl) 10.D-100) ;
  346. chfl = 'CHAINE' '0' ;
  347. 'SINON' ;
  348. *! sans le 1.D-10, ca ne fonctionne pas
  349. *! qd on entre pile poil une puissance de 10
  350. lfl = LOG10 ('ABS' fl) ;
  351. * lfl = '+' (LOG10 ('ABS' fl)) 1.D-10 ;
  352. slfl = 'SIGNE' ('ENTIER' lfl) ;
  353. 'SI' ('EGA' slfl 1) ;
  354. elfl = 'ENTIER' lfl ;
  355. 'SINON' ;
  356. elfl = '-' ('ENTIER' lfl) 1 ;
  357. 'FINSI' ;
  358. man = '/' fl ('**' 10.D0 elfl) ;
  359. *
  360. * Une verrue pour des histoires de précision...
  361. *
  362. 'SI' ('EGA' man 10.D0 ('**' 10.D0 ('*' vir -1.D0))) ;
  363. man = '/' man 10.D0 ;
  364. elfl = '+' elfl 1 ;
  365. 'FINSI' ;
  366. *
  367. sman = 'SIGNE' man ;
  368. 'SI' ('EGA' sman 1) ;
  369. fman = 'CHAINE' '(F' ('+' vir 2) '.0' vir ')' ;
  370. 'SINON' ;
  371. fman = 'CHAINE' '(F' ('+' vir 3) '.0' vir ')' ;
  372. 'FINSI' ;
  373. 'SI' ('NEG' vir 0) ;
  374. 'SI' ('NEG' elfl 0) ;
  375. chfl = 'CHAINE' 'FORMAT' fman man 'E' elfl ;
  376. 'SINON' ;
  377. chfl = 'CHAINE' 'FORMAT' fman man ;
  378. 'FINSI' ;
  379. 'SINON' ;
  380. man2 = 'ENTIER' ('+' man ('*' 0.5D0 sman)) ;
  381. 'SI' ('NEG' elfl 0) ;
  382. chfl = 'CHAINE' man2 'E' elfl ;
  383. 'SINON' ;
  384. chfl = 'CHAINE' man2 ;
  385. 'FINSI' ;
  386. 'FINSI' ;
  387. 'FINSI' ;
  388. 'RESPRO' chfl ;
  389. *
  390. * End of procedure file FORMAR
  391. *
  392. 'FINPROC' ;
  393. *ENDPROCEDUR formar
  394. *BEGINPROCEDUR getcoo
  395. ************************************************************************
  396. * NOM : GETCOO
  397. * DESCRIPTION :
  398. * Renvoie les coordonnées des points dans un champ type déplacement
  399. *
  400. *
  401. *
  402. * LANGAGE : GIBIANE-CAST3M
  403. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  404. * mél : gounand@semt2.smts.cea.fr
  405. **********************************************************************
  406. * VERSION : v1, 22/04/2011, version initiale
  407. * HISTORIQUE : v1, 22/04/2011, création
  408. * HISTORIQUE :
  409. * HISTORIQUE :
  410. ************************************************************************
  411. *
  412. *
  413. 'DEBPROC' GETCOO ;
  414. 'ARGUMENT' mail*'MAILLAGE' ;
  415. 'ARGUMENT' incop*'LISTMOTS' ;
  416. *
  417. dim = 'VALEUR' 'DIME' ;
  418. 'REPETER' iidim dim ;
  419. idim= &iidim ;
  420. icoo = 'NOMC' ('EXTRAIRE' incop idim)
  421. ('COORDONNEE' idim mail) ;
  422. 'SI' ('EGA' idim 1) ;
  423. vcoo = icoo ;
  424. 'SINON' ;
  425. vcoo = 'ET' vcoo icoo ;
  426. 'FINSI' ;
  427. 'FIN' iidim ;
  428. 'RESPRO' vcoo ;
  429. *
  430. * End of procedure file GETCOO
  431. *
  432. 'FINPROC' ;
  433. *ENDPROCEDUR getcoo
  434. *BEGINPROCEDUR ggravi
  435. ************************************************************************
  436. * NOM : GGRAVI
  437. * DESCRIPTION : Calcul de la force associée au potentiel gravitaire
  438. * (\rho g z si g vertical)
  439. *
  440. *
  441. *
  442. * LANGAGE : GIBIANE-CAST3M
  443. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  444. * mél : gounand@semt2.smts.cea.fr
  445. **********************************************************************
  446. * VERSION : v1, 22/04/2011
  447. * HISTORIQUE : v1, 22/04/2011, création
  448. * HISTORIQUE :
  449. * HISTORIQUE :
  450. ************************************************************************
  451. *
  452. *
  453. 'DEBPROC' GGRAVI ;
  454. 'ARGUMENT' _surf*'MAILLAGE' ;
  455. 'ARGUMENT' tdisc*'TABLE' ;
  456. 'ARGUMENT' coef*'FLOTTANT' ;
  457. 'ARGUMENT' ang*'FLOTTANT' ;
  458. *
  459. vdim = 'VALEUR' 'DIME' ;
  460. pgrax = '*' ('COORDONNEE' 1 _surf) ('*' +1. ('SIN' ang)) ;
  461. pgraz = '*' ('COORDONNEE' vdim _surf) ('*' -1. ('COS' ang)) ;
  462. DISCG = TDISC . 'GEOM' . 'DISC' ;
  463. fpgrax = GNOR _surf tdisc 'NPRI' discg 'CPRI' pgrax 'NDUA' 'XN' ;
  464. fpgraz = GNOR _surf tdisc 'NPRI' discg 'CPRI' pgraz 'NDUA' 'XN' ;
  465. fpgra = '+' fpgrax fpgraz ;
  466. fpgra = '*' fpgra ('*' -1. coef) ;
  467. 'RESPRO' fpgra ;
  468. *
  469. * End of procedure file GGRAVI
  470. *
  471. 'FINPROC' ;
  472. *ENDPROCEDUR ggravi
  473. *BEGINPROCEDUR gkgravi
  474. ************************************************************************
  475. * NOM : GKGRAVI
  476. * DESCRIPTION : Calcul de la matrice tangente de la force
  477. * associée au potentiel gravitaire (calculée par GGRAVI)
  478. * en fonction des déplacements des points de la surface.
  479. *
  480. *
  481. *
  482. * LANGAGE : GIBIANE-CAST3M
  483. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  484. * mél : gounand@semt2.smts.cea.fr
  485. **********************************************************************
  486. * VERSION : v1, 22/04/2011
  487. * HISTORIQUE : v1, 22/04/2011, création
  488. * HISTORIQUE :
  489. * HISTORIQUE :
  490. ************************************************************************
  491. *
  492. *
  493. 'DEBPROC' GKGRAVI ;
  494. 'ARGUMENT' _surf*'MAILLAGE' ;
  495. 'ARGUMENT' tdisc*'TABLE' ;
  496. 'ARGUMENT' ijaco*'ENTIER' ;
  497. *'SI' ('NON' ('EXISTE' ijaco)) ;
  498. * ijaco = 0 ;
  499. *'FINSI' ;
  500. 'ARGUMENT' coef*'FLOTTANT' ;
  501. 'ARGUMENT' ang*'FLOTTANT' ;
  502. *
  503. vdim = 'VALEUR' 'DIME' ;
  504. pgrax = '*' ('COORDONNEE' 1 _surf) ('*' +1. ('SIN' ang)) ;
  505. pgraz = '*' ('COORDONNEE' vdim _surf) ('*' -1. ('COS' ang)) ;
  506. *pgra = '*' ('-' ('COORDONNEE' vdim _surf) H) -1. ;
  507. *pgra = '*' ('COORDONNEE' vdim _surf) -1. ;
  508. DISCG = TDISC . 'GEOM' . 'DISC' ;
  509. NOMDEP = @STBL (TDISC . 'XN' . 'NOMINC') ;
  510. *fpgra = GNOR _surf tdisc 'NPRI' discg 'CPRI' pgra 'NDUA' 'XN' ;
  511. k1x = GNOR _surf tdisc 'NPRI' discg 'NDUA' 'XN' ;
  512. k1x = '*' k1x ('*' +1. ('SIN' ang)) ;
  513. k1x = 'CHANGER' 'INCO' k1x ('MOTS' 'SCAL')
  514. ('MOTS' ('EXTRAIRE' NOMDEP 1)) NOMDEP NOMDEP ;
  515. k2x = GNORKTAN _surf tdisc 'NPRI' 'XN'
  516. 'NCOF' discg 'CCOF' pgrax 'NDUA' 'XN' ;
  517. k1z = GNOR _surf tdisc 'NPRI' discg 'NDUA' 'XN' ;
  518. k1z = '*' k1z ('*' -1. ('COS' ang)) ;
  519. k1z = 'CHANGER' 'INCO' k1z ('MOTS' 'SCAL')
  520. ('MOTS' ('EXTRAIRE' NOMDEP vdim)) NOMDEP NOMDEP ;
  521. k2z = GNORKTAN _surf tdisc 'NPRI' 'XN'
  522. 'NCOF' discg 'CCOF' pgraz 'NDUA' 'XN' ;
  523. 'SI' ('EGA' ijaco 0) ;
  524. ktgra = k1x 'ET' k1z 'ET' k2x 'ET' k2z ;
  525. 'FINSI' ;
  526. 'SI' ('EGA' ijaco 1) ;
  527. ktgra = k1x 'ET' k1z ;
  528. 'FINSI' ;
  529. 'SI' ('EGA' ijaco 2) ;
  530. ktgra = k2x 'ET' k2z ;
  531. 'FINSI' ;
  532. ktgra = '*' ktgra coef ;
  533. 'RESPRO' ktgra ;
  534. *
  535. * End of procedure file GKGRAVI
  536. *
  537. 'FINPROC' ;
  538. *ENDPROCEDUR gkgravi
  539. *BEGINPROCEDUR gkvol
  540. ************************************************************************
  541. * NOM : GKVOL
  542. * DESCRIPTION : Matrice tangente associée à la variation du volume
  543. * contenu dans une surface (calculé par GVOL)
  544. * en fonction des déplacements des points de la surface.
  545. *
  546. *
  547. * LANGAGE : GIBIANE-CAST3M
  548. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  549. * mél : gounand@semt2.smts.cea.fr
  550. **********************************************************************
  551. * VERSION : v1, 22/04/2011, version initiale
  552. * HISTORIQUE : v1, 22/04/2011, création
  553. * HISTORIQUE :
  554. * HISTORIQUE :
  555. ************************************************************************
  556. *
  557. *
  558. 'DEBPROC' GKVOL ;
  559. 'ARGUMENT' _surf*'MAILLAGE' ;
  560. 'ARGUMENT' tdisc*'TABLE' ;
  561. 'ARGUMENT' ijaco/'ENTIER' ;
  562. 'SI' ('NON' ('EXISTE' ijaco)) ;
  563. ijaco = 0 ;
  564. 'FINSI' ;
  565. * Vecteur position et calcul du volume
  566. NOMVIT = @STBL (TDISC . 'XN' . 'NOMINC') ;
  567. DISCG = TDISC . 'GEOM' . 'DISC' ;
  568. vdim = 'VALEUR' 'DIME' ;
  569. 'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'AXIS')) ;
  570. fdim = 3 ;
  571. 'SINON' ;
  572. fdim = vdim ;
  573. 'FINSI' ;
  574. vpos = GETCOO _surf nomvit ;
  575. kvol1 = GNOR _surf tdisc 'NPRI' ('CHAINE' discg 'V')
  576. 'NDUA' 'XN' 'FDUA' ('PROG' vdim * 1.) ;
  577. kvol2 = GNORKTAN _surf tdisc 'NPRI' ('CHAINE' discg 'V')
  578. 'NCOF' ('CHAINE' discg 'V') 'CCOF' vpos
  579. 'NDUA' 'XN' 'FDUA' ('PROG' vdim * 1.) ;
  580. 'SI' ('EGA' ijaco 0) ;
  581. kvol = '/' ('+' kvol1 kvol2) fdim ;
  582. 'FINSI' ;
  583. 'SI' ('EGA' ijaco 1) ;
  584. kvol = '/' kvol1 fdim ;
  585. 'FINSI' ;
  586. 'SI' ('EGA' ijaco 2) ;
  587. kvol = '/' kvol2 fdim ;
  588. 'FINSI' ;
  589. 'RESPRO' kvol ;
  590. *
  591. * End of procedure file GKVOL
  592. *
  593. 'FINPROC' ;
  594. *ENDPROCEDUR gkvol
  595. *BEGINPROCEDUR gmass2
  596. ************************************************************************
  597. * NOM : GMASS2
  598. * DESCRIPTION : Une matrice de masse
  599. *
  600. *
  601. *
  602. * LANGAGE : GIBIANE-CAST3M
  603. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  604. * mél : gounand@semt2.smts.cea.fr
  605. **********************************************************************
  606. * VERSION : v2, 14/03/2006, mise à jour NLIN évolué
  607. * VERSION : v1, 13/05/2004, version initiale
  608. * HISTORIQUE : v1, 13/05/2004, création
  609. * HISTORIQUE :
  610. * HISTORIQUE :
  611. ************************************************************************
  612. * Prière de PRENDRE LE TEMPS de compléter les commentaires
  613. * en cas de modification de ce sous-programme afin de faciliter
  614. * la maintenance !
  615. ************************************************************************
  616. *
  617. *
  618. 'DEBPROC' GMASS2 ;
  619. 'ARGUMENT' _mt*'MAILLAGE' ;
  620. 'ARGUMENT' _smt/'MAILLAGE' ;
  621. 'ARGUMENT' tdisc*'TABLE' ;
  622. *
  623. * Lectures
  624. *
  625. debug = FAUX ;
  626. lmotcle = 'MOTS' 'NPRI' 'FPRI' 'CPRI' 'NDUA' 'FDUA' 'CDUA'
  627. 'NCOF' 'FCOF' 'CCOF' ;
  628. * Il faut initialiser valt et valq, sinon on peut capturer ceux de
  629. * la procédure appelante
  630. valt = 'valt' ; valq = 'valq' ;
  631. 'REPETER' imotcle ;
  632. 'ARGUMENT' motcle/'MOT' ;
  633. 'SI' ('NON' ('EXISTE' motcle)) ; 'QUITTER' imotcle ; 'FINSI' ;
  634. 'SI' ('NON' ('EXISTE' lmotcle motcle)) ;
  635. cherr = 'CHAINE' 'Keyword ' motcle ' unknown.' ; 'ERREUR' cherr ;
  636. 'FINSI' ;
  637. 'SI' ('EGA' motcle 'NPRI') ; 'ARGUMENT' nomt*'MOT' ; 'FINSI' ;
  638. 'SI' ('EGA' motcle 'NDUA') ; 'ARGUMENT' nomq*'MOT' ; 'FINSI' ;
  639. 'SI' ('EGA' motcle 'NCOF') ; 'ARGUMENT' nomo*'MOT' ; 'FINSI' ;
  640. tst1 = 'EGA' motcle 'FPRI' ; tst2 = 'EGA' motcle 'FDUA' ;
  641. tst = tst1 'OU' tst2 ;
  642. 'SI' tst ;
  643. 'SI' tst1 ; tt = TDISC . nomt ; 'FINSI' ;
  644. 'SI' tst2 ; tt = TDISC . nomq ; 'FINSI' ;
  645. isvec = ('>' ('DIME' (tt. 'NOMINC')) 1) ;
  646. 'SI' isvec ; 'ARGUMENT' valv*'LISTREEL' ; 'SINON' ;
  647. 'ARGUMENT' valv*'FLOTTANT' ;
  648. 'FINSI' ;
  649. 'SI' tst1 ; valt = valv ; 'FINSI' ;
  650. 'SI' tst2 ; valq = valv ; 'FINSI' ;
  651. 'FINSI' ;
  652. 'SI' ('EGA' motcle 'FCOF') ; 'ARGUMENT' valo*'FLOTTANT' ; 'FINSI' ;
  653. 'SI' ('EGA' motcle 'CPRI') ; 'ARGUMENT' valt*'CHPOINT' ; 'FINSI' ;
  654. 'SI' ('EGA' motcle 'CDUA') ; 'ARGUMENT' valq*'CHPOINT' ; 'FINSI' ;
  655. 'SI' ('EGA' motcle 'CCOF') ; 'ARGUMENT' valo*'CHPOINT' ; 'FINSI' ;
  656. 'FIN' imotcle ;
  657. *
  658. * Tests
  659. *
  660. discg = TDISC . 'GEOM' . 'DISC' ;
  661. 'SI' ('EXISTE' tdisc 'methgau') ;
  662. methgau = tdisc . 'methgau' . 'mass' ;
  663. 'SINON' ;
  664. methgau = 'GAU7' ;
  665. 'FINSI' ;
  666. tnomt = TDISC . nomt ;
  667. lvalt = 'NEG' ('TYPE' valt) 'MOT' ;
  668. tnomq = TDISC . nomq ;
  669. lvalq = 'NEG' ('TYPE' valq) 'MOT' ;
  670. * Scalaire ou vecteur
  671. ninct = 'DIME' (tnomt . 'NOMINC') ;
  672. nincq = 'DIME' (tnomq . 'NOMINC') ;
  673. 'SI' ('NEG' ninct nincq) ;
  674. cherr = 'CHAINE'
  675. 'les primales et duales nont pas le meme nombre de composantes' ;
  676. 'ERREUR' cherr ;
  677. 'FINSI' ;
  678. ninc = ninct ;
  679. *
  680. lcof = 'EXISTE' TDISC nomo ;
  681. 'SI' lcof ; ncof = 1 ; tcof = TDISC . nomo ;
  682. 'SINON' ; ncof = 0 ;
  683. 'FINSI' ;
  684. *
  685. 'SI' debug ;
  686. 'SI' lcof ; 'MESSAGE' 'Un coef a ete detecte' ;
  687. 'SINON' ; 'MESSAGE' 'pas de coef detecte' ;
  688. 'FINSI' ;
  689. 'FINSI' ;
  690. *
  691. vdim = 'VALEUR' 'DIME' ;
  692. vmod = 'VALEUR' 'MODE' ;
  693. idim = 0 ;
  694. 'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'PLANDEFO')) ;
  695. idim = 2 ;
  696. iaxi = FAUX ;
  697. 'FINSI' ;
  698. 'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'AXIS')) ;
  699. idim = 2 ;
  700. iaxi = VRAI ;
  701. 'FINSI' ;
  702. 'SI' ('ET' ('EGA' vdim 3) ('EGA' vmod 'TRID')) ;
  703. idim = 3 ;
  704. iaxi = FAUX ;
  705. 'FINSI' ;
  706. 'SI' ('EGA' vdim 1) ;
  707. idim = 1 ;
  708. iaxi = FAUX ;
  709. 'FINSI' ;
  710. * 'MESSAGE' ('CHAINE' 'iaxi=' iaxi );
  711. 'SI' ('EGA' idim 0) ;
  712. 'ERREUR' ('CHAINE' 'vmod=' vmod ' et vdim=' vdim ' non prevu') ;
  713. 'FINSI' ;
  714. 'SI' iaxi ;
  715. dprmt = '*' ('COORDONNEE' 1 _mt) ('*' PI 2.D0) ;
  716. 'FINSI' ;
  717. *
  718. * Optimisation possible : construire la matrice par blocs
  719. * qd valt et valq ne sont pas donnés
  720. *
  721. numop = ninc ; numder = idim ; numvar = ninc ;
  722. numdat = ncof ; numcof = ncof ;
  723. A = ININLIN numop numvar numdat numcof numder ;
  724. 'SI' lcof ;
  725. A . 'DAT' . 1 . 'NOMDDL' = tcof . 'NOMINC' . 1 ;
  726. A . 'DAT' . 1 . 'DISC' = tcof . 'DISC' ;
  727. A . 'DAT' . 1 . 'VALEUR' = valo ;
  728. A . 'COF' . 1 . 'COMPOR' = 'IDEN' ;
  729. A . 'COF' . 1 . 'LDAT' = 'LECT' 1 ;
  730. ll = 'LECT' 1 ;
  731. 'SINON' ;
  732. ll = 'LECT' ;
  733. 'FINSI' ;
  734. lvt = 'EGA' ('TYPE' valt) 'LISTREEL' ;
  735. 'REPETER' iiinc ninc ;
  736. iinc = &iiinc ;
  737. A . 'VAR' . iinc . 'NOMDDL' = tnomt . 'NOMINC' . iinc ;
  738. A . 'VAR' . iinc . 'DISC' = tnomt . 'DISC' ;
  739. 'SI' lvalt ;
  740. 'SI' lvt ;
  741. A . 'VAR' . iinc . 'VALEUR' = 'EXTRAIRE' valt iinc ;
  742. 'SINON' ;
  743. A . 'VAR' . iinc . 'VALEUR' = valt ;
  744. 'FINSI' ;
  745. 'FINSI' ;
  746. A . iinc . iinc . 0 = ll ;
  747. 'FIN' iiinc ;
  748. *
  749. 'SI' iaxi ;
  750. numdat = 1 ;
  751. numcof = 1 ;
  752. 'SINON' ;
  753. numdat = 0 ;
  754. numcof = 0 ;
  755. 'FINSI' ;
  756. B = ININLIN numop numvar numdat numcof numder ;
  757. 'SI' iaxi ;
  758. B . 'DAT' . 1 . 'NOMDDL' = 'MOTS' 'SCAL' ;
  759. B . 'DAT' . 1 . 'DISC' = discg ;
  760. B . 'DAT' . 1 . 'VALEUR' = dprmt ;
  761. B . 'COF' . 1 . 'COMPOR' = 'IDEN' ;
  762. B . 'COF' . 1 . 'LDAT' = 'LECT' 1 ;
  763. ll = 'LECT' 1 ;
  764. 'SINON' ;
  765. ll = 'LECT' ;
  766. 'FINSI' ;
  767. lvq = 'EGA' ('TYPE' valq) 'LISTREEL' ;
  768. 'REPETER' iiinc ninc ;
  769. iinc = &iiinc ;
  770. B . 'VAR' . iinc . 'NOMDDL' = tnomq . 'NOMINC' . iinc ;
  771. B . 'VAR' . iinc . 'DISC' = tnomq . 'DISC' ;
  772. 'SI' lvalq ;
  773. 'SI' lvq ;
  774. B . 'VAR' . iinc . 'VALEUR' = 'EXTRAIRE' valq iinc ;
  775. 'SINON' ;
  776. B . 'VAR' . iinc . 'VALEUR' = valq ;
  777. 'FINSI' ;
  778. 'FINSI' ;
  779. B . iinc . iinc . 0 = ll ;
  780. 'FIN' iiinc ;
  781. *
  782. 'SI' ('EXISTE' _smt) ;
  783. mgmass2 = 'NLIN' discg _mt _smt A B methgau ;
  784. 'SINON' ;
  785. mgmass2 = NLINP discg _mt A B methgau ;
  786. 'FINSI' ;
  787. *
  788. 'RESPRO' mgmass2 ;
  789. 'FINPROC' ;
  790. *
  791. * End of procedure file GMASS2
  792. *
  793. *ENDPROCEDUR gmass2
  794. *BEGINPROCEDUR gnorktan
  795. ************************************************************************
  796. * NOM : GNORKTAN
  797. * DESCRIPTION : Matrice tangente associée à la variation de la normale
  798. * à une surface (calculée par GNOR)
  799. * en fonction des déplacements des points de la surface.
  800. *
  801. *
  802. *
  803. * LANGAGE : GIBIANE-CAST3M
  804. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  805. * mél : gounand@semt2.smts.cea.fr
  806. **********************************************************************
  807. * VERSION : v1, 22/04/2011, version initiale
  808. * HISTORIQUE : v1, 22/04/2011, création
  809. * HISTORIQUE :
  810. * HISTORIQUE :
  811. ************************************************************************
  812. *
  813. *
  814. 'DEBPROC' GNORKTAN ;
  815. 'ARGUMENT' _mt*'MAILLAGE' ;
  816. 'ARGUMENT' tdisc*'TABLE' ;
  817. *
  818. * Lectures
  819. *
  820. dim = 'VALEUR' 'DIME' ;
  821. mdim = DEADUTIL 'DIMM' _mt ;
  822. 'SI' ('NEG' mdim ('-' dim 1)) ;
  823. 'ERREUR' 'Dim. maillage .neq. dim. espace - 1' ;
  824. 'FINSI' ;
  825. loi = 'CHAINE' 'VNOJ' ;
  826. debug = FAUX ;
  827. lmotcle = 'MOTS' 'NPRI' 'FPRI' 'CPRI' 'NDUA' 'FDUA' 'CDUA'
  828. 'NCOF' 'FCOF' 'CCOF' ;
  829. * Il faut initialiser valt et valq, sinon on peut capturer ceux de
  830. * la procédure appelante
  831. valt = 'valt' ; valq = 'valq' ;
  832. 'REPETER' imotcle ;
  833. 'ARGUMENT' motcle/'MOT' ;
  834. 'SI' ('NON' ('EXISTE' motcle)) ; 'QUITTER' imotcle ; 'FINSI' ;
  835. 'SI' ('NON' ('EXISTE' lmotcle motcle)) ;
  836. cherr = 'CHAINE' 'Keyword ' motcle ' unknown.' ; 'ERREUR' cherr ;
  837. 'FINSI' ;
  838. 'SI' ('EGA' motcle 'NPRI') ; 'ARGUMENT' nomt*'MOT' ; 'FINSI' ;
  839. 'SI' ('EGA' motcle 'NDUA') ; 'ARGUMENT' nomq*'MOT' ; 'FINSI' ;
  840. 'SI' ('EGA' motcle 'NCOF') ; 'ARGUMENT' nomo*'MOT' ; 'FINSI' ;
  841. tst1 = 'EGA' motcle 'FPRI' ; tst2 = 'EGA' motcle 'FDUA' ;
  842. tst = tst1 'OU' tst2 ;
  843. 'SI' tst ;
  844. 'SI' tst1 ; tt = TDISC . nomt ; 'FINSI' ;
  845. 'SI' tst2 ; tt = TDISC . nomq ; 'FINSI' ;
  846. isvec = ('>' ('DIME' (tt. 'NOMINC')) 1) ;
  847. 'SI' isvec ; 'ARGUMENT' valv*'LISTREEL' ; 'SINON' ;
  848. 'ARGUMENT' valv*'FLOTTANT' ;
  849. 'FINSI' ;
  850. 'SI' tst1 ; valt = valv ; 'FINSI' ;
  851. 'SI' tst2 ; valq = valv ; 'FINSI' ;
  852. 'FINSI' ;
  853. 'SI' ('EGA' motcle 'FCOF') ; 'ARGUMENT' valo*'FLOTTANT' ; 'FINSI' ;
  854. 'SI' ('EGA' motcle 'CPRI') ; 'ARGUMENT' valt*'CHPOINT' ; 'FINSI' ;
  855. 'SI' ('EGA' motcle 'CDUA') ; 'ARGUMENT' valq*'CHPOINT' ; 'FINSI' ;
  856. 'SI' ('EGA' motcle 'CCOF') ; 'ARGUMENT' valo*'CHPOINT' ; 'FINSI' ;
  857. 'FIN' imotcle ;
  858. *
  859. * Tests
  860. *
  861. discg = TDISC . 'GEOM' . 'DISC' ;
  862. 'SI' ('EXISTE' tdisc 'methgau') ;
  863. methgau = tdisc . 'methgau' . 'mass' ;
  864. 'SINON' ;
  865. methgau = 'GAU7' ;
  866. 'FINSI' ;
  867. tnomt = TDISC . nomt ;
  868. lvalt = 'NEG' ('TYPE' valt) 'MOT' ;
  869. tnomq = TDISC . nomq ;
  870. lvalq = 'NEG' ('TYPE' valq) 'MOT' ;
  871. * Scalaire ou vecteur
  872. ninct = 'DIME' (tnomt . 'NOMINC') ;
  873. nincq = 'DIME' (tnomq . 'NOMINC') ;
  874. 'SI' ('NEG' ninct dim) ;
  875. cherr = 'CHAINE'
  876. 'la primale doit etre un vecteur' ;
  877. 'ERREUR' cherr ;
  878. 'FINSI' ;
  879. 'SI' ('NEG' nincq dim) ;
  880. cherr = 'CHAINE'
  881. 'la duale doit etre un vecteur' ;
  882. 'ERREUR' cherr ;
  883. 'FINSI' ;
  884. ninc = dim ;
  885. *
  886. lcof = 'EXISTE' TDISC nomo ;
  887. 'SI' lcof ; tcof = TDISC . nomo ;
  888. ncof = 'DIME' (tcof . 'NOMINC') ;
  889. 'SINON' ; ncof = 0 ;
  890. 'FINSI' ;
  891. *
  892. 'SI' debug ;
  893. 'SI' lcof ; 'MESSAGE' 'Un coef a ete detecte' ;
  894. 'SINON' ; 'MESSAGE' 'pas de coef detecte' ;
  895. 'FINSI' ;
  896. 'FINSI' ;
  897. *
  898. vdim = 'VALEUR' 'DIME' ;
  899. vmod = 'VALEUR' 'MODE' ;
  900. idim = 0 ;
  901. 'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'PLANDEFO')) ;
  902. idim = 2 ;
  903. iaxi = FAUX ;
  904. 'FINSI' ;
  905. 'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'AXIS')) ;
  906. idim = 2 ;
  907. iaxi = VRAI ;
  908. 'FINSI' ;
  909. 'SI' ('ET' ('EGA' vdim 3) ('EGA' vmod 'TRID')) ;
  910. idim = 3 ;
  911. iaxi = FAUX ;
  912. 'FINSI' ;
  913. 'SI' ('EGA' vdim 1) ;
  914. idim = 1 ;
  915. iaxi = FAUX ;
  916. 'FINSI' ;
  917. * 'MESSAGE' ('CHAINE' 'iaxi=' iaxi );
  918. 'SI' ('EGA' idim 0) ;
  919. 'ERREUR' ('CHAINE' 'vmod=' vmod ' et vdim=' vdim ' non prevu') ;
  920. 'FINSI' ;
  921. 'SI' iaxi ;
  922. deupi = '*' PI 2.D0 ;
  923. dprmt = '*' ('COORDONNEE' 1 _mt) deupi ;
  924. 'FINSI' ;
  925. *
  926. * Optimisation possible : construire la matrice par blocs
  927. * qd valt et valq ne sont pas donnés
  928. *
  929. numop = idim '*' idim '*' idim ;
  930. 'SI' iaxi ;
  931. numop = numop '+' idim ;
  932. 'FINSI' ;
  933. numder = idim ; numvar = ninct ;
  934. numdat = ncof ; numcof = ncof ;
  935. A = ININLIN numop numvar numdat numcof numder ;
  936. 'SI' lcof ;
  937. lvo = 'EGA' ('TYPE' valo) 'LISTREEL' ;
  938. 'REPETER' iicof ncof ;
  939. icof = &iicof ;
  940. A . 'DAT' . icof . 'NOMDDL' = tcof . 'NOMINC' . icof ;
  941. A . 'DAT' . icof . 'DISC' = tcof . 'DISC' ;
  942. 'SI' lvo ;
  943. A . 'DAT' . icof . 'VALEUR' = 'EXTRAIRE' valo icof ;
  944. 'SINON' ;
  945. A . 'DAT' . icof . 'VALEUR' = valo ;
  946. 'FINSI' ;
  947. A . 'COF' . icof . 'COMPOR' = 'IDEN' ;
  948. A . 'COF' . icof . 'LDAT' = 'LECT' icof ;
  949. 'FIN' iicof ;
  950. 'SINON' ;
  951. ll = 'LECT' ;
  952. 'FINSI' ;
  953. lvt = 'EGA' ('TYPE' valt) 'LISTREEL' ;
  954. iop = 0 ;
  955. 'REPETER' iiinct ninct ;
  956. iinct = &iiinct ;
  957. A . 'VAR' . iinct . 'NOMDDL' = tnomt . 'NOMINC' . iinct ;
  958. A . 'VAR' . iinct . 'DISC' = tnomt . 'DISC' ;
  959. 'SI' lvalt ;
  960. 'SI' lvt ;
  961. A . 'VAR' . iinct . 'VALEUR' = 'EXTRAIRE' valt iinct ;
  962. 'SINON' ;
  963. A . 'VAR' . iinct . 'VALEUR' = valt ;
  964. 'FINSI' ;
  965. 'FINSI' ;
  966. 'REPETER' iiincq nincq ;
  967. 'REPETER' iiider numder ;
  968. iop = '+' iop 1 ;
  969. 'SI' lcof ;
  970. icof = 'MINIMUM' ('LECT' &iiincq ncof) ;
  971. A . iop . iinct . &iiider = 'LECT' icof ;
  972. 'SINON' ;
  973. A . iop . iinct . &iiider = ll ;
  974. 'FINSI' ;
  975. 'FIN' iiider ;
  976. 'FIN' iiincq ;
  977. 'FIN' iiinct ;
  978. 'SI' iaxi ;
  979. 'REPETER' iiincq nincq ;
  980. iop = '+' iop 1 ;
  981. 'SI' lcof ;
  982. icof = 'MINIMUM' ('LECT' &iiincq ncof) ;
  983. A . iop . 1 . 0 = 'LECT' icof ;
  984. 'SINON' ;
  985. A . iop . 1 . 0 = ll ;
  986. 'FINSI' ;
  987. 'FIN' iiincq ;
  988. 'FINSI' ;
  989. *
  990. * 'SI' iaxi ;
  991. * numdat = 1 ;
  992. * numcof = dim '+' 1 ;
  993. * 'SINON' ;
  994. numdat = 0 ;
  995. numcof = idim '*' idim '*' idim ;
  996. * 'FINSI' ;
  997. 'SI' iaxi ;
  998. numdat = '+' numdat 2 ;
  999. numcof = '+' numcof ('+' idim 2) ;
  1000. 'FINSI' ;
  1001. numvar = nincq ;
  1002. B = ININLIN numop numvar numdat numcof numder ;
  1003. *
  1004. lvq = 'EGA' ('TYPE' valq) 'LISTREEL' ;
  1005. 'REPETER' iiinc nincq ;
  1006. iinc = &iiinc ;
  1007. B . 'VAR' . iinc . 'NOMDDL' = tnomq . 'NOMINC' . iinc ;
  1008. B . 'VAR' . iinc . 'DISC' = tnomq . 'DISC' ;
  1009. 'SI' lvalq ;
  1010. 'SI' lvq ;
  1011. B . 'VAR' . iinc . 'VALEUR' = 'EXTRAIRE' valq iinc ;
  1012. 'SINON' ;
  1013. B . 'VAR' . iinc . 'VALEUR' = valq ;
  1014. 'FINSI' ;
  1015. 'FINSI' ;
  1016. 'FIN' iiinc ;
  1017. idat = 0 ;
  1018. icof = 0 ;
  1019. 'SI' iaxi ;
  1020. 'REPETER' iiidim idim ;
  1021. icof = '+' icof 1 ;
  1022. B . 'COF' . icof . 'COMPOR' = 'CHAINE' 'VNOR' &iiidim ;
  1023. B . 'COF' . icof . 'LDAT' = 'LECT' ;
  1024. 'FIN' iiidim ;
  1025. idat = '+' idat 1 ;
  1026. icof = '+' icof 1 ;
  1027. B . 'DAT' . idat . 'NOMDDL' = 'MOTS' 'SCAL' ;
  1028. B . 'DAT' . idat . 'DISC' = discg ;
  1029. B . 'DAT' . idat . 'VALEUR' = dprmt ;
  1030. B . 'COF' . icof . 'COMPOR' = 'IDEN' ;
  1031. B . 'COF' . icof . 'LDAT' = 'LECT' idat ;
  1032. ll = 'LECT' icof ;
  1033. idat = '+' idat 1 ;
  1034. icof = '+' icof 1 ;
  1035. B . 'DAT' . idat . 'NOMDDL' = 'MOTS' 'SCAL' ;
  1036. B . 'DAT' . idat . 'DISC' = 'CSTE' ;
  1037. B . 'DAT' . idat . 'VALEUR' = deupi ;
  1038. B . 'COF' . icof . 'COMPOR' = 'IDEN' ;
  1039. B . 'COF' . icof . 'LDAT' = 'LECT' idat ;
  1040. ll2 = 'LECT' icof ;
  1041. 'SINON' ;
  1042. ll = 'LECT' ;
  1043. 'FINSI' ;
  1044. *
  1045. iop = 0 ;
  1046. 'REPETER' iiinct ninct ;
  1047. 'REPETER' iiincq nincq ;
  1048. 'REPETER' iiider numder ;
  1049. iop = '+' iop 1 ;
  1050. icof = '+' icof 1 ;
  1051. lcomp = 'CHAINE' loi &iiincq &iiinct &iiider ;
  1052. * lcomp = 'CHAINE' loi &iiinct &iiincq &iiider ;
  1053. B . 'COF' . icof . 'COMPOR' = lcomp ;
  1054. B . 'COF' . icof . 'LDAT' = 'LECT' ;
  1055. B . iop . &iiincq . 0 = ('LECT' icof) 'ET' ll ;
  1056. 'FIN' iiider ;
  1057. 'FIN' iiincq ;
  1058. 'FIN' iiinct ;
  1059. 'SI' iaxi ;
  1060. 'REPETER' iiincq nincq ;
  1061. iincq = &iiincq ;
  1062. iop = '+' iop 1 ;
  1063. B . iop . iincq . 0 = ('LECT' iincq) 'ET' ll2 ;
  1064. 'FIN' iiincq ;
  1065. 'FINSI' ;
  1066. *
  1067. * mgnorkt = NLIN discg _mt A B 'CRES' methgau ;
  1068. mgnorkt = NLIN discg _mt A B methgau ;
  1069. *
  1070. 'RESPRO' mgnorkt ;
  1071. 'FINPROC' ;
  1072. *
  1073. * End of procedure file GNORKTAN
  1074. *
  1075. *ENDPROCEDUR gnorktan
  1076. *BEGINPROCEDUR gnor
  1077. ************************************************************************
  1078. * NOM : GNOR
  1079. * DESCRIPTION : Calcule le champ de normales à une surface.
  1080. * Peut servir à calculer une pression, un potentiel
  1081. * lié à la gravité, un volume contenu dans une surface.
  1082. * Attention à l'orientation de la surface !
  1083. *
  1084. * Computes a field of normal to a surface.
  1085. * Also useful to compute a pressure field,
  1086. * a gravity potential field, a volume enclosed
  1087. * by a surface.
  1088. * WARNING : The orientation of the surface matters !
  1089. *
  1090. *
  1091. *
  1092. *
  1093. * LANGAGE : GIBIANE-CAST3M
  1094. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  1095. * mél : gounand@semt2.smts.cea.fr
  1096. **********************************************************************
  1097. * VERSION : v1, 22/04/2011
  1098. * HISTORIQUE : v1, 22/04/2011, création
  1099. * HISTORIQUE :
  1100. * HISTORIQUE :
  1101. ************************************************************************
  1102. *
  1103. *
  1104. 'DEBPROC' GNOR ;
  1105. 'ARGUMENT' _mt*'MAILLAGE' ;
  1106. 'ARGUMENT' tdisc*'TABLE' ;
  1107. *
  1108. * Lectures
  1109. *
  1110. dim = 'VALEUR' 'DIME' ;
  1111. debug = FAUX ;
  1112. lmotcle = 'MOTS' 'NPRI' 'FPRI' 'CPRI' 'NDUA' 'FDUA' 'CDUA'
  1113. 'NCOF' 'FCOF' 'CCOF' ;
  1114. * Il faut initialiser valt et valq, sinon on peut capturer ceux de
  1115. * la procédure appelante
  1116. valt = 'valt' ; valq = 'valq' ;
  1117. 'REPETER' imotcle ;
  1118. 'ARGUMENT' motcle/'MOT' ;
  1119. 'SI' ('NON' ('EXISTE' motcle)) ; 'QUITTER' imotcle ; 'FINSI' ;
  1120. 'SI' ('NON' ('EXISTE' lmotcle motcle)) ;
  1121. cherr = 'CHAINE' 'Keyword ' motcle ' unknown.' ; 'ERREUR' cherr ;
  1122. 'FINSI' ;
  1123. 'SI' ('EGA' motcle 'NPRI') ; 'ARGUMENT' nomt*'MOT' ; 'FINSI' ;
  1124. 'SI' ('EGA' motcle 'NDUA') ; 'ARGUMENT' nomq*'MOT' ; 'FINSI' ;
  1125. 'SI' ('EGA' motcle 'NCOF') ; 'ARGUMENT' nomo*'MOT' ; 'FINSI' ;
  1126. tst1 = 'EGA' motcle 'FPRI' ; tst2 = 'EGA' motcle 'FDUA' ;
  1127. tst3 = 'EGA' motcle 'FCOF' ;
  1128. tst = tst1 'OU' tst2 'OU' tst3 ;
  1129. 'SI' tst ;
  1130. 'SI' tst1 ; tt = TDISC . nomt ; 'FINSI' ;
  1131. 'SI' tst2 ; tt = TDISC . nomq ; 'FINSI' ;
  1132. 'SI' tst3 ; tt = TDISC . nomo ; 'FINSI' ;
  1133. isvec = ('>' ('DIME' (tt. 'NOMINC')) 1) ;
  1134. 'SI' isvec ; 'ARGUMENT' valv*'LISTREEL' ; 'SINON' ;
  1135. 'ARGUMENT' valv*'FLOTTANT' ;
  1136. 'FINSI' ;
  1137. 'SI' tst1 ; valt = valv ; 'FINSI' ;
  1138. 'SI' tst2 ; valq = valv ; 'FINSI' ;
  1139. 'SI' tst3 ; valo = valv ; 'FINSI' ;
  1140. 'FINSI' ;
  1141. 'SI' ('EGA' motcle 'CPRI') ; 'ARGUMENT' valt*'CHPOINT' ; 'FINSI' ;
  1142. 'SI' ('EGA' motcle 'CDUA') ; 'ARGUMENT' valq*'CHPOINT' ; 'FINSI' ;
  1143. 'SI' ('EGA' motcle 'CCOF') ; 'ARGUMENT' valo*'CHPOINT' ; 'FINSI' ;
  1144. 'FIN' imotcle ;
  1145. *
  1146. * Tests
  1147. *
  1148. discg = TDISC . 'GEOM' . 'DISC' ;
  1149. 'SI' ('EXISTE' tdisc 'methgau') ;
  1150. methgau = tdisc . 'methgau' . 'mass' ;
  1151. 'SINON' ;
  1152. methgau = 'GAU7' ;
  1153. 'FINSI' ;
  1154. tnomt = TDISC . nomt ;
  1155. lvalt = 'NEG' ('TYPE' valt) 'MOT' ;
  1156. tnomq = TDISC . nomq ;
  1157. lvalq = 'NEG' ('TYPE' valq) 'MOT' ;
  1158. * Scalaire ou vecteur
  1159. ninct = 'DIME' (tnomt . 'NOMINC') ;
  1160. nincq = 'DIME' (tnomq . 'NOMINC') ;
  1161. 'SI' ('ET' ('NEG' ninct 1) ('NEG' ninct dim)) ;
  1162. cherr = 'CHAINE'
  1163. 'la primale doit etre un scalaire ou un vecteur' ;
  1164. 'ERREUR' cherr ;
  1165. 'FINSI' ;
  1166. 'SI' ('NEG' nincq dim) ;
  1167. cherr = 'CHAINE'
  1168. 'la duale doit etre un vecteur' ;
  1169. 'ERREUR' cherr ;
  1170. 'FINSI' ;
  1171. *ninc = ninct ;
  1172. *
  1173. lcof = 'EXISTE' TDISC nomo ;
  1174. 'SI' lcof ; tcof = TDISC . nomo ;
  1175. ncof = 'DIME' (tcof . 'NOMINC') ;
  1176. 'SINON' ; ncof = 0 ;
  1177. 'FINSI' ;
  1178. *
  1179. 'SI' debug ;
  1180. 'SI' lcof ; 'MESSAGE' 'Un coef a ete detecte' ;
  1181. 'SINON' ; 'MESSAGE' 'pas de coef detecte' ;
  1182. 'FINSI' ;
  1183. 'FINSI' ;
  1184. *
  1185. vdim = 'VALEUR' 'DIME' ;
  1186. vmod = 'VALEUR' 'MODE' ;
  1187. idim = 0 ;
  1188. 'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'PLANDEFO')) ;
  1189. idim = 2 ;
  1190. iaxi = FAUX ;
  1191. 'FINSI' ;
  1192. 'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'AXIS')) ;
  1193. idim = 2 ;
  1194. iaxi = VRAI ;
  1195. 'FINSI' ;
  1196. 'SI' ('ET' ('EGA' vdim 3) ('EGA' vmod 'TRID')) ;
  1197. idim = 3 ;
  1198. iaxi = FAUX ;
  1199. 'FINSI' ;
  1200. 'SI' ('EGA' vdim 1) ;
  1201. idim = 1 ;
  1202. iaxi = FAUX ;
  1203. 'FINSI' ;
  1204. * 'MESSAGE' ('CHAINE' 'iaxi=' iaxi );
  1205. 'SI' ('EGA' idim 0) ;
  1206. 'ERREUR' ('CHAINE' 'vmod=' vmod ' et vdim=' vdim ' non prevu') ;
  1207. 'FINSI' ;
  1208. 'SI' iaxi ;
  1209. dprmt = '*' ('COORDONNEE' 1 _mt) ('*' PI 2.D0) ;
  1210. 'FINSI' ;
  1211. *
  1212. * Optimisation possible : construire la matrice par blocs
  1213. * qd valt et valq ne sont pas donnés
  1214. *
  1215. numop = nincq ; numder = idim ; numvar = ninct ;
  1216. numdat = ncof ; numcof = ncof ;
  1217. A = ININLIN numop numvar numdat numcof numder ;
  1218. 'SI' lcof ;
  1219. lvo = 'EGA' ('TYPE' valo) 'LISTREEL' ;
  1220. 'REPETER' iicof ncof ;
  1221. icof = &iicof ;
  1222. A . 'DAT' . icof . 'NOMDDL' = tcof . 'NOMINC' . icof ;
  1223. A . 'DAT' . icof . 'DISC' = tcof . 'DISC' ;
  1224. 'SI' lvo ;
  1225. A . 'DAT' . icof . 'VALEUR' = 'EXTRAIRE' valo icof ;
  1226. 'SINON' ;
  1227. A . 'DAT' . icof . 'VALEUR' = valo ;
  1228. 'FINSI' ;
  1229. A . 'COF' . icof . 'COMPOR' = 'IDEN' ;
  1230. A . 'COF' . icof . 'LDAT' = 'LECT' icof ;
  1231. 'FIN' iicof ;
  1232. 'SINON' ;
  1233. ll = 'LECT' ;
  1234. 'FINSI' ;
  1235. lvt = 'EGA' ('TYPE' valt) 'LISTREEL' ;
  1236. 'REPETER' iiincq nincq ;
  1237. iincq = &iiincq ;
  1238. iinct = 'MINIMUM' ('LECT' iincq ninct) ;
  1239. A . 'VAR' . iinct . 'NOMDDL' = tnomt . 'NOMINC' . iinct ;
  1240. A . 'VAR' . iinct . 'DISC' = tnomt . 'DISC' ;
  1241. 'SI' lvalt ;
  1242. 'SI' lvt ;
  1243. A . 'VAR' . iinct . 'VALEUR' = 'EXTRAIRE' valt iinct ;
  1244. 'SINON' ;
  1245. A . 'VAR' . iinct . 'VALEUR' = valt ;
  1246. 'FINSI' ;
  1247. 'FINSI' ;
  1248. 'SI' lcof ;
  1249. icof = 'MINIMUM' ('LECT' iincq ncof) ;
  1250. A . iincq . iinct . 0 = 'LECT' icof ;
  1251. 'SINON' ;
  1252. A . iincq . iinct . 0 = ll ;
  1253. 'FINSI' ;
  1254. 'FIN' iiincq ;
  1255. *
  1256. 'SI' iaxi ;
  1257. numdat = 1 ;
  1258. numcof = dim '+' 1 ;
  1259. 'SINON' ;
  1260. numdat = 0 ;
  1261. numcof = dim ;
  1262. 'FINSI' ;
  1263. numvar = nincq ;
  1264. B = ININLIN numop numvar numdat numcof numder ;
  1265. icof = 0 ;
  1266. 'REPETER' iiidim idim ;
  1267. icof = '+' icof 1 ;
  1268. B . 'COF' . icof . 'COMPOR' = 'CHAINE' 'VNOR' &iiidim ;
  1269. B . 'COF' . icof . 'LDAT' = 'LECT' ;
  1270. 'FIN' iiidim ;
  1271. *
  1272. 'SI' iaxi ;
  1273. icof = '+' icof 1 ;
  1274. B . 'DAT' . 1 . 'NOMDDL' = 'MOTS' 'SCAL' ;
  1275. B . 'DAT' . 1 . 'DISC' = discg ;
  1276. B . 'DAT' . 1 . 'VALEUR' = dprmt ;
  1277. B . 'COF' . icof . 'COMPOR' = 'IDEN' ;
  1278. B . 'COF' . icof . 'LDAT' = 'LECT' 1 ;
  1279. ll = 'LECT' icof ;
  1280. 'SINON' ;
  1281. ll = 'LECT' ;
  1282. 'FINSI' ;
  1283. lvq = 'EGA' ('TYPE' valq) 'LISTREEL' ;
  1284. 'REPETER' iiincq nincq ;
  1285. iincq = &iiincq ;
  1286. B . 'VAR' . iincq . 'NOMDDL' = tnomq . 'NOMINC' . iincq ;
  1287. B . 'VAR' . iincq . 'DISC' = tnomq . 'DISC' ;
  1288. 'SI' lvalq ;
  1289. 'SI' lvq ;
  1290. B . 'VAR' . iincq . 'VALEUR' = 'EXTRAIRE' valq iincq ;
  1291. 'SINON' ;
  1292. B . 'VAR' . iincq . 'VALEUR' = valq ;
  1293. 'FINSI' ;
  1294. 'FINSI' ;
  1295. B . iincq . iincq . 0 = ('LECT' iincq) 'ET' ll ;
  1296. 'FIN' iiincq ;
  1297. *
  1298. mgnor = NLIN discg _mt A B methgau ;
  1299. *
  1300. 'RESPRO' mgnor ;
  1301. 'FINPROC' ;
  1302. *
  1303. * End of procedure file GNOR
  1304. *
  1305. *ENDPROCEDUR gnor
  1306. *BEGINPROCEDUR gvol
  1307. ************************************************************************
  1308. * NOM : GVOL
  1309. * DESCRIPTION :
  1310. * Calcule le volume compris dans une surface fermée
  1311. * La normale doit être vers l'intérieur pour que le volume soit positif
  1312. *
  1313. *
  1314. * LANGAGE : GIBIANE-CAST3M
  1315. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  1316. * mél : gounand@semt2.smts.cea.fr
  1317. **********************************************************************
  1318. * VERSION : v1, 22/04/2011, version initiale
  1319. * HISTORIQUE : v1, 22/04/2011, création
  1320. * HISTORIQUE :
  1321. * HISTORIQUE :
  1322. ************************************************************************
  1323. *
  1324. *
  1325. 'DEBPROC' GVOL ;
  1326. 'ARGUMENT' _surf*'MAILLAGE' ;
  1327. 'ARGUMENT' tdisc*'TABLE' ;
  1328. 'ARGUMENT' dbg/'LOGIQUE' ;
  1329. *
  1330. 'SI' ('NON' ('EXISTE' dbg)) ;
  1331. dbg = FAUX ;
  1332. 'FINSI' ;
  1333. *
  1334. * Vecteur position et calcul du volume
  1335. NOMVIT = @STBL (TDISC . 'XN' . 'NOMINC') ;
  1336. DISCG = TDISC . 'GEOM' . 'DISC' ;
  1337. vdim = 'VALEUR' 'DIME' ;
  1338. 'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'AXIS')) ;
  1339. fdim = 3 ;
  1340. 'SINON' ;
  1341. fdim = vdim ;
  1342. 'FINSI' ;
  1343. vposc = GETCOO _surf nomvit ;
  1344. * 'SI' iaxi ;
  1345. * rs zs = 'COORDONNEE' _surf ;
  1346. * nr = 'EXTRAIRE' nomvit 1 ;
  1347. * nz = 'EXTRAIRE' nomvit 2 ;
  1348. * vposc =
  1349. * 'FINSI' ;
  1350. * fvol = GNOR _surf tdisc 'NPRI' ('CHAINE' discg 'V') 'CPRI' vpos
  1351. * 'NDUA' 'CSTEV' ;
  1352. * rfvol = 'RESULT' fvol ;
  1353. * volx = 'MAXIMUM' ('EXCO' 'UX' rfvol) ;
  1354. * voly = 'MAXIMUM' ('EXCO' 'UY' rfvol) ;
  1355. * vol = '/' ('+' volx voly) vdim ;
  1356. fvolc = GNOR _surf tdisc 'NPRI' discg
  1357. 'NCOF' (chai discg 'V')
  1358. 'CCOF' vposc
  1359. 'NDUA' (chai discg 'V')
  1360. 'FDUA' ('PROG' vdim * 1.) ;
  1361. volc = '/' ('MAXIMUM' ('RESULT' fvolc))
  1362. fdim ;
  1363. vol = volc '*' -1. ;
  1364. 'RESPRO' vol ;
  1365. *
  1366. * End of procedure file GVOL
  1367. *
  1368. 'FINPROC' ;
  1369. *ENDPROCEDUR gvol
  1370. *BEGINPROCEDUR log10
  1371. ************************************************************************
  1372. * NOM : LOG10
  1373. * DESCRIPTION : Log_10
  1374. *
  1375. *
  1376. *
  1377. * LANGAGE : GIBIANE-CAST3M
  1378. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  1379. * mél : gounand@semt2.smts.cea.fr
  1380. **********************************************************************
  1381. * VERSION : v1, 18/02/2003, version initiale
  1382. * HISTORIQUE : v1, 18/02/2003, création
  1383. * HISTORIQUE :
  1384. * HISTORIQUE :
  1385. ************************************************************************
  1386. * Prière de PRENDRE LE TEMPS de compléter les commentaires
  1387. * en cas de modification de ce sous-programme afin de faciliter
  1388. * la maintenance !
  1389. ************************************************************************
  1390. *
  1391. *
  1392. 'DEBPROC' LOG10 ;
  1393. 'REPETER' bouc ;
  1394. ok = FAUX ;
  1395. 'ARGUMENT' fl/'FLOTTANT' ;
  1396. 'ARGUMENT' lr/'LISTREEL' ;
  1397. 'ARGUMENT' cp/'CHPOINT ' ;
  1398. 'ARGUMENT' cm/'MCHAML ' ;
  1399. 'SI' ('EXISTE' fl) ;
  1400. ok = VRAI ;
  1401. 'RESPRO' ('/' ('LOG' fl) ('LOG' 10.D0)) ;
  1402. 'FINSI' ;
  1403. 'SI' ('EXISTE' lr) ;
  1404. ok = VRAI ;
  1405. 'RESPRO' ('/' ('LOG' lr) ('LOG' 10.D0)) ;
  1406. 'FINSI' ;
  1407. 'SI' ('EXISTE' cp) ;
  1408. ok = VRAI ;
  1409. 'RESPRO' ('/' ('LOG' cp) ('LOG' 10.D0)) ;
  1410. 'FINSI' ;
  1411. 'SI' ('EXISTE' cm) ;
  1412. ok = VRAI ;
  1413. 'RESPRO' ('/' ('LOG' cm) ('LOG' 10.D0)) ;
  1414. 'FINSI' ;
  1415. 'SI' ('NON' ok) ;
  1416. 'QUITTER' bouc ;
  1417. 'FINSI' ;
  1418. 'FIN' bouc ;
  1419. *
  1420. * End of procedure file LOG10
  1421. *
  1422. 'FINPROC' ;
  1423. *ENDPROCEDUR log10
  1424. *BEGINPROCEDUR modulo
  1425. ************************************************************************
  1426. * NOM : MODULO
  1427. * DESCRIPTION : Calcule un entier modulo un autre...
  1428. *
  1429. *
  1430. *
  1431. * LANGAGE : GIBIANE-CAST3M
  1432. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  1433. * mél : gounand@semt2.smts.cea.fr
  1434. **********************************************************************
  1435. * VERSION : v1, 15/10/2002, version initiale
  1436. * HISTORIQUE : v1, 15/10/2002, création
  1437. * HISTORIQUE :
  1438. * HISTORIQUE :
  1439. ************************************************************************
  1440. * Prière de PRENDRE LE TEMPS de compléter les commentaires
  1441. * en cas de modification de ce sous-programme afin de faciliter
  1442. * la maintenance !
  1443. ************************************************************************
  1444. *
  1445. *
  1446. 'DEBPROC' MODULO ;
  1447. 'ARGUMENT' i*'ENTIER' j*'ENTIER' ;
  1448. 'SI' ('EGA' j 0) ;
  1449. 'MESSAGE' 'Impossible de faire modulo 0' ;
  1450. 'ERREUR' 5 ;
  1451. 'SINON' ;
  1452. k=i '/' j ;
  1453. mod=i '-' ( k '*'j ) ;
  1454. 'RESPRO' mod ;
  1455. 'FINSI' ;
  1456. *
  1457. * End of procedure file MODULO
  1458. *
  1459. 'FINPROC' ;
  1460. *ENDPROCEDUR modulo
  1461. *BEGINPROCEDUR projsysc
  1462. ************************************************************************
  1463. * NOM : PROJSYSC
  1464. * DESCRIPTION : Calcul matrice et second membre projetés suivant
  1465. * un champ de directions données
  1466. *
  1467. * Project a linear system with respect to a given
  1468. * vector field
  1469. *
  1470. *
  1471. *
  1472. * LANGAGE : GIBIANE-CAST3M
  1473. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  1474. * mél : gounand@semt2.smts.cea.fr
  1475. **********************************************************************
  1476. * VERSION : v1, 22/04/2011, version initiale
  1477. * HISTORIQUE : v1, 22/04/2011, création
  1478. * HISTORIQUE :
  1479. * HISTORIQUE :
  1480. ************************************************************************
  1481. *
  1482. *
  1483. 'DEBPROC' PROJSYSC ;
  1484. 'ARGUMENT' tdisc*'TABLE' ;
  1485. 'ARGUMENT' vnorn*'CHPOINT' ;
  1486. 'ARGUMENT' ktgra*'RIGIDITE' ;
  1487. 'ARGUMENT' fpgra*'CHPOINT' ;
  1488. 'ARGUMENT' kvol/'CHPOINT' ;
  1489. lcnt = 'EXISTE' kvol ;
  1490. 'SI' lcnt ;
  1491. 'ARGUMENT' dvol*'FLOTTANT' ;
  1492. 'FINSI' ;
  1493. vdim = 'VALEUR' 'DIME' ;
  1494. NOMVIT = @STBL (TDISC . 'XN' . 'NOMINC') ;
  1495. fpgran = 'PSCAL' fpgra vnorn nomvit nomvit ;
  1496. * Condensation de la matrice
  1497. * tknor = 'TABLE' 'ESCLAVE' ;
  1498. * 'REPETER' idim vdim ;
  1499. * lnomi = 'EXTRAIRE' nomvit ('LECT' &idim) ;
  1500. * tknor . &idim = 'KOPS' 'MATDIAGO' lnomi
  1501. * ('EXCO' lnomi vnorn lnomi) ;
  1502. * 'FIN' idim ;
  1503. * knord = 'ET' tknor ;
  1504. 'SI' lmatrik ;
  1505. knord = 'KOPS' 'MATDIAGO' vnorn 'MATRIK' ;
  1506. 'SINON' ;
  1507. knord = 'KOPS' 'MATDIAGO' vnorn ;
  1508. knord = 'CHANGER' 'INCO' knord nomvit nomvit
  1509. nomfor nomvit ;
  1510. 'FINSI' ;
  1511. 'SI' ('EGA' vdim 2) ;
  1512. nomscal = 'MOTS' 'SCAL' 'SCAL' ;
  1513. 'SINON' ;
  1514. nomscal = 'MOTS' 'SCAL' 'SCAL' 'SCAL' ;
  1515. 'FINSI' ;
  1516. knor = 'CHANGER' 'INCO' knord nomvit nomscal nomvit nomvit ;
  1517. knort = 'CHANGER' 'INCO' knord nomvit nomvit nomvit nomscal ;
  1518. 'SI' lmatrik ;
  1519. ktgrak = 'KOPS' 'RIMA' ktgra ;
  1520. 'SINON' ;
  1521. ktgrak = ktgra ;
  1522. 'FINSI' ;
  1523. ktg1 = 'KOPS' 'CMCT' ktgrak knort ;
  1524. ktg2 = 'KOPS' 'TRANSPOS' ktg1 ;
  1525. ktg3 = 'KOPS' 'CMCT' knort ktg2 ;
  1526. ktot = ktg3 ;
  1527. ftot = fpgran ;
  1528. 'SI' lcnt ;
  1529. * 'MESSAGE' 'Une contrainte dans projsysc' ;
  1530. ktvol = 'PSCAL' kvol vnorn nomvit nomvit ;
  1531. ktv = rela ('NOMC' 'T' ktvol) ;
  1532. smbvol = 'DEPIMPOSE' ktv dvol ;
  1533. 'SI' lmatrik ;
  1534. ktv = 'KOPS' 'RIMA' ktv ;
  1535. 'FINSI' ;
  1536. lpr1 = 'MOTS' 'T' 'LX' ; l2 = 'MOTS' 'SCAL' 'LX' ;
  1537. ldu1 = 'MOTS' 'Q' 'FLX' ;
  1538. ktv = 'CHANGER' 'INCO' ktv lpr1 l2 ldu1 l2 'MULT' ;
  1539. smbvol = 'NOMC' ldu1 l2 smbvol ;
  1540. ktot = ktot 'ET' ktv ;
  1541. ftot = ftot '+' smbvol ;
  1542. 'SINON' ;
  1543. * 'MESSAGE' 'Pas de contrainte dans projsysc' ;
  1544. 'FINSI' ;
  1545. 'RESPRO' ktot ftot ;
  1546. *
  1547. * End of procedure file PROJSYSC
  1548. *
  1549. 'FINPROC' ;
  1550. *ENDPROCEDUR projsysc
  1551. *BEGINPROCEDUR quafme
  1552. ************************************************************************
  1553. * NOM : QUAFME
  1554. * DESCRIPTION :
  1555. *
  1556. *
  1557. *
  1558. * LANGAGE : GIBIANE-CAST3M
  1559. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  1560. * mél : gounand@semt2.smts.cea.fr
  1561. **********************************************************************
  1562. * VERSION : v1, 01/12/2004, version initiale
  1563. * HISTORIQUE : v1, 01/12/2004, création
  1564. * HISTORIQUE :
  1565. * HISTORIQUE :
  1566. ************************************************************************
  1567. * Prière de PRENDRE LE TEMPS de compléter les commentaires
  1568. * en cas de modification de ce sous-programme afin de faciliter
  1569. * la maintenance !
  1570. ************************************************************************
  1571. *
  1572. *
  1573. 'DEBPROC' QUAFME ;
  1574. 'REPETER' bcl ;
  1575. 'ARGUMENT' mquad/'MAILLAGE' ;
  1576. 'SI' ('EXISTE' mquad) ;
  1577. mquaf = 'CHANGER' mquad 'QUAF' ;
  1578. * mlin = 'CHANGER' mquad 'LINEAIRE' ;
  1579. 'RESPRO' mquaf ;
  1580. 'SINON' ;
  1581. 'QUITTER' bcl ;
  1582. 'FINSI' ;
  1583. 'FIN' bcl ;
  1584. 'FINPROC' ;
  1585. *
  1586. * End of procedure file QUAFME
  1587. *
  1588. *ENDPROCEDUR quafme
  1589. *BEGINPROCEDUR trvec
  1590. ************************************************************************
  1591. * NOM : TRVEC
  1592. * DESCRIPTION : Trace des champs de vecteurs.
  1593. * Utile pour tracer des bilans de forces
  1594. *
  1595. * Display vector fields.
  1596. * Useful for visualization of force balance.
  1597. *
  1598. *
  1599. *
  1600. * LANGAGE : GIBIANE-CAST3M
  1601. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  1602. * mél : gounand@semt2.smts.cea.fr
  1603. **********************************************************************
  1604. * VERSION : v1, 22/04/2011, version initiale
  1605. * HISTORIQUE : v1, 22/04/2011, création
  1606. * HISTORIQUE :
  1607. * HISTORIQUE :
  1608. ************************************************************************
  1609. *
  1610. *
  1611. 'DEBPROC' TRVEC ;
  1612. *'ARGUMENT' tdisc*'TABLE' ;
  1613. *'ARGUMENT' motdom*'MOT' ;
  1614. 'ARGUMENT' tdom*'MAILLAGE' ;
  1615. *
  1616. tvec = 'TABLE' ; ttit = 'TABLE' ;
  1617. i = 0 ;
  1618. lcoul = 'MOTS' 'JAUN' 'ROUG' 'BLAN' 'TURQ' 'VERT' 'OLIV'
  1619. 'AZUR' 'ORAN' 'VIOL' 'GRIS' 'OCEA' ;
  1620. *
  1621. 'REPETER' livec ;
  1622. 'SI' ('EGA' i 0) ;
  1623. 'ARGUMENT' ccvec*'CHPOINT' ;
  1624. 'SINON' ;
  1625. 'ARGUMENT' ccvec/'CHPOINT' ;
  1626. 'FINSI' ;
  1627. 'SI' ('EXISTE' ccvec) ;
  1628. 'ARGUMENT' ttvec*'MOT' ;
  1629. 'SINON' ;
  1630. 'QUITTER' livec ;
  1631. 'FINSI' ;
  1632. i = '+' i 1 ;
  1633. * 'MESSAGE' ('CHAINE' 'i=' i) ;
  1634. * 'LISTE' ccvec ;
  1635. * 'LISTE' tvec ;
  1636. tvec . i = ccvec ;
  1637. ttit . i = ttvec ;
  1638. 'FIN' livec ;
  1639. 'ARGUMENT' echv/'FLOTTANT' ;
  1640. 'ARGUMENT' lnclk/'LOGIQUE' ;
  1641. 'SI' ('NON' ('EXISTE' lnclk)) ;
  1642. lnclk = faux ;
  1643. 'FINSI' ;
  1644. *
  1645. lmax = 'PROG' ;
  1646. 'REPETER' ii i ;
  1647. mm = 'MAXIMUM' (tvec . &ii) 'ABS' ;
  1648. lmax = 'ET' lmax ('PROG' mm) ;
  1649. 'FIN' ii ;
  1650. mm = '+' ('MAXIMUM' lmax) 1.D-60 ;
  1651. 'SI' ('NON' ('EXISTE' echv)) ;
  1652. vmodo = 'EGA' ('VALEUR' 'MODE') 'AXIS' ;
  1653. 'SI' vmodo ; 'OPTI' 'MODE' 'PLAN' ; 'FINSI' ;
  1654. ctail = gmass2 ('CHANGER' tdom 'QUAF') tdisc
  1655. 'NPRI' 'CSTE' 'FPRI' 1. 'NDUA' 'CSTE' 'FDUA' 1. ;
  1656. * 'LISTE' ctail ;
  1657. 'SI' vmodo ; 'OPTI' 'MODE' 'AXIS' ; 'FINSI' ;
  1658. vdim = 'VALEUR' 'DIME' ;
  1659. dimm = DEADUTIL 'DIMM' tdom ;
  1660. * ctail = '**' ctail ('/' 1. ('-' vdim 1)) ;
  1661. ctail = '**' ctail ('/' 1. dimm) ;
  1662. tail = '**' ('*' ('MAXIMUM' ctail) ('MINIMUM' ctail)) 0.5 ;
  1663. *'LISTE' tail ;
  1664. *'LISTE' mm ;
  1665. echv = '/' ('*' 0.9 tail) mm ;
  1666. 'FINSI' ;
  1667. tit = 'CHAINE' 'Max. =' (formar mm 2) ;
  1668. *'MESSAGE' ('CHAINE' 'mm=' mm) ;
  1669. 'REPETER' ii i ;
  1670. cou = EXMOMOD lcoul &ii ;
  1671. tvec . &ii = 'VECT' (tvec . &ii) echv 'DEPL' cou ;
  1672. tit = 'CHAINE' tit ' ' cou '=' (ttit . &ii) ;
  1673. 'FIN' ii ;
  1674. 'SI' lnclk ;
  1675. 'TRACER' (@stbl tvec) tdom 'TITR' tit 'NCLK' ;
  1676. 'SINON' ;
  1677. 'TRACER' (@stbl tvec) tdom 'TITR' tit ;
  1678. 'FINSI' ;
  1679. *
  1680. * End of procedure file TRVEC
  1681. *
  1682. 'FINPROC' ;
  1683. *ENDPROCEDUR trvec
  1684. *BEGINPROCEDUR tsurfonc
  1685. ************************************************************************
  1686. * NOM : TSURFONC
  1687. * DESCRIPTION : La fonctionnelle à minimiser pour la tension
  1688. * de surface
  1689. *
  1690. *
  1691. * LANGAGE : GIBIANE-CAST3M
  1692. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  1693. * mél : gounand@semt2.smts.cea.fr
  1694. **********************************************************************
  1695. * VERSION : v1, 17/01/2007, version initiale
  1696. * HISTORIQUE : v1, 17/01/2007, création
  1697. * HISTORIQUE :
  1698. * HISTORIQUE :
  1699. ************************************************************************
  1700. * Prière de PRENDRE LE TEMPS de compléter les commentaires
  1701. * en cas de modification de ce sous-programme afin de faciliter
  1702. * la maintenance !
  1703. ************************************************************************
  1704. *
  1705. *
  1706. 'DEBPROC' TSURFONC ;
  1707. 'ARGUMENT' _mt*'MAILLAGE' ;
  1708. 'ARGUMENT' gdisc*'MOT' ;
  1709. 'ARGUMENT' methgau*'MOT' ;
  1710. *
  1711. idim = 'VALEUR' 'DIME' ;
  1712. vdim = DEADUTIL 'DIMM' _mt ;
  1713. *
  1714. 'ARGUMENT' coef/'FLOTTANT' ;
  1715. 'SI' ('NON' ('EXISTE' coef)) ;
  1716. 'ARGUMENT' coef2/'CHPOINT ' ;
  1717. 'SI' ('NON' ('EXISTE' coef2)) ;
  1718. 'ERREUR' 'Il faut donner un coef FLOTTANT ou CHPOINT' ;
  1719. 'SINON' ;
  1720. coef = coef2 ;
  1721. 'ARGUMENT' discc*'MOT ' ;
  1722. 'FINSI' ;
  1723. 'SINON' ;
  1724. discc = 'CSTE' ;
  1725. 'FINSI' ;
  1726. *'ARGUMENT' met/'CHPOINT' ;
  1727. *lmet = 'EXISTE' met ;
  1728. *'SI' lmet ;
  1729. * debloi = 'CHAINE' 'ADD' ;
  1730. * ncmet = '/' ('*' idim ('+' idim 1)) 2 ;
  1731. * 'ARGUMENT' metdisc*'MOT' ;
  1732. ** metdisc = gdisc ;
  1733. *'SINON' ;
  1734. * debloi = 'CHAINE' 'ADC' ;
  1735. * ncmet = 0 ;
  1736. *'FINSI' ;
  1737. *loi = 'CHAINE' debloi 'F' ;
  1738. 'ARGUMENT' optelem/'MOT' ;
  1739. *
  1740. 'SI' ('EXISTE' optelem) ;
  1741. 'SI' ('EGA' optelem 'ELEM') ;
  1742. lelem = VRAI ;
  1743. 'SINON' ;
  1744. cherr = 'CHAINE' 'Option ' optelem ' inconnue' ;
  1745. 'ERREUR' cherr ;
  1746. 'FINSI' ;
  1747. 'SINON' ;
  1748. lelem = FAUX ;
  1749. 'FINSI' ;
  1750. *
  1751. * Calcul de la fonctionnelle
  1752. *
  1753. numop = 1 ;
  1754. numvar = 1 ;
  1755. numder = vdim ;
  1756. numdat = 0 ;
  1757. numcof = 0 ;
  1758. *
  1759. A = ININLIN numop numvar numdat numcof numder ;
  1760. A . 'VAR' . 1 . 'NOMDDL' = 'MOTS' 'DUMM' ;
  1761. A . 'VAR' . 1 . 'DISC' = 'CSTE' ;
  1762. A . 'VAR' . 1 . 'VALEUR' = 1.D0 ;
  1763. *
  1764. numvar = 1 ;
  1765. numdat = 1 ;
  1766. numcof = 1 ;
  1767. B = ININLIN numop numvar numdat numcof numder ;
  1768. B . 'VAR' . 1 . 'NOMDDL' = 'MOTS' 'DUMM' ;
  1769. B . 'VAR' . 1 . 'DISC' = 'CSTE' ;
  1770. B . 'VAR' . 1 . 'VALEUR' = 1.D0 ;
  1771. B . 'DAT' . 1 . 'NOMDDL' = 'MOTS' 'SCAL' ;
  1772. B . 'DAT' . 1 . 'DISC' = discc ;
  1773. B . 'DAT' . 1 . 'VALEUR' = coef ;
  1774. *
  1775. B . 'COF' . 1 . 'COMPOR' = 'TSUF' ;
  1776. B . 'COF' . 1 . 'LDAT' = 'LECT' 1 ;
  1777. *
  1778. A . 1 . 1 . 0 = 'LECT' ;
  1779. B . 1 . 1 . 0 = 'LECT' 1 ;
  1780. *
  1781. vfonc = 'NLIN' gdisc _mt A B 'EREF' methgau ;
  1782. *
  1783. 'SI' ('NON' lelem) ;
  1784. vfonc = 'MAXIMUM' ('RESULT' vfonc) ;
  1785. 'FINSI' ;
  1786. *
  1787. 'RESPRO' vfonc ;
  1788. *
  1789. * End of procedure file TSURFONC
  1790. *
  1791. 'FINPROC' ;
  1792. *ENDPROCEDUR tsurfonc
  1793. *BEGINPROCEDUR tsurktan
  1794. ************************************************************************
  1795. * NOM : TSURKTAN
  1796. * DESCRIPTION : La matrice tangente pour la tension de surface
  1797. *
  1798. *
  1799. *
  1800. * LANGAGE : GIBIANE-CAST3M
  1801. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  1802. * mél : gounand@semt2.smts.cea.fr
  1803. **********************************************************************
  1804. * VERSION : v1, 17/01/2007, version initiale
  1805. * HISTORIQUE : v1, 17/01/2007, création
  1806. * HISTORIQUE :
  1807. * HISTORIQUE :
  1808. ************************************************************************
  1809. * Prière de PRENDRE LE TEMPS de compléter les commentaires
  1810. * en cas de modification de ce sous-programme afin de faciliter
  1811. * la maintenance !
  1812. ************************************************************************
  1813. *
  1814. *
  1815. 'DEBPROC' TSURKTAN ;
  1816. 'ARGUMENT' _mt*'MAILLAGE' ;
  1817. 'ARGUMENT' gdisc*'MOT' ;
  1818. 'ARGUMENT' methgau*'MOT' ;
  1819. 'ARGUMENT' dppri*'LISTMOTS' ;
  1820. 'ARGUMENT' dpdua*'LISTMOTS' ;
  1821. *
  1822. dpdis = gdisc ;
  1823. *
  1824. idim = 'VALEUR' 'DIME' ;
  1825. vdim = DEADUTIL 'DIMM' _mt ;
  1826. laxi = DEADUTIL 'AXI?' ;
  1827. lsph = DEADUTIL 'SPH?' ;
  1828. *
  1829. loi = 'CHAINE' 'TSUJ' ;
  1830. loij = 'CHAINE' 'TSU' ;
  1831. *
  1832. 'ARGUMENT' coef/'FLOTTANT' ;
  1833. 'SI' ('NON' ('EXISTE' coef)) ;
  1834. 'ARGUMENT' coef2/'CHPOINT ' ;
  1835. 'SI' ('NON' ('EXISTE' coef2)) ;
  1836. 'ERREUR' 'Il faut donner un coef FLOTTANT ou CHPOINT' ;
  1837. 'SINON' ;
  1838. coef = coef2 ;
  1839. 'ARGUMENT' discc*'MOT ' ;
  1840. 'FINSI' ;
  1841. 'SINON' ;
  1842. discc = 'CSTE' ;
  1843. 'FINSI' ;
  1844. *
  1845. 'ARGUMENT' jaco/'ENTIER' ;
  1846. 'SI' ('NON' ('EXISTE' jaco)) ;
  1847. jaco = 1 ;
  1848. dir1 = VRAI ;
  1849. 'FINSI' ;
  1850. 'SI' ('OU' ('EGA' jaco 2) ('EGA' jaco 3)) ;
  1851. 'ARGUMENT' idir/'ENTIER' ;
  1852. 'SI' ('EXISTE' idir) ;
  1853. ldir = 'LECT' idir ;
  1854. 'SI' ('EGA' idir 1) ;
  1855. dir1 = VRAI ;
  1856. 'FINSI' ;
  1857. 'SINON' ;
  1858. 'ARGUMENT' ldir/'LISTENTI' ;
  1859. 'SI' ('NON' ('EXISTE' ldir)) ;
  1860. ldir = 'LECT' 1 'PAS' 1 idim ;
  1861. dir1 = VRAI ;
  1862. 'FINSI' ;
  1863. 'FINSI' ;
  1864. 'FINSI' ;
  1865. 'ARGUMENT' lterm/'LISTENTI' ;
  1866. llterm = 'EXISTE' lterm ;
  1867. 'SI' llterm ;
  1868. dlterm = 'DIME' lterm ;
  1869. 'SINON' ;
  1870. dlterm = 1 ;
  1871. 'FINSI' ;
  1872. *
  1873. * Calcul du jacobien complet (jaco = 1)
  1874. *
  1875. 'SI' ('EGA' jaco 1) ;
  1876. numop = '*' ('**' vdim 2) ('**' idim 2) ;
  1877. numop = '*' numop dlterm ;
  1878. 'SI' ('OU' laxi lsph) ;
  1879. numop = '+' numop ('*' (vdim '*' idim) 2) ;
  1880. 'FINSI' ;
  1881. 'SI' lsph ;
  1882. numop = '+' numop 1 ;
  1883. 'FINSI' ;
  1884. numder = vdim ;
  1885. numvar = idim ;
  1886. numdat = 0 ;
  1887. numcof = 0 ;
  1888. *
  1889. A = ININLIN numop numvar numdat numcof numder ;
  1890. numdat = 1 ;
  1891. numcof = numop ;
  1892. B = ININLIN numop numvar numdat numcof numder ;
  1893. 'REPETER' ivar numvar ;
  1894. A . 'VAR' . &ivar . 'NOMDDL' = 'MOTS' ('EXTRAIRE' dppri &ivar) ;
  1895. A . 'VAR' . &ivar . 'DISC' = dpdis ;
  1896. B . 'VAR' . &ivar . 'NOMDDL' = 'MOTS' ('EXTRAIRE' dpdua &ivar) ;
  1897. B . 'VAR' . &ivar . 'DISC' = dpdis ;
  1898. 'FIN' ivar ;
  1899. iop = 0 ;
  1900. 'REPETER' h dlterm ;
  1901. 'REPETER' i idim ;
  1902. 'REPETER' j vdim ;
  1903. 'REPETER' k idim ;
  1904. 'REPETER' l vdim ;
  1905. iop = iop '+' 1 ;
  1906. A . iop . &i . &j = 'LECT' ;
  1907. 'SI' llterm ;
  1908. nl = 'EXTRAIRE' lterm &h ;
  1909. nomloi = 'CHAINE' loij nl &i &j &k &l ;
  1910. 'SINON' ;
  1911. nomloi = 'CHAINE' loi &i &j &k &l ;
  1912. 'FINSI' ;
  1913. * 'MESSAGE' ('CHAINE' 'Nomloi=' nomloi) ;
  1914. B . 'COF' . iop . 'COMPOR' = nomloi ;
  1915. B . 'COF' . iop . 'LDAT' = 'LECT' 1 ;
  1916. B . iop . &k . &l = 'LECT' iop ;
  1917. 'FIN' l ;
  1918. 'FIN' k ;
  1919. 'FIN' j ;
  1920. 'FIN' i ;
  1921. 'FIN' h ;
  1922. 'SI' ('OU' laxi lsph) ;
  1923. 'REPETER' i idim ;
  1924. 'REPETER' j vdim ;
  1925. iop = iop '+' 1 ;
  1926. A . iop . &i . &j = 'LECT' ;
  1927. nomloi = 'CHAINE' loi &i &j '10' ;
  1928. * 'MESSAGE' ('CHAINE' 'Nomloi=' nomloi) ;
  1929. B . 'COF' . iop . 'COMPOR' = nomloi ;
  1930. B . 'COF' . iop . 'LDAT' = 'LECT' 1 ;
  1931. B . iop . 1 . 0 = 'LECT' iop ;
  1932. 'FIN' j ;
  1933. 'FIN' i ;
  1934. 'REPETER' k idim ;
  1935. 'REPETER' l vdim ;
  1936. iop = iop '+' 1 ;
  1937. A . iop . 1 . 0 = 'LECT' ;
  1938. nomloi = 'CHAINE' loi '10' &k &l ;
  1939. * 'MESSAGE' ('CHAINE' 'Nomloi=' nomloi) ;
  1940. B . 'COF' . iop . 'COMPOR' = nomloi ;
  1941. B . 'COF' . iop . 'LDAT' = 'LECT' 1 ;
  1942. B . iop . &k . &l = 'LECT' iop ;
  1943. 'FIN' l ;
  1944. 'FIN' k ;
  1945. 'FINSI' ;
  1946. 'SI' lsph ;
  1947. iop = iop '+' 1 ;
  1948. A . iop . 1 . 0 = 'LECT' ;
  1949. nomloi = 'CHAINE' loi '1010' ;
  1950. * 'MESSAGE' ('CHAINE' 'Nomloi=' nomloi) ;
  1951. B . 'COF' . iop . 'COMPOR' = nomloi ;
  1952. B . 'COF' . iop . 'LDAT' = 'LECT' 1 ;
  1953. B . iop . 1 . 0 = 'LECT' iop ;
  1954. 'FINSI' ;
  1955. 'FINSI' ;
  1956. 'SI' ('EGA' jaco 2) ;
  1957. nldir = 'DIME' ldir ;
  1958. numop = '*' nldir ('**' vdim 2) ;
  1959. 'SI' dir1 ;
  1960. 'SI' ('OU' laxi lsph) ;
  1961. numop = '+' numop ('*' vdim 2) ;
  1962. 'FINSI' ;
  1963. 'SI' lsph ;
  1964. numop = '+' numop 1 ;
  1965. 'FINSI' ;
  1966. 'FINSI' ;
  1967. *
  1968. numder = vdim ;
  1969. numvar = idim ;
  1970. numdat = 0 ;
  1971. numcof = 0 ;
  1972. *
  1973. A = ININLIN numop numvar numdat numcof numder ;
  1974. numdat = 1 ;
  1975. numcof = numop ;
  1976. B = ININLIN numop numvar numdat numcof numder ;
  1977. 'REPETER' ivar numvar ;
  1978. A . 'VAR' . &ivar . 'NOMDDL' = 'MOTS' ('EXTRAIRE' dppri &ivar) ;
  1979. A . 'VAR' . &ivar . 'DISC' = dpdis ;
  1980. B . 'VAR' . &ivar . 'NOMDDL' = 'MOTS' ('EXTRAIRE' dpdua &ivar) ;
  1981. B . 'VAR' . &ivar . 'DISC' = dpdis ;
  1982. 'FIN' ivar ;
  1983. iop = 0 ;
  1984. 'REPETER' i nldir ;
  1985. idir = 'EXTRAIRE' ldir &i ;
  1986. 'REPETER' j vdim ;
  1987. 'REPETER' l vdim ;
  1988. iop = iop '+' 1 ;
  1989. A . iop . idir . &j = 'LECT' ;
  1990. nomloi = 'CHAINE' loi idir &j idir &l ;
  1991. * 'MESSAGE' ('CHAINE' 'Nomloi=' nomloi) ;
  1992. B . 'COF' . iop . 'COMPOR' = nomloi ;
  1993. B . 'COF' . iop . 'LDAT' = 'LECT' 1 ;
  1994. B . iop . idir . &l = 'LECT' iop ;
  1995. 'FIN' l ;
  1996. 'FIN' j ;
  1997. 'FIN' i ;
  1998. 'SI' dir1 ;
  1999. 'SI' ('OU' laxi lsph) ;
  2000. 'REPETER' j vdim ;
  2001. iop = iop '+' 1 ;
  2002. A . iop . 1 . &j = 'LECT' ;
  2003. nomloi = 'CHAINE' loi '1' &j '10' ;
  2004. * 'MESSAGE' ('CHAINE' 'Nomloi=' nomloi) ;
  2005. B . 'COF' . iop . 'COMPOR' = nomloi ;
  2006. B . 'COF' . iop . 'LDAT' = 'LECT' 1 ;
  2007. B . iop . 1 . 0 = 'LECT' iop ;
  2008. 'FIN' j ;
  2009. 'REPETER' l vdim ;
  2010. iop = iop '+' 1 ;
  2011. A . iop . 1 . 0 = 'LECT' ;
  2012. nomloi = 'CHAINE' loi '101' &l ;
  2013. * 'MESSAGE' ('CHAINE' 'Nomloi=' nomloi) ;
  2014. B . 'COF' . iop . 'COMPOR' = nomloi ;
  2015. B . 'COF' . iop . 'LDAT' = 'LECT' 1 ;
  2016. B . iop . 1 . &l = 'LECT' iop ;
  2017. 'FIN' l ;
  2018. 'FINSI' ;
  2019. 'SI' lsph ;
  2020. iop = iop '+' 1 ;
  2021. A . iop . 1 . 0 = 'LECT' ;
  2022. nomloi = 'CHAINE' loi '1010' ;
  2023. * 'MESSAGE' ('CHAINE' 'Nomloi=' nomloi) ;
  2024. B . 'COF' . iop . 'COMPOR' = nomloi ;
  2025. B . 'COF' . iop . 'LDAT' = 'LECT' 1 ;
  2026. B . iop . 1 . 0 = 'LECT' iop ;
  2027. 'FINSI' ;
  2028. 'FINSI' ;
  2029. 'FINSI' ;
  2030. *
  2031. 'SI' ('EGA' jaco 3) ;
  2032. nldir = 'DIME' ldir ;
  2033. * numop = '**' vdim 2 ;
  2034. numop = '*' nldir vdim ;
  2035. 'SI' ('ET' dir1 lsph) ;
  2036. numop = '+' numop 1 ;
  2037. 'FINSI' ;
  2038. *
  2039. numder = vdim ;
  2040. numvar = idim ;
  2041. numdat = 0 ;
  2042. numcof = 0 ;
  2043. *
  2044. A = ININLIN numop numvar numdat numcof numder ;
  2045. numdat = 1 ;
  2046. numcof = numop ;
  2047. B = ININLIN numop numvar numdat numcof numder ;
  2048. 'REPETER' ivar numvar ;
  2049. A . 'VAR' . &ivar . 'NOMDDL' = 'MOTS' ('EXTRAIRE' dppri &ivar) ;
  2050. A . 'VAR' . &ivar . 'DISC' = dpdis ;
  2051. B . 'VAR' . &ivar . 'NOMDDL' = 'MOTS' ('EXTRAIRE' dpdua &ivar) ;
  2052. B . 'VAR' . &ivar . 'DISC' = dpdis ;
  2053. 'FIN' ivar ;
  2054. iop = 0 ;
  2055. 'REPETER' i nldir ;
  2056. idir = 'EXTRAIRE' ldir &i ;
  2057. 'REPETER' j vdim ;
  2058. iop = iop '+' 1 ;
  2059. A . iop . idir . &j = 'LECT' ;
  2060. nomloi = 'CHAINE' loi idir &j idir &j ;
  2061. * 'MESSAGE' ('CHAINE' 'Nomloi=' nomloi) ;
  2062. B . 'COF' . iop . 'COMPOR' = nomloi ;
  2063. B . 'COF' . iop . 'LDAT' = 'LECT' 1 ;
  2064. B . iop . idir . &j = 'LECT' iop ;
  2065. 'FIN' j ;
  2066. 'FIN' i ;
  2067. 'SI' ('ET' dir1 lsph) ;
  2068. iop = iop '+' 1 ;
  2069. A . iop . 1 . 0 = 'LECT' ;
  2070. nomloi = 'CHAINE' loi '1010' ;
  2071. * 'MESSAGE' ('CHAINE' 'Nomloi=' nomloi) ;
  2072. B . 'COF' . iop . 'COMPOR' = nomloi ;
  2073. B . 'COF' . iop . 'LDAT' = 'LECT' 1 ;
  2074. B . iop . 1 . 0 = 'LECT' iop ;
  2075. 'FINSI' ;
  2076. 'FINSI' ;
  2077. *
  2078. * Partie commune
  2079. *
  2080. B . 'DAT' . 1 . 'NOMDDL' = 'MOTS' 'SCAL' ;
  2081. B . 'DAT' . 1 . 'DISC' = discc ;
  2082. B . 'DAT' . 1 . 'VALEUR' = coef ;
  2083. *
  2084. jac = 'NLIN' gdisc _mt A B 'EREF' methgau ;
  2085. *
  2086. 'RESPRO' jac ;
  2087. *
  2088. * End of procedure file TSURKTAN
  2089. *
  2090. 'FINPROC' ;
  2091. *ENDPROCEDUR tsurktan
  2092. *BEGINPROCEDUR tsurresi
  2093. ************************************************************************
  2094. * NOM : TSURRESI
  2095. * DESCRIPTION : Le résidu à annuler pour la tension de surface
  2096. *
  2097. *
  2098. *
  2099. * LANGAGE : GIBIANE-CAST3M
  2100. * AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  2101. * mél : gounand@semt2.smts.cea.fr
  2102. **********************************************************************
  2103. * VERSION : v1, 17/01/2007, version initiale
  2104. * HISTORIQUE : v1, 17/01/2007, création
  2105. * HISTORIQUE :
  2106. * HISTORIQUE :
  2107. ************************************************************************
  2108. * Prière de PRENDRE LE TEMPS de compléter les commentaires
  2109. * en cas de modification de ce sous-programme afin de faciliter
  2110. * la maintenance !
  2111. ************************************************************************
  2112. *
  2113. *
  2114. 'DEBPROC' TSURRESI ;
  2115. 'ARGUMENT' _mt*'MAILLAGE' ;
  2116. 'ARGUMENT' gdisc*'MOT' ;
  2117. 'ARGUMENT' methgau*'MOT' ;
  2118. 'ARGUMENT' dpdua*'LISTMOTS' ;
  2119. *
  2120. dpdis = gdisc ;
  2121. *
  2122. idim = 'VALEUR' 'DIME' ;
  2123. vdim = DEADUTIL 'DIMM' _mt ;
  2124. laxi = DEADUTIL 'AXI?' ;
  2125. lsph = DEADUTIL 'SPH?' ;
  2126. *
  2127. loi = 'CHAINE' 'TSUR' ;
  2128. *
  2129. 'ARGUMENT' coef/'FLOTTANT' ;
  2130. 'SI' ('NON' ('EXISTE' coef)) ;
  2131. 'ARGUMENT' coef2/'CHPOINT ' ;
  2132. 'SI' ('NON' ('EXISTE' coef2)) ;
  2133. 'ERREUR' 'Il faut donner un coef FLOTTANT ou CHPOINT' ;
  2134. 'SINON' ;
  2135. coef = coef2 ;
  2136. 'ARGUMENT' discc*'MOT ' ;
  2137. 'FINSI' ;
  2138. 'SINON' ;
  2139. discc = 'CSTE' ;
  2140. 'FINSI' ;
  2141. *
  2142. dir1 = FAUX ;
  2143. 'ARGUMENT' idir/'ENTIER' ;
  2144. 'SI' ('EXISTE' idir) ;
  2145. ldir = 'LECT' idir ;
  2146. 'SI' ('EGA' idir 1) ;
  2147. dir1 = VRAI ;
  2148. 'FINSI' ;
  2149. 'SINON' ;
  2150. 'ARGUMENT' ldir/'LISTENTI' ;
  2151. 'SI' ('NON' ('EXISTE' ldir)) ;
  2152. ldir = 'LECT' 1 'PAS' 1 idim ;
  2153. dir1 = VRAI ;
  2154. 'FINSI' ;
  2155. 'FINSI' ;
  2156. *
  2157. * Calcul du résidu
  2158. *
  2159. nldir = 'DIME' ldir ;
  2160. *
  2161. numop = '*' nldir vdim ;
  2162. term1 = ('OU' laxi lsph) 'ET' dir1 ;
  2163. 'SI' term1 ;
  2164. numop = '+' numop 1 ;
  2165. 'FINSI' ;
  2166. numder = vdim ;
  2167. numvar = 1 ;
  2168. numdat = 0 ;
  2169. numcof = 0 ;
  2170. *
  2171. A = ININLIN numop numvar numdat numcof numder ;
  2172. A . 'VAR' . 1 . 'NOMDDL' = 'MOTS' 'DUMM' ;
  2173. A . 'VAR' . 1 . 'DISC' = 'CSTE' ;
  2174. A . 'VAR' . 1 . 'VALEUR' = 1.D0 ;
  2175. *
  2176. numvar = idim ;
  2177. numdat = 1 ;
  2178. numcof = numop ;
  2179. B = ININLIN numop numvar numdat numcof numder ;
  2180. 'REPETER' ivar numvar ;
  2181. B . 'VAR' . &ivar . 'NOMDDL' = 'MOTS' ('EXTRAIRE' dpdua &ivar) ;
  2182. B . 'VAR' . &ivar . 'DISC' = dpdis ;
  2183. 'FIN' ivar ;
  2184. *
  2185. B . 'DAT' . 1 . 'NOMDDL' = 'MOTS' 'SCAL' ;
  2186. B . 'DAT' . 1 . 'DISC' = discc ;
  2187. B . 'DAT' . 1 . 'VALEUR' = coef ;
  2188. *
  2189. iop = 0 ;
  2190. 'REPETER' k nldir ;
  2191. idir = 'EXTRAIRE' ldir &k ;
  2192. 'REPETER' l vdim ;
  2193. iop = '+' iop 1 ;
  2194. A . iop . 1 . 0 = 'LECT' ;
  2195. nomloi = 'CHAINE' loi idir &l ;
  2196. * 'MESSAGE' ('CHAINE' 'Nomloi=' nomloi) ;
  2197. B . 'COF' . iop . 'COMPOR' = nomloi ;
  2198. B . 'COF' . iop . 'LDAT' = 'LECT' 1 ;
  2199. B . iop . idir . &l = 'LECT' iop ;
  2200. 'FIN' l ;
  2201. 'FIN' k ;
  2202. * 'LISTE' A ; 'LISTE' iop ;
  2203. 'SI' term1 ;
  2204. iop = '+' iop 1 ;
  2205. A . iop . 1 . 0 = 'LECT' ;
  2206. nomloi = 'CHAINE' loi '10' ;
  2207. B . 'COF' . iop . 'COMPOR' = nomloi ;
  2208. B . 'COF' . iop . 'LDAT' = 'LECT' 1 ;
  2209. B . iop . 1 . 0 = 'LECT' iop ;
  2210. 'FINSI' ;
  2211. *
  2212. res = 'NLIN' gdisc _mt A B 'EREF' methgau ;
  2213. *
  2214. 'RESPRO' res ;
  2215. *
  2216. * End of procedure file TSURRESI
  2217. *
  2218. 'FINPROC' ;
  2219. *ENDPROCEDUR tsurresi
  2220. **
  2221. * COPTCLIM Copy a table with boundary conditions
  2222. **
  2223. 'DEBPROC' COPTCLIM ;
  2224. 'ARGUMENT' otclim*'TABLE' ;
  2225. tclim = 'TABLE' ;
  2226. iotclim = index otclim ;
  2227. 'REPETER' ii ('DIME' iotclim) ;
  2228. indx = iotclim . &ii ;
  2229. tclim . indx = otclim . indx ;
  2230. 'FIN' ii ;
  2231. 'RESPRO' tclim ;
  2232. 'FINPROC' ;
  2233.  
  2234. **
  2235. * DIRDEP
  2236. * Procédure créant le champ de direction selon lequel
  2237. * les points se déplacent.
  2238. * idir = 0 : les points se déplacent suivant la normale à la surface
  2239. * idir = 1 : les points se déplacent suivant la direction passant par
  2240. * 0. 0.
  2241. *
  2242. 'DEBPROC' DIRDEP ;
  2243. 'ARGUMENT' _cmt*'MAILLAGE' ;
  2244. 'ARGUMENT' cmt*'MAILLAGE' ;
  2245. 'ARGUMENT' sur*'MAILLAGE' ;
  2246. 'ARGUMENT' tdisc*'TABLE' ;
  2247. 'ARGUMENT' idir/'ENTIER' ;
  2248. 'SI' ('NON' ('EXISTE' idir)) ;
  2249. idir = 0 ;
  2250. 'FINSI' ;
  2251. *
  2252. NOMVIT = @STBL (TDISC . 'XN' . 'NOMINC') ;
  2253. *
  2254. vdim = 'VALEUR' 'DIME' ;
  2255. DISCG = TDISC . 'GEOM' . 'DISC' ;
  2256. 'SI' ('EGA' idir 0) ;
  2257. vnor = GNOR _cmt tdisc 'NPRI' discg 'FPRI' 1. 'NDUA' 'XN' ;
  2258. vnor = '*' vnor -1. ;
  2259. nvnor = '**' ('PSCAL' vnor vnor nomvit nomvit) 0.5 ;
  2260. nvnor = '+' nvnor ('MASQUE' nvnor 'INFERIEUR' 1.D-100) ;
  2261. vnorn = '/' vnor nvnor ;
  2262. 'FINSI' ;
  2263. 'SI' ('EGA' idir 1) ;
  2264. vnor = GETCOO cmt nomvit ;
  2265. nvnor = '**' ('PSCAL' vnor vnor nomvit nomvit) 0.5 ;
  2266. nvnor = '+' nvnor ('MASQUE' nvnor 'INFERIEUR' 1.D-100) ;
  2267. vnorn = '/' vnor nvnor ;
  2268. 'FINSI' ;
  2269. 'SI' ('EGA' idir 2) ;
  2270. vnor1 = GNOR _cmt tdisc 'NPRI' discg 'FPRI' 1. 'NDUA' 'XN' ;
  2271. vnor1 = '*' vnor1 -1. ;
  2272. nvnor1 = '**' ('PSCAL' vnor1 vnor1 nomvit nomvit) 0.5 ;
  2273. nvnor1 = '+' nvnor1 ('MASQUE' nvnor1 'INFERIEUR' 1.D-100) ;
  2274. vnor1n = '/' vnor1 nvnor1 ;
  2275. 'FINSI' ;
  2276. * Cette formule peut poser problème en axi !
  2277. nvnor = '**' ('PSCAL' vnor vnor nomvit nomvit) 0.5 ;
  2278. * 'LISTE' ('MAXIMUM' nvnor) ;
  2279. * 'LISTE' ('MINIMUM' nvnor) ;
  2280. nvnor = '+' nvnor ('MASQUE' nvnor 'INFERIEUR' 1.D-100) ;
  2281. vnorn = '/' vnor nvnor ;
  2282. * Correction de vnorn aux points extrémités
  2283. pcmt = 'CHANGER' 'POI1' cmt ;
  2284. mcorr = pB 'ET' pC ;
  2285. vnorn2 = 'MANUEL' 'CHPO' pB ('EXTRAIRE' nomvit 1) 1. ;
  2286. vnorn3 = 'MANUEL' 'CHPO' pC ('EXTRAIRE' nomvit 2) 1. ;
  2287. pmcorr = 'CHANGER' 'POI1' mcorr ;
  2288. pcmtr = 'DIFF' pcmt pmcorr ;
  2289. vnorn1 = 'REDU' vnorn pcmtr ;
  2290. vnorn = vnorn1 '+' vnorn2 '+' vnorn3 ;
  2291. 'RESPRO' vnorn ;
  2292. 'FINPROC' ;
  2293. **
  2294. ************************************************************************
  2295. *
  2296. *
  2297. * END OF PROCEDURES
  2298. *
  2299. *
  2300. ************************************************************************
  2301. ************************************************************************
  2302. *
  2303. *
  2304. * MAIN : 1) MESH
  2305. * 2) COMPUTATIONAL LOOP
  2306. * 3) TESTs if interact=FAUX ;
  2307. * GUI if interact=VRAI ;
  2308. *
  2309. ************************************************************************
  2310. *
  2311. * Construction du "modèle" (maillage)
  2312. * et des paramètres de départ
  2313. *
  2314. idisc = 'QUAF' ;
  2315. 'SI' complet ;
  2316. nx = 25 ;
  2317. 'SINON' ;
  2318. nx = 5 ;
  2319. 'FINSI' ;
  2320. 'SI' interact ;
  2321. critnewt = 2.D-3 ; nitermax = 25 ; omeg = 0.45 ;
  2322. 'SINON' ;
  2323. critnewt = 1.D-4 ; nitermax = 20 ; omeg = 1.0 ;
  2324. 'FINSI' ;
  2325. methgau = 'GAU7' ;
  2326. idir = 1 ;
  2327. jacoxf = 3 ; jacoxg = 0 ; jacoxv = 0 ;
  2328. jacoxt = 'LECT' 1 2 3 4 ;
  2329. *
  2330. * Création du maillage
  2331. *
  2332. pA = 0. 0. ; pB = 1. 0. ; pC = 0. 1. ;
  2333. bas = 'DROIT' 1 pA pB ; sur = 'CERCLE' nx pB pA pC ;
  2334. gau = 'DROIT' 1 pC pA ;
  2335. cmt = bas 'ET' sur 'ET' gau ;
  2336. tol = 1.D-5 ;
  2337. _bas _sur _gau _cmt = QUAFME bas sur gau cmt ;
  2338. 'ELIMINATION' (_bas 'ET' _sur 'ET' _gau 'ET' _cmt) 1.D-5 ;
  2339. 'SI' ('EGA' idisc 'QUAF') ; cmt = _cmt ; 'FINSI' ;
  2340. *
  2341. ************************************************************************
  2342. *
  2343. * COMPUTATIONAL LOOP
  2344. *
  2345. ************************************************************************
  2346. *
  2347. * Bo : nombre de Bond (gravité / tension de surface)
  2348. * ang : angle de la gravité par rapport à la verticale
  2349. *
  2350. * Structure de la table TCLIM pour les conditions aux limites
  2351. *
  2352. * En ENTREE :
  2353. * cvol = vrai => Contrainte sur le volume
  2354. * volv : volume cible
  2355. * cvol = faux
  2356. * dpv : différence de pression voulue
  2357. * blocb = VRAI => Contrainte sur la position du point pB
  2358. * rbv : rayon voulu du point B
  2359. * blocb = FAUX
  2360. * abv : angle de contact voulu au point B
  2361. * blocc = VRAI => Contrainte sur la position du point pC
  2362. * rcv : rayon voulu du point C
  2363. * blocc = FAUX
  2364. * acv : angle de contact voulu au point C
  2365. * En SORTIE :
  2366. * les mêmes indices sans le v final indiquent
  2367. * les quantités trouvées par le calcul
  2368. *
  2369. *
  2370. 'DEBPROC' calcul ;
  2371. 'ARGUMENT' Bo*'FLOTTANT' ;
  2372. 'ARGUMENT' ang*'FLOTTANT' ;
  2373. 'ARGUMENT' tclim*'TABLE' ;
  2374. *
  2375. vdim = 'VALEUR' 'DIME' ;
  2376. vmod = 'VALEUR' 'MODE' ;
  2377. *
  2378. * Création du "modèle"
  2379. *
  2380. 'SI' ('EGA' vdim 2) ;
  2381. 'SI' ('NEG' vmod 'AXIS') ;
  2382. nomvit = 'MOTS' 'UX' 'UY' ;
  2383. nomfor = 'MOTS' 'FX' 'FY' ;
  2384. 'SINON' ;
  2385. nomvit = 'MOTS' 'UR' 'UZ' ;
  2386. nomfor = 'MOTS' 'FR' 'FZ' ;
  2387. 'FINSI' ;
  2388. 'SINON' ;
  2389. nomvit = 'MOTS' 'UX' 'UY' 'UZ' ;
  2390. nomfor = 'MOTS' 'FX' 'FY' 'FZ' ;
  2391. 'FINSI' ;
  2392. *
  2393. nomsca = 'MOTS' 'SCAL' ;
  2394. *
  2395. TDISC = 'TABLE' ;
  2396. TDISC . 'GEOM' = 'TABLE' ;
  2397. TDISC . 'GEOM' . 'DISC' = idisc ;
  2398. TDISC . 'XN' = 'TABLE' ;
  2399. TDISC . 'XN' . 'DISC' = TDISC . 'GEOM' . 'DISC' ;
  2400. TDISC . 'XN' . 'NOMINC' = 'TABLE' ;
  2401. 'REPETER' idim vdim ;
  2402. TDISC . 'XN' . 'NOMINC' . &idim = 'MOTS' ('EXTRAIRE' nomvit &idim) ;
  2403. 'FIN' idim ;
  2404. TDISC . 'FN' = 'TABLE' ;
  2405. TDISC . 'FN' . 'DISC' = TDISC . 'GEOM' . 'DISC' ;
  2406. TDISC . 'FN' . 'NOMINC' = 'TABLE' ;
  2407. 'REPETER' idim vdim ;
  2408. TDISC . 'FN' . 'NOMINC' . &idim = 'MOTS' ('EXTRAIRE' nomfor &idim) ;
  2409. 'FIN' idim ;
  2410. lmdisc = 'MOTS' 'CSTE' 'LINE' 'QUAF' ;
  2411. 'REPETER' iidisc ('DIME' lmdisc) ;
  2412. mdisc = 'EXTRAIRE' lmdisc &iidisc ;
  2413. tdisc . mdisc = 'TABLE' ;
  2414. tdisc . mdisc . 'DISC' = mdisc ;
  2415. tdisc . mdisc . 'NOMINC' = 'TABLE' ;
  2416. tdisc . mdisc . 'NOMINC' . 1 = nomsca ;
  2417. mdiscv = 'CHAINE' mdisc 'V' ;
  2418. tdisc . mdiscv = 'TABLE' ;
  2419. tdisc . mdiscv . 'DISC' = mdisc ;
  2420. tdisc . mdiscv . 'NOMINC' = 'TABLE' ;
  2421. 'REPETER' idim vdim ;
  2422. TDISC . mdiscv. 'NOMINC' . &idim =
  2423. 'MOTS' ('EXTRAIRE' nomvit &idim) ;
  2424. 'FIN' idim ;
  2425. 'FIN' iidisc ;
  2426. *
  2427. DISCG = TDISC . 'GEOM' . 'DISC' ;
  2428. DISCU = TDISC . 'XN' . 'DISC' ;
  2429. NOMVIT = @STBL (TDISC . 'XN' . 'NOMINC') ;
  2430. NOMFOR = @STBL (TDISC . 'FN' . 'NOMINC') ;
  2431. *
  2432. tdisc . 'cmt' = 'TABLE' ;
  2433. tdisc . 'cmt' .'QUAF' = _cmt ; tdisc . 'cmt' .'LINE' = cmt ;
  2434. tdisc . 'bas' = 'TABLE' ;
  2435. tdisc . 'bas' .'QUAF' = _bas ; tdisc . 'bas' .'LINE' = bas ;
  2436. tdisc . 'sur' = 'TABLE' ;
  2437. tdisc . 'sur' .'QUAF' = _sur ; tdisc . 'sur' .'LINE' = sur ;
  2438. tdisc . 'gau' = 'TABLE' ;
  2439. tdisc . 'gau' .'QUAF' = _gau ; tdisc . 'gau' .'LINE' = gau ;
  2440. *
  2441. cmt = tdisc . 'cmt' . discg ;
  2442. bas = tdisc . 'bas' . discg ;
  2443. gau = tdisc . 'gau' . discg ;
  2444. sur = tdisc . 'sur' . discg ;
  2445. *
  2446. *
  2447. lok = vrai ;
  2448. lquit = faux ;
  2449. *
  2450. * Boucle d'itérations (Newton)
  2451. *
  2452. tclim . 'fini' = 'FORME' ;
  2453. vol = GVOL _cmt tdisc faux ;
  2454. debug = faux ;
  2455. 'SI' debug ;
  2456. 'MESSAGE' ('CHAINE' 'Volume initial = ' (formar vol 2)) ;
  2457. 'FINSI' ;
  2458. 'REPETER' it nitermax ;
  2459. tabres = 'TABLE' ;
  2460. * Gravité
  2461. fpgra = GGRAVI _cmt tdisc 1. ang ;
  2462. ktgra = GKGRAVI _cmt tdisc jacoxg 1. ang ;
  2463. fpgra = '*' fpgra Bo ;
  2464. ktgra = '*' ktgra Bo ;
  2465. tabres = append tabres 'ftot' fpgra ;
  2466. tabres = append tabres 'ktot' ktgra ;
  2467. * Tension de surface
  2468. ftsur = TSURRESI _sur discg methgau nomvit -1. ;
  2469. ktsur = TSURKTAN _sur discg methgau nomvit nomvit +1.
  2470. jacoxt ;
  2471. tabres = append tabres 'ftot' ftsur ;
  2472. tabres = append tabres 'ktot' ktsur ;
  2473. * Contrainte éventuelle sur le volume
  2474. 'SI' (tclim . 'cvol') ;
  2475. volc = GVOL _cmt tdisc ;
  2476. dvol = ('-' (tclim . 'volv') volc) '*' -1. ;
  2477. * 'MESSAGE' ('CHAINE' 'volc=' volc) ;
  2478. 'SI' debug ;
  2479. 'MESSAGE' ('CHAINE' '-volv volc=' dvol) ;
  2480. 'FINSI' ;
  2481. kvol = GKVOL _cmt tdisc jacoxv ;
  2482. tabres = append tabres 'kcnt' kvol ;
  2483. tabres . 'fcnt' = dvol ;
  2484. 'SINON' ;
  2485. fpvol = GNOR _cmt tdisc 'NPRI' 'CSTE'
  2486. 'FPRI' ('*' (tclim . 'dpv') -1.)
  2487. 'NDUA' 'XN' ;
  2488. ktvol = GNORKTAN _cmt tdisc 'NPRI' 'XN' 'NCOF' 'CSTE'
  2489. 'FCOF' ('*' (tclim . 'dpv') +1.)
  2490. 'NDUA' 'XN' ;
  2491. tabres = append tabres 'ftot' fpvol ;
  2492. tabres = append tabres 'ktot' ktvol ;
  2493. 'FINSI' ;
  2494. * Conditions aux limites
  2495. 'SI' (tclim . 'blob') ;
  2496. rbc = 'COORDONNEE' 1 pB ;
  2497. dr = ('-' (tclim . 'rbv') rbc) '*' +1. ;
  2498. 'SI' debug ;
  2499. 'MESSAGE' ('CHAINE' '- rbv rbc=' dr) ;
  2500. 'FINSI' ;
  2501. ccl = 'MANUEL' 'CHPO' pB 1 'SCAL' dr ;
  2502. tabres = append tabres 'ccl' ccl ;
  2503. 'SINON' ;
  2504. * Force de bord
  2505. delta = tclim . 'abv' ;
  2506. sdel = 'SIN' delta ;
  2507. fborb = TSURRESI _bas discg methgau nomvit ('*' -1. sdel) ;
  2508. kborb = TSURKTAN _bas discg methgau nomvit nomvit ('*' +1. sdel) ;
  2509. tabres = append tabres 'ftot' fborb ;
  2510. tabres = append tabres 'ktot' kborb ;
  2511. 'FINSI' ;
  2512. * Conditions aux limites
  2513. 'SI' (tclim . 'bloc') ;
  2514. rcc = 'COORDONNEE' 2 pC ;
  2515. dr = ('-' (tclim . 'rcv') rcc) '*' +1. ;
  2516. 'SI' debug ;
  2517. 'MESSAGE' ('CHAINE' '- rcv rcc=' dr) ;
  2518. 'FINSI' ;
  2519. ccl = 'MANUEL' 'CHPO' pC 1 'SCAL' dr ;
  2520. append tabres 'ccl' ccl ;
  2521. 'SINON' ;
  2522. * Force de bord
  2523. delta = tclim . 'acv' ;
  2524. sdel = 'SIN' delta ;
  2525. fborc = TSURRESI _gau discg methgau nomvit ('*' -1. sdel) ;
  2526. kborc = TSURKTAN _gau discg methgau nomvit nomvit ('*' +1. sdel) ;
  2527. tabres = append tabres 'ftot' fborc ;
  2528. tabres = append tabres 'ktot' kborc ;
  2529. 'FINSI' ;
  2530. * Direction du déplacement des points de la surface
  2531. vnor = GNOR _sur tdisc 'NPRI' discg 'FPRI' 1. 'NDUA' 'XN' ;
  2532. vnorn = DIRDEP _cmt cmt sur tdisc idir ;
  2533. * trvec cmt vnorn 'Vnorn' fpgra 'Fpgra' ftsur 'Ftsur' ;
  2534. * fpvol 'Fpvol' ;
  2535. *
  2536. * Réduction du système sur l'inconnue déplacement normal
  2537. *
  2538. ktg = tabres . 'ktot' ;
  2539. *
  2540. * Convergence avec backtracking repris sur l'algorithme de
  2541. * DEDU ADAP pas parfait...
  2542. *
  2543. fback = 2. ; fvdet = 6.D0 ; nback = 10 ; damp = omeg ;
  2544. det0 = DEADJACO _sur discg methgau ;
  2545. 'OPTI' 'MODE' plan ;
  2546. mtail = GMASS2 _sur tdisc 'NPRI' discg 'NDUA' discg ;
  2547. 'OPTI' 'MODE' vmod ;
  2548. tail = 'KOPS' 'EXTRDIAG' mtail ;
  2549. suri = 'DIFF' ('CHANGER' sur 'POI1') (pB 'ET' pC) ;
  2550. taili = 'REDU' tail suri ;
  2551. dampi = damp ;
  2552. deti = det0 ;
  2553. backok = FAUX ;
  2554. 'REPETER' iback nback ;
  2555. 'SI' ('>' &iback 1) ; dampi = '/' dampi fback ; 'FINSI' ;
  2556. 'SI' debug ;
  2557. ch = 'CHAINE' ' dampi=' dampi ;
  2558. 'MESSAGE' ch ;
  2559. 'FINSI' ;
  2560. ktg2 = '*' ktg ('/' 1. dampi) ;
  2561. * ktg2 = ktg ;
  2562. fpg = tabres . 'ftot' ;
  2563. lcnt = 'EXISTE' tabres 'kcnt' ;
  2564. 'SI' lcnt ;
  2565. kvol = tabres . 'kcnt' ;
  2566. dvol = tabres . 'fcnt' ;
  2567. ktotp ftotp = PROJSYSC tdisc vnorn ktg2 fpg kvol dvol ;
  2568. 'SINON' ;
  2569. ktotp ftotp = PROJSYSC tdisc vnorn ktg2 fpg ;
  2570. 'FINSI' ;
  2571. domb = 'DIFF' ('CHANGER' 'POI1' cmt)
  2572. ('CHANGER' 'POI1' sur) ;
  2573. ccl = 'MANUEL' 'CHPO' domb 1 'SCAL' 0. ;
  2574. tabres = append tabres 'ccl' ccl ;
  2575. 'SI' ('EGA' ('TYPE' ktotp) 'RIGIDITE') ;
  2576. ktotpr = 'KOPS' 'NINCPRDU' ktotp ;
  2577. ftotpr = 'KOPS' 'NINCPRDU' ftotp ;
  2578. 'SINON' ;
  2579. ktotpr = ktotp ;
  2580. ftotpr = ftotp ;
  2581. 'FINSI' ;
  2582. * 'LISTE' (tabres . 'ccl') ;
  2583. sol = 'KRES' ktotpr ftotpr 'CLIM' (tabres . 'ccl') ;
  2584. * 'LDEPE' FAUX ;
  2585. * sol = '*' sol dampi ;
  2586. dep = '*' ('EXCO' 'SCAL' sol) vnorn ;
  2587. deps = 'REDU' dep sur ;
  2588. oldconf = 'FORME' ;
  2589. 'FORME' deps ;
  2590. depsi = 'REDU' dep suri ;
  2591. * depsi = deps ;
  2592. ndepsi = '**' ('PSCAL' depsi depsi nomvit nomvit) 0.5 ;
  2593. ndr = '/' ndepsi taili ;
  2594. mndr = 'MAXIMUM' ndr ;
  2595. xs ys = 'COORDONNEE' sur ;
  2596. detip = DEADJACO _sur discg methgau ;
  2597. 'FORME' oldconf ;
  2598. tyde = 'TYPE' detip ;
  2599. 'SI' ('EGA' tyde 'ENTIER') ;
  2600. 'SI' debug ;
  2601. ch = 'CHAINE' ' Warning : inv. loc. jacobien !' ;
  2602. 'MESSAGE' ch ;
  2603. 'FINSI' ;
  2604. 'SINON' ;
  2605. vardet = ('/' detip deti) ;
  2606. mivd = 'MINIMUM' vardet ;
  2607. mavd = 'MAXIMUM' vardet ;
  2608. mixs = 'MINIMUM' xs ;
  2609. miys = 'MINIMUM' ys ;
  2610. 'SI' debug ;
  2611. 'MESSAGE' ('CHAINE' 'Maxi dep rel = ' mndr ) ;
  2612. 'MESSAGE' ('CHAINE' 'Mini var jaco = ' mivd ) ;
  2613. 'MESSAGE' ('CHAINE' 'Maxi var jaco = ' mavd ) ;
  2614. 'MESSAGE' ('CHAINE' 'Mini xsurf = ' mixs ) ;
  2615. 'MESSAGE' ('CHAINE' 'Mini ysurf = ' miys ) ;
  2616. 'FINSI' ;
  2617. bigvar = 'OU' ('>' mavd fvdet) ('<' mivd ('/' 1.D0 fvdet))
  2618. 'OU' ('>' mndr ('-' fvdet 1.)) 'OU' ('&lt;EG' mixs -1.D-5)
  2619. 'OU' ('&lt;EG' miys -1.D-5) ;
  2620. 'SI' bigvar ;
  2621. 'SI' debug ;
  2622. ch = 'CHAINE'
  2623. ' Warn : trop grande variation du jaco !' ;
  2624. 'MESSAGE' ch ;
  2625. 'FINSI' ;
  2626. 'SINON' ;
  2627. backok = VRAI ;
  2628. 'QUITTER' iback ;
  2629. 'FINSI' ;
  2630. 'FINSI' ;
  2631. 'FIN' iback ;
  2632. *
  2633. * On peut avoir un problème à la première itération si
  2634. * les contraintes ne sont pas vérifiées au départ
  2635. *
  2636. 'SI' ('ET' ('NON' backok) ('NEG' &it 1)) ;
  2637. chinfo1 = 'CHAINE'
  2638. 'Backtracking failed to converge !' ;
  2639. 'MESSAGE' chinfo1 ;
  2640. chinfo2 = 'CHAINE' 'Please check the output displacement' ;
  2641. 'MESSAGE' chinfo2 ;
  2642. * 'SI' ('NEG' &it 1) ;
  2643. * lok = faux ;
  2644. * lquit = vrai ;
  2645. * 'FINSI' ;
  2646. 'FINSI' ;
  2647. *
  2648. freac = (ktotpr '*' sol) '-' ftotpr ;
  2649. odep = dep ;
  2650. *
  2651. * Bilan des forces
  2652. *
  2653. fpgrann = '*' ('PSCAL' fpgra vnorn nomvit nomvit) vnorn ;
  2654. ftsurnn = '*' ('PSCAL' ftsur vnorn nomvit nomvit) vnorn ;
  2655. 'SI' (tclim . 'cvol') ;
  2656. mulag = 'EXCO' 'LX' sol 'LX' ;
  2657. fpvolnn = '*' ('*' ktotp ('*' mulag -1.)) vnorn ;
  2658. 'SINON' ;
  2659. fpvolnn = '*' ('PSCAL' fpvol vnorn nomvit nomvit) vnorn ;
  2660. 'FINSI' ;
  2661. 'SI' (tclim . 'blob') ;
  2662. fborbnn = '*' ('REDU' freac pB) vnorn ;
  2663. 'SINON' ;
  2664. fborbnn = '*' ('PSCAL' fborb vnorn nomvit nomvit) vnorn ;
  2665. 'FINSI' ;
  2666. 'SI' (tclim . 'bloc') ;
  2667. fborcnn = '*' ('REDU' freac pC) vnorn ;
  2668. 'SINON' ;
  2669. fborcnn = '*' ('PSCAL' fborc vnorn nomvit nomvit) vnorn ;
  2670. 'FINSI' ;
  2671. * Bilan des forces normales en surface
  2672. desfnn = fpvolnn '+' fpgrann '+' ftsurnn '+' fborbnn '+' fborcnn ;
  2673. fpvolnn = 'REDU' fpvolnn sur ;
  2674. fpgrann = 'REDU' fpgrann sur ;
  2675. ftsurnn = 'REDU' ftsurnn sur ;
  2676. desfnn = 'REDU' desfnn sur ;
  2677. *
  2678. 'SI' (tclim . 'cvol') ;
  2679. mulag = 'EXCO' 'LX' sol 'LX' ;
  2680. tclim . 'dp' = 'MAXIMUM' mulag ;
  2681. 'SINON' ;
  2682. tclim . 'dp' = tclim . 'dpv' ;
  2683. 'FINSI' ;
  2684. *
  2685. maxdepr = 'MAXIMUM' dep 'ABS' ;
  2686. 'SI' debug ;
  2687. 'MESSAGE' ('CHAINE' 'Maxdepr=' maxdepr) ;
  2688. 'FINSI' ;
  2689. fbornn = fborcnn '+' fborbnn ;
  2690. 'SI' graph ;
  2691. * TRVEC sur fpgrann 'gran' fpvolnn 'voln'
  2692. * ftsurnn 'tsun' fbornn 'born' desfnn 'desn' VRAI ;
  2693. TRVEC cmt dep ('CHAINE' 'Deplacement') 1. VRAI ;
  2694. 'FINSI' ;
  2695. 'SI' ('<' maxdepr critnewt) ;
  2696. lquit = vrai ;
  2697. 'FINSI' ;
  2698. * 'SI' ('>' maxdepr ('+' ('ABS' dx) 0.5)) ;
  2699. * 'SI' ('>' maxdepr 6.) ;
  2700. * lok = faux ;
  2701. * lquit = vrai ;
  2702. * 'FINSI' ;
  2703. * 'OPTION' 'DONN' 5 ;
  2704. * Extension du déplacement sur les points de l'axe et du haut
  2705. * pour les quadratiques !
  2706. *
  2707. 'SI' lok ;
  2708. 'SI' ('EGA' idisc 'QUAF') ;
  2709. 'OPTI' 'MODE' 'PLAN' ;
  2710. nomvit2 = 'MOTS' 'UX' 'UY' ;
  2711. nv21 = 'EXTRAIRE' nomvit2 1 ;
  2712. nv22 = 'EXTRAIRE' nomvit2 2 ;
  2713. nv1 = 'EXTRAIRE' nomvit 1 ;
  2714. nv2 = 'EXTRAIRE' nomvit 2 ;
  2715. bux = 'BLOQUE' nv21 gau ;
  2716. buy = 'BLOQUE' nv22 bas ;
  2717. bs = 'BLOQUE' 'DEPL' sur ;
  2718. fs = 'DEPIMPOSE' bs ('EXCO' nomvit odep nomvit2) ;
  2719. btot = bux 'ET' buy 'ET' bs ;
  2720. ftot = fs ;
  2721. dx = 'DEDU' adap cmt btot ftot ;
  2722. 'OPTI' 'MODE' vmod ;
  2723. dx = 'NOMC' nomvit2 dx nomvit 'NATURE' 'DIFFUS' ;
  2724. odep = 'CHANGER' 'ATTRIBUT' ('REDU' odep sur)
  2725. 'NATURE' 'DIFFUS' ;
  2726. * 'LISTE' odep ;
  2727. * 'LISTE' dx ;
  2728. odep = odep 'ET' dx ;
  2729. 'FINSI' ;
  2730. 'FORME' odep ;
  2731. * Eventuelle régularisation ?
  2732. 'SI' faux ;
  2733. 'OPTI' 'MODE' 'PLAN' ;
  2734. nomvit2 = 'MOTS' 'UX' 'UY' ;
  2735. vnor = GNOR _sur tdisc 'NPRI' discg 'FPRI' +1. 'NDUA' 'XN' ;
  2736. vnor = 'NOMC' nomvit nomvit2 vnor ;
  2737. b1 = 'BLOQUE' 'DEPL' 'DIRECTION' vnor sur ;
  2738. b2 = 'BLOQUE' 'DEPL' (pB 'ET' pC) ;
  2739. dx = deduad2 sur (b1 'ET' b2) 'NITM' 1 ;
  2740. 'OPTI' 'MODE' vmod ;
  2741. dx = 'NOMC' nomvit2 dx nomvit 'NATURE' 'DIFFUS' ;
  2742. depsi = 'REDU' dx suri ;
  2743. ndepsi = '**' ('PSCAL' depsi depsi nomvit nomvit) 0.5 ;
  2744. ndr = '/' ndepsi taili ;
  2745. mndr = 'MAXIMUM' ndr ;
  2746. 'MESSAGE' ('CHAINE' 'Maxi dep reg rel = ' mndr ) ;
  2747. 'SI' ('>' mndr 0.5) ;
  2748. cof = '/' 0.5 mndr ;
  2749. dx = '*' dx cof ;
  2750. 'FINSI' ;
  2751. * trvec sur dx 'dxreg' 1. ;
  2752. 'FORME' dx ;
  2753. 'FINSI' ;
  2754. 'FINSI' ;
  2755. * La normale intégrée est nulle sur l'axe avec des quadratiques...
  2756. vmod = 'VALEUR' 'MODE' ;
  2757. 'OPTI' 'MODE' 'PLAN' ;
  2758. vnor = GNOR _sur tdisc 'NPRI' discg 'FPRI' +1. 'NDUA' 'XN' ;
  2759. 'OPTI' 'MODE' vmod ;
  2760. * trvec cmt vnor 'VNor' ;
  2761. vnbx = 'EXTRAIRE' vnor ('EXTRAIRE' nomvit 1) pB ;
  2762. vnby = 'EXTRAIRE' vnor ('EXTRAIRE' nomvit 2) pB ;
  2763. vnb = vnbx vnby ;
  2764. fb = 'PVEC' vnb ;
  2765. fbux = 'COORDONNEE' 1 fb ;
  2766. fbuy = 'COORDONNEE' 2 fb ;
  2767. * fbux = 'EXTRAIRE' ftsur ('EXTRAIRE' nomvit 1) pB ;
  2768. * fbuy = 'EXTRAIRE' ftsur ('EXTRAIRE' nomvit 2) pB ;
  2769. fb = fbux fbuy ;
  2770. nfb = '**' ('PSCA' fb fb) 0.5 ;
  2771. * 'LISTE' fb ; 'LISTE' nfb ;
  2772. tclim . 'ab' = 'ASIN' ('/' ('*' fbux -1.) nfb) ;
  2773. * tclim . 'ab' = 'ATG' fbux fbuy ;
  2774. vncx = 'EXTRAIRE' vnor ('EXTRAIRE' nomvit 1) pC ;
  2775. vncy = 'EXTRAIRE' vnor ('EXTRAIRE' nomvit 2) pC ;
  2776. vnc = vncx vncy ;
  2777. fc = 'PVEC' vnc ;
  2778. fcux = 'COORDONNEE' 1 fc ;
  2779. fcuy = 'COORDONNEE' 2 fc ;
  2780. * fcux = 'EXTRAIRE' ftsur ('EXTRAIRE' nomvit 1) pC ;
  2781. * fcuy = 'EXTRAIRE' ftsur ('EXTRAIRE' nomvit 2) pC ;
  2782. fc = fcux fcuy ;
  2783. nfc = '**' ('PSCA' fc fc) 0.5 ;
  2784. * 'LISTE' fc ; liste nfc ;
  2785. tclim . 'ac' = 'ASIN' ('/' ('*' fcuy +1.) nfc) ;
  2786. * tclim . 'ac' = 'ATG' fcuy fcux ;
  2787. tclim . 'vol' = GVOL _cmt tdisc faux ;
  2788. tclim . 'rb' = 'COORDONNEE' 1 pB ;
  2789. tclim . 'rc' = 'COORDONNEE' 2 pC ;
  2790. * 'LISTE' tclim ;
  2791. * 'OPTION' 'DONN' 5 ;
  2792. 'SI' lquit ; 'QUITTER' it ; 'FINSI' ;
  2793. 'FIN' it ;
  2794. 'SI' ('NON' lquit) ;
  2795. lok = faux ;
  2796. 'FINSI' ;
  2797. 'SI' (graph 'ET' faux) ;
  2798. * TRVEC cmt dep ('CHAINE' 'Depl omeg=' (formar omeg 1)) 1. ;
  2799. fbornn = fborcnn '+' fborbnn ;
  2800.  
  2801. TRVEC sur fpgrann 'gran' fpvolnn 'voln'
  2802. ftsurnn 'tsun' fbornn 'born' desfnn 'desn' ;
  2803. fpvol = GNOR _cmt tdisc 'NPRI' 'CSTE'
  2804. 'FPRI' ('*' (tclim . 'dp') -1.)
  2805. 'NDUA' 'XN' ;
  2806. desf = fpvol '+' fpgra '+' ftsur ;
  2807. suri = 'DIFF' ('CHANGER' sur 'POI1') (pB 'ET' pC) ;
  2808. fpvoli = 'REDU' fpvol suri ;
  2809. fpgrai = 'REDU' fpgra suri ;
  2810. ftsuri = 'REDU' ftsur suri ;
  2811. desfi = 'REDU' desf suri ;
  2812. TRVEC sur fpgrai 'grai' fpvoli 'voli'
  2813. ftsuri 'tsui' desfi 'desi' ;
  2814. TRVEC sur desfi 'desi' ;
  2815. * 'OPTION' 'DONN' 5 ;
  2816. * Bilan des forces en surface
  2817. * fpgrat = '-' fpgra fpgrann ;
  2818. * ftsurt = '-' ftsur ftsurnn ;
  2819. * desft = fpgrat '+' ftsurt ;
  2820. * psurm = 'DIFF' ('CHANGER' 'POI1' sur) (pB 'ET' pC) ;
  2821. * fpgrat = 'REDU' fpgrat psurm ;
  2822. * ftsurt = 'REDU' ftsurt psurm ;
  2823. * desft = 'REDU' desft psurm ;
  2824. * TRVEC sur fpgrat 'grat'
  2825. * ftsurt 'tsut' desft 'dest' ;
  2826. * 'OPTION' 'DONN' 5 ;
  2827. 'FINSI' ;
  2828.  
  2829. 'RESPRO' lok ;
  2830. 'FINPROC' ;
  2831. ************************************************************************
  2832. *
  2833. * END OF COMPUTATIONAL LOOP
  2834. *
  2835. ************************************************************************
  2836. *
  2837.  
  2838. 'SI' interact ;
  2839. ************************************************************************
  2840. *
  2841. * GUI PART
  2842. *
  2843. ************************************************************************
  2844. tclim = 'TABLE' ;
  2845. *
  2846. 'SI' ('EGA' ('VALEUR' 'MODE') 'AXIS') ;
  2847. tclim . 'cvol' = vrai ; tclim . 'volv' = '*' PI ('/' 2. 3.) ;
  2848. 'SINON' ;
  2849. tclim . 'cvol' = vrai ; tclim . 'volv' = '/' PI 4. ;
  2850. 'FINSI' ;
  2851. tclim . 'blob' = faux ; tclim . 'abv' = 0. ;
  2852. tclim . 'bloc' = faux ; tclim . 'acv' = 0. ;
  2853. *
  2854. Bo = 1.D-6 ; ang = 0. ;
  2855. lok = CALCUL Bo ang tclim ;
  2856. *'LISTE' tclim ;
  2857. 'SI' ('NON' lok) ;
  2858. 'MESSAGE' 'Pb: pas detat initial!' ;
  2859. 'FINSI' ;
  2860. 'REPETER' bouc1 ;
  2861. laxi = 'EGA' ('VALEUR' 'MODE') 'AXIS' ;
  2862. 'SI' laxi ; tit = 'CHAINE' 'Axi' ;
  2863. 'SINON' ; tit = 'CHAINE' 'Plane' ; 'FINSI' ;
  2864. tit = 'CHAINE' tit ' Bo=' (formar Bo 2) ;
  2865. 'SI' ('NON' laxi) ;
  2866. tit = 'CHAINE' tit ' angg=' (formar ang 2) ; 'FINSI' ;
  2867. 'SI' (tclim . 'cvol') ;
  2868. tit = 'CHAINE' tit
  2869. ('CHAINE' ' vol=' (formar (tclim . 'volv') 2)) ;
  2870. 'SINON' ;
  2871. tit = 'CHAINE' tit
  2872. ('CHAINE' ' dpv=' (formar (tclim . 'dpv') 2)) ;
  2873. 'FINSI' ;
  2874. 'SI' (tclim . 'blob') ;
  2875. tit = 'CHAINE' tit
  2876. ('CHAINE' ' rbv=' (formar (tclim . 'rbv') 2)) ;
  2877. 'SINON' ;
  2878. tit = 'CHAINE' tit
  2879. ('CHAINE' ' abv=' (formar (tclim . 'abv') 2)) ;
  2880. 'FINSI' ;
  2881. 'SI' (tclim . 'bloc') ;
  2882. tit = 'CHAINE' tit
  2883. ('CHAINE' ' rcv=' (formar (tclim . 'rcv') 2)) ;
  2884. 'SINON' ;
  2885. tit = 'CHAINE' tit
  2886. ('CHAINE' ' acv=' (formar (tclim . 'acv') 2)) ;
  2887. 'FINSI' ;
  2888. 'TRACER' cmt 'TITR' tit 'NCLK' ;
  2889. omod = 'VALEUR' 'MODE' ;
  2890. oBo = Bo ;
  2891. oang = ang ;
  2892. otclim = coptclim tclim ;
  2893. ofor = 'FORME' ;
  2894. *
  2895. * Menu items
  2896. *
  2897. tmenu = 'TABLE' ; imenu = 0 ;
  2898. tquest = 'TABLE' ;
  2899. imenu = '+' imenu 1 ;
  2900. tmenu . imenu = 'Bond' ;
  2901. tquest . imenu = 'CHAINE' 'Bond number (gravity/surface tension) ? ('
  2902. (formar Bo 2) ')' ;
  2903. 'SI' ('NEG' ('VALEUR' 'MODE') 'AXIS') ;
  2904. imenu = '+' imenu 1 ; tmenu . imenu = 'Angg' ;
  2905. tquest . imenu = 'CHAINE'
  2906. 'Gravity direction/downward in degrees ? ('
  2907. (formar ang 2) ')' ;
  2908. 'FINSI' ;
  2909. imenu = '+' imenu 1 ;
  2910. 'SI' (tclim . 'cvol') ;
  2911. tmenu . imenu = 'Volume' ;
  2912. tquest . imenu = 'CHAINE'
  2913. 'Prescribed volume ? (' (formar (tclim . 'volv') 2) ')' ;
  2914. 'SINON' ;
  2915. tmenu . imenu = 'DeltaP' ;
  2916. tquest . imenu = 'CHAINE'
  2917. 'Prescribed pressure difference ? ('
  2918. (formar (tclim . 'dpv') 2) ')' ;
  2919. 'FINSI' ;
  2920. imenu = '+' imenu 1 ;
  2921. 'SI' (tclim . 'blob') ;
  2922. tmenu . imenu = 'RadB' ;
  2923. tquest . imenu = 'CHAINE'
  2924. 'Prescribed radius at B (lower right) ? ('
  2925. (formar (tclim . 'rbv') 2) ')' ;
  2926. 'SINON' ;
  2927. tmenu . imenu = 'AngB' ;
  2928. tquest . imenu = 'CHAINE'
  2929. 'Prescribed angle at B (lower right) in degrees ? ('
  2930. (formar (tclim . 'abv') 2) ')' ;
  2931. 'FINSI' ;
  2932. 'SI' ('NEG' ('VALEUR' 'MODE') 'AXIS') ;
  2933. imenu = '+' imenu 1 ;
  2934. 'SI' (tclim . 'bloc') ;
  2935. tmenu . imenu = 'RadC' ;
  2936. tquest . imenu = 'CHAINE'
  2937. 'Prescribed radius at C (upper left) ? ('
  2938. (formar (tclim . 'rcv') 2) ')' ;
  2939. 'SINON' ;
  2940. tmenu . imenu = 'AngC' ;
  2941. tquest . imenu = 'CHAINE'
  2942. 'Prescribed angle at C (upper left) in degrees ? ('
  2943. (formar (tclim . 'acv') 2) ')' ;
  2944. 'FINSI' ;
  2945. 'FINSI' ;
  2946. imenu = '+' imenu 1 ; tmenu . imenu = 'Options' ;
  2947. *
  2948. * Menu display
  2949. *
  2950. 'SI' ('EGA' imenu 4) ;
  2951. ret = 'MENU' tit (tmenu . 1) (tmenu . 2) (tmenu . 3) (tmenu . 4) ;
  2952. 'SINON' ;
  2953. ret = 'MENU' tit (tmenu . 1) (tmenu . 2) (tmenu . 3) (tmenu . 4)
  2954. (tmenu . 5) (tmenu . 6) ;
  2955. 'FINSI' ;
  2956. *
  2957. * Menu actions
  2958. *
  2959. 'SI' ('ET' ('NEG' ret 'Quitter') ('NEG' ret 'Options')) ;
  2960. irep = 0 ;
  2961. 'REPETER' ii ('DIME' tmenu) ;
  2962. 'SI' ('EGA' ret (tmenu . &ii)) ;
  2963. irep = &ii ;
  2964. 'FINSI' ;
  2965. 'FIN' ii ;
  2966. 'SI' ('EGA' irep 0) ;
  2967. cherr = 'CHAINE' 'Option ' ret ' unknown' ;
  2968. 'ERREUR' cherr ;
  2969. 'FINSI' ;
  2970. val = sais ('CHAINE' (tquest . irep) ' : ') 'FLOTTANT' ;
  2971. 'SI' ('EGA' ret 'Bond') ; Bo = val ; 'FINSI' ;
  2972. 'SI' ('EGA' ret 'Angg') ; ang = val ; 'FINSI' ;
  2973. 'SI' ('EGA' ret 'Volume') ; tclim . 'volv' = 'ABS' val ; 'FINSI' ;
  2974. 'SI' ('EGA' ret 'DeltaP') ; tclim . 'dpv' = val ; 'FINSI' ;
  2975. 'SI' ('EGA' ret 'RadB') ; tclim . 'rbv' = 'ABS' val ; 'FINSI' ;
  2976. 'SI' ('EGA' ret 'AngB') ; tclim . 'abv' = val ; 'FINSI' ;
  2977. 'SI' ('EGA' ret 'RadC') ; tclim . 'rcv' = 'ABS' val ; 'FINSI' ;
  2978. 'SI' ('EGA' ret 'AngC') ; tclim . 'acv' = val ; 'FINSI' ;
  2979. 'FINSI' ;
  2980. 'SI' ('EGA' ret 'Quitter') ;
  2981. val = SAIS 'Are you sure you want to quit ? (Y/N)' 'MOT' ;
  2982. 'REPETER' bquit ;
  2983. 'SI' ('EGA' val 'Y') ; 'QUITTER' bouc1 ; 'FINSI' ;
  2984. 'SI' ('EGA' val 'N') ; 'QUITTER' bquit ; 'FINSI' ;
  2985. val = SAIS 'Please answer Yes or No : quit ? (Y/N)' 'MOT' ;
  2986. 'FIN' bquit ;
  2987. 'FINSI' ;
  2988. 'SI' ('EGA' ret 'Options') ;
  2989. axi = 'EGA' ('VALEUR' 'MODE') 'AXIS' ;
  2990. volume = tclim . 'cvol' ;
  2991. radiusb = tclim . 'blob' ;
  2992. radiusc = tclim . 'bloc' ;
  2993. titchoi = 'CHAINE' 'Check options : plane/axi vol./press. '
  2994. 'radius/angle@B radius/angle@C' ;
  2995. axi volume radiusb radiusc = 'CHOI' titchoi
  2996. axi volume radiusb radiusc ;
  2997. 'SI' axi ;
  2998. 'OPTI' 'MODE' 'AXIS' ;
  2999. ang = 0 ; radiusc = FAUX ; tclim . 'acv' = 0. ;
  3000. 'SINON' ;
  3001. 'OPTI' 'MODE' 'PLAN' ;
  3002. 'FINSI' ;
  3003. tclim . 'cvol' = volume ;
  3004. tclim . 'blob' = radiusb ;
  3005. tclim . 'bloc' = radiusc ;
  3006. 'SI' volume ;
  3007. 'SI' ('NON' ('EXISTE' tclim 'volv')) ;
  3008. tclim . 'volv' = tclim . 'vol' ;
  3009. 'OUBLIER' tclim 'dpv' ;
  3010. 'FINSI' ;
  3011. 'SINON' ;
  3012. 'SI' ('NON' ('EXISTE' tclim 'dpv')) ;
  3013. tclim . 'dpv' = tclim . 'dp' ;
  3014. 'OUBLIER' tclim 'volv' ;
  3015. 'FINSI' ;
  3016. 'FINSI' ;
  3017. 'SI' radiusb ;
  3018. 'SI' ('NON' ('EXISTE' tclim 'rbv')) ;
  3019. tclim . 'rbv' = tclim . 'rb' ;
  3020. 'OUBLIER' tclim 'abv' ;
  3021. 'FINSI' ;
  3022. 'SINON' ;
  3023. 'SI' ('NON' ('EXISTE' tclim 'abv')) ;
  3024. tclim . 'abv' = tclim . 'ab' ;
  3025. 'OUBLIER' tclim 'rbv' ;
  3026. 'FINSI' ;
  3027. 'FINSI' ;
  3028. 'SI' radiusc ;
  3029. 'SI' ('NON' ('EXISTE' tclim 'rcv')) ;
  3030. tclim . 'rcv' = tclim . 'rc' ;
  3031. 'OUBLIER' tclim 'acv' ;
  3032. 'FINSI' ;
  3033. 'SINON' ;
  3034. 'SI' ('NON' ('EXISTE' tclim 'acv')) ;
  3035. tclim . 'acv' = tclim . 'ac' ;
  3036. 'OUBLIER' tclim 'rcv' ;
  3037. 'FINSI' ;
  3038. 'FINSI' ;
  3039. 'FINSI' ;
  3040. *
  3041. * New shape
  3042. *
  3043. 'SI' ('NEG' ret 'Quitter') ;
  3044. lok = CALCUL Bo ang tclim ;
  3045. * 'LISTE' tclim ;
  3046. 'SI' ('NON' lok) ;
  3047. 'MESSAGE' 'Convergence error !!!' ;
  3048. 'OPTI' 'MODE' omod ;
  3049. Bo = oBo ;
  3050. ang = oang ;
  3051. tclim = otclim ;
  3052. 'FORME' ofor ;
  3053. 'FINSI' ;
  3054. 'FINSI' ;
  3055. 'FIN' bouc1 ;
  3056. ************************************************************************
  3057. *
  3058. * END OF GUI PART
  3059. *
  3060. ************************************************************************
  3061. 'FINSI' ;
  3062.  
  3063. 'SI' ('NON' interact) ;
  3064. ************************************************************************
  3065. *
  3066. * TEST PART
  3067. *
  3068. ************************************************************************
  3069. lpass = VRAI ;
  3070. *
  3071. 'OPTI' 'MODE' 'PLAN' ;
  3072. 'SAUTER' 2 'LIGNE' ; 'OPTI' 'ECHO' 1 ;
  3073. ***
  3074. *** Test 1 Plane circular drop, Prescribed contact angle = 90 degrees.
  3075. *** Prescribed volume = 1.
  3076. *** No gravity
  3077. *** Test on radius R and pressure difference \delta P
  3078. *** 2D Laplace-Young's law tells \delta P = \gamma / R.
  3079. *** Here \gamma = 1 (non dimensional)
  3080. 'OPTI' 'ECHO' 0 ;
  3081. *
  3082. * Computation
  3083. *
  3084. vol = 1. ;
  3085. tclim = 'TABLE' ;
  3086. tclim . 'cvol' = vrai ; tclim . 'volv' = vol ;
  3087. tclim . 'blob' = faux ; tclim . 'abv' = 0. ;
  3088. tclim . 'bloc' = faux ; tclim . 'acv' = 0. ;
  3089. Bo = 1.D-6 ; ang = 0. ;
  3090. lok = CALCUL Bo ang tclim ;
  3091. * list tclim ;
  3092. *
  3093. * Reference values
  3094. *
  3095. * V = \frac{\pi R^2}{4}
  3096. rref = '**' ('/' ('*' vol 4.) PI) 0.5 ;
  3097. dpref = '/' 1. rref ;
  3098. *
  3099. * Tests
  3100. *
  3101. rb = tclim . 'rb' ; rc = tclim . 'rc' ; dp = tclim . 'dp' ;
  3102. AFFVAR rref 'rref' rb 'rb' rc 'rc' ;
  3103. AFFVAR dpref 'dpref' dp 'dp' ;
  3104. errv = 2.D-4 ;
  3105. err1 = errrel rb rref ; tst1 = '<' err1 errv ;
  3106. err2 = errrel rc rref ; tst2 = '<' err2 errv ;
  3107. err3 = errrel dp dpref ; tst3 = '<' err3 errv ;
  3108. tst = lok 'ET' tst1 'ET' tst2 'ET' tst3 ;
  3109. 'MESSAGE' ('CHAINE' 'Test 1 : lok = ' lok) ;
  3110. 'MESSAGE' ' err1 = ' err1 ;
  3111. 'MESSAGE' ' err2 = ' err2 ;
  3112. 'MESSAGE' ' err3 = ' err3 ;
  3113. 'SI' tst ;
  3114. 'MESSAGE' 'Test 1 OK' ;
  3115. 'SINON' ;
  3116. 'MESSAGE' '!!! Test 1 not passed ' ;
  3117. 'FINSI' ;
  3118. *'OPTION' 'DONN' 5 ;
  3119. lpass = lpass 'ET' tst ;
  3120. 'SAUTER' 2 'LIGNE' ; 'OPTI' 'ECHO' 1 ;
  3121. ***
  3122. *** Test 2 Plane circular drop, Prescribed contact angle at B = 70 degrees.
  3123. *** Prescribed volume = 1.
  3124. *** No gravity
  3125. *** Test on position of B (xB, 0) and C (0, yC)
  3126. *** versus analytical solution
  3127. 'OPTI' 'ECHO' 0 ;
  3128. *
  3129. * Computation
  3130. *
  3131. vol = 1. ; alpha = 20. ; alphar = '*' alpha ('/' PI 180.) ;
  3132. tclim = 'TABLE' ;
  3133. tclim . 'cvol' = vrai ; tclim . 'volv' = vol ;
  3134. tclim . 'blob' = faux ; tclim . 'abv' = -1. '*' alpha ;
  3135. tclim . 'bloc' = faux ; tclim . 'acv' = 0. ;
  3136. Bo = 1.D-6 ; ang = 0. ;
  3137. lok = CALCUL Bo ang tclim ;
  3138. *'LISTE' tclim ;
  3139. *
  3140. * Reference values
  3141. *
  3142. * V = \frac{R^2}{4} (\pi - 2\alpha - sin 2\alpha) (\alpha in radians)
  3143. * xB = R \cos \alpha
  3144. * yC = R (1 - \sin \alpha)
  3145. dar = '*' alphar 2. ;
  3146. da = '*' alpha 2. ;
  3147. rref = '**' ('/' ('*' vol 4.) (PI '-' dar '-' ('SIN' da))) 0.5 ;
  3148. xbref = rref '*' ('COS' alpha) ;
  3149. ycref = rref '*' ('-' 1. ('SIN' alpha)) ;
  3150. *
  3151. * Tests
  3152. *
  3153. xb = tclim . 'rb' ; yc = tclim . 'rc' ;
  3154. AFFVAR xbref 'xbref' xb 'xb' ;
  3155. AFFVAR ycref 'ycref' yc 'yc' ;
  3156. errv = 2.D-4 ;
  3157. err1 = errrel xb xbref ; tst1 = '<' err1 errv ;
  3158. err2 = errrel yc ycref ; tst2 = '<' err2 errv ;
  3159. tst = lok 'ET' tst1 'ET' tst2 ;
  3160. 'MESSAGE' ('CHAINE' 'Test 2 : lok = ' lok) ;
  3161. 'MESSAGE' ' err1 = ' err1 ;
  3162. 'MESSAGE' ' err2 = ' err2 ;
  3163. 'SI' tst ;
  3164. 'MESSAGE' 'Test 2 OK' ;
  3165. 'SINON' ;
  3166. 'MESSAGE' '!!! Test 2 not passed ' ;
  3167. 'FINSI' ;
  3168. lpass = lpass 'ET' tst ;
  3169. 'SAUTER' 2 'LIGNE' ; 'OPTI' 'ECHO' 1 ;
  3170. ***
  3171. *** Test 3 Plane drop, Prescribed contact angle at B = 70 degrees.
  3172. *** Prescribed volume = 1.
  3173. *** Bond = -1. (upward gravity)
  3174. *** Test on position of C (0, yC)
  3175. *** versus numerical solution \cite{sumesh}, figure 2
  3176. 'OPTI' 'ECHO' 0 ;
  3177. *
  3178. * Computation
  3179. *
  3180. vol = 1. ; alpha = 20. ; alphar = '*' alpha ('/' PI 180.) ;
  3181. tclim = 'TABLE' ;
  3182. tclim . 'cvol' = vrai ; tclim . 'volv' = vol ;
  3183. tclim . 'blob' = faux ; tclim . 'abv' = -1. '*' alpha ;
  3184. tclim . 'bloc' = faux ; tclim . 'acv' = 0. ;
  3185. Bo = -1.D0 ; ang = 0. ;
  3186. lok = CALCUL Bo ang tclim ;
  3187. *'LISTE' tclim ;
  3188. *
  3189. * Reference values
  3190. *
  3191. ycref = 1.20 ;
  3192. *
  3193. * Tests
  3194. *
  3195. yc = tclim . 'rc' ;
  3196. AFFVAR ycref 'ycref' yc 'yc' ;
  3197. errv = 1.D-2 ;
  3198. err1 = errrel yc ycref ; tst1 = '<' err1 errv ;
  3199. tst = lok 'ET' tst1 ;
  3200. 'MESSAGE' ('CHAINE' 'Test 3 : lok = ' lok) ;
  3201. 'MESSAGE' ' err1 = ' err1 ;
  3202. 'SI' tst ;
  3203. 'MESSAGE' 'Test 3 OK' ;
  3204. 'SINON' ;
  3205. 'MESSAGE' '!!! Test 3 not passed ' ;
  3206. 'FINSI' ;
  3207. lpass = lpass 'ET' tst ;
  3208. 'OPTI' 'MODE' 'AXIS' ;
  3209. 'SAUTER' 2 'LIGNE' ; 'OPTI' 'ECHO' 1 ;
  3210. ***
  3211. *** Test 4 Axisymmetric spherical drop,
  3212. *** Prescribed contact angle = 90 degrees.
  3213. *** Prescribed pressure difference
  3214. *** such that volume= 3 PI / 2
  3215. *** No gravity
  3216. *** Test on radius at B and C and volume V
  3217. *** 3D Laplace-Young's law tells \delta P = 2 \gamma / R.
  3218. *** Here \gamma = 1 (non dimensional)
  3219. 'OPTI' 'ECHO' 0 ;
  3220. *
  3221. * Computation
  3222. *
  3223. volref = '/' ('*' PI 3.) 2. ;
  3224. *
  3225. * Reference values
  3226. *
  3227. * V = \frac{2 \pi R^3}{3}
  3228. rref = '**' ('/' ('*' volref 3.) ('*' 2. PI)) ('/' 1. 3.) ;
  3229. dpref = '/' 2. rref ;
  3230. *
  3231. tclim = 'TABLE' ;
  3232. tclim . 'cvol' = faux ; tclim . 'dpv' = dpref ;
  3233. tclim . 'blob' = faux ; tclim . 'abv' = 0. ;
  3234. tclim . 'bloc' = faux ; tclim . 'acv' = 0. ;
  3235. Bo = 1.D-6 ; ang = 0. ;
  3236. lok = CALCUL Bo ang tclim ;
  3237. *'LISTE' tclim ;
  3238. *
  3239. * Tests
  3240. *
  3241. rb = tclim . 'rb' ; rc = tclim . 'rc' ; vol = tclim . 'vol' ;
  3242. AFFVAR rref 'rref' rb 'rb' rc 'rc' ;
  3243. AFFVAR volref 'volref' vol 'vol' ;
  3244. errv = 2.D-4 ;
  3245. err1 = errrel rb rref ; tst1 = '<' err1 errv ;
  3246. err2 = errrel rc rref ; tst2 = '<' err2 errv ;
  3247. err3 = errrel vol volref ; tst3 = '<' err3 errv ;
  3248. tst = lok 'ET' tst1 'ET' tst2 'ET' tst3 ;
  3249. 'MESSAGE' ('CHAINE' 'Test 4 : lok = ' lok) ;
  3250. 'MESSAGE' ' err1 = ' err1 ;
  3251. 'MESSAGE' ' err2 = ' err2 ;
  3252. 'MESSAGE' ' err3 = ' err3 ;
  3253. 'SI' tst ;
  3254. 'MESSAGE' 'Test 4 OK' ;
  3255. 'SINON' ;
  3256. 'MESSAGE' '!!! Test 4 not passed ' ;
  3257. 'FINSI' ;
  3258. lpass = lpass 'ET' tst ;
  3259. 'SAUTER' 2 'LIGNE' ; 'OPTI' 'ECHO' 1 ;
  3260. ***
  3261. *** Test 5 Axisymmetric spherical drop,
  3262. *** Prescribed contact angle at B = 50 degrees.
  3263. *** Prescribed volume = 3 PI / 2
  3264. *** No gravity
  3265. *** Test on position of B (xB, 0) and C (0, yC)
  3266. *** versus analytical solution
  3267. 'OPTI' 'ECHO' 0 ;
  3268. *
  3269. * Computation
  3270. *
  3271. vol = '/' ('*' PI 3.) 2. ; alpha = 40. ;
  3272. tclim = 'TABLE' ;
  3273. tclim . 'cvol' = vrai ; tclim . 'volv' = vol ;
  3274. tclim . 'blob' = faux ; tclim . 'abv' = -1. '*' alpha ;
  3275. tclim . 'bloc' = faux ; tclim . 'acv' = 0. ;
  3276. Bo = 1.D-6 ; ang = 0. ;
  3277. lok = CALCUL Bo ang tclim ;
  3278. *'LISTE' tclim ;
  3279. *
  3280. * Reference values
  3281. *
  3282. * V = \pi \frac{R^3}{3} (1 - sin \alpha)^2 (2 + sin \alpha)
  3283. * rB = R \cos \alpha
  3284. * zC = R (1 - \sin \alpha)
  3285. rref = '**' ('/' ('*' vol 3.) (PI
  3286. '*' ('**' ('-' 1. ('SIN' alpha)) 2.)
  3287. '*' ('+' 2. ('SIN' alpha)))
  3288. ) ('/' 1. 3.) ;
  3289. rbref = rref '*' ('COS' alpha) ;
  3290. zcref = rref '*' ('-' 1. ('SIN' alpha)) ;
  3291. *
  3292. * Tests
  3293. *
  3294. rb = tclim . 'rb' ; zc = tclim . 'rc' ;
  3295. AFFVAR rbref 'rbref' rb 'rb' ;
  3296. AFFVAR zcref 'zcref' zc 'zc' ;
  3297. errv = 2.D-4 ;
  3298. err1 = errrel rb rbref ; tst1 = '<' err1 errv ;
  3299. err2 = errrel zc zcref ; tst2 = '<' err2 errv ;
  3300. tst = lok 'ET' tst1 'ET' tst2 ;
  3301. 'MESSAGE' ('CHAINE' 'Test 5 : lok = ' lok) ;
  3302. 'MESSAGE' ' err1 = ' err1 ;
  3303. 'MESSAGE' ' err2 = ' err2 ;
  3304. 'SI' tst ;
  3305. 'MESSAGE' 'Test 5 OK' ;
  3306. 'SINON' ;
  3307. 'MESSAGE' '!!! Test 5 not passed ' ;
  3308. 'FINSI' ;
  3309. lpass = lpass 'ET' tst ;
  3310. 'SAUTER' 2 'LIGNE' ; 'OPTI' 'ECHO' 1 ;
  3311. ***
  3312. *** Test 6 Axisymmetric drop, Prescribed contact angle at B = 50 degrees.
  3313. *** Prescribed volume = 3 PI / 2
  3314. *** Bond = -1. (upward gravity)
  3315. *** Test on position of C (0, zC)
  3316. *** versus numerical solution \cite{sumesh}, figure 8
  3317. 'OPTI' 'ECHO' 0 ;
  3318. *
  3319. * Computation
  3320. *
  3321. vol = '/' ('*' PI 3.) 2. ; alpha = 40. ;
  3322. tclim = 'TABLE' ;
  3323. tclim . 'cvol' = vrai ; tclim . 'volv' = vol ;
  3324. tclim . 'blob' = faux ; tclim . 'abv' = -1. '*' alpha ;
  3325. tclim . 'bloc' = faux ; tclim . 'acv' = 0. ;
  3326. Bo = -1.D0 ; ang = 0. ;
  3327. lok = CALCUL Bo ang tclim ;
  3328. *'LISTE' tclim ;
  3329. *
  3330. * Reference values
  3331. *
  3332. *zcref = 0.90 ;
  3333. zcref = 1.094 ;
  3334. *
  3335. * Tests
  3336. *
  3337. zc = tclim . 'rc' ;
  3338. AFFVAR zcref 'zcref' zc 'zc' ;
  3339. errv = 1.D-2 ;
  3340. err1 = errrel zc zcref ; tst1 = '<' err1 errv ;
  3341. tst = lok 'ET' tst1 ;
  3342. 'MESSAGE' ('CHAINE' 'Test 6 : lok = ' lok) ;
  3343. 'MESSAGE' ' err1 = ' err1 ;
  3344. 'SI' tst ;
  3345. 'MESSAGE' 'Test 6 OK' ;
  3346. 'SINON' ;
  3347. 'MESSAGE' '!!! Test 6 not passed ' ;
  3348. 'FINSI' ;
  3349. lpass = lpass 'ET' tst ;
  3350.  
  3351.  
  3352.  
  3353. 'SAUTER' 2 'LIGNE' ;
  3354. 'SI' lpass ;
  3355. 'MESSAGE' 'Tout sest bien passe' ;
  3356. 'SINON' ;
  3357. 'MESSAGE' 'Il y a eu des erreurs' ;
  3358. 'FINSI' ;
  3359. 'SAUTER' 2 'LIGNE' ;
  3360. 'SI' ('NON' lpass) ;
  3361. 'ERREUR' 5 ;
  3362. 'FINSI' ;
  3363. ************************************************************************
  3364. *
  3365. * END OF TEST PART
  3366. *
  3367. ************************************************************************
  3368. 'FINSI' ;
  3369. 'SI' interact ;
  3370. 'OPTION' 'DONN' 5 'ECHO' 1 ;
  3371. 'FINSI' ;
  3372. *
  3373. 'FIN' ;
  3374.  
  3375.  
  3376.  
  3377.  
  3378.  
  3379.  
  3380.  
  3381.  
  3382.  
  3383.  
  3384.  
  3385.  

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