Télécharger cadgsi.eso

Retour à la liste

Numérotation des lignes :

cadgsi
  1. C CADGSI SOURCE CB215821 20/11/25 13:19:03 10792
  2. SUBROUTINE CADGSI
  3. C************************************************************************
  4. C
  5. C OBJET :
  6. C
  7. C CALCUL DE LA MATRICE MASSE DIAGONALE ---> Creation d'un CHAMPOIN
  8. C Pour les SOMMETS D0=NI ( MASSE LUMPE )
  9. C Pour les FACES D0=1/2(Vol1 + Vol2)
  10. C Pour les CENTRES D0=Vol Elt
  11. C
  12. C SYNTAXE :
  13. C
  14. C RES = DGSI OBJ1 <TYPE> <'IMPR'> ;
  15. C
  16. C OBJ1 : Table DOMAINE
  17. C TYPE ; SOMMET , FACE , CENTRE (par defaut SOMMET) MSOMMET CENTREP1
  18. C
  19. C
  20. C************************************************************************
  21. IMPLICIT INTEGER(I-N)
  22. IMPLICIT REAL*8 (A-H,O-Z)
  23. INTEGER INEFMD
  24. -INC CCGEOME
  25. -INC SMCHAML
  26.  
  27. -INC PPARAM
  28. -INC CCOPTIO
  29. -INC SMELEME
  30. POINTEUR MELEMP.MELEME
  31. -INC SMTABLE
  32. POINTEUR MTABD.MTABLE
  33. -INC SMCOORD
  34. -INC SMLENTI
  35. -INC SMCHPOI
  36. -INC SIZFFB
  37. POINTEUR IZF1.IZFFM,IZH2.IZHR
  38. CHARACTER*8 NOM0,CHAI,LISMO(5),TYPE,TYPC,NOM
  39. DATA LISMO/'SOMMET ','FACE ','CENTRE ','MSOMMET ','CENTREP1'/
  40. C ***************************************************************
  41. C----------------------------------------------------------------------
  42. C KPOIN = 0->SOMMET 1-> FACE 2-> CENTRE 3-> CENTREP0 4-> CENTREP1 5-> MSOMMET
  43. C INEFMD : Type formulation INEFMD=1 LINE,=2 MACRO,=3 QUADRATIQUE, INEFMD=4 LINB
  44. C************************************************************************
  45.  
  46.  
  47. INEFMD=0
  48. IKAS=1
  49. IMPR=0
  50. IAXI=0
  51. IF(IFOMOD.EQ.0)IAXI=2
  52.  
  53. CALL LITABS('DOMAINE ',MTABD,1,1,IRET)
  54. IF(IRET.EQ.0)THEN
  55. WRITE(6,*)' On attend une table de soustype DOMAINE'
  56. RETURN
  57. ENDIF
  58.  
  59. 19 CONTINUE
  60. CALL LIRCHA(CHAI,0,IRET)
  61. IF(IRET.EQ.0)GO TO 20
  62. CALL OPTLI(IP,LISMO,CHAI,5)
  63. IF(IP.EQ.0)THEN
  64. WRITE(6,*)' On attend un mot cle parmi SOMMET FACE CENTRE ',
  65. & 'MSOMMET CENTREP1 '
  66. RETURN
  67. ENDIF
  68. IKAS=IP
  69.  
  70. 20 CONTINUE
  71.  
  72. IF(IKAS.EQ.1.OR.IKAS.EQ.4)THEN
  73. C SOMMET et MSOMMET
  74. C SOMMET
  75. IF(IKAS.EQ.1)THEN
  76. TYPE=' '
  77. CALL ACMO(MTABD,'MAILLAGE',TYPE,MELEME)
  78. IF(TYPE.NE.'MAILLAGE')RETURN
  79. TYPE=' '
  80. CALL ACMO(MTABD,'SOMMET',TYPE,MELEMS)
  81. IF(TYPE.NE.'MAILLAGE')RETURN
  82. CALL CRCHPT('SOMMET',MELEMS,1,MCHPOI)
  83. ENDIF
  84. C MSOMMET
  85. IF(IKAS.EQ.4)THEN
  86. c write(6,*)' CADGSI: IKAS=4'
  87. TYPE=' '
  88. CALL ACMO(MTABD,'MMAIL ',TYPE,MELEME)
  89. IF(TYPE.NE.'MAILLAGE')RETURN
  90. TYPE=' '
  91. CALL ACMO(MTABD,'MSOMMET',TYPE,MELEMS)
  92. IF(TYPE.NE.'MAILLAGE')RETURN
  93. CALL CRCHPT('MSOMMET',MELEMS,1,MCHPOI)
  94. ENDIF
  95.  
  96. CALL KRIPAD(MELEMS,MLENTI)
  97.  
  98. C CREATION DE LA DIAGONALE
  99. CALL LICHTM(MCHPOI,MPOVAL,TYPC,IGEOM)
  100. SEGACT MELEME
  101. NBSOUS=LISOUS(/1)
  102. IF(NBSOUS.EQ.0)NBSOUS=1
  103. C
  104. C BOUCLE SUR LES TYPES D'ELEMENTS ET CALCUL
  105. C
  106. DO 1 KSOUS=1,NBSOUS
  107. IF(NBSOUS.EQ.1)IPT1=MELEME
  108. IF(NBSOUS.GT.1)IPT1=LISOUS(KSOUS)
  109. SEGACT IPT1
  110.  
  111. NP=IPT1.NUM(/1)
  112. NEL=IPT1.NUM(/2)
  113. C
  114. NOM0=NOMS(IPT1.ITYPEL)//' '
  115. IF(INEFMD.EQ.1.AND.IKAS.EQ.4)NOM0=NOMS(IPT1.ITYPEL)//'P1P1'
  116. IF(INEFMD.EQ.2.AND.IKAS.EQ.4)NOM0=NOMS(IPT1.ITYPEL)//'MCF1'
  117. IF(INEFMD.EQ.3.AND.IKAS.EQ.4)NOM0=NOMS(IPT1.ITYPEL)//'PFP1'
  118. c write(6,*)'4 NOM0=',NOM0
  119. CALL KALPBG(NOM0,'FONFORM ',IZFFM)
  120. SEGACT IZFFM*MOD
  121. IZHR=KZHR(1)
  122. SEGACT IZHR*MOD
  123. C
  124. NPG=FN(/2)
  125. NES=GR(/1)
  126.  
  127. IF(IMPR.NE.0)THEN
  128. WRITE(6,*)' SUB CADGSI : NES,NP,NPG,IDIM,NEL='
  129. & ,NES,NP,NPG,IDIM,NEL
  130. ENDIF
  131. C
  132. DO 10 K=1,NEL
  133. C
  134. NPGR=0
  135. IF(IAXI.NE.0)NPGR=NPG
  136. C
  137. DO 12 I=1,NP
  138. J=IPT1.NUM(I,K)
  139. DO 13 N=1,IDIM
  140. XYZ(N,I)=XCOOR((J-1)*(IDIM+1) +N)
  141. 13 CONTINUE
  142. 12 CONTINUE
  143.  
  144. CALL CALJBC(FN,GR,PG,XYZ,HR,PGSQ,RPG,NES,IDIM,NP,NPG,IAXI,AIRE)
  145.  
  146. IF(IMPR.NE.0)THEN
  147. WRITE(6,*)' SUB CADGSI : AIRE=',AIRE
  148. WRITE(6,*)' SUB CADGSI : LER '
  149. WRITE(6,1001)(IPT1.NUM(I,K),I=1,NP)
  150. WRITE(6,*)' SUB CADGSI : XYZ '
  151. WRITE(6,1002)((XYZ(N,I),N=1,2),I=1,NP)
  152. ENDIF
  153. C
  154. DO 3 J=1,NP
  155. SJ=0.D0
  156. DO 4 L=1,NPG
  157. SJ=SJ+FN(J,L)*PGSQ(L)
  158. 4 CONTINUE
  159. C SD(J,K)=SJ
  160. JU=LECT(IPT1.NUM(J,K))
  161. C D0(JU)=D0(JU)+SJ
  162. VPOCHA(JU,1)=VPOCHA(JU,1)+SJ
  163. 3 CONTINUE
  164. C
  165. 10 CONTINUE
  166.  
  167. SEGDES IPT1
  168. SEGSUP IZFFM,IZHR
  169. 1 CONTINUE
  170. SEGDES MELEME
  171. SEGSUP MLENTI
  172. SEGDES MPOVAL
  173.  
  174. IF(IMPR.NE.0)THEN
  175. WRITE(6,*)' SUB CADGSI : CALCUL DE LA DIAGONALE'
  176. C WRITE(6,1003)(I,VPOCHA(I,1),I=1,NPT)
  177. WRITE(6,*)' FIN DE CADGSI'
  178. ENDIF
  179.  
  180. C FACE
  181. ELSEIF(IKAS.EQ.2)THEN
  182.  
  183. TYPE=' '
  184. CALL ACMO(MTABD,'FACE',TYPE,MELEMF)
  185. IF(TYPE.NE.'MAILLAGE')RETURN
  186. TYPE=' '
  187. CALL ACMO(MTABD,'CENTRE',TYPE,MELEMC)
  188. IF(TYPE.NE.'MAILLAGE')RETURN
  189. CALL LEKTAB(MTABD,'FACEL',MELEME)
  190. IF(MELEME.EQ.0)RETURN
  191. CALL LEKTAB(MTABD,'XXVOLUM',MCHPO1)
  192. IF(MCHPO1.EQ.0)RETURN
  193.  
  194. CALL KRIPAD(MELEMC,MLENT1)
  195. CALL KRIPAD(MELEMF,MLENT2)
  196.  
  197. C CREATION DE LA DIAGONALE
  198. CALL CRCHPT('FACE',MELEMF,1,MCHPOI)
  199. CALL LICHTM(MCHPOI,MPOVAL,TYPC,IGEOM)
  200. CALL LICHTM(MCHPO1,MPOVA1,TYPC,IGEOM)
  201. SEGACT MELEME
  202. C
  203. C BOUCLE SUR LES TYPES D'ELEMENTS ET CALCUL
  204. C
  205. NBEL=NUM(/2)
  206. DO 2 K=1,NBEL
  207. I1=NUM(1,K)
  208. I2=NUM(2,K)
  209. I3=NUM(3,K)
  210. I1=MLENT1.LECT(I1)
  211. I2=MLENT2.LECT(I2)
  212. I3=MLENT1.LECT(I3)
  213. V=(MPOVA1.VPOCHA(I1,1)+MPOVA1.VPOCHA(I3,1) )*0.5D0
  214. VPOCHA(I2,1)=V
  215. 2 CONTINUE
  216.  
  217. SEGDES MELEME
  218. SEGDES MPOVA1,MPOVAL
  219. SEGSUP MLENT1,MLENT2
  220.  
  221. ELSEIF(IKAS.EQ.3)THEN
  222. C CENTRE
  223. CALL LEKTAB(MTABD,'XXVOLUM',MCHPOI)
  224. IF(MCHPOI.EQ.0)RETURN
  225.  
  226. ELSEIF(IKAS.EQ.5.OR.IKAS.EQ.6)THEN
  227. C CENTREP1 et CENTREP0
  228.  
  229. CALL ACME(MTABD,'INEFMD',INEFMD)
  230.  
  231. IF(INEFMD.EQ.1)THEN
  232. C% Le type d'élément fini %m1:8 ne convient pas.
  233. MOTERR( 1: 8) = 'LINE'
  234. CALL ERREUR(927)
  235. RETURN
  236. ENDIF
  237.  
  238. TYPE = ' '
  239. CALL ACMO(MTABD,'XXPSOML',TYPE,MCHELM)
  240. TYPE = ' '
  241. IF (TYPE.NE.'MCHAML ') THEN
  242. CALL ECROBJ('TABLE',MTABD)
  243. CALL KPSOML
  244. TYPE = 'MCHAML'
  245. CALL LIROBJ(TYPE,MCHELM,1,IRET)
  246. IF (IRET.EQ.0)THEN
  247. NOM='XXPSOML'
  248. GOTO 5000
  249. ENDIF
  250. CALL ECMO(MTABD,'XXPSOML','MCHAML',MCHELM)
  251. ENDIF
  252. SEGACT MCHELM
  253.  
  254. TYPE=' '
  255. IF(INEFMD.EQ.2)THEN
  256. CALL ACMO(MTABD,'MACRO1 ',TYPE,MELEME)
  257. IF(TYPE.NE.'MAILLAGE')RETURN
  258. ELSEIF(INEFMD.EQ.3)THEN
  259. CALL ACMO(MTABD,'MAILLAGE',TYPE,MELEME)
  260. NOM='MAILLAGE'
  261. IF(TYPE.NE.'MAILLAGE')GO TO 5000
  262. ENDIF
  263.  
  264. TYPE=' '
  265. CALL ACMO(MTABD,'CENTREP1',TYPE,MELEMS)
  266. IF(TYPE.NE.'MAILLAGE')THEN
  267. CALL KCTRP1(MTABD,MELEMS,1)
  268. ENDIF
  269. CALL KRIPAD(MELEMS,MLENTI)
  270. CALL ACMO(MTABD,'ELTP1NC ',TYPE,MELEMP)
  271. NOM='ELTP1NC '
  272. IF(TYPE.NE.'MAILLAGE')GO TO 5000
  273. CALL CRCHPT('CENTREP1',MELEMS,1,MCHPOI)
  274.  
  275.  
  276. CALL LICHTM(MCHPOI,MPOVAL,TYPC,IGEOM)
  277. SEGACT MELEME
  278. NBSOUS=LISOUS(/1)
  279. IF(NBSOUS.EQ.0)NBSOUS=1
  280. NUTOEL=0
  281.  
  282. NPTD=VPOCHA(/1)
  283. IES=IDIM
  284. MP10=0
  285.  
  286. DO 11 L=1,NBSOUS
  287. IPT1=MELEME
  288. IF(NBSOUS.NE.1)IPT1=LISOUS(L)
  289. SEGACT IPT1
  290.  
  291. MCHAML=ICHAML(L)
  292. SEGACT MCHAML
  293. MELVAL=IELVAL(1)
  294. SEGACT MELVAL
  295.  
  296. NP =IPT1.NUM(/1)
  297. NBEL=IPT1.NUM(/2)
  298.  
  299. IPT2=MELEMP
  300. IF(NBSOUS.NE.1)IPT2=LISOUS(L)
  301. SEGACT IPT2
  302.  
  303. IF(INEFMD.EQ.2.AND.IKAS.EQ.5)NOM0=NOMS(IPT1.ITYPEL)//'MCP1'
  304. IF(INEFMD.EQ.3.AND.IKAS.EQ.5)NOM0=NOMS(IPT1.ITYPEL)//'PRP1'
  305.  
  306. CALL KALPBG(NOM0,'FONFORM ',IZFFM)
  307.  
  308.  
  309. SEGACT IZFFM*MOD
  310. IZHR=KZHR(1)
  311. IZH2=KZHR(2)
  312. SEGACT IZHR*MOD,IZH2*MOD
  313. NES=GR(/1)
  314. NPG=GR(/3)
  315. IZF1=KTP(1)
  316. SEGACT IZF1*MOD
  317. MP1=IZF1.FN(/1)
  318. NP = IPT1.NUM(/1)
  319.  
  320. DO 21 K=1,NBEL
  321. DO 109 I=1,NP
  322. J=IPT1.NUM(I,K)
  323. JC = (J-1)*(IDIM+1)
  324. DO 110 N=1,IDIM
  325. XYZ(N,I)=XCOOR( JC + N )
  326. 110 CONTINUE
  327. 109 CONTINUE
  328.  
  329. CALL CALJBC(FN,GR,PG,XYZ,HR,PGSQ,RPG,NES,IDIM,NP,NPG,IAXI,AIRE)
  330.  
  331. DO 39 M=1,MP1
  332. M11=LECT(IPT2.NUM(M,K))
  333. M1=M+MP10
  334. c IF(KPOIND.EQ.5)M1=M11
  335.  
  336. U=0.D0
  337. DO 33 LL=1,NPG
  338. U=U+IZF1.FN(M,LL)*PGSQ(LL)
  339. 33 CONTINUE
  340.  
  341. VPOCHA(M1,1)=VPOCHA(M1,1)+U
  342. 39 CONTINUE
  343.  
  344. MP10=MP10+MP1
  345. 21 CONTINUE
  346.  
  347. SEGSUP IZFFM,IZF1,IZHR,IZH2
  348.  
  349.  
  350. SEGDES IPT1,MCHAML,MELVAL,IPT2
  351. NUTOEL=NUTOEL+NBEL
  352.  
  353. 11 CONTINUE
  354.  
  355. SEGSUP MLENTI
  356.  
  357. ENDIF
  358.  
  359.  
  360. SEGDES MTABD,MPOVAL,MCHPOI
  361. CALL ACTOBJ('CHPOINT ',MCHPOI,1)
  362. CALL ECROBJ('CHPOINT ',MCHPOI)
  363. RETURN
  364. 5000 CONTINUE
  365. C Indice %m1:8 : Problème de données détecté dans lektab
  366. IPOINT = 0
  367. MOTERR(1:8) = NOM
  368. CALL ERREUR(792)
  369. RETURN
  370. 1001 FORMAT(20(1X,I5))
  371. 1002 FORMAT(10(1X,1PE11.4))
  372. 1003 FORMAT(6(1X,I7,1X,1PE11.4))
  373. END
  374.  
  375.  
  376.  
  377.  
  378.  
  379.  
  380.  
  381.  
  382.  
  383.  

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