Télécharger sigmax.eso

Retour à la liste

Numérotation des lignes :

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

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