Télécharger cadgsi.eso

Retour à la liste

Numérotation des lignes :

cadgsi
  1. C CADGSI SOURCE CB215821 23/01/25 21:15:05 11573
  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. INTEGER TMTABD(1)
  34. -INC SMCOORD
  35. -INC SMLENTI
  36. -INC SMCHPOI
  37. -INC SIZFFB
  38. POINTEUR IZF1.IZFFM,IZH2.IZHR
  39. CHARACTER*8 NOM0,CHAI,LISMO(5),TYPE,TYPC,NOM
  40. DATA LISMO/'SOMMET ','FACE ','CENTRE ','MSOMMET ','CENTREP1'/
  41. C ***************************************************************
  42. C----------------------------------------------------------------------
  43. C KPOIN = 0->SOMMET 1-> FACE 2-> CENTRE 3-> CENTREP0 4-> CENTREP1 5-> MSOMMET
  44. C INEFMD : Type formulation INEFMD=1 LINE,=2 MACRO,=3 QUADRATIQUE, INEFMD=4 LINB
  45. C************************************************************************
  46.  
  47.  
  48. INEFMD=0
  49. IKAS=1
  50. IMPR=0
  51. IAXI=0
  52.  
  53. SEGACT,MCOORD
  54.  
  55. IF(IFOMOD.EQ.0)IAXI=2
  56.  
  57. CALL LITABS('DOMAINE ',TMTABD,1,1,IRET)
  58. MTABD=TMTABD(1)
  59. IF(IRET.EQ.0)THEN
  60. WRITE(6,*)' On attend une table de soustype DOMAINE'
  61. RETURN
  62. ENDIF
  63.  
  64. 19 CONTINUE
  65. CALL LIRCHA(CHAI,0,IRET)
  66. IF(IRET.EQ.0)GO TO 20
  67. CALL OPTLI(IP,LISMO,CHAI,5)
  68. IF(IP.EQ.0)THEN
  69. WRITE(6,*)' On attend un mot cle parmi SOMMET FACE CENTRE ',
  70. & 'MSOMMET CENTREP1 '
  71. RETURN
  72. ENDIF
  73. IKAS=IP
  74.  
  75. 20 CONTINUE
  76.  
  77. IF(IKAS.EQ.1.OR.IKAS.EQ.4)THEN
  78. C SOMMET et MSOMMET
  79. C SOMMET
  80. IF(IKAS.EQ.1)THEN
  81. TYPE=' '
  82. CALL ACMO(MTABD,'MAILLAGE',TYPE,MELEME)
  83. IF(TYPE.NE.'MAILLAGE')RETURN
  84. TYPE=' '
  85. CALL ACMO(MTABD,'SOMMET',TYPE,MELEMS)
  86. IF(TYPE.NE.'MAILLAGE')RETURN
  87. CALL CRCHPT('SOMMET',MELEMS,1,MCHPOI)
  88. ENDIF
  89. C MSOMMET
  90. IF(IKAS.EQ.4)THEN
  91. c write(6,*)' CADGSI: IKAS=4'
  92. TYPE=' '
  93. CALL ACMO(MTABD,'MMAIL ',TYPE,MELEME)
  94. IF(TYPE.NE.'MAILLAGE')RETURN
  95. TYPE=' '
  96. CALL ACMO(MTABD,'MSOMMET',TYPE,MELEMS)
  97. IF(TYPE.NE.'MAILLAGE')RETURN
  98. CALL CRCHPT('MSOMMET',MELEMS,1,MCHPOI)
  99. ENDIF
  100.  
  101. CALL KRIPAD(MELEMS,MLENTI)
  102.  
  103. C CREATION DE LA DIAGONALE
  104. CALL LICHTM(MCHPOI,MPOVAL,TYPC,IGEOM)
  105. SEGACT MELEME
  106. NBSOUS=LISOUS(/1)
  107. IF(NBSOUS.EQ.0)NBSOUS=1
  108. C
  109. C BOUCLE SUR LES TYPES D'ELEMENTS ET CALCUL
  110. C
  111. DO 1 KSOUS=1,NBSOUS
  112. IF(NBSOUS.EQ.1)IPT1=MELEME
  113. IF(NBSOUS.GT.1)IPT1=LISOUS(KSOUS)
  114. SEGACT IPT1
  115.  
  116. NP=IPT1.NUM(/1)
  117. NEL=IPT1.NUM(/2)
  118. C
  119. NOM0=NOMS(IPT1.ITYPEL)//' '
  120. IF(INEFMD.EQ.1.AND.IKAS.EQ.4)NOM0=NOMS(IPT1.ITYPEL)//'P1P1'
  121. IF(INEFMD.EQ.2.AND.IKAS.EQ.4)NOM0=NOMS(IPT1.ITYPEL)//'MCF1'
  122. IF(INEFMD.EQ.3.AND.IKAS.EQ.4)NOM0=NOMS(IPT1.ITYPEL)//'PFP1'
  123. c write(6,*)'4 NOM0=',NOM0
  124. CALL KALPBG(NOM0,'FONFORM ',IZFFM)
  125. SEGACT IZFFM*MOD
  126. IZHR=KZHR(1)
  127. SEGACT IZHR*MOD
  128. C
  129. NPG=FN(/2)
  130. NES=GR(/1)
  131.  
  132. IF(IMPR.NE.0)THEN
  133. WRITE(6,*)' SUB CADGSI : NES,NP,NPG,IDIM,NEL='
  134. & ,NES,NP,NPG,IDIM,NEL
  135. ENDIF
  136. C
  137. DO 10 K=1,NEL
  138. C
  139. NPGR=0
  140. IF(IAXI.NE.0)NPGR=NPG
  141. C
  142. DO 12 I=1,NP
  143. J=IPT1.NUM(I,K)
  144. DO 13 N=1,IDIM
  145. XYZ(N,I)=XCOOR((J-1)*(IDIM+1) +N)
  146. 13 CONTINUE
  147. 12 CONTINUE
  148.  
  149. CALL CALJBC(FN,GR,PG,XYZ,HR,PGSQ,RPG,NES,IDIM,NP,NPG,IAXI,AIRE)
  150.  
  151. IF(IMPR.NE.0)THEN
  152. WRITE(6,*)' SUB CADGSI : AIRE=',AIRE
  153. WRITE(6,*)' SUB CADGSI : LER '
  154. WRITE(6,1001)(IPT1.NUM(I,K),I=1,NP)
  155. WRITE(6,*)' SUB CADGSI : XYZ '
  156. WRITE(6,1002)((XYZ(N,I),N=1,2),I=1,NP)
  157. ENDIF
  158. C
  159. DO 3 J=1,NP
  160. SJ=0.D0
  161. DO 4 L=1,NPG
  162. SJ=SJ+FN(J,L)*PGSQ(L)
  163. 4 CONTINUE
  164. C SD(J,K)=SJ
  165. JU=LECT(IPT1.NUM(J,K))
  166. C D0(JU)=D0(JU)+SJ
  167. VPOCHA(JU,1)=VPOCHA(JU,1)+SJ
  168. 3 CONTINUE
  169. C
  170. 10 CONTINUE
  171.  
  172. SEGSUP IZFFM,IZHR
  173. 1 CONTINUE
  174. SEGSUP MLENTI
  175.  
  176. IF(IMPR.NE.0)THEN
  177. WRITE(6,*)' SUB CADGSI : CALCUL DE LA DIAGONALE'
  178. C WRITE(6,1003)(I,VPOCHA(I,1),I=1,NPT)
  179. WRITE(6,*)' FIN DE CADGSI'
  180. ENDIF
  181.  
  182. C FACE
  183. ELSEIF(IKAS.EQ.2)THEN
  184.  
  185. TYPE=' '
  186. CALL ACMO(MTABD,'FACE',TYPE,MELEMF)
  187. IF(TYPE.NE.'MAILLAGE')RETURN
  188. TYPE=' '
  189. CALL ACMO(MTABD,'CENTRE',TYPE,MELEMC)
  190. IF(TYPE.NE.'MAILLAGE')RETURN
  191. CALL LEKTAB(MTABD,'FACEL',MELEME)
  192. IF(MELEME.EQ.0)RETURN
  193. CALL LEKTAB(MTABD,'XXVOLUM',MCHPO1)
  194. IF(MCHPO1.EQ.0)RETURN
  195.  
  196. CALL KRIPAD(MELEMC,MLENT1)
  197. CALL KRIPAD(MELEMF,MLENT2)
  198.  
  199. C CREATION DE LA DIAGONALE
  200. CALL CRCHPT('FACE',MELEMF,1,MCHPOI)
  201. CALL LICHTM(MCHPOI,MPOVAL,TYPC,IGEOM)
  202. CALL LICHTM(MCHPO1,MPOVA1,TYPC,IGEOM)
  203. SEGACT MELEME
  204. C
  205. C BOUCLE SUR LES TYPES D'ELEMENTS ET CALCUL
  206. C
  207. NBEL=NUM(/2)
  208. DO 2 K=1,NBEL
  209. I1=NUM(1,K)
  210. I2=NUM(2,K)
  211. I3=NUM(3,K)
  212. I1=MLENT1.LECT(I1)
  213. I2=MLENT2.LECT(I2)
  214. I3=MLENT1.LECT(I3)
  215. V=(MPOVA1.VPOCHA(I1,1)+MPOVA1.VPOCHA(I3,1) )*0.5D0
  216. VPOCHA(I2,1)=V
  217. 2 CONTINUE
  218.  
  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. SEGACT,MCOORD
  321. DO 21 K=1,NBEL
  322. DO 109 I=1,NP
  323. J=IPT1.NUM(I,K)
  324. JC = (J-1)*(IDIM+1)
  325. DO 110 N=1,IDIM
  326. XYZ(N,I)=XCOOR( JC + N )
  327. 110 CONTINUE
  328. 109 CONTINUE
  329.  
  330. CALL CALJBC(FN,GR,PG,XYZ,HR,PGSQ,RPG,NES,IDIM,NP,NPG,IAXI,AIRE)
  331.  
  332. DO 39 M=1,MP1
  333. M11=LECT(IPT2.NUM(M,K))
  334. M1=M+MP10
  335. c IF(KPOIND.EQ.5)M1=M11
  336.  
  337. U=0.D0
  338. DO 33 LL=1,NPG
  339. U=U+IZF1.FN(M,LL)*PGSQ(LL)
  340. 33 CONTINUE
  341.  
  342. VPOCHA(M1,1)=VPOCHA(M1,1)+U
  343. 39 CONTINUE
  344.  
  345. MP10=MP10+MP1
  346. 21 CONTINUE
  347.  
  348. SEGSUP IZFFM,IZF1,IZHR,IZH2
  349.  
  350. NUTOEL=NUTOEL+NBEL
  351.  
  352. 11 CONTINUE
  353.  
  354. SEGSUP MLENTI
  355.  
  356. ENDIF
  357.  
  358.  
  359. SEGDES MTABD
  360. CALL ACTOBJ('CHPOINT ',MCHPOI,1)
  361. CALL ECROBJ('CHPOINT ',MCHPOI)
  362. RETURN
  363.  
  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.  

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