Télécharger epsix.eso

Retour à la liste

Numérotation des lignes :

  1. C EPSIX SOURCE FANDEUR 16/01/07 21:15:20 8756
  2. subroutine EPSIX (IMODEL,IREPS2,IVADEP,IVAEPS,
  3. & UZDPG,RYDPG,RXDPG,IIPDPG,IRETER)
  4. c
  5. C PROCEDURE UTILISEE DANS LE CAS D'ELEMENTS XFEM
  6. c POUR LE CALCUL DE la deformation
  7. C
  8. C*********************************************************
  9. C PARTIE DECLARATIVE
  10. C*********************************************************
  11. C
  12. IMPLICIT REAL*8(A-H,O-Z)
  13. C
  14. PARAMETER (NDDLMAX=30,NBNIMAX=10)
  15. C
  16. PARAMETER (NBENRMAX=5)
  17. DIMENSION MLRE(NBENRMAX+1)
  18. C
  19. -INC CCOPTIO
  20. -INC SMCOORD
  21. -INC SMMODEL
  22. -INC SMCHAML
  23. -INC SMINTE
  24. -INC SMELEME
  25. -INC SMLREEL
  26. C
  27. POINTEUR MCHEX1.MCHELM,MINT1.MINTE,MINT2.MINTE
  28. C
  29. C Segment (type LISTENTI) contenant les informations sur un element
  30. SEGMENT INFO
  31. INTEGER INFELL(JG)
  32. ENDSEGMENT
  33. c
  34. SEGMENT WRK1
  35. REAL*8 XE(3,NBBB)
  36. REAL*8 XDDL(LRE),XSTRS(LHOOK)
  37. ENDSEGMENT
  38. c
  39. SEGMENT WRK2
  40. REAL*8 SHPWRK(6,NBNO),BGENE(LHOOK,LRE)
  41. c REAL*8 LV7WRK(NBENRMA2,2,6,NBNO)
  42. REAL*8 LV7WRK(NBENRMA2,2,6,NBBB)
  43. ENDSEGMENT
  44. C
  45. SEGMENT MPTVAL
  46. INTEGER IPOS(NS),NSOF(NS)
  47. INTEGER IVAL(NCOSOU)
  48. CHARACTER*16 TYVAL(NCOSOU)
  49. ENDSEGMENT
  50. C
  51. SEGMENT MRACC
  52. INTEGER TLREEL(NBENRMA2,NBI)
  53. ENDSEGMENT
  54. C
  55. SEGMENT MTRACE
  56. REAL*8 TRACE(NBPTEL)
  57. ENDSEGMENT
  58. DIMENSION UDPGE(3)
  59. C
  60. C write (*,*) '############################'
  61. C write (*,*) '##### DEBUT DE EPSIX #####'
  62. C write (*,*) '############################'
  63. C
  64. C*********************************************************
  65. c Introduction du point autour duquel se fait le mouvement
  66. c de la section en defo plane generalisee
  67. C*********************************************************
  68. C IIPDPG > 0 (argument) si noeud/point support defini dans le modele
  69. C NDPGE > 0 si besoin de calcul avec point support
  70. IF (IIPDPG.GT.0) THEN
  71. IF (IFOUR.EQ.-3) THEN
  72. NDPGE=3
  73. UDPGE(1)=UZDPG
  74. UDPGE(2)=RYDPG
  75. UDPGE(3)=RXDPG
  76. C SEGACT,MCOORD
  77. IREF=(IIPDPG-1)*(IDIM+1)
  78. XDPGE=XCOOR(IREF+1)
  79. YDPGE=XCOOR(IREF+2)
  80. C* ELSE IF (IFOUR.EQ.11) THEN
  81. C* NDPGE=2
  82. C* UDPGE(1)=UZDPG
  83. C* UDPGE(2)=RXDPG
  84. C* UDPGE(3)=0.D0
  85. C* XDPGE=0.D0
  86. C* YDPGE=0.D0
  87. C* ELSE IF (IFOUR.EQ. 7 .OR. IFOUR.EQ. 8 .OR. IFOUR.EQ. 9 .OR.
  88. C* & IFOUR.EQ.10 .OR. IFOUR.EQ.14) THEN
  89. C* NDPGE=1
  90. C* UDPGE(1)=UZDPG
  91. C* UDPGE(2)=0.D0
  92. C* UDPGE(3)=0.D0
  93. C* XDPGE=0.D0
  94. C* YDPGE=0.D0
  95. else
  96. write(ioimp,*) 'EPSI4 : erreur NDPGE'
  97. call erreur(5)
  98. return
  99. ENDIF
  100. ELSE
  101. NDPGE=0
  102. UDPGE(1)=0.D0
  103. UDPGE(2)=0.D0
  104. UDPGE(3)=0.D0
  105. XDPGE=0.D0
  106. YDPGE=0.D0
  107. ENDIF
  108. C
  109. C
  110. C*********************************************************
  111. c
  112. C RECUPERATION + ACTIVATIONS + VALEURS UTILES
  113. c
  114. C*********************************************************
  115. C SEGACT MMODEL,IMODEL deja activé dans epsi1
  116. c
  117. C++++ Activation au cas ou ++++++++++++++++++++++++++++++
  118. SEGACT MCOORD
  119.  
  120. C++++ Recup + Activation de la geometrie ++++++++++++++++
  121. MELEME= IMAMOD
  122. c SEGACT MELEME deja activé dans epsi1
  123.  
  124.  
  125. C++++ RECUP DES INFOS EF ++++++++++++++++++++++++++++++++
  126. c + OBTENUES PAR ELQUOI DANS RIGI1 PENDANT PHASE 1
  127. C segment INFO deja actif dans RIGI1
  128. c + rigi1 n appelle pas elquoi, c'est modeli qui l'a fait
  129. c mais du coup, on na pas de segment minte
  130. c (car depend de si pt de g pour rigi, pour sigma....)
  131. c c'est + simple de rappeler elquoi ici
  132. MELE = NEFMOD
  133. call elquoi(MELE,0,3,IPINF,IMODEL)
  134. INFO = IPINF
  135. c MELE = INFELL(1)
  136. c NBPGA2= INFELL(2)
  137. c NBPGAU= INFELL(3)
  138. c NBPGAU= INFELL(4)
  139. c ICARA = INFELL(5)
  140. NGAU1 = INFELL(6)
  141. c LW = INFELL(7)
  142. LRE = INFELL(9)
  143. LHOOK = INFELL(10)
  144. MINT1 = INFELL(11)
  145. segact,MINT1
  146. MINT2 = INFELL(12)
  147. if(MINT2.ne.0) segact,MINT2
  148. MFR = INFELL(13)
  149. IELE = INFELL(14)
  150. NDDL = INFELL(15)
  151. NSTRS = INFELL(16)
  152. c write(6,*)'-> EPSIX infell',(infell(iou),iou=1,16)
  153.  
  154. c + AUTRES INFOS
  155. C nbre de noeuds par element
  156. NBNN1 = NUM(/1)
  157. C nbre d elements
  158. NBEL1 = NUM(/2)
  159.  
  160. c REM: pour se passer du dimensionnement du nbre d'enrichissement dans
  161. c elquoi et le realiser localement , on pourrait ecrire:
  162. c LRE = NDDLMAX*NBNN1
  163. c NDDL= NDDLMAX
  164.  
  165. C sous decoupage et points de Gauss de l element geometrique de base
  166. if(MELE.eq.263.or.mele.eq.264) then
  167. NGAU2 = MINT2.POIGAU(/1)
  168. endif
  169. C write(*,*) 'dim de MINT2=6,',(MINT2.SHPTOT(/2)),(MINT2.SHPTOT(/3))
  170. C write(*,*) 'MINT2',(MINT2.QSIGAU(iou),iou=1,NGAU)
  171. c segdes,MINT2
  172.  
  173.  
  174. c nbre maxi de fonction de forme par noeud (fonction std comprise)
  175. c NBNI = NDDL/IDIM inutile!
  176.  
  177.  
  178. C++++ Recup des infos d enrichissement +++++++++++++++++++
  179. c recup du MCHEX1 d'enrichissement
  180. NOBMO1 = IVAMOD(/1)
  181. if(NOBMO1.ne.0) then
  182. do iobmo1=1,NOBMO1
  183. if((TYMODE(iobmo1)).eq.'MCHAML') then
  184. MCHEX1 = IVAMOD(iobmo1)
  185. segact,MCHEX1
  186. if((MCHEX1.TITCHE).eq.'ENRICHIS') then
  187. MCHAM1 = MCHEX1.ICHAML(1)
  188. segact,MCHAM1
  189. segdes,MCHEX1
  190. goto 1000
  191. endif
  192. segdes,MCHEX1
  193. endif
  194. enddo
  195. write(*,*) 'Le modele est vide (absence d enrichissement)'
  196. * return
  197. else
  198. write(*,*) 'Il n y a pas de MCHEML enrichissement dans le Modele'
  199. * return
  200. endif
  201.  
  202. 1000 continue
  203. c niveau d enrichissement(s) du modele (ddls std u exclus)
  204. c NBENR1= 0 si std, 1 si H seul, 2 si H+F1, 3 si H+F1+F2, etc...
  205. if(NOBMO1.ne.0) then
  206. NBENR1= MCHAM1.IELVAL(/1)
  207. else
  208. NBENR1 = 0
  209. endif
  210. c write(*,*) 'niveau d enrichissement(s) du modele',NBENR1
  211. c
  212. C
  213. C
  214. C*********************************************************
  215. C INITIALISATIONS...
  216. C*********************************************************
  217. IRETER = 0
  218. c
  219. c preparation des tables avec:
  220.  
  221. if(NBENR1.ne.0) then
  222. do ienr=1,NBENR1
  223. c -le nombre d'inconnues de chaque sous-zone
  224. c determinee depuis le nombre de fonction de forme
  225. c ienr= 1: U+H(1+1=2), 2: U+H+F1(2+4=6), 3: U+H+F1+F2(6+4=10)
  226. nbniJ = 2 + ((ienr-1)*4)
  227. MLRE(1+ienr) = IDIM*NBNN1*nbniJ
  228. enddo
  229. endif
  230. C Tables + longues car 1er indice correspond au fontion de forme std
  231. MLRE(1) = IDIM*NBNN1*1
  232.  
  233. if(NBENR1.lt.(NBENRMAX+1)) then
  234. do ienr=(NBENR1+1),(NBENRMAX)
  235. MLRE(1+ienr) = 0
  236. enddo
  237. endif
  238. c
  239. c ...DU SEGMENT WRK1
  240. NBENRMA2 = NBENRMAX
  241. NBBB = NBNN1
  242. SEGINI,WRK1
  243.  
  244. c ...DU SEGMENT WRK2
  245. c NBNO = NBNI
  246. NBNO = LRE/IDIM
  247. SEGINI,WRK2
  248. C
  249. c ...DU SEGMENT MRACC
  250. NBENRMA2 = NBENRMAX
  251. NBI = NBNN1
  252. segini,MRACC
  253. C
  254. C du nombre d erreur sur le noms de composantes
  255. NBERR1=0
  256.  
  257.  
  258.  
  259. C*********************************************************
  260. C
  261. C>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> BOUCLE SUR LES ELEMENTS
  262. C
  263. DO 2000 J=1,NBEL1
  264. c write(*,*) '==================================='
  265. c write(*,*) '========= EF',J,' /NBEL1 ========='
  266. C
  267. C
  268. C*********************************************************
  269. C POUR CHAQUE ELEMENT, ...
  270. C*********************************************************
  271. C
  272. C++++ ON RECUPERE LES COORDONNEES DES NOEUDS DE L ELEMENT
  273. CALL DOXE(XCOOR,IDIM,NBNN1,NUM,J,XE)
  274. C
  275. c
  276. C++++ NBENRJ = niveau d enrichissement de l element ++++
  277. C =0 si EF std =1 si U+H =2 si U+H+F1 =3 si U+H+F1+F2
  278. NBENRJ=0
  279. if(NBENR1.ne.0) then
  280. do IENR=1,NBENR1
  281. MELVA1 = MCHAM1.IELVAL(IENR)
  282. segact,MELVA1
  283. do I=1,NBNN1
  284. mlree1 = MELVA1.IELCHE(I,J)
  285. C on en profite pour remplir MRACC table de raccourcis pour cet element
  286. TLREEL(IENR,I) = mlree1
  287. if(mlree1.ne.0) then
  288. NBENRJ = max(NBENRJ,IENR)
  289. C et on active la listreel
  290. segact,mlree1
  291. endif
  292. enddo
  293. segdes,MELVA1
  294. enddo
  295. endif
  296. if(NBENRMA2.gt.NBENR1) then
  297. do IENR2=(NBENR1+1),NBENRMA2
  298. do I=1,NBNN1
  299. TLREEL(IENR2,I) = 0
  300. enddo
  301. enddo
  302. endif
  303. C
  304. c
  305. c++++ quelques indices pour dimensionner au plus juste
  306. c nbre total de ddl de l'élément considéré
  307. NLIGRD = MLRE(1+NBENRJ)
  308. NLIGRP = MLRE(1+NBENRJ)
  309. c nbre fonction de forme=((Ni_std+Ni_enrichi)*nbnoeud)=Ninconnues/idim
  310. NBNI = (MLRE(1+NBENRJ)) / IDIM
  311.  
  312. C write(*,*) 'EF',J,' NBENRJ=',NBENRJ,
  313. C &'donc',NLIGRD,' ddls et ',NBNI,' fonctions de forme'
  314. c
  315. C
  316. C++++ ON CALCULE XDDL ++++
  317. MPTVAL = IVADEP
  318. NCOSOU = IVAL(/1)
  319. C
  320. INITOT = 0
  321. C-->> BOUCLE SUR LES niveaux d'ENRICHISSEMENTS (0:U 1:A 2:BCDE1 3:BCDE2)
  322. DO IENR=0,NBENRJ
  323. *nbre de fonction(s) de ce niveau d'enrichissement (=1 si std ou H, =4 pour F1,2,...)
  324. if(IENR .le. 1) then
  325. NBNIENR = 1
  326. else
  327. NBNIENR = 4
  328. endif
  329. C---->> BOUCLE SUR LES fonctions de forme de ce type d'enrichissement
  330. DO INI=1,NBNIENR
  331. INITOT = INITOT + 1
  332. C------>> BOUCLE SUR LA DIMENSION
  333. DO 2220 KDIM=1,IDIM
  334. ICOMP = (INITOT-1)*IDIM + KDIM
  335.  
  336. c --cas ou on n'a pas trouvé assez de composantes--
  337. if(ICOMP.GT.NCOSOU) GOTO 2221
  338.  
  339. MELVAL = IVAL(ICOMP)
  340. c --cas ou on a pas trouvé cette composante dans cette zone du
  341. c chpoint solution devenu mchaml --
  342. if(MELVAL.eq.0)then
  343. c Avait on besoin de cette composante?
  344. c oui, si c'est une composante obligatoire
  345. if(IENR.eq.0) goto 2991
  346. c oui, si l'un des noeuds est enrichi
  347. do I=1,NBNN1
  348. if(TLREEL(IENR,I).ne.0) goto 2991
  349. enddo
  350. c non, si c'est facultatif et qu'on n'est pas enrichi -> on saute
  351. goto 2220
  352. c ->AVERTISSEMENT puis on saute
  353. 2991 NBERR1=NBERR1+1
  354. if(IIMPI.lt.1) goto 2220
  355. c write(IOIMP,*) 'PB OPERATEUR EPSI :'
  356. write(IOIMP,991) ICOMP,IENR,INI,KDIM
  357. 991 format(2X,'ABSENCE DANS LE CHPOINT DEPLACEMENT DE LA ',I3,
  358. $ 'ieme COMPOSANTE (enrichissement',I3,
  359. $ ', fonction',I3,', direction ',I3,
  360. $ ') NECESSAIRE POUR L UN DES NOEUDS SUIVANTS :')
  361. write(IOIMP,*)' noeuds :',(NUM(iou,J),iou=1,NBNN1)
  362. goto 2220
  363. endif
  364.  
  365. C---------->> BOUCLE SUR LES NOEUDS
  366. DO I=1,NBNN1
  367. IQ = ((INITOT-1)*NBNN1*IDIM) + ((I-1)*IDIM) + KDIM
  368. XDDL(IQ) = VELCHE(I,J)
  369. ENDDO
  370.  
  371. 2220 CONTINUE
  372. ENDDO
  373. ENDDO
  374.  
  375. c --cas normal (toutes les composantes souhaitees etaient presentes)--
  376. GOTO 2223
  377.  
  378. c --cas ou on n'a pas trouvé assez de composantes--
  379. 2221 CONTINUE
  380. if (IIMPI.ge.1) then
  381. WRITE(IOIMP,2222) J,NCOSOU,ICOMP
  382. 2222 FORMAT(2X,'ATTENTION : ELEMENT ',I6,
  383. & ' LE CHAMP DE DEPLACEMENT CONTIENT ',I3,' COMPOSANTES',
  384. & ' PAR NOEUD AU LIEU DE ',I3,' ATTENDUES')
  385. endif
  386. NDDL=NCOSOU*NBNN1
  387. NBENRJ=IENR
  388.  
  389. 2223 CONTINUE
  390. c --cas ou on a une ou des erreurs--
  391. IF (NBERR1.gt.0.and.J.eq.NBEL1) THEN
  392. write(IOIMP,*) 'OPERATEUR GRAD : ABSENCE DANS LE CHPOINT ',
  393. & 'DEPLACEMENT DE CERTAINES INCONNUES ATTENDUES PAR LE MODELE'
  394. ENDIF
  395. c
  396. c rem: il serait probablement interessant au niveau du tems cpu
  397. c d'utiliser moins de pts de Gauss lorsque l element est élastique
  398. c On pourrait par ex. utiliser MINT2 = infell(12) qui contient
  399. c le segment d'integration de l'EF std (QUA4 par ex.)
  400. * if((NBENRJ.eq.0).and.(MINT2.ne.0)) then
  401. * MINTE = MINT2
  402. * NBPGAU= NGAU2
  403. * else
  404. MINTE = MINT1
  405. NBPGAU= NGAU1
  406. * endif
  407. c
  408. c pour les def quadratiques
  409. ISDJC=0
  410. NBPTEL=NBPGAU
  411. IF (IREPS2.EQ.1) SEGINI MTRACE
  412. c
  413. C
  414. C*********************************************************
  415. C
  416. C>>>>>>>>>> BOUCLE SUR LES POINTS DE GAUSS
  417. C
  418. DO 2500 KGAU=1,NBPGAU
  419. C
  420. C*********************************************************
  421. C Initialisation à 0
  422. C*********************************************************
  423.  
  424. c ZERO ne serait-il pas facultatif?
  425. CALL ZERO(SHPWRK,6,NBNO)
  426. XGAU = 0.D0
  427. YGAU = 0.D0
  428. ZGAU = 0.D0
  429. C
  430. i6zz = 3
  431. IF (IDIM.EQ.3) i6zz = 4
  432. C
  433. c do ienr7=1,NBENRMAX
  434. do ienr7=1,NBENRJ
  435. do inod7=1,NBNN1
  436. c do i6=1,6
  437. do i6=1,i6zz
  438. LV7WRK(ienr7,1,i6,inod7) = 0.D0
  439. LV7WRK(ienr7,2,i6,inod7) = 0.D0
  440. enddo
  441. enddo
  442. enddo
  443.  
  444.  
  445. c*********************************************************
  446. c Calcul des fonction de forme std dans repere local
  447. c*********************************************************
  448.  
  449. ccccc BOUCLE SUR LES NOEUDS ccccccccccccccccccccccccccccc
  450. c (et donc sur les Ni std)
  451. DO 2510 I=1,NBNN1
  452.  
  453. C++++ Calcul des Ni std
  454. c (rappel: 1:Ni 2:Ni,qsi 3:Ni,eta avec i=1,4)
  455. SHPWRK(1,I) = SHPTOT(1,I,KGAU)
  456. SHPWRK(2,I) = SHPTOT(2,I,KGAU)
  457. SHPWRK(3,I) = SHPTOT(3,I,KGAU)
  458. IF (IDIM.EQ.3) SHPWRK(4,I) = SHPTOT(4,I,KGAU)
  459.  
  460. C++++ Calcul des coordonnees dans le repere global
  461. XGAU = XGAU + (SHPWRK(1,I)*XE(1,I))
  462. YGAU = YGAU + (SHPWRK(1,I)*XE(2,I))
  463. IF (IDIM.EQ.3) ZGAU = ZGAU + (SHPWRK(1,I)*XE(3,I))
  464.  
  465. 2510 CONTINUE
  466. ccccc fin de BOUCLE SUR LES NOEUDS ccccccccccccccccccccccc
  467.  
  468.  
  469.  
  470. c*********************************************************
  471. c Passage des fonctions de forme std dans repere global
  472. c*********************************************************
  473.  
  474. C++++ CALCUL DES Ni,x Ni,y (i=1,4) + CALCUL DE det(J)
  475. CALL JACOBI(XE,SHPWRK,IDIM,NBNN1,DJAC)
  476. c if(J.eq.1.and.KGAU.eq.1)
  477. c &write(*,*) 'Ni(i=1,4)=',(SHPWRK(1,iou),iou=1,NBNN1)
  478.  
  479. c*********************************************************
  480. c Si on est pas du tout enrichi on peut sauter une grosse partie
  481. c*********************************************************
  482. if(NBENRJ.eq.0) goto 2999
  483.  
  484. c*********************************************************
  485. c Calcul des level set + leurs derivees dans repere global
  486. c*********************************************************
  487.  
  488. ccccc BOUCLE SUR LES enrichissements ccccccccccccccccccc
  489. do 2520 ienr=1,NBENRJ
  490.  
  491. c MELVA1=MCHAM1.IELVAL(IENR)
  492. c segact,MELVA1
  493.  
  494. ccccc BOUCLE SUR LES NOEUDS ccccccccccccccccccccccccccc
  495. DO 2521 I=1,NBNN1
  496.  
  497. C++++ Le I eme noeud est-il ienr-enrichi?
  498. mlree1= TLREEL(IENR,I)
  499. if(mlree1.eq.0) goto 2521
  500.  
  501.  
  502. c Calcul du repere local de fissure(=PSI,PHI)
  503. c (rappel: 1,1:psi 1,2:phi 2,1 psi,x ...etc...)
  504. do 2522 inode=1,NBNN1
  505. c pour le H-enrichissement, on n a pas gardé PSI (inutile)
  506. if(ienr.ne.1) then
  507. c valeur de PSI au inode^ieme noeud
  508. xpsi1 = mlree1.PROG(inode)
  509. c qu on multiplie par la valeur de Ni^std au pt de G considéré
  510. LV7WRK(ienr,1,1,I)= LV7WRK(ienr,1,1,I)
  511. & + (SHPWRK(1,inode)*xpsi1)
  512. LV7WRK(ienr,1,2,I)= LV7WRK(ienr,1,2,I)
  513. & + (SHPWRK(2,inode)*xpsi1)
  514. LV7WRK(ienr,1,3,I)= LV7WRK(ienr,1,3,I)
  515. & + (SHPWRK(3,inode)*xpsi1)
  516. IF (IDIM.EQ.3) LV7WRK(ienr,1,4,I)= LV7WRK(ienr,1,4,I)
  517. & + (SHPWRK(4,inode)*xpsi1)
  518. c valeur de PHI au inode^ieme noeud
  519. xphi1 = mlree1.PROG(NBNN1+inode)
  520. else
  521. xphi1 = mlree1.PROG(inode)
  522. endif
  523. LV7WRK(ienr,2,1,I)= LV7WRK(ienr,2,1,I)
  524. & + (SHPWRK(1,inode)*xphi1)
  525. LV7WRK(ienr,2,2,I)= LV7WRK(ienr,2,2,I)
  526. & + (SHPWRK(2,inode)*xphi1)
  527. LV7WRK(ienr,2,3,I)= LV7WRK(ienr,2,3,I)
  528. & + (SHPWRK(3,inode)*xphi1)
  529. IF (IDIM.EQ.3) LV7WRK(ienr,2,4,I)= LV7WRK(ienr,2,4,I)
  530. & + (SHPWRK(4,inode)*xphi1)
  531. 2522 continue
  532.  
  533. 2521 continue
  534. ccccc fin de BOUCLE SUR LES NOEUDS ccccccccccccccccccccccc
  535.  
  536.  
  537. 2520 CONTINUE
  538. ccccc fin de BOUCLE SUR LES enrichissements cccccccccccccccc
  539.  
  540. c on a construit
  541. C LV7WRK(ienr, PSI/PHI, valeur/deriveeparqsi/pareta, iNOEUD)
  542.  
  543.  
  544. c*********************************************************
  545. c Ajout des fonctions de forme d enrichissement
  546. c + leurs derivees dans repere global
  547. c*********************************************************
  548. CALL SHAPX(XGAU,YGAU,ZGAU,MELE,LV7WRK,NBENRMAX,NBENRJ,
  549. &TLREEL,J,SHPWRK,IRET)
  550.  
  551.  
  552. c retour a la partie commune aux EF enrichi et non enrichi
  553. 2999 continue
  554.  
  555. C*********************************************************
  556. C CALCUL DE B'
  557. C*********************************************************
  558. c ZERO ne serait-il pas facultatif?
  559. c call ZERO(BGENE,LHOOK,NLIGRP)
  560. KB=1
  561. C boucle sur tous les Ni
  562. DO 3001 II=1,NBNI
  563.  
  564. BGENE(1,KB) = SHPWRK(2,II)
  565. BGENE(2,KB+1) = SHPWRK(3,II)
  566. BGENE(4,KB) = SHPWRK(3,II)
  567. BGENE(4,KB+1) = SHPWRK(2,II)
  568.  
  569. IF(IDIM.EQ.3) THEN
  570. BGENE(3,KB+2)=SHPWRK(4,II)
  571. BGENE(5,KB)=SHPWRK(4,II)
  572. BGENE(5,KB+2)=SHPWRK(2,II)
  573. BGENE(6,KB+1)=SHPWRK(4,II)
  574. BGENE(6,KB+2)=SHPWRK(3,II)
  575. ENDIF
  576.  
  577. KB = KB + IDIM
  578.  
  579. 3001 CONTINUE
  580. C
  581. c if(J.eq.5.and.KGAU.eq.1) then
  582. c if(KGAU.eq.1) then
  583. c write(*,*) 'BGENE(1,..)=',(BGENE(1,iou),iou=1,2*NBNI)
  584. c write(*,*) 'BGENE(2,..)=',(BGENE(2,iou),iou=1,2*NBNI)
  585. c write(*,*) 'BGENE(6,..)=',(BGENE(6,iou),iou=1,2*NBNI)
  586. c endif
  587. c endif
  588.  
  589.  
  590. c
  591.  
  592. C*********************************************************
  593. C CALCUL DE B.q'
  594. C*********************************************************
  595.  
  596. C APPEL A LA PROCEDURE DE CALCUL
  597. CALL BST(BGENE,XDDL,(NBNI*IDIM),LHOOK,XSTRS)
  598. c
  599. c cas de la priâe en compte des termes quadratiques
  600. IF (IREPS2.EQ.1)
  601. & CALL BST2(SHPWRK,XDDL,XE,NBNI,IFOUR,XSTRS,TRACE,KGAU,
  602. & XDPGE,YDPGE,UDPGE,NIFOUR)
  603.  
  604. C*********************************************************
  605. C ECRITURE DES DIFFERENTES COMPOSANTES DES DEFORMATIONS
  606. C*********************************************************
  607. MPTVAL = IVAEPS
  608. DO 7000 ICOMP=1,LHOOK
  609. MELVAL = IVAL(ICOMP)
  610. VELCHE(KGAU,J) = XSTRS(ICOMP)
  611. 7000 CONTINUE
  612. C if(KGAU.eq.1) then
  613. C write(*,*) J,KGAU,'EPSI(..)=',(XSTRS(iou),iou=1,LHOOK)
  614. C endif
  615. c
  616. c
  617. 2500 CONTINUE
  618. C FIN DE BOUCLE SUR LES POINTS DE GAUSS <<<<<<<<<<<<<<
  619. C
  620.  
  621. *
  622. c quelques suppressions
  623. c (ici element non-incompressible=> pas besoin de MTRACE (cf epsi2)
  624. IF (IREPS2.EQ.1) THEN
  625. SEGSUP MTRACE
  626. ENDIF
  627.  
  628. c
  629. c quelques desactivations
  630. do IENR=1,NBENR1
  631. do I=1,NBNN1
  632. mlree1 = TLREEL(IENR,I)
  633. if(mlree1.ne.0) segdes,mlree1
  634. enddo
  635. enddo
  636. c
  637. c
  638. 2000 CONTINUE
  639. C FIN DE BOUCLE SUR LES ELEMENTS <<<<<<<<<<<<<<<<<<<<<
  640. c
  641. c
  642.  
  643. C*********************************************************
  644. C SUPPRESSION ET DESACTIVATION DE SEGMENTS
  645. C*********************************************************
  646. C
  647. SEGSUP,WRK1,WRK2
  648. SEGSUP,MRACC
  649.  
  650. segdes,MINT1
  651. if(MINT2.ne.0) segdes,MINT2
  652. C
  653. RETURN
  654. END
  655.  
  656.  
  657.  

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