Télécharger massx.eso

Retour à la liste

Numérotation des lignes :

  1. C MASSX SOURCE BP208322 18/04/12 21:15:08 9802
  2. C RIGIX SOURCE BP208322 08/12/10 21:15:26 6216
  3. subroutine MASSX (ivamat,ivacar,NMATT,CMATE,
  4. & IMODEL,LOCIRI)
  5. c
  6. C PROCEDURE UTILISEE DANS LE CAS D'ELEMENTS XFEM
  7. c POUR LE CALCUL DE RIGIDITE ELEMENTAIRES
  8. C
  9. C
  10. C*********************************************************
  11. C PARTIE DECLARATIVE
  12. C*********************************************************
  13.  
  14. C
  15. IMPLICIT REAL*8(A-H,O-Z)
  16. C
  17. CHARACTER*8 CMATE
  18. CYT PARAMETER (NDDLMAX=20,NBNIMAX=10)
  19. PARAMETER (NDDLMAX=30,NBNIMAX=10)
  20. CHARACTER*4 MOTINC(NDDLMAX),MOTDUA(NDDLMAX)
  21. CYT DATA MOTINC/'UX ','UY ','AX ','AY ',
  22. CYT >'B1X ','B1Y ','C1X ','C1Y ','D1X ','D1Y ','E1X ','E1Y ',
  23. CYT >'B2X ','B2Y ','C2X ','C2Y ','D2X ','D2Y ','E2X ','E2Y '/
  24. CYT DATA MOTDUA/ 'FX ','FY ','FAX ','FAY ',
  25. CYT >'FB1X','FB1Y','FC1X','FC1Y','FD1X','FD1Y','FE1X','FE1Y',
  26. CYT >'FB2X','FB2Y','FC2X','FC2Y','FD2X','FD2Y','FE2X','FE2Y'/
  27.  
  28. DATA MOTINC/'UX ','UY ','UZ ','AX ','AY ','AZ ',
  29. >'B1X ','B1Y ','B1Z ','C1X ','C1Y ','C1Z ','D1X ','D1Y ',
  30. >'D1Z ','E1X ','E1Y ','E1Z ','B2X ','B2Y ','B2Z ','C2X ',
  31. >'C2Y ','C2Z ','D2X ','D2Y ','D2Z ','E2X ','E2Y ','E2Z '/
  32. DATA MOTDUA/'FX ','FY ', 'FZ ','FAX ','FAY ','FAZ ',
  33. >'FB1X','FB1Y','FB1Z','FC1X','FC1Y','FC1Z','FD1X','FD1Y',
  34. >'FD1Z','FE1X','FE1Y','FE1Z','FB2X','FB2Y','FB2Z','FC2X',
  35. >'FC2Y','FC2Z','FD2X','FD2Y','FD2Z','FE2X','FE2Y','FE2Z'/
  36. C
  37. PARAMETER (NBENRMAX=5)
  38. C NBENRMAX deja defini dans rigixr.eso
  39. INTEGER LOCIRI(10,(1+NBENRMAX))
  40. DIMENSION MLRE(NBENRMAX+1)
  41. C
  42. -INC CCOPTIO
  43. -INC SMCOORD
  44. -INC SMMODEL
  45. -INC SMCHAML
  46. -INC SMINTE
  47. -INC SMELEME
  48. -INC SMRIGID
  49. -INC SMLREEL
  50. C
  51. POINTEUR MCHEX1.MCHELM,MINT1.MINTE,MINT2.MINTE
  52. C
  53. C Segment (type LISTENTI) contenant les informations sur un element
  54. SEGMENT INFO
  55. INTEGER INFELL(JG)
  56. ENDSEGMENT
  57.  
  58. SEGMENT WRK1
  59. REAL*8 XE(3,NBBB)
  60. REAL*8 REL(LRE,LRE),RINT(LRE,LRE)
  61. ENDSEGMENT
  62. C
  63. SEGMENT WRK2
  64. REAL*8 SHPWRK(6,NBNO),BGENE(LHOOK,LRE)
  65. c REAL*8 LV7WRK(NBENRMA2,2,6,NBNO)
  66. REAL*8 LV7WRK(NBENRMA2,2,6,NBBB)
  67. ENDSEGMENT
  68. C
  69. SEGMENT MVELCH
  70. REAL*8 VALMAT(NV1)
  71. ENDSEGMENT
  72. C
  73. SEGMENT MPTVAL
  74. INTEGER IPOS(NS),NSOF(NS)
  75. INTEGER IVAL(NCOSOU)
  76. CHARACTER*16 TYVAL(NCOSOU)
  77. ENDSEGMENT
  78. C
  79. SEGMENT MRACC
  80. INTEGER TLREEL(NBENRMA2,NBI)
  81. ENDSEGMENT
  82. C
  83. c write (ioimp,*) '############################'
  84. c write (ioimp,*) '# DEBUT DE MASSX #'
  85. c write (ioimp,*) '############################'
  86. C
  87. C*********************************************************
  88. c
  89. C RECUPERATION + ACTIVATIONS + VALEURS UTILES
  90. c
  91. C*********************************************************
  92. C sont deja actifs dans rigi1.eso :
  93. C SEGACT MMODEL,IMODEL,MELVAL
  94. c IVAMAT
  95. MPTVAL=IVAMAT
  96. c
  97. C++++ Activation au cas ou ++++++++++++++++++++++++++++++
  98. SEGACT MCOORD
  99.  
  100. C++++ Recup + Activation de la geometrie ++++++++++++++++
  101. MELEME= IMAMOD
  102. c SEGACT MELEME deja actif
  103.  
  104. C++++ RECUP DES INFOS EF ++++++++++++++++++++++++++++++++
  105. c + OBTENUES PAR ELQUOI DANS MASSE1
  106. C segment INFO OBTENUES PAR ELQUOI DANS MASSE1 mais supprimé dans MASSE1
  107. MELE = NEFMOD
  108. call elquoi(MELE,0,4,IPINF,IMODEL)
  109. INFO = IPINF
  110. c MELE = INFELL(1)
  111. c NBPGA2= INFELL(2)
  112. c NBPGAU= INFELL(3)
  113. c NBPGAU= INFELL(4)
  114. c ICARA = INFELL(5)
  115. NGAU1 = INFELL(6)
  116. c LW = INFELL(7)
  117. LRE = INFELL(9)
  118. LHOOK = INFELL(10)
  119. MINT1 = INFELL(11)
  120. segact,MINT1
  121. MINT2 = INFELL(12)
  122. if(MINT2.ne.0) segact,MINT2
  123. MFR = INFELL(13)
  124. IELE = INFELL(14)
  125. NDDL = INFELL(15)
  126. NSTRS = INFELL(16)
  127. c write(ioimp,*)'-> MASSX infell',(infell(iou),iou=1,16)
  128.  
  129. c + AUTRES INFOS
  130. C nbre de noeuds par element
  131. NBNN1 = NUM(/1)
  132. C nbre d elements
  133. NBEL1 = NUM(/2)
  134.  
  135. c REM: pour se passer du dimensionnement du nbre d'enrichissement dans
  136. c elquoi et le realiser localement , on pourrait ecrire:
  137. c LRE = NDDLMAX*NBNN1
  138. c NDDL= NDDLMAX
  139.  
  140. C sous decoupage et points de Gauss de l element geometrique de base
  141. CYT if(MELE.eq.263) then
  142. if(MELE.eq.263.or.MELE.eq.264) then
  143. c NGAU2 = 4
  144. NGAU2 = MINT2.POIGAU(/1)
  145. c NBSSEF = (NGAU1/NGAU2)**0.5
  146. endif
  147. c write(ioimp,*) 'dim de MINT2=6,',(MINT2.SHPTOT(/2)),(MINT2.SHPTOT(/3))
  148. c write(ioimp,*) 'MINT2',(MINT2.QSIGAU(iou),iou=1,NGAU)
  149. c segdes,MINT2
  150.  
  151.  
  152. c nbre maxi de fonction de forme par noeud (fonction std comprise)
  153. c NBNI = NDDL/IDIM inutile!
  154.  
  155.  
  156. C++++ Recup des infos d enrichissement +++++++++++++++++++
  157. c recup du MCHEX1 d'enrichissement
  158. NOBMO1 = IVAMOD(/1)
  159. if(NOBMO1.ne.0) then
  160. do iobmo1=1,NOBMO1
  161. if((TYMODE(iobmo1)).eq.'MCHAML') then
  162. MCHEX1 = IVAMOD(iobmo1)
  163. segact,MCHEX1
  164. if((MCHEX1.TITCHE).eq.'ENRICHIS') then
  165. MCHAM1 = MCHEX1.ICHAML(1)
  166. segact,MCHAM1
  167. segdes,MCHEX1
  168. goto 1000
  169. endif
  170. segdes,MCHEX1
  171. endif
  172. enddo
  173. write(ioimp,*) 'Le modele est vide (absence d enrichissement)'
  174. * return
  175. else
  176. write(ioimp,*) 'Aucun MCHAML enrichissement dans le Modele'
  177. * return
  178. endif
  179.  
  180. 1000 continue
  181. c niveau d enrichissement(s) du modele (ddls std u exclus)
  182. c NBENR1= 0 si std, 1 si H seul, 2 si H+F1, 3 si H+F1+F2, etc...
  183. if(NOBMO1.ne.0) then
  184. NBENR1= MCHAM1.IELVAL(/1)
  185. else
  186. NBENR1 = 0
  187. endif
  188. c write(ioimp,*) 'niveau d enrichissement(s) du modele',NBENR1
  189.  
  190.  
  191. C*********************************************************
  192. c
  193. c PHASE 2 : PREPARATION DES OBJETS RESULTATS
  194. c
  195. C*********************************************************
  196. C
  197. C++++ RECHERCHE DES NOMS D'INCONNUES ET DES DUAUX
  198. c MOTINC et MODUA en dur pour l instant
  199.  
  200. C++++ REMPLISSAGE DE DESCR de taille maxi
  201. c (on ajustera en enlevant si besoin)
  202. NLIGRP = LRE
  203. NLIGRD = LRE
  204. SEGINI DESCR
  205. IPDSCR = DESCR
  206. C
  207. KINC = 0
  208. C----->> BOUCLE SUR LES fonctions de Forme
  209. DO IENR=1,NBNIMAX
  210. C------->> BOUCLE SUR LES NOEUDS
  211. DO I=1,NBNN1
  212. C--------->> BOUCLE SUR LA DIMENSION
  213. DO KDIM=1,IDIM
  214. KINC = KINC + 1
  215. LISINC(KINC) = MOTINC((3*(IENR-1))+KDIM)
  216. LISDUA(KINC) = MOTDUA((3*(IENR-1))+KDIM)
  217. CYT LISINC(KINC) = MOTINC((IDIM*(IENR-1))+KDIM)
  218. CYT LISDUA(KINC) = MOTDUA((IDIM*(IENR-1))+KDIM)
  219. NOELEP(KINC) = I
  220. NOELED(KINC) = I
  221. ENDDO
  222. C--------------|
  223. ENDDO
  224. C-------------|
  225. ENDDO
  226. C------------|
  227.  
  228. IF(KINC.NE.LRE) THEN
  229. WRITE(*,*) 'IL Y A UNE ERREUR DANS DESCR'
  230. WRITE(*,*) 'KINC , LRE = ',KINC,LRE
  231. RETURN
  232. ENDIF
  233. C
  234. SEGDES DESCR
  235. C
  236. C*********************************************************
  237. C INITIALISATIONS...
  238. C*********************************************************
  239. C
  240. c preparation des tables avec:
  241.  
  242. if(NBENR1.ne.0) then
  243. do ienr=1,NBENR1
  244. c -le nombre d'inconnues de chaque sous-zone
  245. c determinee depuis le nombre de fonction de forme
  246. c ienr= 1: U+H(1+1=2), 2: U+H+F1(2+4=6), 3: U+H+F1+F2(6+4=10)
  247. nbniJ = 2 + ((ienr-1)*4)
  248. MLRE(1+ienr) = IDIM*NBNN1*nbniJ
  249. c -LOCIRI
  250. LOCIRI(5,1+ienr)= NIFOUR
  251. enddo
  252. endif
  253. C Tables + longues car 1er indice correspond au fontion de forme std
  254. MLRE(1) = IDIM*NBNN1*1
  255. LOCIRI(5,1)= NIFOUR
  256.  
  257. if(NBENR1.lt.(NBENRMAX+1)) then
  258. do ienr=(NBENR1+1),(NBENRMAX)
  259. MLRE(1+ienr) = 0
  260. enddo
  261. endif
  262. c
  263. c ...DU SEGMENT WRK1
  264. NBENRMA2 = NBENRMAX
  265. NBBB = NBNN1
  266. SEGINI,WRK1
  267.  
  268. c ...DU SEGMENT WRK2
  269. c NBNO = NBNI
  270. NBNO = LRE/IDIM
  271. c ici, pour le calcul de la masse, BGENE ne designe pas B mais N => lhook=2
  272. CTY LHOOK=2
  273. LHOOK=IDIM
  274. SEGINI,WRK2
  275.  
  276. c ...DU SEGMENT MVELCH
  277. NV1 = NMATT
  278. SEGINI,MVELCH
  279. C
  280. c ...DU SEGMENT MRACC
  281. NBENRMA2 = NBENRMAX
  282. NBI = NBNN1
  283. segini,MRACC
  284. C
  285. C
  286. C*********************************************************
  287. C
  288. C>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> BOUCLE SUR LES ELEMENTS
  289. C
  290. DO 2000 J=1,NBEL1
  291. c write(ioimp,*) '========= EF',J,' /NBEL1 ========='
  292. C
  293. C
  294. C*********************************************************
  295. C POUR CHAQUE ELEMENT, ...
  296. C*********************************************************
  297. C
  298. C++++ ON RECUPERE LES COORDONNEES DES NOEUDS DE L ELEMENT
  299. CALL DOXE(XCOOR,IDIM,NBNN1,NUM,J,XE)
  300. C
  301. c
  302. c CALL ZEROI(TLREEL,NBENRMA2,NBI)
  303. C++++ NBENRJ = niveau d enrichissement de l element ++++
  304. C =0 si EF std =1 si U+H =2 si U+H+F1 =3 si U+H+F1+F2
  305. NBENRJ=0
  306. if(NBENR1.ne.0) then
  307. do IENR=1,NBENR1
  308. MELVA1 = MCHAM1.IELVAL(IENR)
  309. segact,MELVA1
  310. do I=1,NBNN1
  311. mlree1 = MELVA1.IELCHE(I,J)
  312. C on en profite pour remplir MRACC table de raccourcis pour cet element
  313. TLREEL(IENR,I) = mlree1
  314. if(mlree1.ne.0) then
  315. NBENRJ = max(NBENRJ,IENR)
  316. C et on active la listreel
  317. segact,mlree1
  318. endif
  319. enddo
  320. segdes,MELVA1
  321. enddo
  322. endif
  323. if(NBENRMA2.gt.NBENR1) then
  324. do IENR2=(NBENR1+1),NBENRMA2
  325. do I=1,NBNN1
  326. TLREEL(IENR2,I) = 0
  327. enddo
  328. enddo
  329. endif
  330. C
  331.  
  332. c++++ quelques indices pour dimensionner au plus juste
  333. c nbre total de ddl de l'élément considéré
  334. NLIGRD = MLRE(1+NBENRJ)
  335. NLIGRP = MLRE(1+NBENRJ)
  336. c nbre fonction de forme=((Ni_std+Ni_enrichi)*nbnoeud)=Ninconnues/idim
  337. NBNI = (MLRE(1+NBENRJ)) / IDIM
  338.  
  339. c write(ioimp,*) 'EF',J,' NBENRJ=',NBENRJ,
  340. c &'donc',NLIGRD,' ddls et ',NBNI,' fonctions de forme'
  341. c
  342.  
  343. c++++ Selon le niveau d enrichissement,
  344. c
  345. C + On copie cet element a la geo correspondante
  346. IPT1 = LOCIRI(1,1+NBENRJ)
  347. c si elle n existe pas, il faut la creer
  348. if(IPT1.eq.0) then
  349. NBNN =NBNN1
  350. NBELEM=1
  351. NBSOUS=0
  352. NBREF=0
  353. segini,IPT1
  354. IPT1.ITYPEL = ITYPEL
  355. LOCIRI(1,1+NBENRJ)=IPT1
  356. c si elle existe, il faut l agrandir
  357. else
  358. NBNN =NBNN1
  359. NBELEM = (IPT1.NUM(/2)) + 1
  360. NBSOUS=0
  361. NBREF=0
  362. segadj,IPT1
  363. endif
  364. c copie en cours
  365. do I=1,NBNN1
  366. IPT1.NUM(I,NBELEM) = NUM(I,J)
  367. enddo
  368.  
  369. C + On crée le descr qui va bien si il n existe pas
  370. DES1 = LOCIRI(3,1+NBENRJ)
  371. if(DES1.eq.0) then
  372. segini,DES1=DESCR
  373. c write(ioimp,*)' on ne garde que les',MLRE(1+NBENRJ),
  374. c &' 1eres inconnues'
  375. c NLIGRP = MLRE(1+NBENRJ) fait + haut (ligne 313)
  376. c NLIGRD = MLRE(1+NBENRJ)
  377. segadj,DES1
  378. LOCIRI(3,1+NBENRJ) = DES1
  379. segdes,DES1
  380. endif
  381.  
  382. C + On en profite pour créer IMATRI si il n existe pas
  383. xMATR1 = LOCIRI(4,1+NBENRJ)
  384. if(xMATR1.eq.0) then
  385. NELRIG = 1
  386. segini,xMATR1
  387. LOCIRI(4,1+NBENRJ) = xMATR1
  388. c si il existe, il faut l agrandir
  389. else
  390. NELRIG = xmatr1.RE(/3) + 1
  391. segadj,xMATR1
  392. endif
  393.  
  394.  
  395. C++++ ON MET A ZERO LA SOMME QUE L ON VA CALCULER
  396. c CALL ZERO(RINT,LRE,NLIGRP)
  397. * CALL ZERO(RINT,NLIGRD,NLIGRP)
  398. c |-> impossible car ZERO considere qu il sagit de la dimension de SHPWRK
  399. CALL ZEROX(RINT,NLIGRD,NLIGRP,LRE)
  400.  
  401. c
  402. c rem: il serait probablement interessant au niveau du tems cpu
  403. c d'utiliser moins de pts de Gauss lorsque l element n est pas enrichi
  404. c car cela n'apprte rien de plus
  405. c On pourrait par ex. utiliser MINT2 = infell(12) qui contiendrait
  406. c le segment d'integration de l'EF std (QUA4 par ex.)
  407. if((NBENRJ.eq.0).and.(MINT2.ne.0)) then
  408. MINTE = MINT2
  409. NBPGAU= NGAU2
  410. else
  411. MINTE = MINT1
  412. NBPGAU= NGAU1
  413. endif
  414. C
  415. C*********************************************************
  416. C
  417. C>>>>>>>>>> BOUCLE SUR LES POINTS DE GAUSS
  418. C
  419. DO 2500 KGAU=1,NBPGAU
  420. C
  421. C*********************************************************
  422. C Initialisation à 0
  423. C*********************************************************
  424.  
  425. c ZERO ne serait-il pas facultatif?
  426. * CALL ZERO(SHPWRK,6,NBNO)
  427. CALL ZERO(SHPWRK,6,NBNI)
  428. C
  429. i6zz = 3
  430. IF (IDIM.EQ.3) i6zz = 4
  431. C
  432. c do ienr7=1,NBENRMAX
  433. do ienr7=1,NBENRJ
  434. do inod7=1,NBNN1
  435. c do i6=1,6
  436. CYT do i6=1,3
  437. do i6=1,i6zz
  438. LV7WRK(ienr7,1,i6,inod7) = 0
  439. LV7WRK(ienr7,2,i6,inod7) = 0
  440. enddo
  441. enddo
  442. enddo
  443.  
  444. c*********************************************************
  445. c Calcul des fonction de forme std dans repere local
  446. c*********************************************************
  447.  
  448. ccccc BOUCLE SUR LES NOEUDS ccccccccccccccccccccccccccccc
  449. c (et donc sur les Ni std)
  450. DO 2510 I=1,NBNN1
  451.  
  452. C++++ Calcul des Ni std
  453. c (rappel: 1:Ni 2:Ni,qsi 3:Ni,eta avec i=1,4)
  454. SHPWRK(1,I) = SHPTOT(1,I,KGAU)
  455. SHPWRK(2,I) = SHPTOT(2,I,KGAU)
  456. SHPWRK(3,I) = SHPTOT(3,I,KGAU)
  457. IF (IDIM.EQ.3) SHPWRK(4,I) = SHPTOT(4,I,KGAU)
  458.  
  459. 2510 CONTINUE
  460. ccccc fin de BOUCLE SUR LES NOEUDS ccccccccccccccccccccccc
  461.  
  462.  
  463.  
  464. c*********************************************************
  465. c Passage des fonctions de forme std dans repere global
  466. c*********************************************************
  467. c
  468. C++++ CALCUL DES Ni,x Ni,y (i=1,4) + CALCUL DE det(J)
  469. CALL JACOBI(XE,SHPWRK,IDIM,NBNN1,DJAC)
  470. c
  471. C PRISE EN COMPTE DU POIDS D'INTEGRATION
  472. DJAC = ABS(DJAC) * POIGAU(KGAU)
  473.  
  474. c CAS DES CONTRAINTES PLANES
  475. C recuperation de l'epaisseur: DIM3 donnée facultative du materiau
  476. IF (IFOUR.EQ.-2) THEN
  477. MPTVAL=IVACAR
  478. IF (IVACAR.NE.0) THEN
  479. MELVAL=IVAL(1)
  480. IF (MELVAL.NE.0) THEN
  481. IGMN=MIN(KGAU,VELCHE(/1))
  482. IBMN=MIN(J,VELCHE(/2))
  483. DIM3=VELCHE(IGMN,IBMN)
  484. ELSE
  485. DIM3=1.D0
  486. ENDIF
  487. ENDIF
  488. DJAC=DJAC*DIM3
  489. MPTVAL=IVAMAT
  490. ENDIF
  491.  
  492.  
  493. c*********************************************************
  494. c Si on est pas du tout enrichi on peut sauter une grosse partie
  495. c*********************************************************
  496. if(NBENRJ.eq.0) goto 2999
  497.  
  498. c*********************************************************
  499. c Calcul des level set + leurs derivees dans repere global
  500. c*********************************************************
  501.  
  502. ccccc BOUCLE SUR LES enrichissements ccccccccccccccccccc
  503. do 2520 ienr=1,NBENRJ
  504.  
  505. c MELVA1=MCHAM1.IELVAL(IENR)
  506. c segact,MELVA1
  507.  
  508. ccccc BOUCLE SUR LES NOEUDS ccccccccccccccccccccccccccc
  509. DO 2521 I=1,NBNN1
  510.  
  511. C++++ Le I eme noeud est-il ienr-enrichi?
  512. mlree1= TLREEL(IENR,I)
  513. if(mlree1.eq.0) goto 2521
  514.  
  515.  
  516. c Calcul du repere local de fissure(=PSI,PHI)
  517. c (rappel: 1,1:psi 1,2:phi 2,1 psi,x ...etc...)
  518. do 2522 inode=1,NBNN1
  519. c pour le H-enrichissement, on n a pas gardé PSI (inutile)
  520. if(ienr.ne.1) then
  521. c valeur de PSI au inode^ieme noeud
  522. xpsi1 = mlree1.PROG(inode)
  523. c qu on multiplie par la valeur de Ni^std au pt de G considéré
  524. LV7WRK(ienr,1,1,I)= LV7WRK(ienr,1,1,I)
  525. & + (SHPWRK(1,inode)*xpsi1)
  526. LV7WRK(ienr,1,2,I)= LV7WRK(ienr,1,2,I)
  527. & + (SHPWRK(2,inode)*xpsi1)
  528. LV7WRK(ienr,1,3,I)= LV7WRK(ienr,1,3,I)
  529. & + (SHPWRK(3,inode)*xpsi1)
  530. IF (IDIM.EQ.3) LV7WRK(ienr,1,4,I)= LV7WRK(ienr,1,4,I)
  531. & + (SHPWRK(4,inode)*xpsi1)
  532. c valeur de PHI au inode^ieme noeud
  533. xphi1 = mlree1.PROG(NBNN1+inode)
  534. else
  535. xphi1 = mlree1.PROG(inode)
  536. endif
  537. LV7WRK(ienr,2,1,I)= LV7WRK(ienr,2,1,I)
  538. & + (SHPWRK(1,inode)*xphi1)
  539. LV7WRK(ienr,2,2,I)= LV7WRK(ienr,2,2,I)
  540. & + (SHPWRK(2,inode)*xphi1)
  541. LV7WRK(ienr,2,3,I)= LV7WRK(ienr,2,3,I)
  542. & + (SHPWRK(3,inode)*xphi1)
  543. IF (IDIM.EQ.3) LV7WRK(ienr,2,4,I)= LV7WRK(ienr,2,4,I)
  544. & + (SHPWRK(4,inode)*xphi1)
  545. 2522 continue
  546.  
  547. 2521 continue
  548. ccccc fin de BOUCLE SUR LES NOEUDS ccccccccccccccccccccccc
  549.  
  550.  
  551. 2520 CONTINUE
  552. ccccc fin de BOUCLE SUR LES enrichissements cccccccccccccccc
  553.  
  554. c on a construit
  555. C LV7WRK(ienr, PSI/PHI, valeur/deriveeparqsi/pareta, iNOEUD)
  556.  
  557.  
  558. c*********************************************************
  559. c Ajout des fonctions de forme d enrichissement
  560. c + leurs derivees dans repere global
  561. c*********************************************************
  562. CYT CALL SHAPX(XGAU,YGAU,MELE,LV7WRK,NBENRMAX,NBENRJ,
  563. CALL SHAPX(MELE,LV7WRK,NBENRMAX,NBENRJ,TLREEL,SHPWRK,IRET)
  564.  
  565.  
  566. c retour a la partie commune aux EF enrichi et non enrichi
  567. 2999 continue
  568.  
  569. C*********************************************************
  570. C CALCUL DE [N]
  571. C*********************************************************
  572. c ZERO ne serait-il pas facultatif?
  573. c call ZERO(BGENE,LHOOK,NLIGRP)
  574. KB=1
  575. C boucle sur tous les Ni
  576. DO 3001 II=1,NBNI
  577.  
  578. BGENE(1,KB) = SHPWRK(1,II)
  579. BGENE(2,KB+1) = SHPWRK(1,II)
  580.  
  581. IF(IDIM.EQ.3) BGENE(3,KB+2) = SHPWRK(1,II)
  582.  
  583. KB = KB + IDIM
  584. CTY KB = KB + 2
  585.  
  586. 3001 CONTINUE
  587. C
  588. c if(J.eq.1.and.KGAU.eq.1) then
  589. c write(ioimp,*) 'BGENE(1,..)=',(BGENE(1,iou),iou=1,2*NBNI)
  590. c write(ioimp,*) 'BGENE(2,..)=',(BGENE(2,iou),iou=1,2*NBNI)
  591. c endif
  592.  
  593. C*********************************************************
  594. C CALCUL DE rho
  595. C*********************************************************
  596. c
  597. C on remplit le segment MVELCH => valmat(1)=RHO
  598. IF (IVAL(1).NE.0) THEN
  599. MELVAL=IVAL(1)
  600. IBMN=MIN(J ,VELCHE(/2))
  601. IGMN=MIN(KGAU,VELCHE(/1))
  602. VALMAT(1)=VELCHE(IGMN,IBMN)
  603. ELSE
  604. VALMAT(1)=0.D0
  605. ENDIF
  606.  
  607. C PRISE EN COMPTE DE RHO
  608. DJAC=DJAC*VALMAT(1)
  609. c if(J.eq.1) write(ioimp,*) 'DJAC=',DJAC
  610. c
  611. c
  612. C*********************************************************
  613. C CALCUL DE rho.Nt.N
  614. C*********************************************************
  615.  
  616. CALL ZEROX(REL,NLIGRD,NLIGRP,LRE)
  617.  
  618. CALL NTNST(BGENE,DJAC,LRE,LHOOK,REL)
  619. c pour gagner du temps cpu il faudrait qqchose du type:
  620. c CALL NTNSTX(BGENE,DJAC,DDHOOK,NLIGRP,2,LRE,REL)
  621. c if(J.eq.1) then
  622. c write(ioimp,*) 'REL',(REL(1,iou),iou=1,8)
  623. c write(ioimp,*) 'REL',(REL(2,iou),iou=1,8)
  624. c write(ioimp,*) 'REL',(REL(3,iou),iou=1,8)
  625. c write(ioimp,*) 'REL',(REL(4,iou),iou=1,8)
  626. c endif
  627.  
  628.  
  629. C*********************************************************
  630. C CUMUL DANS RINT(II,JJ)
  631. C*********************************************************
  632. DO 4201 II=1,NLIGRD
  633. DO 4202 JJ=1,NLIGRP
  634. RINT(II,JJ) = RINT(II,JJ) + REL(II,JJ)
  635. 4202 CONTINUE
  636. 4201 CONTINUE
  637. c
  638. c
  639. 2500 CONTINUE
  640. C FIN DE BOUCLE SUR LES POINTS DE GAUSS <<<<<<<<<<<<<<
  641. C
  642.  
  643.  
  644. C*********************************************************
  645. c write(ioimp,*) 'C STOCKAGE DANS XMATRI de RE(II,JJ)'
  646. C et de XMATRI dans IMATRI
  647. C*********************************************************
  648. C
  649. * SEGINI,XMATRI
  650. * IMATR1.IMATTT(NELRIG) = XMATRI
  651. c
  652. DO 6001 II=1,NLIGRD
  653.  
  654. DO 6002 JJ=1,NLIGRP
  655.  
  656. xmatr1.RE(II,JJ,nelrig) = RINT(II,JJ)
  657.  
  658. c si NON enrichi, on met tout a 0 sauf la diago
  659. c if(mlree1.eq.0) RE(II,JJ) = 0. inutile
  660.  
  661. 6002 CONTINUE
  662. c if(J.eq.1) write(ioimp,*) 'RE_',II,'=',(RE(II,jou),jou=1,NLIGRP)
  663.  
  664. 6001 CONTINUE
  665. c
  666. c quelques desactivations
  667. CYT SEGDES,XMATRI
  668. do IENR=1,NBENR1
  669. do I=1,NBNN1
  670. mlree1 = TLREEL(IENR,I)
  671. if(mlree1.ne.0) segdes,mlree1
  672. enddo
  673. enddo
  674. c
  675. c
  676. c
  677. 2000 CONTINUE
  678. C FIN DE BOUCLE SUR LES ELEMENTS <<<<<<<<<<<<<<<<<<<<<
  679. c
  680. c
  681.  
  682. C*********************************************************
  683. C SUPPRESSION ET DESACTIVATION DE SEGMENTS
  684. C*********************************************************
  685. C
  686. SEGSUP,WRK1,WRK2,MVELCH
  687. SEGSUP,MRACC
  688.  
  689. segdes,MELEME
  690. segdes,MINT1
  691. if(MINT2.ne.0) segdes,MINT2
  692. c desactivation des maillages et segment imatri
  693. do ienr=1,(1+NBENR1)
  694. IPT1 = LOCIRI(1,ienr)
  695. if(IPT1.ne.0) segdes,IPT1
  696. xMATR1 = LOCIRI(4,ienr)
  697. if(xMATR1.ne.0) segdes,xMATR1
  698. enddo
  699.  
  700. c SEGDES dans RIGI1 = IMODEL
  701. C
  702. RETURN
  703. END
  704.  
  705.  
  706.  
  707.  
  708.  

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