Télécharger matp.eso

Retour à la liste

Numérotation des lignes :

  1. C MATP SOURCE CB215821 19/08/20 21:19:35 10287
  2. SUBROUTINE MATP
  3. C-----------------------------------------------------------------------
  4. C Cette subroutine permet de creer les matrices elementaires du système
  5. C matriciel en trace de charge dans le cadre de la résolution des
  6. C équations de DARCY par une méthode d'éléments finis mixtes hybride.
  7. C On obtient ce systeme en exprimant les inconnues de charge et de
  8. C flux en fonction des traces de charge.
  9. C -1 t -1
  10. C M1 = M2 * ( COEF * DIV * DIV * M2 - Id )
  11. C
  12. C avec M2 : Matrice massse hybride
  13. C DIV : Matrice uniligne representant la divergence
  14. C COEF : Scalaire valant 1/d dans le cas permanant, THETA*BETA/d
  15. C dans le cas transitoire.
  16. C -1 t
  17. C d : Scalaire valant ( DIV * M2 * DIV )
  18. C
  19. C THETA : Scalaire de la THETA-méthode | Données
  20. C BETA : Scalaire d*dt/(Ck|K| + THETA*d*dt) | utiles
  21. C dt : Scalaire pas de temps | uniquement
  22. C Ck : Scalaire coefficient d'emmagasinement | en
  23. C |K| : Scalaire aire de l'élément K | transitoire
  24. C
  25. C---------------------------
  26. C Phrase d'appel (GIBIANE) :
  27. C---------------------------
  28. C
  29. C MRIGI1 = MATP MMODEL MRIGI2 (TABLE2) ;
  30. C
  31. C------------------------
  32. C Operandes et resultat :
  33. C------------------------
  34. C
  35. C MRIGI1 : Matrices elementaires en trace de charge.
  36. C MMODEL : Objet modele specifiant la formulation.
  37. C MRIGI2 : Matrices masses hybrides elementaires creees par MHYB.
  38. C TABLE2 : TABLE DARCY_TRANSITOIRE contenant les infos transitoires.
  39. C
  40. C-----------------------------------------------------------------------
  41. C
  42. C Langage : ESOPE + FORTRAN77
  43. C
  44. C Auteurs : 09/93 F.DABBENE - Cas permanent
  45. C 09/94 X.NOUVELLON - Extension au cas transitoire
  46. C
  47. C-----------------------------------------------------------------------
  48. IMPLICIT INTEGER(I-N)
  49. IMPLICIT REAL*8 (A-H,O-Z)
  50. *
  51.  
  52. -INC PPARAM
  53. -INC CCOPTIO
  54. -INC SMELEME
  55. -INC SMRIGID
  56. -INC SMMODEL
  57. -INC SMTABLE
  58. -INC SMCHAML
  59. *
  60. SEGMENT IPMAHY
  61. INTEGER MAHYBR(NSOUS)
  62. ENDSEGMENT
  63. *
  64. LOGICAL LOGRE,LOGIN
  65. CHARACTER*6 CHAR6
  66. CHARACTER*8 TAPIND,TYPOBJ,CHARIN,CHARRE,LETYPE
  67. *
  68. * Initialisations
  69. *
  70. IVALIN = 0
  71. XVALIN = 0.D0
  72. LOGIN = .TRUE.
  73. IOBIN = 0
  74. TAPIND = 'MOT '
  75. *
  76. * Lecture du MMODEL
  77. *
  78. CALL LIROBJ('MMODEL',IPMODE,1,IRET1)
  79. IF (IERR.NE.0) RETURN
  80. MMODEL = IPMODE
  81. *
  82. * Lecture de la TABLE domaine et récupération des infos
  83. *
  84. IPTABL = 0
  85. CALL LEKMOD(IPMODE,IPTABL,IRET)
  86. IF (IERR.NE.0) RETURN
  87. CHARIN = 'MAILLAGE'
  88. TYPOBJ = 'MAILLAGE'
  89. CALL LEKTAB(IPTABL,CHARIN,IOBRE)
  90. IF (IERR.NE.0) RETURN
  91. IPGEOM = IOBRE
  92. CALL LEKTAB(IPTABL,'ELTFA',IOBRE)
  93. IF (IERR.NE.0) RETURN
  94. IELTFA = IOBRE
  95. *
  96. * Lecture de RIGIDITE
  97. *
  98. CALL LIROBJ('RIGIDITE',IPRIGI,1,IRET1)
  99. IF (IERR.NE.0) RETURN
  100. MRIGID = IPRIGI
  101. *
  102. * Lecture éventuelle de la TABLE DARCY_TRANSITOIRE
  103. *
  104. IPTABL = 0
  105. CALL LIRTAB('DARCY_TRANSITOIRE',IPTABL,0,IRET)
  106. IF (IERR.NE.0) RETURN
  107. IF (IRET.EQ.1) THEN
  108. TYPOBJ = 'FLOTTANT'
  109. CALL ACCTAB(IPTABL,TAPIND,IVALIN,XVALIN,'THETA',LOGIN,IOBIN,
  110. S TYPOBJ,IVALRE,XVALRE,CHARRE,LOGRE,IOBRE)
  111. IF (IERR.NE.0) RETURN
  112. THETA = XVALRE
  113. THEMIN = -1.D-12
  114. THEMAX = 1.D0 - THEMIN
  115. IF ((THETA.LT.THEMIN).OR.(THETA.GT.THEMAX)) THEN
  116. REAERR(1) = REAL(THETA)
  117. REAERR(2) = REAL(0.D0)
  118. REAERR(3) = REAL(1.D0)
  119. MOTERR(1:8) = ' THETA '
  120. CALL ERREUR(42)
  121. RETURN
  122. ENDIF
  123. IF (THETA.LT.0.D0) THETA=0.D0
  124. IF (THETA.GT.1.D0) THETA=1.D0
  125. CALL ACCTAB(IPTABL,TAPIND,IVALIN,XVALIN,'PAS',LOGIN,IOBIN,
  126. S TYPOBJ,IVALRE,XVALRE,CHARRE,LOGRE,IOBRE)
  127. IF (IERR.NE.0) RETURN
  128. DELTAT = XVALRE
  129. IF (DELTAT.LT.0.D0) THEN
  130. REAERR(1) = REAL(DELTAT)
  131. REAERR(2) = REAL(0.D0)
  132. MOTERR(1:8) = ' DELTAT '
  133. CALL ERREUR(41)
  134. RETURN
  135. ENDIF
  136. TYPOBJ = 'MCHAML '
  137. CALL ACCTAB(IPTABL,TAPIND,IVALIN,XVALIN,'SURF',LOGIN,IOBIN,
  138. S TYPOBJ,IVALRE,XVALRE,CHARRE,LOGRE,IOBRE)
  139. CALL ACTOBJ(TYPOBJ,IOBRE,1)
  140. IF (IERR.NE.0) RETURN
  141. IPCK = IOBRE
  142. ELSE
  143. THETA = 0.D0
  144. DELTAT = 0.D0
  145. IPCK = 0
  146. ENDIF
  147. *
  148. * Test de l'ojet MODELE
  149. * Récupération des zones élémentaires ou DARCY est défini
  150. *
  151. SEGACT MMODEL
  152. NSOUS = KMODEL(/1)
  153. SEGINI IPMAHY
  154. IDARCY = 0
  155. IPT2=IPGEOM
  156. IPT1=IPGEOM
  157. DO 10 ISOUS=1,NSOUS
  158. IF(NSOUS.GT.1)THEN
  159. SEGACT IPT2
  160. IPT1=IPT2.LISOUS(ISOUS)
  161. ENDIF
  162. IMODEL = KMODEL(ISOUS)
  163. SEGACT IMODEL
  164. LETYPE = FORMOD(1)
  165. IF (LETYPE.EQ.'DARCY') THEN
  166. IDARCY = IDARCY + 1
  167. MAHYBR(ISOUS) = IPT1
  168. ENDIF
  169. 10 CONTINUE
  170. IF (IDARCY.EQ.0) THEN
  171. MOTERR = LETYPE
  172. CALL ERREUR(193)
  173. GOTO 100
  174. ENDIF
  175. *
  176. * Recuperation des pointeurs ELTFA pour les zones ou DARCY est defini
  177. *
  178. MELEME = IELTFA
  179. SEGACT MELEME
  180. LZONES = LISOUS(/1)
  181. IF (LZONES.EQ.0) LZONES = 1
  182. IPT1 = IPGEOM
  183. SEGACT IPT1
  184. DO 30 ISOUS=1,NSOUS
  185. IMAMEL = MAHYBR(ISOUS)
  186. IF (IMAMEL.NE.0) THEN
  187. DO 20 ISZ=1,LZONES
  188. IF (LZONES.EQ.1) THEN
  189. IPT2 = IPT1
  190. IPT3 = MELEME
  191. ELSE
  192. IPT2 = IPT1.LISOUS(ISZ)
  193. IPT3 = LISOUS(ISZ)
  194. ENDIF
  195. IF (IPT2.EQ.IMAMEL) THEN
  196. MAHYBR(ISOUS) = IPT3
  197. GOTO 30
  198. ENDIF
  199. 20 CONTINUE
  200. IF (IMAMEL.EQ.MAHYBR(ISOUS)) THEN
  201. MOTERR(1:8) = ' MODELE '
  202. MOTERR(9:16) = ' ELTFA '
  203. INTERR(1) = ISOUS
  204. CALL ERREUR(698)
  205. GOTO 100
  206. ENDIF
  207. ENDIF
  208. 30 CONTINUE
  209. *
  210. * Test du sous-type de la matrice de rigidité récupérée
  211. * Controle des MELEME de l'objet RIGIDITE % au ELTFA liée au MMODEL
  212. *
  213. SEGACT MRIGID
  214. LETYPE = MTYMAT
  215. IF (LETYPE.NE.'DARCY') THEN
  216. MOTERR(1:8) = 'RIGIDITE'
  217. MOTERR(9:16) = 'DARCY'
  218. CALL ERREUR(79)
  219. SEGDES MRIGID
  220. GOTO 100
  221. ENDIF
  222. DO 40 ISOUS=1,NSOUS
  223. IF (MAHYBR(ISOUS).NE.0) THEN
  224. IPTEST = IRIGEL(1,ISOUS)
  225. IF (IPTEST.NE.MAHYBR(ISOUS)) THEN
  226. MOTERR(1:8) = 'DARCY'
  227. MOTERR(9:16) = ' ELTFA '
  228. INTERR(1) = ISOUS
  229. CALL ERREUR(698)
  230. SEGDES MRIGID
  231. GOTO 100
  232. ENDIF
  233. ENDIF
  234. 40 CONTINUE
  235. SEGDES MRIGID
  236. *
  237. * Vérification du support du MCHAML % au MMODEL
  238. * Controle du nombre de composantes reelles
  239. *
  240. IF (IPCK.NE.0) THEN
  241. ISUP = 2
  242. ICOND = 1
  243. CALL QUESUP(IPMODE,IPCK,ISUP,ICOND,IRET,IRET2)
  244. IF (IRET.GT.1) GOTO 100
  245. MCHELM = IPCK
  246. SEGACT MCHELM
  247. DO 50 ISOUS=1,NSOUS
  248. IF (MAHYBR(ISOUS).NE.0) THEN
  249. MCHAML = ICHAML(ISOUS)
  250. SEGACT MCHAML
  251. N2 = TYPCHE(/2)
  252. IF (N2.GT.1) THEN
  253. CALL ERREUR(320)
  254. GOTO 100
  255. ENDIF
  256. CHAR6 = TYPCHE(1)(1:6)
  257. IF (CHAR6.NE.'REAL*8') THEN
  258. MOTERR(1:8) = NOMCHE(1)
  259. CALL ERREUR(679)
  260. GOTO 100
  261. ENDIF
  262. ENDIF
  263. 50 CONTINUE
  264. ENDIF
  265. SEGDES IPMAHY
  266. *
  267. * Construction de IPRIG1
  268. *
  269. CALL MATP1(NSOUS,IPGEOM,IPMAHY,IPRIGI,THETA,DELTAT,IPCK,IPRIG1)
  270. CALL ECROBJ('RIGIDITE',IPRIG1)
  271. *
  272. * Ménage
  273. *
  274. 100 CONTINUE
  275. SEGSUP IPMAHY
  276. END
  277.  
  278.  
  279.  

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