Télécharger dynamic.procedur

Retour à la liste

Numérotation des lignes :

  1. * DYNAMIC PROCEDUR CB215821 17/11/14 21:15:11 9612
  2. * ======================================================================
  3. * CALCUL DE LA REPONSE DYNAMIQUE PAR INTEGRATION TEMPORELLE
  4. * (SCHEMA DE NEWMARK OU HHT OU ALPHA-GENERALISE)
  5. * ======================================================================
  6. * Creation : ???, 1988 ???
  7. * Refonte : Benoit Prabel, 23/07/2015
  8. * Modifications : JCARDO 15/02/2017 => post-traitement tableur + VTK
  9. * ======================================================================
  10. DEBP DYNAMIC ETAB*'TABLE' ;
  11.  
  12. IECHO = VALE 'ECHO' ;
  13.  
  14. * +-----------------------------------------------------------------+
  15. * | |
  16. * | L E C T U R E D E S P A R A M E T R E S |
  17. * | |
  18. * +-----------------------------------------------------------------+
  19.  
  20. * ================
  21. * OPTIONS INTERNES
  22. * ================
  23.  
  24. KMENAG = VRAI ;
  25.  
  26.  
  27. * ===================================================
  28. * PARAMETRES DE L'ALGORITHME D'INTEGRATION TEMPORELLE
  29. * ===================================================
  30.  
  31. SI ((EXIS ETAB 'ALPHA_F') ET (EXIS ETAB 'RHO_INF')) ;
  32. MESS 'Veuillez fournir 1 seul parametre parmi :' ;
  33. MESS ' ALPHA_F => algorithme HHT ' ;
  34. MESS ' RHO_INF => algorithme alpha-generalise' ;
  35. ERRE 21 ;
  36. FINS ;
  37.  
  38. * => Newmark acceleration moyenne (IALGO = 0)
  39. IALGO = 0 ;
  40. ALPHAM = 0. ;
  41. ALPHAF = 0. ;
  42.  
  43. * => HHT (IALGO = 1)
  44. SI (EXIS ETAB 'ALPHA_F') ;
  45. ALPHAF = ETAB.'ALPHA_F' ;
  46. SI (IECHO >EG 0) ;
  47. MESS 'ALGO = HHT (ALPHA_F=' ALPHAF ')' ;
  48. FINS ;
  49. IALGO = 1 ;
  50. FINS ;
  51.  
  52. * => alpha-generalise (IALGO = 2)
  53. SI (EXIS ETAB 'RHO_INF') ;
  54. RHOINF = ETAB . 'RHO_INF' ;
  55. ALPHAM = ((2.*RHOINF) - 1.) / (1.+RHOINF) ;
  56. ALPHAF = RHOINF / (1.+RHOINF) ;
  57. SI (IECHO >EG 0) ;
  58. MESS 'ALGO = alpha-generalise (ALPHA_M=' ALPHAM
  59. ' & ALPHA_F=' ALPHAF ')' ;
  60. FINS ;
  61. IALGO = 2 ;
  62. FINS ;
  63.  
  64.  
  65. * ====================
  66. * CONDITIONS INITIALES
  67. * ====================
  68.  
  69. U0 = ETAB.'DEPL' ;
  70. V0 = ETAB.'VITE' ;
  71.  
  72.  
  73. * =======================
  74. * MATRICES ET CHARGEMENTS
  75. * =======================
  76.  
  77. * Matrices et chargements constants
  78. K = ETAB.'RIGI' ;
  79. M = ETAB.'MASS' ;
  80. KAMOR = EXIS ETAB 'AMOR' ;
  81. SI KAMOR ;
  82. C = EXTR (ETAB.'AMOR') 'RIGI' 'NOMU' ;
  83. SI (IECHO >EG 0) ;
  84. MESS 'Presence d une matrice d amortissement' ;
  85. FINS ;
  86. FINS ;
  87. CHAR1 = ETAB.'CHAR' ;
  88.  
  89. * Matrices et chargements "utilisateur" via la procedure CHARMECA
  90. KCHARM = FAUX ;
  91. SI (EXIS ETAB 'CHARMECA') ;
  92. KCHARM = ETAB.'CHARMECA' ;
  93. FINS ;
  94.  
  95.  
  96. * ====================================
  97. * APPEL A VITEUNIL POUR LES CONTACTS ?
  98. * ====================================
  99.  
  100. SI (EXIS ETAB 'VITEUNIL') ;
  101. KVITUN = ETAB.'VITEUNIL' ;
  102. SINON ;
  103. KVITUN = (EGA IALGO 0) ;
  104. FINS ;
  105.  
  106. SI (KVITUN ET (NEG IALGO 0)) ;
  107. MESS 'VITEUNIL compatible uniquement avec le schema de Newmark '
  108. 'acceleration moyenne' ;
  109. ERRE 21 ;
  110. FINS ;
  111.  
  112.  
  113. * ====================
  114. * PARAMETRES TEMPORELS
  115. * ====================
  116.  
  117. * Nouvelle syntaxe :
  118. SI (EXIS ETAB 'TEMPS_CALCULES') ;
  119. TPROG = ETAB.'TEMPS_CALCULES' ;
  120.  
  121. SI ((MINI ((ENLE TPROG 1) - (ENLE TPROG (DIME TPROG)))) <EG 0.) ;
  122. MESS 'La liste des TEMPS_CALCULES doit etre strictement '
  123. 'croissante !' ;
  124. ERRE 21 ;
  125. FINS ;
  126.  
  127. * Ancienne syntaxe : on avertit et on traduit
  128. SINON ; SI (EXIS ETAB 'FREQ') ;
  129. MESS '!!! SYNTAXE BIENTOT OSBOLETE : CONSULTEZ LA NOTICE !!!' ;
  130. MESS '!!! (LE CALCUL SE POURSUIT MALGRE TOUT...) !!!' ;
  131. DTFREQ = 0.25 / ETAB.'FREQ' ;
  132. SI (EXIS ETAB 'DEBU') ;
  133. T0 = ETAB.'DEBU' ;
  134. SINON ;
  135. T0 = 0. ;
  136. FINS ;
  137.  
  138. TSORT = ORDO ETAB.'INST' ;
  139. ETAB.'TEMPS_SAUVES' = TSORT ;
  140.  
  141. TFIN = DTFREQ * (ENTI 'SUPERIEUR' (MAXI TSORT / DTFREQ)) ;
  142. TPROG = PROG T0 PAS DTFREQ TFIN ;
  143. ETAB.'TEMPS_CALCULES' = TPROG ;
  144.  
  145. SINON ;
  146. MESS 'Liste des TEMPS_CALCULES absente de la table !' ;
  147. ERRE 21 ;
  148. FINS ; FINS ;
  149.  
  150. NPAS = DIME TPROG ;
  151.  
  152.  
  153. * ========================
  154. * PARAMETRES DE SAUVEGARDE
  155. * ========================
  156.  
  157. IGIB = 0 ;
  158. ICSV = 0 ;
  159. IVTK = 0 ;
  160.  
  161.  
  162. * SAUVEGARDE VTK
  163. * --------------
  164.  
  165. * Instants de sauvegarde
  166. SI (EXIS ETAB 'PAS_SAUVES_VTK') ;
  167. SI (EGA (TYPE ETAB.'PAS_SAUVES_VTK') 'ENTIER') ;
  168. IPASVTK = ETAB.'PAS_SAUVES_VTK' ;
  169. SINON ;
  170. SI (EGA (ETAB.'PAS_SAUVES_VTK') 'TOUS') ;
  171. IPASVTK = 1 ;
  172. SINON ; SI (EGA (ETAB.'PAS_SAUVES_VTK') 'FINAL') ;
  173. IPASVTK = NPAS ;
  174. SINON ;
  175. MESS 'PAS_SAUVES_VTK est de type incorrect !' ;
  176. ERRE 21 ;
  177. FINS ; FINS ;
  178. FINS ;
  179. IVTK = 1 ;
  180. SINON ; SI (EXIS ETAB 'TEMPS_SAUVES_VTK') ;
  181. SI (NEG (TYPE ETAB.'TEMPS_SAUVES_VTK') 'LISTREEL') ;
  182. MESS 'TEMPS_SAUVES_VTK doit etre de type LISTREEL !' ;
  183. ERRE 21 ;
  184. FINS ;
  185. IVTK = 2 ;
  186. FINS ; FINS ;
  187.  
  188. * Emplacement de sauvegarde
  189. SI (EXIS ETAB 'FICHIER_VTK') ;
  190. FICVTK = ETAB.'FICHIER_VTK' ;
  191. SINON ;
  192. FICVTK = 'DYNAMIC' ;
  193. FINS ;
  194.  
  195. * Liste des maillages a sauvegarder
  196. SI (IVTK > 0) ;
  197. SI (EXIS ETAB 'MAILLAGE_VTK') ;
  198. MAIVTK = ETAB.'MAILLAGE_VTK' ;
  199. SINON ;
  200. MESS 'L indice MAILLAGE_VTK est obligatoire !' ;
  201. ERRE 21 ;
  202. FINS ;
  203. FINS ;
  204.  
  205.  
  206. * SAUVEGARDE CSV (TABLEUR)
  207. * ------------------------
  208.  
  209. * Instants de sauvegarde
  210. SI (EXIS ETAB 'PAS_SAUVES_CSV') ;
  211. SI (EGA (TYPE ETAB.'PAS_SAUVES_CSV') 'ENTIER') ;
  212. IPASCSV = ETAB.'PAS_SAUVES_CSV' ;
  213. SINON ;
  214. SI (EGA (ETAB.'PAS_SAUVES_CSV') 'TOUS') ;
  215. IPASCSV = 1 ;
  216. SINON ; SI (EGA (ETAB.'PAS_SAUVES_CSV') 'FINAL') ;
  217. IPASCSV = NPAS ;
  218. SINON ;
  219. MESS 'PAS_SAUVES_CSV est de type incorrect !' ;
  220. ERRE 21 ;
  221. FINS ; FINS ;
  222. FINS ;
  223. ICSV = 1 ;
  224. SINON ; SI (EXIS ETAB 'TEMPS_SAUVES_CSV') ;
  225. SI (NEG (TYPE ETAB.'TEMPS_SAUVES_CSV') 'LISTREEL') ;
  226. MESS 'TEMPS_SAUVES_CSV doit etre de type LISTREEL !' ;
  227. ERRE 21 ;
  228. FINS ;
  229. ICSV = 2 ;
  230. FINS ; FINS ;
  231.  
  232. * Emplacement de sauvegarde (eventuellement a la fin du calcul)
  233. SI (EXIS ETAB 'FICHIER_CSV') ;
  234. KFICCSV = VRAI ;
  235. FICCSV = ETAB.'FICHIER_CSV' ;
  236. SINON ;
  237. KFICCSV = FAUX ;
  238. FINS ;
  239.  
  240. * Liste des noeuds et composantes a sauvegarder
  241. SI (ICSV > 0) ;
  242. SI (EXIS ETAB 'MAILLAGE_CSV') ;
  243. MAICSV = CHAN 'POI1' ETAB.'MAILLAGE_CSV' ;
  244. NELCSV = NBEL MAICSV ;
  245. SI (NELCSV EGA 0) ;
  246. MESS 'Le maillage MAILLAGE_CSV est vide !' ;
  247. ERRE 21 ;
  248. FINS ;
  249. SINON ;
  250. MESS 'L indice MAILLAGE_CSV est obligatoire !' ;
  251. ERRE 21 ;
  252. FINS ;
  253.  
  254. KLISCO = FAUX ;
  255. SI (EXIS ETAB 'COMPOSANTES_CSV') ;
  256. LISCO = ETAB.'COMPOSANTES_CSV' ;
  257. KLISCO = VRAI ;
  258. SI ((DIME LISCO) EGA 0) ;
  259. MESS 'La liste COMPOSANTES_CSV est vide !' ;
  260. ERRE 21 ;
  261. FINS ;
  262. FINS ;
  263. FINS ;
  264.  
  265.  
  266. * SAUVEGARDE GIBIANE
  267. * ------------------
  268.  
  269. * Instants de sauvegarde
  270. SI (EXIS ETAB 'PAS_SAUVES') ;
  271. SI (EGA (TYPE ETAB.'PAS_SAUVES') 'ENTIER') ;
  272. IPASGIB = ETAB.'PAS_SAUVES' ;
  273. SINON ;
  274. SI (EGA (ETAB.'PAS_SAUVES') 'TOUS') ;
  275. IPASGIB = 1 ;
  276. SINON ; SI (EGA (ETAB.'PAS_SAUVES') 'FINAL') ;
  277. IPASGIB = NPAS ;
  278. SINON ;
  279. MESS 'PAS_SAUVES est de type incorrect !' ;
  280. ERRE 21 ;
  281. FINS ; FINS ;
  282. FINS ;
  283. IGIB = 1 ;
  284. SINON ; SI (EXIS ETAB 'TEMPS_SAUVES') ;
  285. SI (NEG (TYPE ETAB.'TEMPS_SAUVES') 'LISTREEL') ;
  286. MESS 'TEMPS_SAUVES doit etre de type LISTREEL !' ;
  287. ERRE 21 ;
  288. FINS ;
  289. IGIB = 2 ;
  290. SINON ;
  291. IPASGIB = 4 ;
  292. IGIB = 1 ;
  293. FINS ; FINS ;
  294.  
  295. * Sauvegarde incrementale ?
  296. * (appel a "SAUV MUET ETAB" => ecriture sur disque)
  297. KINCR = FAUX ;
  298. SI (EXIS ETAB 'SAUV') ;
  299. KINCR = ETAB.'SAUV' ;
  300. SI (KINCR ET (IECHO >EG 0)) ;
  301. MESS 'SAUVegarde incrementale activee' ;
  302. FINS ;
  303. FINS ;
  304.  
  305. * Fantomisation des resultats au fur et a mesure du calcul
  306. KECON = FAUX ;
  307. SI (EXIS ETAB 'ECON') ;
  308. KECON = ETAB.'ECON' ;
  309. FINS ;
  310.  
  311. * Restriction de la geometrie a sauvegarder
  312. KMAIGIB = EXIS ETAB 'MAILLAGE_SAUVE' ;
  313. SI KMAIGIB ;
  314. MAIGIB = ETAB.'MAILLAGE_SAUVE' ;
  315. FINS ;
  316.  
  317.  
  318.  
  319. * +-----------------------------------------------------------------+
  320. * | |
  321. * | I N I T I A L I S A T I O N S |
  322. * | |
  323. * +-----------------------------------------------------------------+
  324.  
  325. * ===================================================
  326. * PARAMETRES DE L'ALGORITHME D'INTEGRATION TEMPORELLE
  327. * ===================================================
  328.  
  329. ALPHA = ALPHAF - ALPHAM ;
  330. GAMMA = 0.5 + ALPHA ;
  331. BETA = 0.25 * ((1.+ALPHA)**2) ;
  332. 1MAF = 1. - ALPHAF ;
  333. 1MAM = 1. - ALPHAM ;
  334. GAMBET = GAMMA / BETA ;
  335. SI (IECHO >EG 0) ;
  336. MESS (CHAI 'Calcul avec GAMMA=' GAMMA ' BETA=' BETA
  337. ' (GAMMA/BETA=' GAMBET ')') ;
  338. FINS ;
  339.  
  340.  
  341. * ====================
  342. * PARAMETRES TEMPORELS
  343. * ====================
  344.  
  345. IPAS = 1 ;
  346. T0 = EXTR TPROG IPAS ;
  347. DT0 = 0. ;
  348. CHACHA = CHAI 'PAS #' IPAS ' TEMPS' T0 ;
  349.  
  350.  
  351. * ===================
  352. * SYSTEME D'EQUATIONS
  353. * ===================
  354.  
  355. OPER = VIDE 'RIGIDITE' ;
  356. F0 = TIRE CHAR1 T0 ;
  357.  
  358.  
  359. * =====================================================
  360. * ACCELERATION INITIALE A0 (NOUVEAUX SCHEMAS TEMPORELS)
  361. * =====================================================
  362.  
  363. SI (EXIS ETAB 'ACCE') ;
  364. A0 = ETAB.'ACCE' ;
  365. SINON ;
  366. K0 = K ;
  367. F2ND = (F0 ENLE 'FLX') - (K * U0) ;
  368. SI KAMOR ;
  369. F2ND = F2ND - (C * V0) ;
  370. FINS ;
  371. SI KCHARM ;
  372. TCHAR = CHARMECA WTAB ;
  373. KCHAF = EXIS TCHAR 'ADDI_SECOND' ;
  374. SI KCHAF ;
  375. F2ND = F2ND + TCHAR.'ADDI_SECOND' ;
  376. FINS ;
  377. KCHAK = EXIS TCHAR 'ADDI_MATRICE' ;
  378. SI KCHAK ;
  379. K0 = K0 ET TCHAR.'ADDI_MATRICE' ;
  380. FINS ;
  381. FINS ;
  382. FREAC0 = REAC K0 U0 ;
  383. F2ND = F2ND + FREAC0 ;
  384. * (1.*M) POUR EVITER PB HORODATAGE SI REUTILISATION DE M PLUS TARD
  385. A0 = RESO (1.*M) F2ND ;
  386. FINS ;
  387.  
  388.  
  389. * ===================
  390. * TABLES DE RESULTATS
  391. * ===================
  392.  
  393. INDVTK = 0 ;
  394. INDCSV = 0 ;
  395. INDGIB = 0 ;
  396. NFANTO = 0 ;
  397. EPOCH0 = DATE 'EPOCH';
  398.  
  399. * SAUVEGARDE VTK
  400. * --------------
  401. SI (IVTK > 0) ;
  402. SI (IVTK EGA 1) ;
  403. KVTK = VRAI ;
  404. SINON ;
  405. KVTK = DANS ETAB.'TEMPS_SAUVES_VTK' T0 ;
  406. FINS ;
  407.  
  408. SI KVTK ;
  409. INDVTK = INDVTK + 1 ;
  410. OPTI 'SORT' FICVTK ;
  411. SORT 'VTK' MAIVTK U0 'DEPL' V0 'VITE' 'TEMP' T0 ;
  412. CHACHA = CHAI CHACHA ' => VTK #' INDVTK ;
  413. FINS ;
  414. FINS ;
  415.  
  416. * SAUVEGARDE CSV
  417. * --------------
  418. SI (ICSV > 0) ;
  419.  
  420. * Initialisation de la table TABCSV
  421. TABCSV = TABL ;
  422. ETAB.'RESULTATS_CSV' = TABCSV ;
  423. TABCSV.'TEMPS' = PROG ;
  424. LCOCSV = MOTS ;
  425. LPOCSV = LECT ;
  426.  
  427. * Titres de colonnes et initialisation des LISTREEL
  428. CHPO1 = REDU (CHPO 'UNIFORME' 0. K) MAICSV ;
  429. REPE I1 NELCSV ;
  430. POIN1 = MAICSV POIN &I1 ;
  431. COMP1 = EXTR (REDU CHPO1 POIN1) 'COMP' ;
  432. SI KLISCO ;
  433. COMP1 = COMP1 SAUF (COMP1 SAUF LISCO) ;
  434. FINS ;
  435. NCO1 = DIME COMP1 ;
  436. SI (NCO1 EGA 0) ;
  437. ITER I1 ;
  438. FINS ;
  439. LPOCSV = LPOCSV ET (LECT NCO1*&I1) ;
  440. REPE I2 NCO1 ;
  441. CHA1 = EXTR COMP1 &I2 ;
  442. REPE I3 4 ;
  443. J3 = 5 - &I3 ;
  444. SI (NEG (EXTR CHA1 J3 J3) ' ') ;
  445. CHA1 = EXTR CHA1 1 J3 ;
  446. QUIT I3 ;
  447. FINS ;
  448. FIN I3 ;
  449. LCOCSV = LCOCSV ET CHA1 ;
  450. NOM1 = CHAI 'DEPL_' CHA1 '_' &I1 ;
  451. NOM2 = CHAI 'VITE_' CHA1 '_' &I1 ;
  452. TABCSV.NOM1 = PROG ;
  453. TABCSV.NOM2 = PROG ;
  454. FIN I2 ;
  455. FIN I1 ;
  456. _TABCSV = INDE TABCSV ;
  457. NBCSV = ((DIME TABCSV) - 1) / 2 ;
  458.  
  459. * Garder l'instant initial ?
  460. SI (ICSV EGA 1) ;
  461. KCSV = VRAI ;
  462. SINON ;
  463. KCSV = DANS ETAB.'TEMPS_SAUVES_CSV' T0 ;
  464. FINS ;
  465.  
  466. SI KCSV ;
  467. INDCSV = INDCSV + 1 ;
  468. TABCSV.'TEMPS' = TABCSV.'TEMPS' ET T0 ;
  469. REPE I1 NBCSV ;
  470. COMP1 = EXTR LCOCSV &I1 ;
  471. POIN1 = MAICSV POIN (EXTR LPOCSV &I1) ;
  472. X1 = EXTR U0 'VALE' COMP1 POIN1 'NOID' ;
  473. X2 = EXTR V0 'VALE' COMP1 POIN1 'NOID' ;
  474. K1 = 2*&I1 ;
  475. K2 = K1 + 1 ;
  476. MOT1 = MOT _TABCSV.K1 ;
  477. MOT2 = MOT _TABCSV.K2 ;
  478. TABCSV.MOT1 = TABCSV.MOT1 ET X1 ;
  479. TABCSV.MOT2 = TABCSV.MOT2 ET X2 ;
  480. FIN I1 ;
  481. CHACHA = CHAI CHACHA ' => CSV #' INDCSV ;
  482. FINS ;
  483. FINS ;
  484.  
  485. * SAUVEGARDE GIBIANE
  486. * ------------------
  487. SI (IGIB > 0) ;
  488.  
  489. * Initialisation de la table STAB (indicee par INDGIB)
  490. STAB = TABL ;
  491. ETAB.'RESULTATS' = STAB ;
  492.  
  493. * Garder l'instant initial ?
  494. SI (EGA IGIB 1) ;
  495. KGIB = VRAI ;
  496. SINON ;
  497. KGIB = DANS ETAB.'TEMPS_SAUVES' T0 ;
  498. FINS ;
  499.  
  500. SI (KGIB OU (EGA IPAS NPAS)) ;
  501. * (on force la sauvegarde si dernier pas)
  502. INDGIB = INDGIB + 1 ;
  503. STN = TABL ;
  504. STAB.INDGIB = STN ;
  505. STN.'TEMP' = T0 ;
  506. SI KMAIGIB ;
  507. STN.'DEPL' = REDU U0 MAIGIB ;
  508. STN.'VITE' = REDU V0 MAIGIB ;
  509. SINON ;
  510. STN.'DEPL' = U0 ;
  511. STN.'VITE' = V0 ;
  512. FINS ;
  513. CHACHA = CHAI CHACHA ' => GIBI #' INDGIB ;
  514. FINS ;
  515. FINS ;
  516.  
  517.  
  518. * =======================
  519. * TABLE DE CALCUL INTERNE
  520. * =======================
  521.  
  522. WTAB = TABL ;
  523.  
  524. * On transmet a WTAB la table de PARAmetres utilisateurs si elle existe
  525. SI (EXIS ETAB 'PARA') ;
  526. WTAB.'PARA' = ETAB.'PARA' ;
  527. FINS ;
  528. WTAB.'TEMP' = T0 ;
  529. WTAB.'DEPL' = U0 ;
  530. WTAB.'VITE' = V0 ;
  531.  
  532.  
  533.  
  534.  
  535. * +-----------------------------------------------------------------+
  536. * | |
  537. * | B O U C L E S U R L E S P A S D E T E M P S |
  538. * | |
  539. * +-----------------------------------------------------------------+
  540.  
  541.  
  542. * AFFICHAGE DU MESSAGE D'INFORMATION POUR L'ETAT INITIAL
  543. MESS CHACHA ;
  544.  
  545. REPE BPAS (NPAS - 1) ;
  546. IPAS = IPAS + 1 ;
  547.  
  548. * ===========================
  549. * PREPARATION DU PAS DE TEMPS
  550. * ===========================
  551.  
  552. * PARAMETRES TEMPORELS
  553. * --------------------
  554. T1 = EXTR TPROG IPAS ;
  555. DT = T1 - T0 ;
  556. DT2 = DT ** 2 ;
  557. CHACHA = CHAI 'PAS #' IPAS ' TEMPS' T1 ;
  558.  
  559. * QUELQUES CONSTANTES
  560. * -------------------
  561. XM = 1. / (BETA*DT2) ;
  562. XC = GAMMA / (BETA*DT) ;
  563.  
  564. * WTAB AVEC U_N ET T_N+1
  565. WTAB.'TEMP' = T1 ;
  566. WTAB.'DEPL' = U0 ;
  567. WTAB.'VITE' = V0 ;
  568. * WTAB.'ACCE' = A0 ;
  569.  
  570. * PROCEDURE "UTILISATEUR" (CHARMECA)
  571. * ----------------------------------
  572. KCHAF = FAUX ;
  573. KCHAK = FAUX ;
  574. KCHAFNL = FAUX ;
  575. KCHAKNL = FAUX ;
  576. KCHACNL = FAUX ;
  577. SI KCHARM ;
  578. TCHAR = CHARMECA WTAB ;
  579.  
  580. * Indices destines a l'ajout de relations (mult. de Lagrange)
  581. * entre deplacements
  582. KCHAF = EXIS TCHAR 'ADDI_SECOND' ;
  583. SI KCHAF ;
  584. FCHAR1 = TCHAR.'ADDI_SECOND' ;
  585. FINS ;
  586. KCHAK = EXIS TCHAR 'ADDI_MATRICE' ;
  587. SI KCHAK ;
  588. KCHAR1 = TCHAR.'ADDI_MATRICE' ;
  589. FINS ;
  590.  
  591. * Indices destines aux autres comportements non lineaires
  592. KCHAFNL = EXIS TCHAR 'ADDI_FNL' ;
  593. SI KCHAFNL ;
  594. FNL0 = TCHAR.'ADDI_FNL' ;
  595. FINS ;
  596. KCHAKNL = EXIS TCHAR 'ADDI_KNL' ;
  597. SI KCHAKNL ;
  598. KNL1 = TCHAR.'ADDI_KNL' ;
  599. FINS ;
  600. KCHACNL = EXIS TCHAR 'ADDI_CNL' ;
  601. SI KCHACNL ;
  602. CNL1 = TCHAR.'ADDI_CNL' ;
  603. FINS ;
  604. FINS ;
  605.  
  606. * CALCUL DU SECOND MEMBRE
  607. * -----------------------
  608. F1 = TIRE CHAR1 T1 ;
  609. * (aucune difference entre ADDI_SECOND et ADDI_FNL si ce n'est
  610. * l'affichage dans CHACHA)
  611. SI KCHAF ;
  612. F1 = F1 + FCHAR1 ;
  613. CHACHA = CHAI CHACHA ' +FLX' ;
  614. FINS ;
  615. SI KCHAFNL ;
  616. F1 = F1 + FNL0 ;
  617. CHACHA = CHAI CHACHA ' +FNL' ;
  618. FINS ;
  619. FEXT1 = (1MAF*(F1 ENLE 'FLX')) + (ALPHAF*(F0 ENLE 'FLX'))
  620. + (EXCO 'FLX' F1 'NOID' 'FLX') ;
  621.  
  622. * ANCIEN SCHEMA => COMPATIBLE AVEC VITEUNIL
  623. SI KVITUN ;
  624. * SI (EGA IALGO 0) ;
  625. K1 = 'EXTRAI' K 'RIGI' 'NOMU' ;
  626. K2 = 'EXTRAI' K 'RIGI' 'MULT' ;
  627. MAICO = EXTR K 'MAIL' 'UNIL' ;
  628. A1PRED = COLI V0 (1MAM / (BETA*DT)) ;
  629. FA = M * A1PRED ;
  630. FB = K1 * (2. * (U0 - (EXCO 'LX' U0 'NOID' 'LX'))) ;
  631. SI ((NBEL MAICO) > 0) ;
  632. U0PRED = U0 - (REDU U0 MAICO) ;
  633. SINON ;
  634. U0PRED = U0 ;
  635. FINS ;
  636. FG = K2 * U0PRED ;
  637. F2ND = FEXT1 + (F0 ENLE 'FLX') - FB + FA - FG ;
  638.  
  639. * SYNTAXE ADAPTEE AU HHT ET ALPHA-GEN (MAIS NON COMPATIBLE AVEC
  640. * RELATIONS UNILATERALES ET VITEUNIL)
  641. SINON ;
  642. A1PRED = COLI V0 ( 1MAM / (BETA*DT) )
  643. A0 ( (1MAM * 0.5 / BETA) - 1. ) ;
  644. V1PRED = COLI V0 ( (1MAF * GAMBET) - 1. )
  645. A0 ( 1MAF * DT * ((0.5*GAMBET) - 1.) ) ;
  646. U0PRED = ENLE U0 'LX' ;
  647. F2ND = FEXT1 - (K * (U0PRED)) + (M * A1PRED) ;
  648. SI KAMOR ;
  649. F2ND = F2ND + (C * V1PRED) ;
  650. FINS ;
  651. SI KCHACNL ;
  652. F2ND = F2ND + (CNL1 * (V1PRED + V0)) ;
  653. FINS ;
  654. * 2nd membre complet = Fext_n+1 - Fint_n + Famor_n + M*(4/dt*V_n+A_n)
  655.  
  656. FINS ;
  657.  
  658. * ON RECALCULE L'OPERATEUR...
  659. * ...SI LE PAS DE TEMPS CHANGE OU SI UNE MATRICE CHANGE
  660. * -----------------------------------------------------
  661. SI ((DT NEG DT0) OU KCHAK) ;
  662. SI (KMENAG) ;
  663. DETR OPER ;
  664. FINS ;
  665. OPER = (1MAF * K) ET ((1MAM * XM) * M) ;
  666. SI KAMOR ;
  667. OPER = OPER ET ((1MAF * XC) * C) ;
  668. FINS ;
  669. SI KCHAK ;
  670. OPER = OPER ET KCHAR1 ;
  671. CHACHA = CHAI CHACHA ' +[LX]' ;
  672. FINS ;
  673. SI KCHAKNL ;
  674. OPER = OPER ET (1MAF * KNL1) ;
  675. CHACHA = CHAI CHACHA ' +[KNL]' ;
  676. FINS ;
  677. SI KCHACNL ;
  678. OPER = OPER ET ((1MAF * XC) * CNL1) ;
  679. CHACHA = CHAI CHACHA ' +[CNL]' ;
  680. FINS ;
  681. FINS ;
  682.  
  683.  
  684. * ======================
  685. * CALCUL DU PAS DE TEMPS
  686. * ======================
  687.  
  688. * RESOLUTION DU SYSTEME LINEAIRE
  689. DU = RESO OPER F2ND ;
  690.  
  691. * MISE A JOUR DES VARIABLES AU TEMPS N+1
  692. * U1 = U0 + DU ;
  693. * A priori inutile d'enlever les LX, mais histoire d'etre sur...
  694. U1 = (ENLE U0 'LX') + DU ;
  695. DU1 = ENLE DU 'LX' ;
  696. V1 = COLI DU1 XC V0 (1. - GAMBET) A0 (DT*(1. - (0.5*GAMBET))) ;
  697. A1 = COLI DU1 XM V0 (-1./(BETA*DT)) A0 (1. - (0.5/BETA)) ;
  698.  
  699. * * Calcul du Residu (a n+1 ou avec les alphas pour NL?)
  700. * * rem : inutile en lineaire
  701. * * K*U1 contient deja les reactions
  702. * Res1 = F1 - (K * U1) - (M * A1);
  703. * MaxRes1 = MAXI Res1 'ABS' 'SANS' (mots 'FLX');
  704. * chacha = chai chacha ' max(R)=' MaxRes1;
  705. *
  706. * * Calcul de l'energie (a n+1)
  707. * Wmeca1 = 0.5 * (XTMX K (U1 enle LX) + (XTMX M V1)) ;
  708. * chacha = chai chacha ' Wmeca=' FORMAT '(F8.1)' Wmeca1 ;
  709.  
  710. * CORRECTION DES VITESSES (SEULEMENT POUR NEWMARK ACC. MOY.)
  711. SI KVITUN ;
  712.  
  713. * Eventuelle mise a jour de l'indicateur de presence de relations
  714. * unilaterales
  715. SI (KCHAK OU (IPAS <EG 2)) ;
  716. MTA = EXTR OPER 'CONTACT' ;
  717. KCONTA = ((DIME MTA) > 0) ;
  718. FINS ;
  719.  
  720. SI KCONTA ;
  721. TBID = TABL ;
  722. V1 = VITEUNIL OPER M V1 DU U0 DT F2ND TBID ;
  723. SI (EXIS TBID 'RATE_VITEUNIL') ;
  724. MESS 'Erreur fatale dans VITEUNIL' ;
  725. ERRE 5 ;
  726. FINS ;
  727. FINS ;
  728. FINS ;
  729.  
  730.  
  731. * ==========================
  732. * SAUVEGARDE DU PAS DE TEMPS
  733. * ==========================
  734.  
  735. * SAUVEGARDE VTK
  736. * --------------
  737. SI (IVTK > 0) ;
  738. SI (IVTK EGA 1) ;
  739. KVTK = MULT (IPAS - 1) IPASVTK ;
  740. SINON ;
  741. KVTK = DANS ETAB.'TEMPS_SAUVES_VTK' T1 ;
  742. FINS ;
  743.  
  744. SI KVTK ;
  745. INDVTK = INDVTK + 1 ;
  746. OPTI 'SORT' FICVTK ;
  747. SORT 'VTK' MAIVTK U1 'DEPL' V1 'VITE' 'TEMP' T1 ;
  748. CHACHA = CHAI CHACHA ' => VTK #' INDVTK ;
  749. FINS ;
  750. FINS ;
  751.  
  752. * SAUVEGARDE CSV
  753. * --------------
  754. SI (ICSV > 0) ;
  755. SI (ICSV EGA 1) ;
  756. KCSV = MULT (IPAS - 1) IPASCSV ;
  757. SINON ;
  758. KCSV = DANS ETAB.'TEMPS_SAUVES_CSV' T1 ;
  759. FINS ;
  760.  
  761. SI KCSV ;
  762. INDCSV = INDCSV + 1 ;
  763. TABCSV.'TEMPS' = TABCSV.'TEMPS' ET T1 ;
  764. REPE I1 NBCSV ;
  765. COMP1 = EXTR LCOCSV &I1 ;
  766. POIN1 = MAICSV POIN (EXTR LPOCSV &I1) ;
  767. X1 = EXTR U1 'VALE' COMP1 POIN1 'NOID' ;
  768. X2 = EXTR V1 'VALE' COMP1 POIN1 'NOID' ;
  769. K1 = 2*&I1 ;
  770. K2 = K1 + 1 ;
  771. MOT1 = MOT _TABCSV.K1 ;
  772. MOT2 = MOT _TABCSV.K2 ;
  773. TABCSV.MOT1 = TABCSV.MOT1 ET X1 ;
  774. TABCSV.MOT2 = TABCSV.MOT2 ET X2 ;
  775. FIN I1 ;
  776. CHACHA = CHAI CHACHA ' => CSV #' INDCSV ;
  777. FINS ;
  778.  
  779. SI ((EGA IPAS NPAS) ET KFICCSV) ;
  780. OPTI 'SORT' FICCSV ;
  781. SORT 'EXCE' TABCSV ;
  782. FINS ;
  783. FINS ;
  784.  
  785. * SAUVEGARDE GIBIANE
  786. * ------------------
  787. SI (IGIB > 0) ;
  788. SI (IGIB EGA 1) ;
  789. KGIB = MULT (IPAS - 1) IPASGIB ;
  790. SINON ;
  791. KGIB = DANS ETAB.'TEMPS_SAUVES' T1 ;
  792. FINS ;
  793.  
  794. EPOCH1='DATE' 'EPOCH';
  795. SI ((IPAS 'EGA' NPAS) 'OU' KGIB);
  796. * (on force la sauvegarde si dernier pas)
  797. INDGIB = INDGIB + 1 ;
  798. STN = TABL ;
  799. STAB.INDGIB = STN ;
  800. STN.'TEMP' = T1 ;
  801. SI KMAIGIB ;
  802. STN.'DEPL' = REDU U1 MAIGIB ;
  803. STN.'VITE' = REDU V1 MAIGIB ;
  804. SINON ;
  805. STN.'DEPL' = U1 ;
  806. STN.'VITE' = V1 ;
  807. FINS ;
  808. CHACHA = CHAI CHACHA ' => GIBI #' INDGIB ;
  809.  
  810. * Ecriture incrementale sur disque dur
  811. SI (KINCR 'ET'
  812. (((EPOCH1 - EPOCH0) '>EG' 300.D0) 'OU'
  813. ((EPOCH0 EGA 0.D0) 'ET' (EPOCH1 EGA 0.D0)) 'OU'
  814. (IPAS 'EGA' NPAS))) ;
  815. * (Il est necessaire de sauver aussi OPER : RIGIDITE)
  816. ETAB.'OPER' = OPER ;
  817. SAUV 'MUET' ETAB ;
  818. EPOCH0='DATE' 'EPOCH';
  819.  
  820. * Fantomisation des resultats sauves precedemment
  821. * (attention : on a encore besoin de ceux en cours pour le
  822. * prochain pas)
  823. SI (KECON) ;
  824. NFANTOL=(INDGIB - 1) - NFANTO ;
  825. SI(NFANTOL > 0);
  826. REPE SURI NFANTOL;
  827. JJ=NFANTO + &SURI;
  828. SI(JJ <EG 0);
  829. ITER SURI;
  830. FINS;
  831. 'FANT' (STAB.JJ) 'DEPL' ;
  832. 'FANT' (STAB.JJ) 'VITE' ;
  833. FIN SURI;
  834. NFANTO=NFANTO+NFANTOL;
  835. FINS ;
  836. FINS ;
  837. FINS ;
  838. FINS ;
  839. FINS ;
  840.  
  841.  
  842. * ===============================
  843. * PASSAGE AU PAS DE TEMPS SUIVANT
  844. * ===============================
  845.  
  846. * AFFICHAGE DU MESSAGE D'INFORMATION DU PAS DE TEMPS COURANT
  847. MESS CHACHA ;
  848.  
  849. T0 = T1 ;
  850. DT0 = DT ;
  851. * U0 = U1 ENLE 'LX' ;
  852. U0 = U1 ;
  853. V0 = V1 ;
  854. DETR A0 ;
  855. A0 = A1 ;
  856. DETR F0 ;
  857. F0 = F1 ;
  858. SI (KMENAG) ;
  859. DETR A1PRED ;
  860. DETR DU ;
  861. FINS ;
  862.  
  863. FIN BPAS ;
  864.  
  865.  
  866. * Sauvegarde de la derniere acceleration (pour reprise eventuelle)
  867. STN.'ACCE' = A1 ;
  868.  
  869. FINP STAB ;
  870.  
  871.  
  872.  

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