Télécharger hhosig.eso

Retour à la liste

Numérotation des lignes :

hhosig
  1. C HHOSIG SOURCE OF166741 24/06/19 21:15:10 11942
  2. C HHOSIG SOURCE FANDEUR
  3.  
  4. *----------------------------------------------------------------------*
  5. * CALCUL DES CONTRAINTES POUR LA FORMULATION 'HHO' *
  6. *----------------------------------------------------------------------*
  7. * ENTREES : *
  8. * MATE Numero du materiau *
  9. * IVAMAT Pointeur sur un segment MPTVAL pour le materiau *
  10. * NVMAT Nombre de composantes obligatoires materiau *
  11. * *
  12. * SORTIES : *
  13. *----------------------------------------------------------------------*
  14.  
  15. SUBROUTINE HHOSIG(imoHHO, ichDEP,nmoDEP,
  16. & IIPDPG,UZDPG,RYDPG,RXDPG,
  17. & MATE,IVAMAT,NVMAT, IPMINT, NBPTEL,
  18. & IVASTR,NCSTR, iret)
  19.  
  20. IMPLICIT INTEGER(I-N)
  21. IMPLICIT REAL*8(A-H,O-Z)
  22.  
  23. -INC PPARAM
  24. -INC CCOPTIO
  25. -INC CCREEL
  26.  
  27. -INC CCHHOPA
  28. -INC CCHHOPR
  29.  
  30. -INC SMCOORD
  31. -INC SMCHAML
  32. -INC SMELEME
  33. -INC SMINTE
  34. -INC SMLENTI
  35. POINTEUR mlent4.mlenti
  36. -INC SMLMOTS
  37. -INC SMLREEL
  38. POINTEUR mlrdef.mlreel, mlrdec.mlreel, mlrmbh.mlreel
  39. -INC SMMODEL
  40.  
  41. SEGMENT MPTVAL
  42. INTEGER IPOS(NS),NSOF(NS)
  43. INTEGER IVAL(NCOSOU)
  44. CHARACTER*16 TYVAL(NCOSOU)
  45. ENDSEGMENT
  46.  
  47. SEGMENT MWKMAT
  48. REAL*8 VALMAT(NTMAT), HOELAS(lhook,lhook)
  49. c** REAL*8 XLOC(3,3),XGLOB(3,3),TXR(IDIM,IDIM)
  50. ENDSEGMENT
  51.  
  52. SEGMENT MWKHHO
  53. INTEGER TABINT(NBINT)
  54. REAL*8 TABFLO(NBFLO)
  55. ENDSEGMENT
  56.  
  57. DIMENSION UDPGE(3)
  58.  
  59. iret = 0
  60.  
  61. imodel = imoHHO
  62. c* segact,imodel <- actif en entree/sortie
  63.  
  64. C- Premieres verifications :
  65. CALL HHONOB(imoHHO, nobHHO, iret)
  66. IF (nobHHO.LE.0) THEN
  67. write(ioimp,*) 'HHOSIG: IMODEL incorrect (not HHO)'
  68. iret = 5
  69. return
  70. END IF
  71.  
  72. C- Introduction du point autour duquel se fait le mouvement
  73. C de la section en defo plane generalisee
  74. C IIPDPG = numero du noeud/point support si defini pour le modele
  75. C NDPGE > 0 si prise en compte du point support
  76. IF (IIPDPG.GT.0) THEN
  77. write(ioimp,*) 'HHOSIG: 2D PLAN GENE mode not implemented!'
  78. iret = 5
  79. return
  80. IF (IFOUR.EQ.-3) THEN
  81. NDPGE = 3
  82. UDPGE(1) = UZDPG
  83. UDPGE(2) = RYDPG
  84. UDPGE(3) = RXDPG
  85. C SEGACT,MCOORD
  86. IREF = (IIPDPG-1)*(IDIM+1)
  87. XDPGE = mcoord.XCOOR(IREF+1)
  88. YDPGE = mcoord.XCOOR(IREF+2)
  89. ELSE IF (IFOUR.EQ.11) THEN
  90. NDPGE = 2
  91. UDPGE(1) = UZDPG
  92. UDPGE(2) = RXDPG
  93. UDPGE(3) = XZero
  94. XDPGE = XZero
  95. YDPGE = XZero
  96. ELSE IF (IFOUR.EQ. 7 .OR. IFOUR.EQ. 8 .OR. IFOUR.EQ. 9 .OR.
  97. & IFOUR.EQ.10 .OR. IFOUR.EQ.14) THEN
  98. NDPGE = 1
  99. UDPGE(1) = UZDPG
  100. UDPGE(2) = XZero
  101. UDPGE(3) = XZero
  102. XDPGE = XZero
  103. YDPGE = XZero
  104. else
  105. write(ioimp,*) 'HHOSIG: ERREUR NDPGE'
  106. iret = 5
  107. return
  108. END IF
  109. ELSE
  110. NDPGE = 0
  111. UDPGE(1) = XZero
  112. UDPGE(2) = XZero
  113. UDPGE(3) = XZero
  114. XDPGE = XZero
  115. YDPGE = XZero
  116. END IF
  117.  
  118. C- Recuperation des donnees du modele en entree
  119. c* MELE = imodel.NEFMOD
  120. c* MFR = imodel.infele(13)
  121. meleme = imodel.IMAMOD
  122. NBNOE = meleme.NUM(/1)
  123. NBELT = meleme.NUM(/2)
  124.  
  125. mlenti = imodel.IVAMOD(nobHHO+1)
  126. c* segact,mlenti
  127. mlent3 = imodel.IVAMOD(nobHHO+3)
  128. c* segact,mlent3
  129. mlent4 = imodel.IVAMOD(nobHHO+4)
  130. c* segact,mlent4
  131.  
  132. c* n_o_face = mlenti.lect(2)
  133. n_d_face = mlenti.lect(3)
  134. c* n_o_cell = mlenti.lect(4)
  135. n_d_cell = mlenti.lect(5)
  136. nb_faces = mlenti.lect(7)
  137. NBPGAU = mlenti.lect(8)
  138. idifo = mlenti.lect(9)
  139. NBDDL = mlenti.lect(11)
  140. NDEPF = idifo * n_d_face
  141. NDEPC = idifo * n_d_cell
  142. lhook = 9
  143.  
  144. nbel3 = mlent3.lect(/1) / 2
  145. nbfa3 = 2 * nb_faces
  146. nbel4 = mlent4.lect(/1) / 2
  147.  
  148. IF (mlenti.lect(6).NE.NBNOE) THEN
  149. write(ioimp,*) 'HHOSIG: Bizarre nb_vertices'
  150. END IF
  151. IF (NBDDL .NE. imodel.INFELE(9)) then
  152. write(ioimp,*) 'HHOSIG: Bizarre NBDDL'
  153. END IF
  154. c NBPGAU =? (NBPTEL = imodel.INFELE(4))
  155. IF (NBPGAU .NE. imodel.INFELE(4)) then
  156. write(ioimp,*) 'HHOSRIG: Bizarre nb.p.gau(1)'
  157. END IF
  158. c NBPGAU =? minte.POIGAU(/1)
  159. minte = IPMINT
  160. c* SEGACT minte <- actif en E/S
  161. if (NBPGAU .NE. minte.POIGAU(/1)) then
  162. write(ioimp,*) 'HHOSIG: Bizarre nb.p.gau (2)'
  163. end if
  164. c-dbg write(ioimp,*) 'HHOSIG:',nb_faces,NBDDL,NBPGAU
  165. if (nbel3.NE.(NBELT*nb_faces)) then
  166. write(ioimp,*) 'HHOSIG: Bizarre nbel3'
  167. end if
  168. if (nbel4.NE.NBELT) then
  169. write(ioimp,*) 'HHOSIG: Bizarre nbel4'
  170. end if
  171.  
  172. ivid = 1
  173. C- Deplacements des faces et des cellules :
  174. nomid = nmoDEP
  175. c* segact,nomid
  176. JGN = nomid.lesobl(/1)
  177. c-dbg write(ioimp,*) 'HHOSIG=',NDEPC,NDEPF,lesobl(/2),lesfac(/2)
  178. nfac = 0
  179. C Deplacements des cellules - Points supports des cellules
  180. JGM = NDEPC + nfac
  181. SEGINI,mlmots
  182. DO i = 1, NDEPC
  183. mlmots.MOTS(i)(1:JGN) = nomid.lesobl(i)(1:JGN)
  184. END DO
  185. c* DO i = 1, nfac
  186. c* mlmots.MOTS(NDEP+i)(1:JGN) = nomid.lesfac(i)(1:JGN)
  187. c* END DO
  188. CALL EXTR23(ichDEP,mlmots,MPCHHO,mlrDEC,ivid)
  189. SEGACT,mlrDEC
  190. c* IF (IERR.NE.0) THEN
  191. c* iret = 21
  192. c* return
  193. c* END IF
  194.  
  195. c-dbg write(ioimp,*) 'mlmots DEPC ',(mots(i),i=1,mots(/2))
  196. c-dbg write(ioimp,*) 'U.CELL',mlrdec.prog(/1),NCEHHO,NDEPC
  197. c-dbg write(ioimp,*) ' ',(mlrdec.prog(i),i=1,mlrdec.prog(/1))
  198. c-dbg write(ioimp,*)
  199.  
  200. C Deplacements des faces - Points supports des faces
  201. JGM = NDEPF + nfac
  202. SEGADJ,mlmots
  203. DO i = 1, NDEPF
  204. mlmots.MOTS(i)(1:JGN) = nomid.lesobl(NDEPC+i)(1:JGN)
  205. END DO
  206. c* DO i = 1, nfac
  207. c* mlmots.MOTS(NDEP+i)(1:JGN) = nomid.lesfac(i)(1:JGN)
  208. c* END DO
  209. CALL EXTR23(ichDEP,mlmots,MPFHHO,mlrDEF,ivid)
  210. SEGACT,mlrDEF
  211. c* IF (IERR.NE.0) THEN
  212. c* iret = 21
  213. c* return
  214. c* END IF
  215.  
  216. c-dbg write(ioimp,*) 'MLMOTS DEPF ',(mots(i),i=1,mots(/2))
  217. c-dbg write(ioimp,*) 'U.FACE',mlrDEF.prog(/1),NFAHHO,NDEPF
  218. c-dbg write(ioimp,*) ' ',(mlrDEF.prog(i),i=1,mlrDEF.prog(/1))
  219. c-dbg write(ioimp,*)
  220.  
  221. C- Verification des proprietes materielles : Champ CONSTANT par ELEMENT !
  222. mptval = IVAMAT
  223. NTMAT = mptval.IVAL(/1)
  224. NUMAT = NTMAT - 4
  225.  
  226. C- Tableau de stockage des proprietes materielles
  227. MWKMAT = 0
  228. SEGINI,MWKMAT
  229. IWKMAT = MWKMAT
  230.  
  231. isz = 0
  232. iunif = 0
  233. DO is = 1, NTMAT
  234. IF (is.LE.(NVMAT-4) .OR. is.GT.NVMAT) THEN
  235. melval = mptval.IVAL(is)
  236. IGMN = melval.VELCHE(/1)
  237. IEMN = melval.VELCHE(/2)
  238. IF (IGMN.NE.1) isz = isz + 1
  239. IF (IGMN.EQ.1 .AND. IEMN.EQ.1) iunif = iunif + 1
  240. VALMAT(is) = melval.VELCHE(1,1)
  241. END IF
  242. END DO
  243. IF (isz.NE.0) THEN
  244. write(ioimp,*) 'HHOSIG: Material properties must be constant ',
  245. & 'by element (for the moment) !'
  246. iret = 21
  247. RETURN
  248. END IF
  249. * Le champ de proprietes est completement uniforme :
  250. IF (iunif.NE.NUMAT) iunif = 0
  251.  
  252. C- Materiau isotrope et proprietes uniformes : -> calcul de HOELAS
  253. IF (MATE.EQ.1 .AND. iunif.NE.0) THEN
  254. CALL DOHMAS(VALMAT,'ISOTROPE',IFOUR,lhook,1,HOELAS,iretc)
  255. r_z = HOELAS(4,4) * 2.D0
  256. DO i = 4, lhook
  257. HOELAS(i,i) = r_z
  258. END DO
  259. END IF
  260. C - Materiaux orthotrope et anisotrope
  261. IF (MATE.EQ.2.OR.MATE.EQ.3) THEN
  262. C** AFAIRE !!!
  263. END IF
  264.  
  265. C- Les donnees liees a la formulation HHO :
  266. IVMBHO = mptval.IVAL(NVMAT-1)
  267. melval = IVMBHO
  268. IGMB = melval.IELCHE(/1)
  269. IEMB = melval.IELCHE(/2)
  270. c-dbg write(ioimp,*) 'IVMBHO',melval,igmb,iemb,tyval(nvmat-1)
  271. mlrmbh = melval.IELCHE(1,1)
  272. c* segact,mlrmbh
  273. c* write(ioimp,*) 'HHOSIG MBHHO SIZE:',mlrmbh.prog(/1),
  274. c* & NBDDL,9*NBDDL,mlenti.lect(14)
  275. IF ((mlrmbh.prog(/1).NE.(9*NBDDL)) .OR.
  276. & (mlrmbh.prog(/1).NE.mlenti.lect(14))) THEN
  277. write(ioimp,*) 'HHOSIG: BHHO matrix size incorrect'
  278. iret = 21
  279. RETURN
  280. END IF
  281.  
  282. C- Indices et tableau de travail
  283. ir_coo = 0
  284. c* si besoin des coordonnees ir_eps = ir_coo + (IDIM*NBNOE)
  285. ir_eps = ir_coo + 0
  286. ir_sig = ir_eps + lhook
  287. ir_uce = ir_sig + lhook
  288. ir_ufa = ir_uce + NDEPC
  289. ir_uge = ir_ufa + (NDEPF*nb_faces)
  290. ir_fin = ir_uge + NDPGE
  291.  
  292. NBINT = 1
  293. NBFLO = ir_fin
  294. SEGINI,MWKHHO
  295.  
  296. IF (NDPGE.GT.0) THEN
  297. DO ic = 1, NDPGE
  298. TABFLO(ir_uge+ic) = UDPGE(ic)
  299. END DO
  300. END IF
  301. c* si besoin des coordonnees SEGACT,MCOORD*NOMOD
  302.  
  303. C-------------------------
  304. C Boucle sur les cellules du modele
  305. C-------------------------
  306. DO IEL = 1,NBELT
  307.  
  308. C- Recuperation des coordonnees des noeuds de l element IEL
  309. c* CALL HHOCOO(meleme.num,NBNOE, mcoord.xcoor, IEL,
  310. c* & TABFLO(ir_coo+1), iret)
  311. c* IF (iret.NE.0) RETURN
  312.  
  313. C- Valeurs des inconnues primales pour l'element IEL (cell+faces)
  314. in1 = IEL * 2
  315. je = mlent4.lect(in1-1)
  316. ip = mlent4.lect(in1)
  317. if (ip.le.0) write(ioimp,*) 'HHOSIG ICEL Bizarre...',iel,je,ip
  318. jp = ip + NBCHHO(je-1)
  319. DO ic = 1, NDEPC
  320. jc = NCEHHO * (ic - 1)
  321. TABFLO(ir_uce + ic) = mlrDEC.prog(jp + jc)
  322. END DO
  323. c-dbg write(ioimp,*) 'HHOSIG ICEL',iel,je,ip,jp,ir_uce
  324. c-dbg write(ioimp,*) (TABFLO(ir_uce+ic),ic=1,ndepc)
  325.  
  326. in1 = (IEL-1) * nbfa3
  327. ir_kc = ir_ufa
  328. DO j1 = 1, nb_faces
  329. in2 = in1 + (2 * j1)
  330. je = mlent3.lect(in2-1)
  331. ip = ABS(mlent3.lect(in2))
  332. if (ip.eq.0) write(ioimp,*) 'HHOSIG IFAE Bizarre...',iel,j1,je,ip
  333. jp = ip + NBFHHO(je-1)
  334. DO ic = 1, NDEPF
  335. jc = NFAHHO * (ic - 1)
  336. TABFLO(ir_kc + ic) = mlrDEF.prog(jp + jc)
  337. END DO
  338. c-dbg write(ioimp,*) 'HHOSIG IFAE',iel,j1,je,ip,jp,ir_kc
  339. c-dbg write(ioimp,*) (TABFLO(ir_kc+ic),ic=1,ndepf)
  340. ir_kc = ir_kc + NDEPF
  341. END DO
  342.  
  343. c-dbg write(ioimp,*) 'HHOSIG ...',ir_uce,ir_kc,NBDDL
  344.  
  345. C- Recuperation des proprietes constantes par element :
  346. IF (iunif.EQ.0) THEN
  347. DO is = 1, NVMAT
  348. IF (is.LE.(NVMAT-4) .OR. is.GT.NVMAT) THEN
  349. melval = mptval.IVAL(is)
  350. VALMAT(is) = melval.VELCHE(1,IEL)
  351. END IF
  352. END DO
  353. END IF
  354. C- Materiau isotrope
  355. IF (MATE.EQ.1 .AND. iunif.EQ.0) THEN
  356. CALL DOHMAS(VALMAT,'ISOTROPE',IFOUR,lhook,1,HOELAS,iretc)
  357. r_z = HOELAS(4,4) * 2.D0
  358. DO i = 4, lhook
  359. HOELAS(i,i) = r_z
  360. END DO
  361. ENDIF
  362. C - Calcul des axes locaux dans les cas orthotrope et anisotrope
  363. IF (MATE.EQ.2.OR.MATE.EQ.3) THEN
  364. ENDIF
  365.  
  366. JEMB = MIN(IEL,IEMB)
  367. C-- -- -- -- -- -- -- -- --
  368. C - Boucle sur les points de Gauss
  369. C-- -- -- -- -- -- -- -- --
  370. DO IGAU = 1, NBPGAU
  371.  
  372. C-AFAIRE C -- Recuperation des proprietes materielles (IGAU)
  373. C-AFAIRE DO is = 1, NVMAT
  374. C-AFAIRE IF (is.LE.(NVMAT-4) .OR. is.GT.NVMAT) THEN
  375. C-AFAIRE melval = mptval.IVAL(is)
  376. C-AFAIRE JEMN = MIN(IEL ,melval.VELCHE(/2))
  377. C-AFAIRE JGMN = MIN(IGAU,melval.VELCHE(/1))
  378. C-AFAIRE VALMAT(is) = melval.VELCHE(JGMN,JEMN)
  379. C-AFAIRE END IF
  380. C-AFAIRE END DO
  381. C-AFAIREC- Materiau isotrope
  382. C-AFAIRE IF (MATE.EQ.1) THEN
  383. C-AFAIRE CALL DOHMAS(VALMAT,'ISOTROPE',IFOUR,lhook,1,HOELAS,iretc)
  384. C-AFAIRE r_z = HOELAS(4,4) * 2.D0
  385. C-AFAIRE DO i = 4, lhook
  386. C-AFAIRE HOELAS(i,i) = r_z
  387. C-AFAIRE END DO
  388. C-AFAIRE END IF
  389.  
  390. melval = IVMBHO
  391. JGMB = MIN(IGAU,IGMB)
  392. mlrmbh = melval.IELCHE(JGMB,JEMB)
  393. c* segact,mlrmbh
  394. c* !! matrice BHHO stockee colonne par colonne : lhook*NBDDL
  395. DO ic = 1, lhook
  396. r_z = XZero
  397. DO jc = 1, NBDDL
  398. jnc = lhook * (jc-1)
  399. r_z = r_z + mlrmbh.prog(ic + jnc) * TABFLO(ir_uce + jc)
  400. END DO
  401. TABFLO(ir_eps + ic) = r_z
  402. END DO
  403. DO ic = 1, lhook
  404. r_z = XZero
  405. DO jc = 1, lhook
  406. r_z = r_z + HOELAS(ic,jc) * TABFLO(ir_eps + jc)
  407. END DO
  408. TABFLO(ir_sig + ic) = r_z
  409. END DO
  410.  
  411. C -- Remplissage du segment contenant les "contraintes"
  412. c* A voir selon ordre - permutation :
  413. mptval = IVASTR
  414. DO ic = 1, NCSTR
  415. melval = mptval.IVAL(ic)
  416. melval.velche(IGAU,IEL) = TABFLO(ir_sig + ic)
  417. END DO
  418.  
  419. C-- -- -- -- -- -- -- -- --
  420. END DO
  421. C-- -- -- -- -- -- -- -- --
  422. C-------------------------
  423. END DO
  424. C-------------------------
  425.  
  426. IF (MATE.EQ.2.OR.MATE.EQ.3) THEN
  427. ENDIF
  428. SEGSUP,MWKMAT,MWKHHO
  429. SEGSUP,mlrDEC,mlrDEF
  430. SEGSUP,mlmots
  431.  
  432. c* RETURN
  433. END
  434.  
  435.  
  436.  
  437.  

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