Télécharger epsi2.eso

Retour à la liste

Numérotation des lignes :

  1. C EPSI2 SOURCE BP208322 16/11/18 21:16:46 9177
  2.  
  3. SUBROUTINE EPSI2(IPMAIL,IPMINT,MELE,IELE,
  4. & IVADEP,NBPTEL,LRE,NSTRS,LHOOK,
  5. & MFR,NDEP,IPORE,IREPS2,NBPGAU,IVAEPS,UZDPG,RYDPG,RXDPG,IIPDPG,
  6. & IDERI,ivamat,ivade2,mate,nmatt,cmate,ngra,noer,kerr)
  7.  
  8. C---------------------------------------------------------------------*
  9. C
  10. C calcul des deformations
  11. C
  12. C massif, poreux, joints poreux, incompressibles
  13. C---------------------------------------------------------------------*
  14. C *
  15. C entrees : *
  16. C ________ *
  17. C *
  18. C ipmail pointeur sur un segment meleme *
  19. C ipmint pointeur sur un segment minte *
  20. C mele numero de l'element fini *
  21. C iele numero geometrique de l'element *
  22. C nbpgau nombre de point d'integration pour la rigidite *
  23. C ivadep pointeur sur le chamelem de deplacements *
  24. C nbptel nombre de points par element *
  25. C lre nombre de ddl dans la matrice de rigidite *
  26. C nstrs nombre de composante de contraintes/deformations *
  27. C pour une matrice de hooke *
  28. C lhook dimension de la matrice de hooke *
  29. C mfr numero de la formulation de l'element fini *
  30. C ndep nombre de composantes de deplacements *
  31. C ipore nombre de fonctions de forme *
  32. C iresp2 flag pour indiquer si on veut les contraintes *
  33. C de piola-kirchhoff *
  34. C uzdpg = deformation au point nsdpge support de la *
  35. C rydpf = deformation plane generalisee *
  36. C rxdpg = *
  37. C *
  38. C sorties : *
  39. C ________ *
  40. C *
  41. C ivaeps pointeur sur un segment mptval contenant les *
  42. C les melvals de deformations *
  43. C---------------------------------------------------------------------*
  44. C Pour MEMOIRE : Si MELE element incompressible alors MFR = 31
  45. C---------------------------------------------------------------------*
  46.  
  47. IMPLICIT INTEGER(I-N)
  48. IMPLICIT REAL*8(A-H,O-Z)
  49.  
  50. -INC CCOPTIO
  51. -INC CCREEL
  52. -INC CCHAMP
  53. -INC CCGEOME
  54. -INC SMCHAML
  55. -INC SMCHPOI
  56. -INC SMELEME
  57. -INC SMCOORD
  58. -INC SMMODEL
  59. -INC SMINTE
  60. -INC SMLREEL
  61.  
  62. SEGMENT WRK1
  63. REAL*8 DDHOOK(NSTRS,NSTRS),XDDL(LRE),XSTRS(NSTRS)
  64. REAL*8 XE(3,NBBB),DDHOMU(NSTRS,NSTRS)
  65. REAL*8 SHPWRK(6,NBNO),BGENE(LHOOK,LRE)
  66. REAL*8 XE1(3,NBBB),XE2(3,NBBB),xstrs2(NSTRS)
  67. REAL*8 xjac(3,3),valmat(20)
  68. ENDSEGMENT
  69. SEGMENT WRK2
  70. REAL*8 BGR(NGRA,LRE),BB(2,NGRA),gradi(ngra),R(ngra),u(ngra)
  71. REAL*8 TENS(9),tentra(9),xddls2(lre)
  72. ENDSEGMENT
  73. C
  74. SEGMENT WRK3
  75. REAL*8 BPSS(3,3),XEL(3,NBBB)
  76. REAL*8 XNTH(LPP,LPP),XNTB(LPP,LPP),XNTT(LPP)
  77. ENDSEGMENT
  78. C
  79. SEGMENT WRK5
  80. REAL*8 XGENE(NSTN,LRN)
  81. ENDSEGMENT
  82. C
  83. SEGMENT NOTYPE
  84. CHARACTER*16 TYPE(NBTYPE)
  85. ENDSEGMENT
  86. C
  87. SEGMENT MPTVAL
  88. INTEGER IPOS(NS),NSOF(NS)
  89. INTEGER IVAL(NCOSOU)
  90. CHARACTER*16 TYVAL(NCOSOU)
  91. ENDSEGMENT
  92. C
  93. SEGMENT MTRACE
  94. REAL*8 TRACE(NBPTEL)
  95. ENDSEGMENT
  96.  
  97. CHARACTER*(NCONCH) CONM
  98. PARAMETER (NINF=3)
  99. DIMENSION A(4,60),BBX(3,60),UDPGE(3)
  100. INTEGER INFOS(NINF)
  101. CHARACTER*8 CMATE
  102.  
  103.  
  104. DIMENSION IN(6),JN(6),ITAB(3,3)
  105. C
  106. DATA IN/1,2,3,1,1,2/
  107. DATA JN/1,2,3,2,3,3/
  108. C
  109. DATA ITAB(1,1),ITAB(1,2),ITAB(1,3)/1,4,5/
  110. DATA ITAB(2,1),ITAB(2,2),ITAB(2,3)/4,2,6/
  111. DATA ITAB(3,1),ITAB(3,2),ITAB(3,3)/5,6,3/
  112. C
  113.  
  114. kerr=0
  115. C Introduction du point autour duquel se fait le mouvement
  116. C de la section en defo plane generalisee
  117. C IIPDPG = numero du noeud/point support si defini pour le modele
  118. C NDPGE > 0 si prise en compte du point support
  119. IF (IIPDPG.GT.0) THEN
  120. IF (IFOUR.EQ.-3) THEN
  121. NDPGE=3
  122. UDPGE(1)=UZDPG
  123. UDPGE(2)=RYDPG
  124. UDPGE(3)=RXDPG
  125. C SEGACT,MCOORD
  126. IREF=(IIPDPG-1)*(IDIM+1)
  127. XDPGE=XCOOR(IREF+1)
  128. YDPGE=XCOOR(IREF+2)
  129. ELSE IF (IFOUR.EQ.11) THEN
  130. NDPGE=2
  131. UDPGE(1)=UZDPG
  132. UDPGE(2)=RXDPG
  133. UDPGE(3)=XZero
  134. XDPGE=XZero
  135. YDPGE=XZero
  136. ELSE IF (IFOUR.EQ. 7 .OR. IFOUR.EQ. 8 .OR. IFOUR.EQ. 9 .OR.
  137. & IFOUR.EQ.10 .OR. IFOUR.EQ.14) THEN
  138. NDPGE=1
  139. UDPGE(1)=UZDPG
  140. UDPGE(2)=XZero
  141. UDPGE(3)=XZero
  142. XDPGE=XZero
  143. YDPGE=XZero
  144. else
  145. write(ioimp,*) 'EPSI2 : ERREUR NDPGE'
  146. call erreur(5)
  147. return
  148. ENDIF
  149. ELSE
  150. NDPGE=0
  151. UDPGE(1)=UZDPG
  152. UDPGE(2)=XZero
  153. UDPGE(3)=XZero
  154. XDPGE=XZero
  155. YDPGE=XZero
  156. ENDIF
  157. C
  158.  
  159. MELEME=IPMAIL
  160. NBNN=NUM(/1)
  161. NBELEM=NUM(/2)
  162. C
  163. NHRM=NIFOUR
  164. MINTE=IPMINT
  165. NBBB=NBNN
  166. C
  167. C_______________________________________________________________________
  168. C
  169. C numero des etiquettes :
  170. C etiquettes de 1 a 98 pour traitement specifique a l element
  171. C dans la zone specifique a chaque element commencant par :
  172. C 5 continue
  173. C element 5 etiquettes 1005 2005 3005 4005 ...
  174. C 44 continue
  175. C element 44 etiquettes 1044 2044 3044 4044 ...
  176. C_______________________________________________________________________
  177. C
  178. GOTO (99,99,99, 4,99, 4,99, 4,99, 4,99,99,99, 4, 4, 4, 4,99,99,99,
  179. 1 99,99, 4, 4, 4, 4,99,99,99,99,99,99,99,99,99,99,99,99,99,99,
  180. 2 99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,
  181. 3 99,99,99,99,99,99,99,99, 4, 4, 4, 4, 4,99,99, 4,99,99,79,79,
  182. 4 79,79,79,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,
  183. 5 99,99,99,99,99,99,99,80,80,80, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4,
  184. 6 4, 4),MELE
  185. C
  186. IF (MELE.EQ.183.OR.MELE.EQ.184.OR.
  187. . MELE.EQ.193.OR.MELE.EQ.194) GOTO 4
  188. IF (MELE.GE.173.AND.MELE.LE.182) GOTO 173
  189. IF (MELE.GE.185.AND.MELE.LE.190)GOTO 185
  190. C
  191. GOTO 99
  192.  
  193. C_______________________________________________________________________
  194. C
  195. C elements massifs et elements incompressibles MECANIQUE
  196. C_______________________________________________________________________
  197. C
  198. 4 CONTINUE
  199. IF (MFR.EQ.71 .OR. MFR.EQ.73) GOTO 97173
  200. IF( ideri.le.2.or.ideri.eq.5 ) then
  201. C ideri le 2 est pour lineaire et quadratique et =5 pour utlisateur
  202.  
  203. C Elements massifs en FORMULATION 'MECANIQUE'
  204. C -------------------------------------------
  205. NBNO=NBNN
  206. NDDD=NDEP-NDPGE
  207. C
  208. C Donnees liees a l'element de reference
  209. C
  210. SEGACT MINTE
  211. C
  212. SEGINI,WRK1
  213. IF (Ideri.eq.2) SEGINI,MTRACE
  214. C
  215. C boucle sur les elements
  216. C
  217. DO 3004 IB=1,NBELEM
  218. C
  219. C on cherche les deplacements
  220. C
  221. MPTVAL=IVADEP
  222. IE=1
  223. DO 4004 IGAU=1,NBNN
  224. DO 4004 ICOMP=1,NDDD
  225. MELVAL=IVAL(ICOMP)
  226. IGMN=MIN(IGAU,VELCHE(/1))
  227. IBMN=MIN(IB ,VELCHE(/2))
  228. XDDL(IE)=VELCHE(IGMN,IBMN)
  229. IE=IE+1
  230. 4004 CONTINUE
  231. IF (NDPGE.GT.0) THEN
  232. DO i=1,NDPGE
  233. XDDL(IE)=UDPGE(i)
  234. IE=IE+1
  235. ENDDO
  236. ENDIF
  237. C
  238. C on cherche les coordonnees des noeuds de l element ib
  239. C
  240. CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  241. if( ideri.eq.5) then
  242. C on se emet à mi-pas
  243. do iou=1,idim
  244. do iyu=1,nbnn
  245. XE(iou,iyu)= xe(iou,iyu) + xddl( iou+ (iyu-1)*nddd)/2.D0
  246. enddo
  247. enddo
  248. endif
  249. C
  250. C Calcul des coeff de modification de b-barre (elts incompres)
  251. C
  252. IF (MFR.EQ.31) THEN
  253. CALL BBCALC(XE,MELE,NBNN,IDIM,NBPGAU,
  254. 1 POIGAU,QSIGAU,ETAGAU,DZEGAU,NSTRS,
  255. 2 LRE,IFOUR,A,BBx,NHRM,SHPTOT,SHPWRK,XDPGE,YDPGE)
  256. ENDIF
  257. C
  258. C boucle sur les points de gauss
  259. C
  260. ISDJC=0
  261. C
  262. DO 5004 IGAU=1,NBPTEL
  263. C
  264. CALL BMATST(IGAU,NBPGAU,POIGAU,QSIGAU,ETAGAU,DZEGAU,
  265. 1 MELE,MFR,NBNN,LRE,IFOUR,NSTRS,NHRM,1.D0,XE,
  266. 2 SHPTOT,SHPWRK,BGENE,DJAC,XDPGE,YDPGE)
  267.  
  268. IF (DJAC.EQ.0.D0) THEN
  269. INTERR(1)=IB
  270. CALL ERREUR(259)
  271. GOTO 9904
  272. ENDIF
  273. IF (DJAC.LT.0.D0) ISDJC=ISDJC+1
  274.  
  275. C En cas d'elements incompressibles : BGENE selon la methode B-BARRE
  276. IF (MFR.EQ.31) THEN
  277. CALL BBAR(IGAU,NBPGAU,POIGAU,QSIGAU,ETAGAU,DZEGAU,
  278. & MELE,NBNN,LRE,IFOUR,NSTRS,XE,DJAC,A,BBX,BGENE)
  279. ENDIF
  280. C
  281. CALL BST(BGENE,XDDL,LRE,NSTRS,XSTRS)
  282. C
  283. C calcul des eps 2
  284. C
  285. IF (Ideri.eq.2)
  286. & CALL BST2(SHPWRK,XDDL,XE,NBNO,IFOUR,XSTRS,TRACE,
  287. & IGAU,XDPGE,YDPGE,UDPGE,NHRM)
  288. C
  289. C remplissage du segment contenant les deformations
  290. C
  291. MPTVAL=IVAEPS
  292. DO 7004 ICOMP=1,NSTRS
  293. MELVAL=IVAL(ICOMP)
  294. IBMN=MIN(IB,VELCHE(/2))
  295. VELCHE(IGAU,IBMN)=XSTRS(ICOMP)
  296. 7004 CONTINUE
  297. C
  298. 5004 CONTINUE
  299. C
  300. C fin de la boucle sur les points de gauss
  301. C
  302. IF (ISDJC.NE.0.AND.ISDJC.NE.NBPGAU) THEN
  303. if (noer.eq.1) then
  304. kerr=195
  305. else
  306. INTERR(1)=IB
  307. CALL ERREUR(195)
  308. endif
  309. GOTO 9904
  310. ENDIF
  311. C
  312. C correction sur la partie quadratique de la deformation dans le cas
  313. C des elements incompressibles
  314. C
  315. IF (Ideri.eq.2) THEN
  316. IF (MFR.EQ.31) THEN
  317. CALL BBST2(TRACE,NBPTEL,IFOUR,MELE,POIGAU,QSIGAU,
  318. & ETAGAU,DZEGAU,SHPTOT,NBNO,SHPWRK,XE)
  319. L=2
  320. IF (IDIM.EQ.3 .OR. IFOUR.EQ.0) L=3
  321. DO 5005 ICOMP=1,L
  322. MELVAL=IVAL(ICOMP)
  323. IBMN=MIN(IB ,VELCHE(/2))
  324. DO 5006 IGAU=1,NBPTEL
  325. VELCHE(IGAU,IBMN)=VELCHE(IGAU,IBMN)+TRACE(IGAU)
  326. 5006 CONTINUE
  327. 5005 CONTINUE
  328. ENDIF
  329. ENDIF
  330.  
  331. 3004 CONTINUE
  332.  
  333. C
  334. C fin de la boucle sur les elements
  335. C
  336. 9904 CONTINUE
  337. SEGDES MINTE
  338. SEGSUP WRK1
  339. IF (IREPS2.EQ.1) SEGSUP MTRACE
  340.  
  341. C cas de la dérivée de Truesdell
  342. elseif( ideri.eq.3) then
  343. NBNO=NBNN
  344. NDDD=NDEP-NDPGE
  345. SEGACT MINTE
  346. SEGINI,WRK1
  347. C IF (IREPS2.EQ.1) SEGINI,MTRACE
  348. C
  349. C boucle sur les elements
  350. C
  351. DO 34 IB=1,NBELEM
  352. C
  353. C on cherche les deplacements
  354. C
  355. MPTVAL=IVADEP
  356. IE=1
  357. DO 44 IGAU=1,NBNN
  358. DO 44 ICOMP=1,NDDD
  359. MELVAL=IVAL(ICOMP)
  360. IGMN=MIN(IGAU,VELCHE(/1))
  361. IBMN=MIN(IB ,VELCHE(/2))
  362. XDDL(IE)=VELCHE(IGMN,IBMN)
  363. IE=IE+1
  364. 44 CONTINUE
  365. IF (NDPGE.GT.0) THEN
  366. DO i=1,NDPGE
  367. XDDL(IE)=UDPGE(i)
  368. IE=IE+1
  369. ENDDO
  370. ENDIF
  371. C on cherche le max des variations des champs pour savoir s'il faut
  372. C appeler hookis plusieurs fois
  373. C
  374. nbgmat=0
  375. nelmat=0
  376. mptval=ivamat
  377. segact mptval
  378. DO IM=1,NMATT
  379. IF (IVAL(IM).NE.0) THEN
  380. MELVAL=IVAL(IM)
  381. nelmat=Max(nelmat ,VELCHE(/2))
  382. nbgmat=Max(nbgmat,VELCHE(/1))
  383. ENDIF
  384. ENDDO
  385. segdes mptval
  386. C
  387. C on cherche les coordonnees des noeuds de l element ib
  388. C
  389. CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  390. C on ajoute aux coordonnées la moitié du deplacements pour faire
  391. C la configuration à mi-pas
  392. do iou=1,idim
  393. do iyu=1,nbnn
  394. XE(iou,iyu)= xe(iou,iyu) + xddl( iou+ (iyu-1)*nddd)/2.D0
  395. enddo
  396. enddo
  397. C
  398. C boucle sur les points de gauss
  399. C
  400. ISDJC=0
  401. C
  402.  
  403. DO 54 IGAU=1,NBPTEL
  404. CALL BMATST(IGAU,NBPGAU,POIGAU,QSIGAU,ETAGAU,DZEGAU,
  405. 1 MELE,MFR,NBNN,LRE,IFOUR,NSTRS,NHRM,1.D0,XE,
  406. 2 SHPTOT,SHPWRK,BGENE,DJAC,XDPGE,YDPGE)
  407.  
  408. C
  409. IF (DJAC.EQ.0.D0) THEN
  410. INTERR(1)=IB
  411. CALL ERREUR(259)
  412. GOTO 994
  413. ELSE IF (DJAC.LT.0.D0) THEN
  414. ISDJC=ISDJC+1
  415. ENDIF
  416. C
  417. C on cherche les matrices de Hooke
  418. C
  419. mptval=ivamat
  420. segact mptval
  421. DO 94 IM=1,NMATT
  422. IF (IVAL(IM).NE.0) THEN
  423. MELVAL=IVAL(IM)
  424. IBMN=MIN(IB ,VELCHE(/2))
  425. IGMN=MIN(IGAU,VELCHE(/1))
  426.  
  427. VALMAT(IM)=VELCHE(IGMN,IBMN)
  428. ELSE
  429. VALMAT(IM)=0.D0
  430. ENDIF
  431. 94 CONTINUE
  432. segdes mptval
  433. kcas=2
  434. if(nbgmat+nelmat.gt.2 . or . ib+igau.eq.2) then
  435. CALL HOOKIS(VALMAT,VALCAR,VAR,MFR,IB,IGAU,EXCEN,EPAIST,
  436. + INAT,MELE,NPINT,IFOUR,KCAS,NBGMAT,Nelmat,
  437. + S,SECT,LHOOK,DDHOMU,DDHOOK,
  438. + COBMA,XMOB,IRETOU)
  439. endif
  440. do iou=1,nstrs
  441. do iyu=1,nstrs
  442. ddhomu(iyu,iou)=ddhook(iyu,iou)
  443. enddo
  444. enddo
  445.  
  446. CALL DBST(BGENE,DDHomu,XDDL,LRE,NSTRS,XSTRS)
  447. C xstrs contient la contrainte on va faire pica xstrs zdep05
  448. DO 220 INO=1,NBNN
  449. DO 220 ID=1,IDIM
  450. XE1(ID,INO)=XE(ID,INO)
  451. XE2(ID,INO)=XE(ID,INO)-xddl( id+ (ino-1)*nddd)/2.D0
  452. 220 CONTINUE
  453. DO iou=1,3
  454. DO IYU=1,3
  455. XJAC(iou,iyu)=0.D0
  456. enddo
  457. enddo
  458. CALL ZERO(XJAC,3,3)
  459. CALL HPRIME(XE1,NBNN,IDIM,SHPtot,IGAU,SHpwrk,DJAC)
  460. C
  461. C CALCUL DE LA MATRICE F
  462. C
  463. DO 140 ID=1,NBNN
  464. DO 140 IE=1,IDIM
  465. DO 140 IF=1,IDIM
  466. XJAC(IE,IF)=SHpwrk(IF+1,ID)*XE2(IE,ID)+XJAC(IE,IF)
  467. 140 CONTINUE
  468. IF(IDIM.EQ.2) THEN
  469. XJAC(3,3)=1.D0
  470. IF(IFOUR.EQ.0) THEN
  471. C
  472. CC CAS AXISYMETRIQUE
  473. C
  474. R1=0.D0
  475. R2=0.D0
  476. DO 150 ID=1,NBNN
  477. R1=R1+SHpwrk(1,ID)*XE1(1,ID)
  478. R2=R2+SHpwrk(1,ID)*XE2(1,ID)
  479. 150 CONTINUE
  480. XJAC(3,3)=R2/(R1+1.D-20)
  481. ENDIF
  482. ENDIF
  483. CC CALCUL DE DETERMINANT DE F
  484. C
  485. IF(IDIM.EQ.2) THEN
  486. DETF=XJAC(1,1)*XJAC(2,2)-XJAC(1,2)*XJAC(2,1)
  487. DETF = DETF * XJAC (3,3)
  488. ENDIF
  489. IF(IDIM.EQ.3) THEN
  490. DETF=XJAC(1,1)*(XJAC(2,2)*XJAC(3,3)-XJAC(3,2)*XJAC(2,3))
  491. DETF=DETF-XJAC(2,1)*(XJAC(1,2)*XJAC(3,3)-XJAC(3,2)*XJAC(1,3))
  492. DETF=DETF+XJAC(3,1)*(XJAC(1,2)*XJAC(2,3)-XJAC(1,3)*XJAC(2,2))
  493. ENDIF
  494. DETF=1.D0/(DETF+1.D-20)
  495. C
  496. C CALCUL DES CONTRAINTES DE CAUCHY
  497. C
  498. DO 160 ID=1,NSTRS
  499. IND=IN(ID)
  500. JND=JN(ID)
  501. xstrs2(ID)=0.D0
  502. DO 170 IE=1,IDIM
  503. DO 170 IF=1,IDIM
  504. ICO=ITAB(IE,IF)
  505. XsTRS2(ID)=XsTRS(ICO)*XJAC(IND,IE)*XJAC(JND,IF)*DETF
  506. 1 +xstrs2(ID)
  507. 170 CONTINUE
  508. 160 CONTINUE
  509.  
  510. C
  511. C PEGON : ON NE FAIT PAS LA TRANSFORMATION SUR LA 3-EME COMPOSANTE
  512. C
  513. IF(IDIM.EQ.2) THEN
  514. xstrs2(3)=xstrs(3)*XJAC(3,3)*XJAC(3,3)*DETF
  515. ENDIF
  516. C fin du caldul de capi dans dans xstrs2 la contrainte de kirchoff
  517. C on continu en calculant epsi sur config initiale
  518. DO 221 INO=1,NBNN
  519. DO 221 ID=1,IDIM
  520. XE(ID,INO)=XE2(ID,INO)+xddl( id+ (ino-1)*nddd)/2.D0
  521. 221 CONTINUE
  522. C inversion loi de hook
  523. CALL INVALM(DDHOMU,LHOOK,LHOOK,KERRE,0.D0)
  524. DO 6000 I=1,LHOOK
  525. xstrs(I)=0.D0
  526. DO 60001 J=1,LHOOK
  527. xstrs(I)=xstrs(I)+DDHOMU(I,J)*xstrs2(J)
  528. 60001 CONTINUE
  529. 6000 CONTINUE
  530. C
  531. C remplissage du segment contenant les deformations
  532. C
  533. MPTVAL=IVAEPS
  534. DO 74 ICOMP=1,NSTRS
  535. MELVAL=IVAL(ICOMP)
  536. IBMN=MIN(IB,VELCHE(/2))
  537. VELCHE(IGAU,IBMN)=XSTRS(ICOMP)
  538. 74 CONTINUE
  539. 54 continue
  540.  
  541. IF (ISDJC.NE.0.AND.ISDJC.NE.NBPGAU) THEN
  542. if (noer.eq.1) then
  543. kerr=195
  544. else
  545. INTERR(1)=IB
  546. CALL ERREUR(195)
  547. GOTO 994
  548. endif
  549. ENDIF
  550. 34 CONTINUE
  551. 994 CONTINUE
  552. SEGSUP WRK1
  553. segdes minte
  554. C fin du truesdell
  555. C debut du jaumann
  556. elseif(ideri.eq.4) then
  557. C==========================
  558. NBNO=NBNN
  559. C* NDDD=NDEP
  560. C* IF (IFOUR.EQ.-3) NDDD=NDEP-3
  561. NDDD=NDEP-NDPGE
  562. C
  563. C Donnees liees a l'element de reference
  564. C
  565. SEGACT MINTE
  566. C
  567. SEGINI,WRK1
  568. SEGINI,MTRACE
  569. segini wrk2
  570.  
  571. C boucle sur les elements
  572. C
  573. DO 394 IB=1,NBELEM
  574. C
  575. C on cherche les deplacements
  576. C
  577. MPTVAL=IVADEP
  578. IE=1
  579. DO 494 IGAU=1,NBNN
  580. DO 494 ICOMP=1,NDDD
  581. MELVAL=IVAL(ICOMP)
  582. IGMN=MIN(IGAU,VELCHE(/1))
  583. IBMN=MIN(IB ,VELCHE(/2))
  584. XDDL(IE)=VELCHE(IGMN,IBMN)
  585. IE=IE+1
  586. 494 CONTINUE
  587. IF (NDPGE.GT.0) THEN
  588. DO i=1,NDPGE
  589. XDDL(IE)=UDPGE(i)
  590. IE=IE+1
  591. ENDDO
  592. ENDIF
  593. C
  594. C on cherche les coordonnees des noeuds de l element ib
  595. C
  596. CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  597. C
  598. C on se met sur la config à mi pas (XE) xe1 est la config initiale
  599. C
  600. do iou=1,idim
  601. do iyu=1,nbnn
  602. XE1(iou,iyu)= xe(iou,iyu)
  603. XE(iou,iyu)= xe(iou,iyu) + xddl( iou+ (iyu-1)*nddd)/2.D0
  604. enddo
  605. enddo
  606. C
  607. C Calcul des coeff de modification de b-barre (elts incompres)
  608. C
  609. IF (MFR.EQ.31) THEN
  610. CALL BBCALC(XE,MELE,NBNN,IDIM,NBPGAU,
  611. 1 POIGAU,QSIGAU,ETAGAU,DZEGAU,NSTRS,
  612. 2 LRE,IFOUR,A,BBX,NHRM,SHPTOT,SHPWRK,XDPGE,YDPGE)
  613. ENDIF
  614. C
  615. C boucle sur les points de gauss
  616. C
  617. ISDJC=0
  618. C
  619. DO 594 IGAU=1,NBPTEL
  620. C
  621. CALL BMATST(IGAU,NBPGAU,POIGAU,QSIGAU,ETAGAU,DZEGAU,
  622. 1 MELE,MFR,NBNN,LRE,IFOUR,NSTRS,NHRM,1.D0,XE,
  623. 2 SHPTOT,SHPWRK,BGENE,DJAC,XDPGE,YDPGE)
  624.  
  625. IF (DJAC.EQ.0.D0) THEN
  626. INTERR(1)=IB
  627. CALL ERREUR(259)
  628. GOTO 9964
  629. ENDIF
  630. IF (DJAC.LT.0.D0) ISDJC=ISDJC+1
  631.  
  632. C En cas d'elements incompressibles : BGENE selon la methode B-BARRE
  633. IF (MFR.EQ.31) THEN
  634. CALL BBAR(IGAU,NBPGAU,POIGAU,QSIGAU,ETAGAU,DZEGAU,
  635. & MELE,NBNN,LRE,IFOUR,NSTRS,XE,DJAC,A,BBX,BGENE)
  636. ENDIF
  637. C
  638. CALL BST(BGENE,XDDL,LRE,NSTRS,XSTRS)
  639. C dans xstrs on a les deformations II sur config mi pas
  640. C on va calculer grad u/2 puis decomposition polaire puis rtens
  641. C on retravaille sur config initiale
  642. xxzero=XZero
  643. iipdpg=0
  644. CALL BGRMAS(iGau,NOELE,NBNO,LRE,IFOUR,NGRA,NIFOUR,XE1,
  645. . xxzero,SHPTOT,SHPWRK,BB,BGR,DJAC,IIPDPG)
  646. do iou=1,lre
  647. xddls2(iou)= 0.5D0 * xddl(iou)
  648. enddo
  649. CALL BGRDEP(BGR,NGRA,XDDLs2,LRE,GRADI)
  650. C on ajoute l'identité au gradient
  651. if(idim.eQ.2) then
  652. gradi(1)=gradi(1)+1.D0
  653. gradi(4)=gradi(4)+1.D0
  654. ELSE IF(IDIM.EQ.3) THEN
  655. gradi(1)=gradi(1)+1.D0
  656. gradi(5)=gradi(5)+1.D0
  657. gradi(9)=gradi(9)+1.D0
  658. ENDIF
  659.  
  660. CALL POLA2(gradi,R,U,IDIM)
  661. C fait le rtens Rt.A.R on utilise u pour mettre Rt
  662. C et on met le tenseur dans le tableau tens
  663. C attention vu le stockage R est en fait Rt
  664. if(idim.eq.2) then
  665. U(1)=r(1)
  666. u(2)=r(3)
  667. U(3)=R(2)
  668. u(4)=R(4)
  669. tens(1)=xstrs(1)
  670. tens(2)=xstrs(4)/2.d0
  671. tens(3)=xstrs(4)/2.D0
  672. tens(4)=xstrs(2)
  673.  
  674. elseif(idim.eq.3) then
  675. U(1)=r(1)
  676. u(2)=r(4)
  677. U(3)=R(7)
  678. u(4)=R(2)
  679. u(5)=r(5)
  680. u(6)=r(8)
  681. u(7)=r(3)
  682. u(8)=r(6)
  683. u(9)=r(9)
  684. tens(1)=xstrs(1)
  685. tens(2)=xstrs(4)/2.D0
  686. tens(3)=xstrs(5)/2.D0
  687. tens(4)=tens(2)
  688. tens(5)=xstrs(2)
  689. tens(6)=xstrs(6)/2.D0
  690. tens(7)=tens(3)
  691. tens(8)=tens(6)
  692. tens(9)=xstrs(3)
  693. else
  694. write(6,*)' idim est ni 2 ni 3 stop'
  695. stop
  696. endif
  697.  
  698. CALL MULMAT(tentra,tens,U,IDIM,IDIM,IDIM)
  699. CALL MULMAT(tens,R,Tentra,IDIM,IDIM,IDIM)
  700. C tens contient le nouveau tenseur on va remplir xstrs
  701. C en 2 D epzz ne change pas
  702. if(idim.eq.2) then
  703. xstrs(1)=tens(1)
  704. xstrs(2)=tens(4)
  705. xstrs(4)=tens(2)*2.D0
  706. else
  707. xstrs(1)=tens(1)
  708. xstrs(2)=tens(5)
  709. xstrs(3)=tens(9)
  710. xstrs(4)=tens(2)*2.D0
  711. xstrs(5)=tens(3)*2.D0
  712. xstrs(6)=tens(6)*2.D0
  713. endif
  714.  
  715. C
  716. C remplissage du segment contenant les deformations
  717. C
  718. MPTVAL=IVAEPS
  719. DO 794 ICOMP=1,NSTRS
  720. MELVAL=IVAL(ICOMP)
  721. IBMN=MIN(IB,VELCHE(/2))
  722. VELCHE(IGAU,IBMN)=XSTRS(ICOMP)
  723. 794 CONTINUE
  724. C
  725. 594 CONTINUE
  726. C
  727. C fin de la boucle sur les points de gauss
  728. C
  729. IF (ISDJC.NE.0.AND.ISDJC.NE.NBPGAU) THEN
  730. INTERR(1)=IB
  731. if (noer.eq.1) then
  732. kerr=195
  733. else
  734. CALL ERREUR(195)
  735. GOTO 9964
  736. endif
  737. ENDIF
  738.  
  739. 394 CONTINUE
  740.  
  741. C
  742. C fin de la boucle sur les elements
  743. C
  744. 9964 CONTINUE
  745. SEGDES MINTE
  746. SEGSUP WRK1,wrk2
  747. SEGSUP MTRACE
  748. endif
  749. C
  750. GOTO 510
  751.  
  752. C Elements massifs en FORMULATIONs 'ELECTROSTATIQUE' et 'DIFFUSION'
  753. C -----------------------------------------------------------------
  754. 97173 CONTINUE
  755. SEGACT,MINTE
  756. NBNO = NBNN
  757. NDDD = NDEP
  758. SEGINI,WRK1
  759. C-------------------------
  760. C Boucle sur les elements
  761. C-------------------------
  762. DO IEL = 1, NBELEM
  763. C - Recuperation des coordonnees des noeuds de l element IEL
  764. CALL DOXE(XCOOR,IDIM,NBNN,NUM,IEL,XE)
  765. C - Recuperation des inconnues primales aux noeuds de l element IEL
  766. MPTVAL = IVADEP
  767. IE = 1
  768. DO IGAU = 1, NBNN
  769. DO ICOMP = 1, NDDD
  770. MELVAL = IVAL(ICOMP)
  771. IGMN = MIN(IGAU,VELCHE(/1))
  772. IEMN = MIN(IEL ,VELCHE(/2))
  773. XDDL(IE) = VELCHE(IGMN,IEMN)
  774. IE = IE+1
  775. ENDDO
  776. ENDDO
  777. C-- -- -- -- -- -- -- -- --
  778. C - Boucle sur les points de Gauss
  779. C-- -- -- -- -- -- -- -- --
  780. ISDJC=0
  781. DO IGAU = 1, NBPTEL
  782. C -- Calcul de la matrice B et du jacobien au point de Gauss IGAU
  783. IF (MFR.EQ.71) THEN
  784. CALL BELEC(XE,SHPTOT(1,1,IGAU),NBNN,LHOOK,-1,
  785. & SHPWRK,BGENE,DJAC)
  786. ELSE IF (MFR.EQ.73) THEN
  787. CALL BDIFF(XE,SHPTOT(1,1,IGAU),NBNN,LHOOK,-1,
  788. & SHPWRK,BGENE,DJAC)
  789. ENDIF
  790. IF (DJAC.EQ.0.D0) THEN
  791. INTERR(1) = IEL
  792. CALL ERREUR(259)
  793. GOTO 98173
  794. ENDIF
  795. IF (DJAC.LT.0.D0) ISDJC = ISDJC+1
  796. CALL BST(BGENE,XDDL,LRE,NSTRS,XSTRS)
  797. C -- Remplissage du segment contenant les "deformations" = -grad(.)
  798. MPTVAL = IVAEPS
  799. DO ICOMP = 1, NSTRS
  800. MELVAL = IVAL(ICOMP)
  801. IEMN = MIN(IEL,VELCHE(/2))
  802. VELCHE(IGAU,IEMN) = XSTRS(ICOMP)
  803. ENDDO
  804. C-- -- -- -- -- -- -- -- --
  805. ENDDO
  806. C-- -- -- -- -- -- -- -- --
  807. IF (ISDJC.NE.0 .AND. ISDJC.NE.NBPGAU) THEN
  808. INTERR(1) = IEL
  809. if (noer.eq.1) then
  810. kerr=195
  811. else
  812. CALL ERREUR(195)
  813. GOTO 98173
  814. endif
  815. ENDIF
  816. C-------------------------
  817. ENDDO
  818. C-------------------------
  819. 98173 CONTINUE
  820. SEGDES,MINTE
  821. SEGSUP,WRK1
  822. GOTO 510
  823.  
  824. C_______________________________________________________________________
  825. C
  826. C milieux poreux
  827. C_______________________________________________________________________
  828. C
  829. 79 CONTINUE
  830. C
  831. C pour ces elements nbbb = nombre de noeuds
  832. C nbno = nombre de fonctions de forme
  833. C
  834. NBNO=IPORE
  835. NSTN=1
  836. LRN=NBNO-NBBB
  837. LRB=LRE-LRN
  838. C
  839. SEGINI WRK1,WRK5
  840. C Initialisation de MTRACE necessaire mais inutilise pour ces elements
  841. IF (IREPS2.EQ.1) SEGINI MTRACE
  842. C
  843. DO 3079 IB=1,NBELEM
  844. C
  845. C on cherche les coordonnees des noeuds de l element ib
  846. C
  847. CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  848. C
  849. C on cherche les deplacements
  850. C
  851. MPTVAL=IVADEP
  852. IE=1
  853. DO 4079 IGAU=1,NBNN
  854. DO 4079 ICOMP=1,NDEP-1
  855. MELVAL=IVAL(ICOMP)
  856. IGMN=MIN(IGAU,VELCHE(/1))
  857. IBMN=MIN(IB ,VELCHE(/2))
  858. XDDL(IE)=VELCHE(IGMN,IBMN)
  859. IE=IE+1
  860. 4079 CONTINUE
  861. C
  862. C puis les pressions
  863. C
  864. MELVAL=IVAL(NDEP)
  865. DO 4179 IGAU=1,LRN
  866. IGAUSO=IBSOM(NSPOS(IELE)+IGAU-1)
  867. IGMN=MIN(IGAUSO,VELCHE(/1))
  868. IBMN=MIN(IB ,VELCHE(/2))
  869. XDDL(IE)=VELCHE(IGMN,IBMN)
  870. IE=IE+1
  871. 4179 CONTINUE
  872. C
  873. C boucle sur les points de gauss
  874. C
  875. ISDJC=0
  876. C
  877. DO 5079 IGAU=1,NBPTEL
  878. C
  879. CALL BNPORE(IGAU,NBNO,NBBB,LRE,IFOUR,LHOOK,NSTN,NHRM,
  880. & 1.D0,XE,SHPTOT,SHPWRK,BGENE,XGENE,DJAC,1)
  881. C
  882. IF (DJAC.EQ.0.D0) THEN
  883. INTERR(1)=IB
  884. CALL ERREUR(259)
  885. GOTO 9979
  886. ENDIF
  887. IF (DJAC.LT.0.D0) ISDJC=ISDJC+1
  888. C
  889. CALL BST(BGENE,XDDL,LRE,LHOOK,XSTRS)
  890. C
  891. C calcul des eps 2
  892. C
  893. IF (IREPS2.EQ.1)
  894. & CALL BST2(SHPWRK,XDDL,XE,NBNN,IFOUR,XSTRS,TRACE,
  895. & IGAU,XDPGE,YDPGE,UDPGE,NHRM)
  896. C
  897. C calcul de la pression
  898. C
  899. XP=0.D0
  900. DO 4279 ID=1,LRN
  901. XP=XP+XGENE(1,ID)*XDDL(LRB+ID)
  902. 4279 CONTINUE
  903. XSTRS(NSTRS)=XP
  904. C
  905. C remplissage du segment contenant les deformations
  906. C
  907. MPTVAL=IVAEPS
  908. DO 7079 ICOMP=1,NSTRS
  909. MELVAL=IVAL(ICOMP)
  910. IGMN=MIN(IGAU,VELCHE(/1))
  911. IBMN=MIN(IB ,VELCHE(/2))
  912. VELCHE(IGMN,IBMN)=XSTRS(ICOMP)
  913. 7079 CONTINUE
  914. C
  915. 5079 CONTINUE
  916. C
  917. C fin de la boucle sur les points de gauss
  918. C
  919. IF (ISDJC.NE.0.AND.ISDJC.NE.NBPGAU) THEN
  920. INTERR(1)=IB
  921. if (noer.eq.1) then
  922. kerr=195
  923. else
  924. CALL ERREUR(195)
  925. GOTO 9979
  926. endif
  927. ENDIF
  928. C
  929. 3079 CONTINUE
  930. C
  931. C fin de la boucle sur les elements
  932. C
  933. 9979 CONTINUE
  934. SEGSUP WRK1,WRK5
  935. IF (IREPS2.EQ.1) SEGSUP MTRACE
  936. C
  937. GOTO 510
  938.  
  939. C_______________________________________________________________________
  940. C
  941. C milieux poreux - SUITE
  942. C_______________________________________________________________________
  943. C
  944. 173 CONTINUE
  945. C
  946. C pour ces elements nbbb = nombre de noeuds
  947. C nbno = nombre de fonctions de forme
  948. C
  949. IF (MELE.GE.173.AND.MELE.LE.177) THEN
  950. IDECAP = 2
  951. ELSE IF (MELE.GE.178.AND.MELE.LE.182) THEN
  952. IDECAP = 3
  953. ENDIF
  954. C
  955. NBNO=IPORE
  956. NSTN=IDECAP
  957. NSTB=4
  958. IF(IFOUR.EQ.1.OR.IFOUR.EQ.2) NSTB=6
  959. C
  960. LPP=NBNO-NBBB
  961. LRN=IDECAP*LPP
  962. LRB=LRE-LRN
  963. C
  964. SEGINI WRK1,WRK5
  965. C Initialise de MTRACE necessaire mais inutilise pour cet element
  966. IF (IREPS2.EQ.1) SEGINI MTRACE
  967. C
  968. DO 3173 IB=1,NBELEM
  969. C
  970. C on cherche les coordonnees des noeuds de l element ib
  971. C
  972. CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  973. C
  974. C on cherche les deplacements
  975. C
  976. MPTVAL=IVADEP
  977. IE=1
  978. DO 4173 IGAU=1,NBNN
  979. DO 4173 ICOMP=1,NDEP-IDECAP
  980. MELVAL=IVAL(ICOMP)
  981. IGMN=MIN(IGAU,VELCHE(/1))
  982. IBMN=MIN(IB ,VELCHE(/2))
  983. XDDL(IE)=VELCHE(IGMN,IBMN)
  984. IE=IE+1
  985. 4173 CONTINUE
  986. C
  987. C puis les pressions
  988. C
  989. DO 4473 IPR = 1,IDECAP
  990. MELVAL=IVAL(NDEP-IDECAP+IPR)
  991. DO 4273 IGAU=1,LPP
  992. IGAUSO=IBSOM(NSPOS(IELE)+IGAU-1)
  993. IGMN=MIN(IGAUSO,VELCHE(/1))
  994. IBMN=MIN(IB ,VELCHE(/2))
  995. XDDL(IE)=VELCHE(IGMN,IBMN)
  996. IE=IE+1
  997. 4273 CONTINUE
  998. 4473 CONTINUE
  999. C
  1000. C boucle sur les points de gauss
  1001. C
  1002. ISDJC=0
  1003. C
  1004. DO 5173 IGAU=1,NBPTEL
  1005. C
  1006. CALL BNQORE(IGAU,NBNO,NBBB,LRE,IFOUR,NSTB,NSTN,NHRM,
  1007. & 1.D0,XE,SHPTOT,SHPWRK,BGENE,XGENE,DJAC,IDECAP,LHOOK,1)
  1008. C
  1009. IF (DJAC.EQ.0.D0) THEN
  1010. INTERR(1)=IB
  1011. CALL ERREUR(259)
  1012. GOTO 9173
  1013. ENDIF
  1014. IF (DJAC.LT.0.D0) ISDJC=ISDJC+1
  1015. C
  1016. CALL BST(BGENE,XDDL,LRE,LHOOK,XSTRS)
  1017. C
  1018. C calcul des eps 2
  1019. C
  1020. IF (IREPS2.EQ.1)
  1021. & CALL BST2(SHPWRK,XDDL,XE,NBNN,IFOUR,XSTRS,TRACE,
  1022. & IGAU,XDPGE,YDPGE,UDPGE,NHRM)
  1023. C
  1024. C calcul des pressions
  1025. C
  1026. IE=LRB
  1027. DO 4673 IPR=1,IDECAP
  1028. XP=0.D0
  1029. IPR1=(IPR-1)*LPP
  1030. DO 4373 ID=1,LPP
  1031. IE=IE+1
  1032. XP=XP+XGENE(IPR,ID+IPR1)*XDDL(IE)
  1033. 4373 CONTINUE
  1034. XSTRS(NSTRS-IDECAP+IPR)=XP
  1035. 4673 CONTINUE
  1036. C
  1037. C remplissage du segment contenant les deformations
  1038. C
  1039. MPTVAL=IVAEPS
  1040. DO 7173 ICOMP=1,NSTRS
  1041. MELVAL=IVAL(ICOMP)
  1042. IGMN=MIN(IGAU,VELCHE(/1))
  1043. IBMN=MIN(IB ,VELCHE(/2))
  1044. VELCHE(IGMN,IBMN)=XSTRS(ICOMP)
  1045. 7173 CONTINUE
  1046. C
  1047. 5173 CONTINUE
  1048. C
  1049. C fin de la boucle sur les points de gauss
  1050. C
  1051. IF (ISDJC.NE.0.AND.ISDJC.NE.NBPGAU) THEN
  1052. INTERR(1)=IB
  1053. if (noer.eq.1) then
  1054. kerr=195
  1055. else
  1056. CALL ERREUR(195)
  1057. GOTO 9173
  1058. endif
  1059. ENDIF
  1060. C
  1061. 3173 CONTINUE
  1062. C
  1063. C fin de la boucle sur les elements
  1064. C
  1065. 9173 CONTINUE
  1066. SEGSUP WRK1,WRK5
  1067. IF (IREPS2.EQ.1) SEGSUP MTRACE
  1068. C
  1069. GOTO 510
  1070.  
  1071. C_______________________________________________________________________
  1072. C
  1073. C joints poreux
  1074. C_______________________________________________________________________
  1075. C
  1076. 80 CONTINUE
  1077. C
  1078. C pour ces elements nbbb = nombre de noeuds
  1079. C nbno = nombre de fonctions de forme
  1080. C
  1081. NBNO=IPORE
  1082. NSTN=1
  1083. LRN=(NBNO-NBBB)*3/2
  1084. LPP = LRN
  1085. LRB=LRE-LRN
  1086. NFAC=(3*NBBB-NBNO)/2
  1087. C
  1088. SEGINI WRK1,WRK3,WRK5
  1089. C
  1090. DO 3080 IB=1,NBELEM
  1091. C
  1092. C on cherche d'abord les deplacements
  1093. C
  1094. MPTVAL=IVADEP
  1095. IE=1
  1096. DO 4180 IGAU=1,NFAC
  1097. DO 4280 ICOMP=1,NDEP-1
  1098. MELVAL=IVAL(ICOMP)
  1099. IGMN=MIN(IGAU,VELCHE(/1))
  1100. IBMN=MIN(IB ,VELCHE(/2))
  1101. XDDL(IE)=VELCHE(IGMN,IBMN)
  1102. IE=IE+1
  1103. 4280 CONTINUE
  1104. 4180 CONTINUE
  1105. C
  1106. C puis les pressions
  1107. C
  1108. MELVAL=IVAL(NDEP)
  1109. DO 4080 IGAU=1,NBNN
  1110. DO 4190 INSOM=1,NBSOM(IELE)
  1111. IF (IGAU.EQ.IBSOM(NSPOS(IELE)+INSOM-1)) GO TO 4191
  1112. 4190 CONTINUE
  1113. IF (IGAU.GT.NFAC) GO TO 4191
  1114. GO TO 4080
  1115. 4191 CONTINUE
  1116. IBMN=MIN(IB ,VELCHE(/2))
  1117. IGMN=MIN(IGAU,VELCHE(/1))
  1118. XDDL(IE)=VELCHE(IGMN,IBMN)
  1119. IE=IE+1
  1120. 4080 CONTINUE
  1121. C
  1122. C on cherche les coordonnees des noeuds de l element ib
  1123. C
  1124. CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  1125. C
  1126. C calcul des exes locaux et des coordonnees locales
  1127. C
  1128. CALL JOPLOC(XE,SHPTOT,NBBB,NBNO,IFOUR,XEL,BPSS)
  1129. C
  1130. CALL INTDEL(XNTH,XNTB,XNTT,LRN,MELE)
  1131. C
  1132. C boucle sur les points de gauss
  1133. C
  1134. ISDJC=0
  1135. C
  1136. DO 5080 IGAU=1,NBPTEL
  1137. C
  1138. CALL BNPORJ(IGAU,NBNO,NBBB,LRE,IFOUR,LHOOK,NSTN,XE,XEL,
  1139. & SHPTOT,SHPWRK,BPSS,BGENE,XGENE,DJAC,1)
  1140. C
  1141. IF (DJAC.EQ.0.D0) THEN
  1142. INTERR(1)=IB
  1143. CALL ERREUR(259)
  1144. GOTO 9980
  1145. ENDIF
  1146. IF (DJAC.LT.0.D0) ISDJC=ISDJC+1
  1147. C
  1148. CALL BST(BGENE,XDDL,LRB,LHOOK,XSTRS)
  1149.  
  1150. C
  1151. C calcul de la pression
  1152. C
  1153. XP=0.D0
  1154. DO 4480 ID=1,LRN
  1155. XP=XP+XNTT(ID)*XGENE(1,ID)*XDDL(LRB+ID)
  1156. 4480 CONTINUE
  1157. XSTRS(NSTRS)=XP
  1158. C
  1159. C remplissage du segment contenant les deformations
  1160. C
  1161. MPTVAL=IVAEPS
  1162. DO 7080 ICOMP=1,NSTRS
  1163. MELVAL=IVAL(ICOMP)
  1164. IGMN=MIN(IGAU,VELCHE(/1))
  1165. IBMN=MIN(IB ,VELCHE(/2))
  1166. VELCHE(IGMN,IBMN)=XSTRS(ICOMP)
  1167. 7080 CONTINUE
  1168. C
  1169. 5080 CONTINUE
  1170. C
  1171. IF (ISDJC.NE.0.AND.ISDJC.NE.NBPGAU) THEN
  1172. INTERR(1)=IB
  1173. if (noer.eq.1) then
  1174. kerr=195
  1175. else
  1176. CALL ERREUR(195)
  1177. GOTO 9980
  1178. endif
  1179. ENDIF
  1180. C
  1181. 3080 CONTINUE
  1182. C
  1183. 9980 CONTINUE
  1184. SEGSUP WRK1,WRK3,WRK5
  1185. C
  1186. GOTO 510
  1187.  
  1188.  
  1189. C_______________________________________________________________________
  1190. C
  1191. C joints poreux - SUITE
  1192. C_______________________________________________________________________
  1193. C
  1194. 185 CONTINUE
  1195. C
  1196. C pour ces elements nbbb = nombre de noeuds
  1197. C nbno = nombre de fonctions de forme
  1198. C
  1199. IF (MELE.GE.185.AND.MELE.LE.187) THEN
  1200. IDECAP = 2
  1201. ELSE IF (MELE.GE.188.AND.MELE.LE.190) THEN
  1202. IDECAP = 3
  1203. ENDIF
  1204.  
  1205. C
  1206. NBNO=IPORE
  1207. NSTN=IDECAP
  1208. NSTB=2
  1209. IF(IFOUR.EQ.1.OR.IFOUR.EQ.2) NSTB=3
  1210.  
  1211. C
  1212. LPP=(NBNO-NBBB)*3/2
  1213. LRN=IDECAP*LPP
  1214. LRB=LRE-LRN
  1215.  
  1216. NFAC=(3*NBBB-NBNO)/2
  1217. C
  1218. SEGINI WRK1,WRK3,WRK5
  1219. C
  1220. DO 3185 IB=1,NBELEM
  1221. C
  1222. C on cherche d'abord les deplacements
  1223. C
  1224. MPTVAL=IVADEP
  1225. IE=1
  1226. DO 4185 IGAU=1,NFAC
  1227. DO 4285 ICOMP=1,NDEP-IDECAP
  1228. MELVAL=IVAL(ICOMP)
  1229. IGMN=MIN(IGAU,VELCHE(/1))
  1230. IBMN=MIN(IB ,VELCHE(/2))
  1231. XDDL(IE)=VELCHE(IGMN,IBMN)
  1232. IE=IE+1
  1233. 4285 CONTINUE
  1234. 4185 CONTINUE
  1235. C
  1236. C puis les pressions
  1237. C
  1238. DO 4785 IPR=1,IDECAP
  1239. MELVAL=IVAL(NDEP-IDECAP+IPR)
  1240. DO 4085 IGAU=1,NBNN
  1241. DO 4195 INSOM=1,NBSOM(IELE)
  1242. IF (IGAU.EQ.IBSOM(NSPOS(IELE)+INSOM-1)) GO TO 4891
  1243. 4195 CONTINUE
  1244. IF (IGAU.GT.NFAC) GO TO 4891
  1245. GO TO 4085
  1246. 4891 CONTINUE
  1247. IBMN=MIN(IB ,VELCHE(/2))
  1248. IGMN=MIN(IGAU,VELCHE(/1))
  1249. XDDL(IE)=VELCHE(IGMN,IBMN)
  1250. IE=IE+1
  1251. 4085 CONTINUE
  1252. 4785 CONTINUE
  1253. C
  1254. C on cherche les coordonnees des noeuds de l element ib
  1255. C
  1256. CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  1257. C
  1258. C calcul des exes locaux et des coordonnees locales
  1259. C
  1260. CALL JOPLOC(XE,SHPTOT,NBBB,NBNO,IFOUR,XEL,BPSS)
  1261. C
  1262. CALL INTDEL(XNTH,XNTB,XNTT,LPP,MELE)
  1263. C
  1264. C boucle sur les points de gauss
  1265. C
  1266. ISDJC=0
  1267. C
  1268. DO 5185 IGAU=1,NBPTEL
  1269. C
  1270. CALL BNPQRJ(IGAU,NBNO,NBBB,LRE,IFOUR,LHOOK,NSTN,XE,XEL,
  1271. & SHPTOT,SHPWRK,BPSS,BGENE,XGENE,DJAC,IDECAP,NSTB,1)
  1272. C
  1273. IF (DJAC.EQ.0.D0) THEN
  1274. INTERR(1)=IB
  1275. CALL ERREUR(259)
  1276. GOTO 9985
  1277. ENDIF
  1278. IF (DJAC.LT.0.D0) ISDJC=ISDJC+1
  1279. C
  1280. CALL BST(BGENE,XDDL,LRB,LHOOK,XSTRS)
  1281. C
  1282. C calcul de la pression
  1283. C
  1284. IE=LRB
  1285. DO 4985 IPR=1,IDECAP
  1286. XP=0.D0
  1287. IPR1=(IPR-1)*LPP
  1288. DO 4485 ID=1,LPP
  1289. IE=IE+1
  1290. XP=XP+XNTT(ID)*XGENE(IPR,ID+IPR1)*XDDL(IE)
  1291. 4485 CONTINUE
  1292. XSTRS(NSTRS-IDECAP+IPR)=XP
  1293. 4985 CONTINUE
  1294.  
  1295. C
  1296. C remplissage du segment contenant les deformations
  1297. C
  1298. MPTVAL=IVAEPS
  1299. DO 7185 ICOMP=1,NSTRS
  1300. MELVAL=IVAL(ICOMP)
  1301. IGMN=MIN(IGAU,VELCHE(/1))
  1302. IBMN=MIN(IB ,VELCHE(/2))
  1303. VELCHE(IGMN,IBMN)=XSTRS(ICOMP)
  1304. 7185 CONTINUE
  1305. C
  1306. 5185 CONTINUE
  1307. C
  1308. IF (ISDJC.NE.0.AND.ISDJC.NE.NBPGAU) THEN
  1309. INTERR(1)=IB
  1310. if (noer.eq.1) then
  1311. kerr=195
  1312. else
  1313. CALL ERREUR(195)
  1314. GOTO 9985
  1315. endif
  1316. ENDIF
  1317. C
  1318. 3185 CONTINUE
  1319. C
  1320. 9985 CONTINUE
  1321. SEGSUP WRK1,WRK3,WRK5
  1322. C
  1323. GOTO 510
  1324.  
  1325. C____________________________________________________________________
  1326. 99 CONTINUE
  1327. MOTERR(1:4)=NOMTP(MELE)
  1328. MOTERR(9:12)='EPSI'
  1329. CALL ERREUR(86)
  1330.  
  1331. 510 CONTINUE
  1332. RETURN
  1333. END
  1334.  
  1335.  
  1336.  
  1337.  
  1338.  
  1339.  
  1340.  
  1341.  
  1342.  
  1343.  
  1344.  
  1345.  

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