Télécharger epsi2.eso

Retour à la liste

Numérotation des lignes :

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

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