Télécharger bsigmx.eso

Retour à la liste

Numérotation des lignes :

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

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