Numérotation des lignes :

1. C MASSX SOURCE BP208322 13/10/02 21:15:20 7830
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
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)
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
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.
429. XGAU = 0.
430. YGAU = 0.
431. ZGAU = 0.
432. C
433. i6zz = 3
434. IF (IDIM.EQ.3) i6zz = 4
435. C
436. c do ienr7=1,NBENRMAX
437. do ienr7=1,NBENRJ
438. do inod7=1,NBNN1
439. c do i6=1,6
440. CYT do i6=1,3
441. do i6=1,i6zz
442. LV7WRK(ienr7,1,i6,inod7) = 0
443. LV7WRK(ienr7,2,i6,inod7) = 0
444. enddo
445. enddo
446. enddo
447.
448. c*********************************************************
449. c Calcul des fonction de forme std dans repere local
450. c*********************************************************
451.
452. ccccc BOUCLE SUR LES NOEUDS ccccccccccccccccccccccccccccc
453. c (et donc sur les Ni std)
454. DO 2510 I=1,NBNN1
455.
456. C++++ Calcul des Ni std
457. c (rappel: 1:Ni 2:Ni,qsi 3:Ni,eta avec i=1,4)
458. SHPWRK(1,I) = SHPTOT(1,I,KGAU)
459. SHPWRK(2,I) = SHPTOT(2,I,KGAU)
460. SHPWRK(3,I) = SHPTOT(3,I,KGAU)
461. IF (IDIM.EQ.3) SHPWRK(4,I) = SHPTOT(4,I,KGAU)
462.
463. C++++ Calcul des coordonnees dans le repere global
464. XGAU = XGAU + (SHPWRK(1,I)*XE(1,I))
465. YGAU = YGAU + (SHPWRK(1,I)*XE(2,I))
466. IF (IDIM.EQ.3) ZGAU = ZGAU + (SHPWRK(1,I)*XE(3,I))
467.
468. 2510 CONTINUE
469. ccccc fin de BOUCLE SUR LES NOEUDS ccccccccccccccccccccccc
470.
471.
472.
473. c*********************************************************
474. c Passage des fonctions de forme std dans repere global
475. c*********************************************************
476. c
477. C++++ CALCUL DES Ni,x Ni,y (i=1,4) + CALCUL DE det(J)
478. CALL JACOBI(XE,SHPWRK,IDIM,NBNN1,DJAC)
479. c
480. C PRISE EN COMPTE DU POIDS D'INTEGRATION
481. DJAC = ABS(DJAC) * POIGAU(KGAU)
482.
483. c CAS DES CONTRAINTES PLANES
484. C recuperation de l'epaisseur: DIM3 donnée facultative du materiau
485. IF (IFOUR.EQ.-2) THEN
486. MPTVAL=IVACAR
487. IF (IVACAR.NE.0) THEN
488. MELVAL=IVAL(1)
489. IF (MELVAL.NE.0) THEN
490. IGMN=MIN(KGAU,VELCHE(/1))
491. IBMN=MIN(J,VELCHE(/2))
492. DIM3=VELCHE(IGMN,IBMN)
493. ELSE
494. DIM3=1.D0
495. ENDIF
496. ENDIF
497. DJAC=DJAC*DIM3
498. MPTVAL=IVAMAT
499. ENDIF
500.
501.
502. c*********************************************************
503. c Si on est pas du tout enrichi on peut sauter une grosse partie
504. c*********************************************************
505. if(NBENRJ.eq.0) goto 2999
506.
507. c*********************************************************
508. c Calcul des level set + leurs derivees dans repere global
509. c*********************************************************
510.
511. ccccc BOUCLE SUR LES enrichissements ccccccccccccccccccc
512. do 2520 ienr=1,NBENRJ
513.
514. c MELVA1=MCHAM1.IELVAL(IENR)
515. c segact,MELVA1
516.
517. ccccc BOUCLE SUR LES NOEUDS ccccccccccccccccccccccccccc
518. DO 2521 I=1,NBNN1
519.
520. C++++ Le I eme noeud est-il ienr-enrichi?
521. mlree1= TLREEL(IENR,I)
522. if(mlree1.eq.0) goto 2521
523.
524.
525. c Calcul du repere local de fissure(=PSI,PHI)
526. c (rappel: 1,1:psi 1,2:phi 2,1 psi,x ...etc...)
527. do 2522 inode=1,NBNN1
528. c pour le H-enrichissement, on n a pas gardé PSI (inutile)
529. if(ienr.ne.1) then
530. c valeur de PSI au inode^ieme noeud
531. xpsi1 = mlree1.PROG(inode)
532. c qu on multiplie par la valeur de Ni^std au pt de G considéré
533. LV7WRK(ienr,1,1,I)= LV7WRK(ienr,1,1,I)
534. & + (SHPWRK(1,inode)*xpsi1)
535. LV7WRK(ienr,1,2,I)= LV7WRK(ienr,1,2,I)
536. & + (SHPWRK(2,inode)*xpsi1)
537. LV7WRK(ienr,1,3,I)= LV7WRK(ienr,1,3,I)
538. & + (SHPWRK(3,inode)*xpsi1)
539. IF (IDIM.EQ.3) LV7WRK(ienr,1,4,I)= LV7WRK(ienr,1,4,I)
540. & + (SHPWRK(4,inode)*xpsi1)
541. c valeur de PHI au inode^ieme noeud
542. xphi1 = mlree1.PROG(NBNN1+inode)
543. else
544. xphi1 = mlree1.PROG(inode)
545. endif
546. LV7WRK(ienr,2,1,I)= LV7WRK(ienr,2,1,I)
547. & + (SHPWRK(1,inode)*xphi1)
548. LV7WRK(ienr,2,2,I)= LV7WRK(ienr,2,2,I)
549. & + (SHPWRK(2,inode)*xphi1)
550. LV7WRK(ienr,2,3,I)= LV7WRK(ienr,2,3,I)
551. & + (SHPWRK(3,inode)*xphi1)
552. IF (IDIM.EQ.3) LV7WRK(ienr,2,4,I)= LV7WRK(ienr,2,4,I)
553. & + (SHPWRK(4,inode)*xphi1)
554. 2522 continue
555.
556. 2521 continue
557. ccccc fin de BOUCLE SUR LES NOEUDS ccccccccccccccccccccccc
558.
559.
560. 2520 CONTINUE
561. ccccc fin de BOUCLE SUR LES enrichissements cccccccccccccccc
562.
563. c on a construit
564. C LV7WRK(ienr, PSI/PHI, valeur/deriveeparqsi/pareta, iNOEUD)
565.
566.
567. c*********************************************************
568. c Ajout des fonctions de forme d enrichissement
569. c + leurs derivees dans repere global
570. c*********************************************************
571. CYT CALL SHAPX(XGAU,YGAU,MELE,LV7WRK,NBENRMAX,NBENRJ,
572. CALL SHAPX(XGAU,YGAU,ZGAU,MELE,LV7WRK,NBENRMAX,NBENRJ,
573. &TLREEL,J,SHPWRK,IRET)
574.
575.
576. c retour a la partie commune aux EF enrichi et non enrichi
577. 2999 continue
578.
579. C*********************************************************
580. C CALCUL DE [N]
581. C*********************************************************
582. c ZERO ne serait-il pas facultatif?
583. c call ZERO(BGENE,LHOOK,NLIGRP)
584. KB=1
585. C boucle sur tous les Ni
586. DO 3001 II=1,NBNI
587.
588. BGENE(1,KB) = SHPWRK(1,II)
589. BGENE(2,KB+1) = SHPWRK(1,II)
590.
591. IF(IDIM.EQ.3) BGENE(3,KB+2) = SHPWRK(1,II)
592.
593. KB = KB + IDIM
594. CTY KB = KB + 2
595.
596. 3001 CONTINUE
597. C
598. c if(J.eq.1.and.KGAU.eq.1) then
599. c write(ioimp,*) 'BGENE(1,..)=',(BGENE(1,iou),iou=1,2*NBNI)
600. c write(ioimp,*) 'BGENE(2,..)=',(BGENE(2,iou),iou=1,2*NBNI)
601. c endif
602.
603. C*********************************************************
604. C CALCUL DE rho
605. C*********************************************************
606. c
607. C on remplit le segment MVELCH => valmat(1)=RHO
608. IF (IVAL(1).NE.0) THEN
609. MELVAL=IVAL(1)
610. IBMN=MIN(J ,VELCHE(/2))
611. IGMN=MIN(KGAU,VELCHE(/1))
612. VALMAT(1)=VELCHE(IGMN,IBMN)
613. ELSE
614. VALMAT(1)=0.D0
615. ENDIF
616.
617. C PRISE EN COMPTE DE RHO
618. DJAC=DJAC*VALMAT(1)
619. c if(J.eq.1) write(ioimp,*) 'DJAC=',DJAC
620. c
621. c
622. C*********************************************************
623. C CALCUL DE rho.Nt.N
624. C*********************************************************
625.
626. CALL ZEROX(REL,NLIGRD,NLIGRP,LRE)
627.
628. CALL NTNST(BGENE,DJAC,LRE,LHOOK,REL)
629. c pour gagner du temps cpu il faudrait qqchose du type:
630. c CALL NTNSTX(BGENE,DJAC,DDHOOK,NLIGRP,2,LRE,REL)
631. c if(J.eq.1) then
632. c write(ioimp,*) 'REL',(REL(1,iou),iou=1,8)
633. c write(ioimp,*) 'REL',(REL(2,iou),iou=1,8)
634. c write(ioimp,*) 'REL',(REL(3,iou),iou=1,8)
635. c write(ioimp,*) 'REL',(REL(4,iou),iou=1,8)
636. c endif
637.
638.
639. C*********************************************************
640. C CUMUL DANS RINT(II,JJ)
641. C*********************************************************
642. DO 4201 II=1,NLIGRD
643. DO 4202 JJ=1,NLIGRP
644. RINT(II,JJ) = RINT(II,JJ) + REL(II,JJ)
645. 4202 CONTINUE
646. 4201 CONTINUE
647. c
648. c
649. 2500 CONTINUE
650. C FIN DE BOUCLE SUR LES POINTS DE GAUSS &lt;&lt;&lt;&lt;&lt;&lt;&lt;&lt;&lt;&lt;&lt;&lt;&lt;&lt;
651. C
652.
653.
654. C*********************************************************
655. c write(ioimp,*) 'C STOCKAGE DANS XMATRI de RE(II,JJ)'
656. C et de XMATRI dans IMATRI
657. C*********************************************************
658. C
659. * SEGINI,XMATRI
660. * IMATR1.IMATTT(NELRIG) = XMATRI
661. c
662. DO 6001 II=1,NLIGRD
663.
664. DO 6002 JJ=1,NLIGRP
665.
666. xmatr1.RE(II,JJ,nelrig) = RINT(II,JJ)
667.
668. c si NON enrichi, on met tout a 0 sauf la diago
669. c if(mlree1.eq.0) RE(II,JJ) = 0. inutile
670.
671. 6002 CONTINUE
672. c if(J.eq.1) write(ioimp,*) 'RE_',II,'=',(RE(II,jou),jou=1,NLIGRP)
673.
674. 6001 CONTINUE
675. c
676. c quelques desactivations
677. CYT SEGDES,XMATRI
678. do IENR=1,NBENR1
679. do I=1,NBNN1
680. mlree1 = TLREEL(IENR,I)
681. if(mlree1.ne.0) segdes,mlree1
682. enddo
683. enddo
684. c
685. c
686. c
687. 2000 CONTINUE
688. C FIN DE BOUCLE SUR LES ELEMENTS &lt;&lt;&lt;&lt;&lt;&lt;&lt;&lt;&lt;&lt;&lt;&lt;&lt;&lt;&lt;&lt;&lt;&lt;&lt;&lt;&lt;
689. c
690. c
691.
692. C*********************************************************
693. C SUPPRESSION ET DESACTIVATION DE SEGMENTS
694. C*********************************************************
695. C
696. SEGSUP,WRK1,WRK2,MVELCH
697. SEGSUP,MRACC
698.
699. segdes,MELEME
700. segdes,MINT1
701. if(MINT2.ne.0) segdes,MINT2
702. c desactivation des maillages et segment imatri
703. do ienr=1,(1+NBENR1)
704. IPT1 = LOCIRI(1,ienr)
705. if(IPT1.ne.0) segdes,IPT1
706. xMATR1 = LOCIRI(4,ienr)
707. if(xMATR1.ne.0) segdes,xMATR1
708. enddo
709.
710. c SEGDES dans RIGI1 = IMODEL
711. C
712. RETURN
713. END
714.
715.
716.
717.

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