Télécharger bsigm3.eso

Retour à la liste

Numérotation des lignes :

  1. C BSIGM3 SOURCE BP208322 15/06/22 21:15:33 8543
  2. SUBROUTINE BSIGM3(IPMAIL,LRE,NSTRS,LW,IVACAR,NCARR,IVECT,
  3. & MELE,CMATE,IVASTR,ISOUS,NBPGAU,NBPTEL,IPMINT,NFOR,IVAFOR
  4. & ,ADPG,BDPG,CDPG,IIPDPG,IVAMAT,NMATT,MFR,dcmate)
  5. *----------------------------------------------------------------------
  6. * _______________________________ *
  7. * | | *
  8. * | calcul des forces aux noeuds| *
  9. * |______________________________| *
  10. * *
  11. * poutre,tuyau,linespring,tuyau fissure,barre,cerce,tuyo *
  12. * poutre de timoschenko,shb8,joi1,zone_cohesive
  13. * *
  14. *---------------------------------------------------------------------*
  15. * *
  16. * entrees : *
  17. * ________ *
  18. * *
  19. * ipmail pointeur sur un segment meleme *
  20. * lre nombre de ddl dans la matrice de rigidite *
  21. * nstrs nombre de composante de contraintes/deformations *
  22. * lw dimension du tableau de travail de l'element *
  23. * ivacar pointeur sur les chamelems de caracteristiques géo. *
  24. * ncarr nombre de caracteristiques geometriques *
  25. * ivamat pointeur sur les chamelems de caracteristiques mat. *
  26. * nmatt nombre de caracteristiques matériau *
  27. * ivect flag indiquant si on a entree les axes locaux *
  28. * mele numero de l'element fini *
  29. * ivastr pointeur sur un segment mptval contenant les *
  30. * les melvals de contraints *
  31. * isous numero de la sous-zone *
  32. * nbpgau nombre de points d'integration pour les contraintes *
  33. * nbptel nombre de points par element *
  34. * ipmint pointeur sur un segment minte *
  35. * nfor nombre de composantes de forces *
  36. * *
  37. * sorties : *
  38. * ________ *
  39. * *
  40. * ivafor pointeur sur un segment mptval contenant les *
  41. * les melvals de forces *
  42.  
  43. * *
  44. *---------------------------------------------------------------------*
  45. IMPLICIT INTEGER(I-N)
  46. IMPLICIT REAL*8(A-H,O-Z)
  47. *
  48. *
  49. -INC CCOPTIO
  50. -INC CCHAMP
  51. -INC CCREEL
  52. -INC SMCHAML
  53. -INC SMCHPOI
  54. -INC SMELEME
  55. -INC SMCOORD
  56. -INC SMMODEL
  57. -INC SMINTE
  58. -INC SMLMOTS
  59. c
  60.  
  61. SEGMENT WRK1
  62. REAL*8 XFORC(LRE), XSTRS(NSTRS), XE(3,NBBB)
  63. ENDSEGMENT
  64. *
  65. SEGMENT WRK2
  66. REAL*8 SHPWRK(6,NBNN), BGENE(NSTRS,LRE)
  67. ENDSEGMENT
  68. *
  69. SEGMENT WRK3
  70. REAL*8 WORK(LW)
  71. ENDSEGMENT
  72. *
  73. SEGMENT WRK4
  74. REAL*8 BPSS(3,3), XEL(3,NBBB), XFOLO(LRE)
  75. ENDSEGMENT
  76. *
  77. SEGMENT WRK5
  78. REAL*8 XGENE(NSTN,LRN)
  79. ENDSEGMENT
  80. * segment pour shb8
  81. SEGMENT WRK7
  82. REAL*8 PROPEL(24),D(3)
  83. REAL*8 rel(24,24),work1(30)
  84. ENDSEGMENT
  85. *
  86. SEGMENT MPTVAL
  87. INTEGER IPOS(NS) ,NSOF(NS)
  88. INTEGER IVAL(NCOSOU)
  89. CHARACTER*16 TYVAL(NCOSOU)
  90. ENDSEGMENT
  91. POINTEUR MPTVA1.MPTVAL,MPTVA2.MPTVAL
  92. *
  93. CHARACTER*4 lesinc(7),lesdua(7)
  94. DATA lesinc/'UX','UY','UZ','RX','RY','RZ','UR'/
  95. DATA lesdua/'FX','FY','FZ','MX','MY','MZ','FR'/
  96. CHARACTER*8 CMATE
  97. logical dcmat2,dcmate
  98.  
  99. dcmat2 = .false.
  100. MELEME=IPMAIL
  101. SEGACT MELEME
  102. NBNN=NUM(/1)
  103. NBELEM=NUM(/2)
  104. *
  105. * introduction des coord du point autour duquel se fait le
  106. * mouvement de la section en defo plane generalisee
  107. * et initialisation des forces au noeud support de la defo
  108. * plane generalisee
  109. *
  110. ADPG=0.D0
  111. BDPG=0.D0
  112. CDPG=0.D0
  113. IF (IFOUR.EQ.-3)THEN
  114. IP=IIPDPG
  115. SEGACT MCOORD
  116. IREF=(IP-1)*(IDIM+1)
  117. XDPGE=XCOOR(IREF+1)
  118. YDPGE=XCOOR(IREF+2)
  119. ELSE
  120. XDPGE=0.D0
  121. YDPGE=0.D0
  122. ENDIF
  123. *
  124. NHRM=NIFOUR
  125. MINTE=IPMINT
  126. *
  127. c_______________________________________________________________________
  128. c_______________________________________________________________________
  129. c
  130. c numero des etiquettes :
  131. c etiquettes de 1 a 98 pour traitement specifique a l element
  132. c dans la zone specifique a chaque element commencant par :
  133. c 5 continue
  134. c element 5 etiquettes 1005 2005 3005 4005 ...
  135. c 44 continue
  136. c element 44 etiquettes 1044 2044 3044 4044 ...
  137. c_______________________________________________________________________
  138. c
  139. IF (MELE.LE.100)
  140. &GOTO (99,2,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,
  141. 1 99,99,99,99,99,99,99,99,29,30,99,99,99,99,99,99,99,99,99,99,
  142. 2 99,29,43,99,45,46,99,99,99,30,99,99,99,99,99,99,99,99,99,99,
  143. 3 99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,
  144. 4 99,99,99,29,99,99,99,99,99,99,99,99,99,99,46,96,99,99,99,99
  145. 5 ),MELE
  146. IF (MELE.LE.200)
  147. &GOTO (99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,
  148. 1 99,99,46,124,125,34,34,34,34,34,34,34,34,34,34,34,34,34,34,
  149. 2 34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,
  150. 3 34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,
  151. 4 34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,
  152. 5 34),MELE-100
  153. IF (MELE.LE.300)
  154. &GOTO (34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,
  155. 1 34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,
  156. 2 34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,34,
  157. 3 260,34,34,34,34,265,266,266,266,99,99,271,272),MELE-200
  158. c
  159. 34 CONTINUE
  160. GOTO 99
  161. C____________________________________________________________________
  162. C
  163. C ELEMENT SEG2 (pour IMPEDANCE)
  164. C____________________________________________________________________
  165. C
  166. 2 CONTINUE
  167. IF (CMATE.EQ.'IMPELAST'.OR.CMATE.EQ.'IMPVOIGT'.OR.
  168. & CMATE.EQ.'IMPREUSS'.OR .cmate.eq.'IMPCOMPL') THEN
  169. * detection impedance "sur mesure"
  170. MPTVA1=IVASTR
  171. MPTVAL=IVAFOR
  172. if (ival(/1).eq.mptva1.ival(/1)*2) dcmat2 = .true.
  173. ENDIF
  174.  
  175. DO 310 IB=1,NBELEM
  176. C ON CHERCHE LES CONTRAINTES -
  177. C
  178. MPTVA1=IVASTR
  179. numstr = mptva1.ival(/1)
  180. C RANGEMENT DANS MELVAL
  181. C
  182. MPTVAL=IVAFOR
  183. DO 910 IGAU=1,NBNN
  184. DO 910 ICOMP=1,NFOR
  185. if (dcmat2) then
  186. if(icomp.le.numstr) then
  187. melva1 = mptva1.ival(icomp)
  188. else
  189. melva1 = mptva1.ival(icomp - numstr)
  190. endif
  191. if (igau.lt.2) then
  192. melval = ival(icomp)
  193. else
  194. melval = ival(icomp + nfor)
  195. endif
  196. else
  197. MELVA1=MPTVA1.IVAL(ICOMP)
  198. MELVAL=IVAL(ICOMP)
  199. endif
  200. IBMN=MIN(IB ,VELCHE(/2))
  201. IBM1=MIN(IB ,MELVA1.VELCHE(/2))
  202. IGM1=MIN(IGAU,MELVA1.VELCHE(/1))
  203. VELCHE(IGAU,IBMN)= MELVA1.VELCHE(IGM1,IBM1)
  204. 910 CONTINUE
  205. 310 CONTINUE
  206. GO TO 510
  207. c_______________________________________________________________________
  208. c_______________________________________________________________________
  209. c
  210. c elements poutre et poutre de timoschenko
  211. c_______________________________________________________________________
  212. c
  213. 29 CONTINUE
  214. if (dcmate) goto 2
  215. NBBB=NBNN
  216. SEGINI WRK1,WRK3
  217. c
  218. NCARR1=NCARR
  219. IF((IVECT.EQ.1).AND.(IFOUR.NE.-2)) NCARR1=NCARR-1
  220. c
  221. DO 3029 IB=1,NBELEM
  222. c
  223. c on cherche les coordonnees des noeuds de l elementib
  224. c
  225. CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  226. c
  227. c il faudrait aussi modifier le vecteur local de la poutre
  228. c
  229. c rangement des caracteristiques dans work
  230. C
  231. MPTVAL=IVACAR
  232. C
  233. DO 6029 IC=1,NCARR1
  234. WORK(IC)=0.D0
  235. MELVAL=IVAL(IC)
  236. IF (MELVAL.NE.0) THEN
  237. IBMN=MIN(IB,VELCHE(/2))
  238. DO 4029 IGAU=1,NBNN
  239. IGMN=MIN(IGAU,VELCHE(/1))
  240. WORK(IC)=WORK(IC)+VELCHE(IGMN,IBMN)
  241. 4029 CONTINUE
  242. WORK(IC)=WORK(IC)/NBNN
  243. ENDIF
  244. 6029 CONTINUE
  245. c
  246. c cas ou on a lu le mot vecteur
  247. C
  248. IF (IFOUR.NE.-2) THEN
  249. C
  250. IF (IVECT.EQ.1) THEN
  251. MELVAL=IVAL(NCARR)
  252. IF (MELVAL.NE.0) THEN
  253. IBMN=MIN(IB,IELCHE(/2))
  254. IP=IELCHE(1,IBMN)
  255. IREF=(IP-1)*(IDIM+1)
  256. DO 6129 IC=1,IDIM
  257. WORK(NCARR+IC-1)=XCOOR(IREF+IC)
  258. 6129 CONTINUE
  259. ELSE
  260. DO 6229 IC=1,IDIM
  261. WORK(NCARR+IC-1)=0.D0
  262. 6229 CONTINUE
  263. ENDIF
  264. ENDIF
  265. C
  266. ENDIF
  267. C
  268. c cas des tuyaux - on calcule les caracteristiques de la poutre
  269. c equivalente
  270. c
  271. IF(MELE.EQ.42) THEN
  272. CISA=WORK(4)
  273. VX =WORK(5)
  274. VY =WORK(6)
  275. VZ =WORK(7)
  276. CALL TUYCAR(WORK,CISA,VX,VY,VZ,KERRE,0)
  277. IF (KERRE.EQ.77) THEN
  278. CALL ERREUR(77)
  279. GOTO 510
  280. ENDIF
  281. ENDIF
  282. c
  283. c on cherche les contraintes -
  284. c
  285. MPTVAL=IVASTR
  286. IE=9
  287. DO 7029 IGAU=1,2
  288. DO 7029 ICOMP=1,NSTRS
  289. IE=IE+1
  290. MELVAL=IVAL(ICOMP)
  291. IBMN=MIN(IB ,VELCHE(/2))
  292. IGMN=MIN(IGAU,VELCHE(/1))
  293. WORK(IE)=VELCHE(IGMN,IBMN)
  294. 7029 CONTINUE
  295. c
  296. c on calcule les forces internes
  297. c
  298. IF(MELE.EQ.84) THEN
  299. IF(CMATE.NE.'SECTION') THEN
  300. IF (IFOUR.EQ.-2.OR.IFOUR.EQ.-1.OR.IFOUR.EQ.-3) THEN
  301. CALL TIMBS2(XFORC,XE,WORK(10),WORK(22),KERRE)
  302. ELSE
  303. CALL TIMBSG(XFORC,WORK(7),XE,WORK(10),WORK(22),KERRE)
  304. ENDIF
  305. ELSE
  306. IF (IFOUR.EQ.-2.OR.IFOUR.EQ.-1.OR.IFOUR.EQ.-3) THEN
  307. CALL TIMBS2(XFORC,XE,WORK(10),WORK(22),KERRE)
  308. ELSE
  309. CALL TIMBSG(XFORC,WORK(1),XE,WORK(10),WORK(22),KERRE)
  310. ENDIF
  311. ENDIF
  312. ELSE
  313. IF (IFOUR.EQ.-2.OR.IFOUR.EQ.-1.OR.IFOUR.EQ.-3) THEN
  314. CALL POUBS2(XFORC,WORK,XE,WORK(10),WORK(22),KERRE)
  315. ELSE
  316. CALL POUBSG(XFORC,WORK,XE,WORK(10),WORK(22),KERRE)
  317. ENDIF
  318. ENDIF
  319. IF(KERRE.NE.0) THEN
  320. INTERR(1)=ISOUS
  321. INTERR(2)=IB
  322. SEGSUP WRK1,WRK3
  323. CALL ERREUR(128)
  324. GO TO 510
  325. ENDIF
  326. c
  327. c rangement dans melval
  328. c
  329. IE=0
  330. MPTVAL=IVAFOR
  331. DO 9029 IGAU=1,NBNN
  332. DO 9029 ICOMP=1,NFOR
  333. IE=IE+1
  334. MELVAL=IVAL(ICOMP)
  335. IBMN=MIN(IB ,VELCHE(/2))
  336. VELCHE(IGAU,IBMN)=XFORC(IE)
  337. 9029 CONTINUE
  338. 3029 CONTINUE
  339. SEGSUP WRK1,WRK3
  340. GO TO 510
  341. c_______________________________________________________________________
  342. c
  343. c elements lisp et lism
  344. c_______________________________________________________________________
  345. c
  346. 30 CONTINUE
  347. NBBB=NBNN
  348. NBNO=SHPTOT(/2)
  349. SEGINI WRK1,WRK3,WRK4
  350. c
  351. DO 3030 IB=1,NBELEM
  352. c
  353. c on cherche les coordonnees des noeuds de l element ib
  354. c
  355. CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  356. c
  357. c
  358. c on cherche les contraintes
  359. c
  360. IE=0
  361. MPTVAL=IVASTR
  362. DO 7030 IGAU=1,NBPGAU
  363. DO 7030 ICOMP=1,NSTRS
  364. IE=IE+1
  365. MELVAL=IVAL(ICOMP)
  366. IGMN=MIN(IGAU,VELCHE(/1))
  367. IBMN=MIN(IB ,VELCHE(/2))
  368. WORK(IE)=VELCHE(IGMN,IBMN)
  369. 7030 CONTINUE
  370. c
  371. c on cherche les caracteristiques
  372. c
  373. MPTVAL=IVACAR
  374. DO 6030 IGAU=1,NBPGAU
  375. DO 6030 ICOMP=1,NCARR
  376. IE=IE+1
  377. MELVAL=IVAL(ICOMP)
  378. IF (MELVAL.NE.0) THEN
  379. IGMN=MIN(IGAU,VELCHE(/1))
  380. IBMN=MIN(IB ,VELCHE(/2))
  381. WORK(IE)=VELCHE(IGMN,IBMN)
  382. ELSE
  383. WORK(IE)=0.D0
  384. ENDIF
  385. 6030 CONTINUE
  386. c
  387. c on calcule b*sigma
  388. c
  389. ICNT=NBPGAU*NSTRS+1
  390. CALL LISPBS(WORK(1),WORK(ICNT),POIGAU,SHPTOT,
  391. 1 NBPGAU,NBNO,XE,XFOLO,BPSS,XFORC)
  392. c
  393. c rangement dans melval
  394. c
  395. IE=0
  396. MPTVAL=IVAFOR
  397. DO 9030 IGAU=1,NBNN
  398. DO 9030 ICOMP=1,6
  399. IE=IE+1
  400. MELVAL=IVAL(ICOMP)
  401. IBMN=MIN(IB ,VELCHE(/2))
  402. VELCHE(IGAU,IBMN)=XFORC(IE)
  403. 9030 CONTINUE
  404. 3030 CONTINUE
  405. SEGSUP WRK1,WRK3,WRK4
  406. GOTO 510
  407. c_______________________________________________________________________
  408. c
  409. c element tuyau fissure
  410. c_______________________________________________________________________
  411. c
  412. 43 CONTINUE
  413. NBBB=NBNN
  414. SEGINI WRK1,WRK3
  415. c
  416. DO 3043 IB=1,NBELEM
  417. c
  418. c on cherche les contraintes
  419. c
  420. IE=0
  421. MPTVAL=IVASTR
  422. DO 4043 IGAU=1,NBPTEL
  423. DO 4043 ICOMP=1,NSTRS
  424. IE=IE+1
  425. MELVAL=IVAL(ICOMP)
  426. IGMN=MIN(IGAU,VELCHE(/1))
  427. IBMN=MIN(IB ,VELCHE(/2))
  428. WORK(IE)=VELCHE(IGMN,IBMN)
  429. 4043 CONTINUE
  430. c
  431. c on cherche les caracteristiques
  432. c
  433. MPTVAL=IVACAR
  434. DO 5043 ICOMP=1,NCARR
  435. MELVAL=IVAL(ICOMP)
  436. IF (MELVAL.NE.0) THEN
  437. IBMN=MIN(IB ,VELCHE(/2))
  438. WORK(ICOMP+8)=VELCHE(1,IBMN)
  439. ELSE
  440. WORK(ICOMP+8)=0.D0
  441. ENDIF
  442. 5043 CONTINUE
  443. c
  444. c on calcule les forces internes
  445. c
  446. CALL TUFIBS(XFORC,WORK,WORK(9),WORK(18),KERRE)
  447. c
  448. c rangement dans melval
  449. c
  450. IE=0
  451. MPTVAL=IVAFOR
  452. DO 6043 IGAU=1,NBNN
  453. DO 6043 ICOMP=1,NFOR
  454. IE=IE+1
  455. MELVAL=IVAL(ICOMP)
  456. IBMN=MIN(IB ,VELCHE(/2))
  457. VELCHE(IGAU,IBMN)=XFORC(IE)
  458. 6043 CONTINUE
  459. 3043 CONTINUE
  460. SEGSUP WRK1,WRK3
  461. GO TO 510
  462. c_______________________________________________________________________
  463. c
  464. c element point (poi1) defos planes generalisees/materiau IMPEDANCE
  465. c_______________________________________________________________________
  466. c
  467. 45 CONTINUE
  468. NBBB=NBNN
  469. *
  470. *
  471. IF ((CMATE.EQ.'IMPELAST').OR.(CMATE.EQ.'IMPVOIGT').OR.
  472. & (CMATE.EQ.'IMPREUSS').OR .cmate.eq.'IMPCOMPL'.or.
  473. & cmate.eq.'MODAL'.or.cmate.eq.'STATIQUE') THEN
  474. mptva1 = ivastr
  475. mptval = ivafor
  476. segact mptva1,mptval
  477. do iv = 1,ival(/1)
  478. melva1 = mptva1.ival(iv)
  479. melval = ival(iv)
  480. segact melva1,melval*mod
  481. DO IB=1,NBELEM
  482.  
  483. DO 4145 IGAU=1,NBNN
  484. IGMN=MIN(IGAU,MELVA1.VELCHE(/1))
  485. IBMN=MIN(IB ,MELVA1.VELCHE(/2))
  486. valstr = MELVA1.VELCHE(IGMN,IBMN)
  487. VELCHE(IGMN,IBMN) = valstr
  488. 4145 CONTINUE
  489. ENDDO
  490.  
  491. enddo
  492.  
  493. GOTO 510
  494. ENDIF
  495. *
  496. IF(MELE.EQ.45.AND.IFOUR.NE.-3) THEN
  497. GO TO 99
  498. ENDIF
  499. *
  500. SEGINI WRK1,WRK3
  501. c
  502. DO 3045 IB=1,NBELEM
  503. c
  504. c on cherche les coordonnees des noeuds de l element ib
  505. c
  506. CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  507. c
  508. c mise a zero des forces internes
  509. c
  510. CALL ZERO(XFORC,1,LRE)
  511. c
  512. c on cherche l'effort
  513. c
  514. MPTVAL=IVASTR
  515. MELVAL=IVAL(1)
  516. IBMN=MIN(IB ,VELCHE(/2))
  517. EFFORT=VELCHE(1,IBMN)
  518. c
  519. c on calcule les forces internes
  520. c
  521. CALL PO1BSG(XE,LRE,XDPGE,YDPGE,EFFORT,XFORC)
  522. ADPG=ADPG+XFORC(3)
  523. BDPG=BDPG+XFORC(4)
  524. CDPG=CDPG+XFORC(5)
  525. c
  526. c rangement dans melval
  527. c
  528. IE=0
  529. MPTVAL=IVAFOR
  530. DO 9045 IGAU=1,NBNN
  531. DO 9045 ICOMP=1,NFOR-3
  532. IE=IE+1
  533. MELVAL=IVAL(ICOMP)
  534. IBMN=MIN(IB ,VELCHE(/2))
  535. VELCHE(IGAU,IBMN)=XFORC(IE)
  536. 9045 CONTINUE
  537. 3045 CONTINUE
  538. c
  539. SEGSUP WRK1,WRK3
  540. GO TO 510
  541. c_______________________________________________________________________
  542. c
  543. c elements barre et cerce
  544. c_______________________________________________________________________
  545. c
  546. 46 CONTINUE
  547. NBBB=NBNN
  548. *
  549. IF(MELE.EQ.95.AND.IFOUR.NE.0.AND.IFOUR.NE.1) THEN
  550. GO TO 99
  551. ENDIF
  552. *
  553. SEGINI WRK1,WRK3
  554. c
  555. DO 3046 IB=1,NBELEM
  556. c
  557. c on cherche les coordonnees des noeuds de l element ib
  558. c
  559. CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  560. C
  561. IF(MELE.EQ.123) THEN
  562. c
  563. c on cherche l'effort
  564. c
  565. IE=0
  566. MPTVAL=IVASTR
  567. DO 4046 IGAU=1,NBPTEL
  568. DO 4046 ICOMP=1,NSTRS
  569. IE=IE+1
  570. MELVAL=IVAL(ICOMP)
  571. IGMN=MIN(IGAU,VELCHE(/1))
  572. IBMN=MIN(IB ,VELCHE(/2))
  573. WORK(IE)=VELCHE(IGMN,IBMN)
  574. 4046 CONTINUE
  575. c
  576. c on calcule les forces internes
  577. c
  578. CALL BARBS3(XFORC,XE,WORK,KERRE,QSIGAU,POIGAU,NBPGAU,IB)
  579. c
  580. ELSE
  581. c
  582. c on cherche l'effort
  583. c
  584. MPTVAL=IVASTR
  585. MELVAL=IVAL(1)
  586. NPPTEL=VELCHE(/1)
  587. IBMN=MIN(IB ,VELCHE(/2))
  588. IF(NPPTEL.EQ.1) THEN
  589. EFFORT=VELCHE(1,IBMN)
  590. ELSE IF(NPPTEL.EQ.2) THEN
  591. EFFOR1=VELCHE(1,IBMN)
  592. EFFOR2=VELCHE(2,IBMN)
  593. EFFORT=0.5D0*(EFFOR1+EFFOR2)
  594. ENDIF
  595. c
  596. c on calcule les forces internes
  597. c
  598. IF(MELE.EQ.46) CALL BARBSG(XFORC,XE,EFFORT,KERRE)
  599. IF(MELE.EQ.95) CALL CERBSG(XFORC,XE,EFFORT,KERRE)
  600. ENDIF
  601. C
  602. IF(KERRE.NE.0) THEN
  603. INTERR(1)=ISOUS
  604. INTERR(2)=IB
  605. SEGSUP WRK1,WRK3
  606. IF(MELE.EQ.46) CALL ERREUR(128)
  607. IF(MELE.EQ.95) CALL ERREUR(128)
  608. GO TO 510
  609. ENDIF
  610. c
  611. c rangement dans melval
  612. c
  613. NFOD=NFOR
  614. IF (IFOUR.EQ.-3.AND.MELE.EQ.46) NFOD=NFOR-3
  615. IE=0
  616. MPTVAL=IVAFOR
  617. DO 9046 IGAU=1,NBNN
  618. DO 9046 ICOMP=1,NFOD
  619. IE=IE+1
  620. MELVAL=IVAL(ICOMP)
  621. IBMN=MIN(IB ,VELCHE(/2))
  622. VELCHE(IGAU,IBMN)=XFORC(IE)
  623. 9046 CONTINUE
  624. 3046 CONTINUE
  625. SEGSUP WRK1,WRK3
  626. GO TO 510
  627. C_______________________________________________________________________
  628. C
  629. C ELEMENT BARRE 3D EXCENTRE (BAEX)
  630. C_______________________________________________________________________
  631. C
  632. 124 CONTINUE
  633. NBBB=NBNN
  634. NSTN=NBNN
  635. LRN =LRE
  636. SEGINI WRK1,WRK3,WRK5
  637. C
  638. C BOUCLE DE CALCUL POUR LES DIFFERENTS ELEMENTS
  639. C
  640. KERRE=0
  641. DO 3124 IB=1,NBELEM
  642. C
  643. C ON RECUPERE LA SECTION DE L'ELEMENT, SES EXCENTREMENTS ET SON
  644. C ORIENTATION. LES CARACTERISTIQUES SONT RANGEES DANS WORK
  645. C SELON L'ORDRE SUIVANT: SECT EXCZ EXCY VX VY VZ
  646. C
  647. MPTVAL=IVACAR
  648. DO IC=1,NCARR
  649. IF(IVAL(IC).NE.0) THEN
  650. MELVAL=IVAL(IC)
  651. IBMN=MIN(IB,VELCHE(/2))
  652. WORK(IC)=VELCHE(1,IBMN)
  653. ELSE
  654. WORK(IC)=0.D0
  655. ENDIF
  656. END DO
  657. SECT=WORK(1)
  658. C
  659. C XGENE STOCKE LA MATRICE DE PASSAGE DE L'ELEMENT EXCENTRE
  660. C
  661. CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  662. CALL MAPAEX(XE,NBNN,WORK,AL,XGENE,LRE,KERRE)
  663. IF(KERRE.NE.0) THEN
  664. INTERR(1)=ISOUS
  665. INTERR(2)=IB
  666. IF(KERRE.EQ.1) CALL ERREUR(128)
  667. ENDIF
  668. C
  669. C MISE A ZERO DES FORCES INTERNES
  670. C
  671. CALL ZERO(XFORC,1,LRE)
  672. C
  673. C ON CHERCHE L'EFFORT
  674. C
  675. MPTVAL=IVASTR
  676. MELVAL=IVAL(1)
  677. IBMN=MIN(IB ,VELCHE(/2))
  678. NPPTEL=VELCHE(/1)
  679. IF(NPPTEL.EQ.1) THEN
  680. EFFORT=VELCHE(1,IBMN)
  681. ELSE IF(NPPTEL.EQ.2) THEN
  682. EFFOR1=VELCHE(1,IBMN)
  683. EFFOR2=VELCHE(2,IBMN)
  684. EFFORT=0.5D0*(EFFOR1+EFFOR2)
  685. ENDIF
  686. CC EFFORT=SECT*EFFORT
  687. C
  688. C ON CALCULE LES FORCES INTERNES
  689. C
  690. CALL BARINT(XFORC,XGENE,EFFORT,LRE)
  691. C
  692. C RANGEMENT DANS MELVAL
  693. C
  694. IE=0
  695. MPTVAL=IVAFOR
  696. DO 9199 IGAU=1,NBNN
  697. DO 9199 ICOMP=1,NFOR
  698. IE=IE+1
  699. MELVAL=IVAL(ICOMP)
  700. IBMN=MIN(IB ,VELCHE(/2))
  701. VELCHE(IGAU,IBMN)=XFORC(IE)
  702. 9199 CONTINUE
  703. 3124 CONTINUE
  704. SEGSUP WRK1,WRK3,WRK5
  705. GO TO 510
  706. c cccccccccccccccccccccccccccccccccccccccccccccccccccccc c
  707. c element coaxial COS2 (3D pour liaison acier-beton) c
  708. c cccccccccccccccccccccccccccccccccccccccccccccccccccccc c
  709. 271 continue
  710. NBBB=NBNN
  711. NSTN=NBNN
  712. LRN =LRE
  713. SEGINI WRK1,WRK3,WRK4,WRK5
  714. do 2713 ib=1,nbelem
  715. C
  716. CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  717. vx1= xe(1,2) - xe(1,1)
  718. vy1= xe(2,2) - xe(2,1)
  719. vz1=0.d0
  720. if(idim.eq.3) vz1= xe(3,2)-xe(3,1)
  721. xl= sqrt( vx1*vx1 + vy1*vy1 + vz1*vz1)
  722. vx1= vx1 / xl
  723. vy1= vy1/ xl
  724. if(idim.eq.3) THEN
  725. vz1=vz1 / xl
  726. ENDIF
  727.  
  728. IE=0
  729. MPTVAL=IVASTR
  730. DO 2711 IGAU=1,2
  731. DO 2711 ICOMP=1,NSTRS
  732. IE=IE+1
  733. MELVAL=IVAL(ICOMP)
  734. IGMN=MIN(IGAU,VELCHE(/1))
  735. IBMN=MIN(IB ,VELCHE(/2))
  736. WORK(IE)=VELCHE(IGMN,IBMN)
  737. 2711 CONTINUE
  738. C
  739. MPTVAL=IVAcar
  740. DO 2712 ICOMP=1,NCARR
  741. MELVAL=IVAL(ICOMP)
  742. IGMN = VELCHE(/1)
  743. IBMN=MIN(IB,VELCHE(/2))
  744. SECA =VELCHE(IGMN,IBMN)
  745. 2712 CONTINUE
  746. diam = sqrt(4.d0*SECA/xpi)
  747. C
  748. C MISE A ZERO DES FORCES INTERNES
  749. C
  750. CALL ZERO(XFORC,1,LRE)
  751. C
  752. C
  753. C FORCES DANS LA DIRECTION TANGENTIELLE
  754. C
  755. ag = 1.d0-0.5773502691896257645d0
  756. agg = ag/(2.d0*(1-ag))
  757. t11 = work(1) + (work(1) -work(1 + idim))*agg
  758. t21 = work(1 + idim) + (work(1 + idim) -work(1))*agg
  759. FO11 = xpi*diam*xl*(t11/2.d0 + (t21 - t11)/8.d0)
  760. FO21 = xpi*diam*xl*(t21/2.d0 + (t11 - t21)/8.d0)
  761. C
  762. C FORCES DANS LA DIRECTION NORMALE
  763. C
  764. t12 = WORK(2) + (work(2) -work(2 + idim))*agg
  765. t22 = WORK(2 + IDIM)+ (work(2 + idim) -work(2))*agg
  766. FO12 = -1.d0*diam*xl*(t12/2.d0 + (t22 - t12)/8.d0)
  767. FO22 = -1.d0*diam*xl*(t22/2.d0 + (t12 - t22)/8.d0)
  768. c write(6,*) 'xstrs 1',ib,t11,t12,t12
  769. c write(6,*) 'xstrs 2',ib,t21,t22,t22
  770. IF (IDIM.EQ.2) THEN
  771. XFORC(1) = (FO11*VX1 + FO12*VY1) + XFORC(1)
  772. XFORC(3) = (FO21*VX1 + FO22*VY1) + XFORC(3)
  773. XFORC(2) = (FO11*VY1 + FO12*VX1) + XFORC(2)
  774. XFORC(4) = (FO21*VY1 + FO22*VX1) + XFORC(4)
  775. XFORC(5) = (-1.d0*(FO11*VX1 + FO12*VY1)) + XFORC(5)
  776. XFORC(7) = (-1.d0*(FO21*VX1 + FO22*VY1)) + XFORC(7)
  777. XFORC(6) = (-1.d0*(FO11*VY1 + FO12*VX1)) + XFORC(6)
  778. XFORC(8) = (-1.d0*(FO21*VY1 + FO22*VX1)) + XFORC(8)
  779. ELSE IF (IDIM.EQ.3) THEN
  780. IF (vy1.EQ.0.0D0.AND.vz1.EQ.0.0D0) THEN
  781. vx22 = 0.0D0
  782. vy22 = 1.0D0
  783. vz22 = 0.0D0
  784. ELSE IF ((vx1.EQ.0.0D0.AND.vy1.EQ.0.0D0).OR.
  785. . (vx1.EQ.0.0D0.AND.vz1.EQ.0.0D0)) THEN
  786. vx22 = 1.0D0
  787. vy22 = 0.0D0
  788. vz22 = 0.0D0
  789. ELSE IF (vy1.NE.0.0D0.AND.vz1.NE.0.0D0) THEN
  790. Vx22 = 0.0D0
  791. Vy22 = -vz1
  792. Vz22 = vy1
  793. LLL = 1
  794. ELSE IF (vx1.NE.0.0D0.AND.vz1.NE.0.0D0) THEN
  795. Vx22 = -vz1
  796. Vy22 = 0.0D0
  797. Vz22 = vx1
  798. LLL = 1
  799. ELSE IF (vy1.NE.0.0D0.AND.vx1.NE.0.0D0) THEN
  800. Vx22 = -vy1
  801. Vy22 = vx1
  802. Vz22 = 0.0D0
  803. LLL = 1
  804. END IF
  805. xl22 = sqrt((vx22*vx22) + (vy22*vy22)+ (vz22*vz22))
  806. vx22 = vx22/xl22
  807. vy22 = vy22/xl22
  808. vz22 = vz22/xl22
  809. vx3 = (vy1*Vz22) - (vz1*Vy22)
  810. vy3 = (vz1*Vx22) - (vx1*vz22)
  811. vz3 = (vx1*vy22) - (vy1*Vx22)
  812. xl3 = sqrt((vx3*vx3) + (vy3*vy3)+ (vz3*vz3))
  813. vx3 = vx3/xl3
  814. vy3 = vy3/xl3
  815. vz3 = vz3/xl3
  816. vx2 = (vy3*vz1) - (vz3*vy1)
  817. vy2 = (vz3*vx1) - (vx3*vz1)
  818. vz2 = (vx3*vy1) - (vy3*vx1)
  819. xl2 = sqrt((vx2*vx2) + (vy2*vy2)+ (vz2*vz2))
  820. vx2 = vx2/xl2
  821. vy2 = vy2/xl2
  822. vz2 = vz2/xl2
  823. F1 = FO11*vx1 + FO12*vx2 + FO12*vx3
  824. F2 = FO11*vy1 + FO12*vy2 + FO12*vy3
  825. F3 = FO11*vz1 + FO12*vz2 + FO12*vz3
  826. F4 = FO21*vx1 + FO22*vx2 + FO22*vx3
  827. F5 = FO21*vy1 + FO22*vy2 + FO22*vy3
  828. F6 = FO21*vz1 + FO22*vz2 + FO22*vz3
  829. XFORC(1) = F1 + XFORC(1)
  830. XFORC(2) = F2 + XFORC(2)
  831. XFORC(3) = F3 + XFORC(3)
  832. XFORC(4) = F4 + XFORC(4)
  833. XFORC(5) = F5 + XFORC(5)
  834. XFORC(6) = F6 + XFORC(6)
  835. XFORC(7) = -1.d0*F4 + XFORC(7)
  836. XFORC(8) = -1.d0*F5 + XFORC(8)
  837. XFORC(9) = -1.d0*F6 + XFORC(9)
  838. XFORC(10) = -1.d0*F1 + XFORC(10)
  839. XFORC(11) = -1.d0*F2 + XFORC(11)
  840. XFORC(12) = -1.d0*F3 + XFORC(12)
  841. ENDIF
  842. c write(6,*) 'xforc cos2', (xforc(io),io = 1,lre)
  843. C
  844. C RANGEMENT DANS MELVAL
  845. C
  846. IE=0
  847. MPTVAL=IVAFOR
  848. C
  849. C NODE=4= NOMBRE DE NOEUDS
  850. C ICOMP=NSTRS= NOMBRE DE COMPOSANTES DE CONTRAINTES PAR NOEUD
  851. C
  852. DO 2714 NODE=1,4
  853. DO 2714 ICOMP=1,NSTRS
  854. IE=IE+1
  855. MELVAL=IVAL(ICOMP)
  856. IBMN=MIN(IB ,VELCHE(/2))
  857. VELCHE(NODE,IBMN)=XFORC(IE)
  858. 2714 CONTINUE
  859.  
  860. 2713 continue
  861.  
  862. SEGSUP WRK1,WRK3,WRK4,WRK5
  863. GO TO 510
  864.  
  865. C_______________________________________________________________________
  866. C
  867. C ELEMENT COAXIAL (COA2)
  868. C_______________________________________________________________________
  869. C
  870. 272 continue
  871. NBNO=NBNN
  872. NBBB=NBNN
  873. SEGINI WRK1,WRK2,WRK3,WRK4
  874. LW = 24
  875. C
  876. DO 2721 IB=1,NBELEM
  877. C
  878. C ON CHERCHE LES COORDONNEES DES NOEUDS DE L ELEMENT IB
  879. C
  880. CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  881. C
  882. C MISE A ZERO DES FORCES INTERNES
  883. C
  884. CALL ZERO(XFORC,1,LRE)
  885. C
  886. C REPERE LOCAL
  887. C
  888. CALL CO2LOC(XE,SHPTOT,NBNN,XEL,BPSS,NOQUAL,IDIM)
  889. C
  890. C BOUCLE SUR LES POINTS DE GAUSS
  891. C
  892. DO 2722 IGAU=1,NBPGAU
  893. C
  894. C CALCUL DE LA MATRICE B ET DU JACOBIEN EN IGAU
  895. C
  896. CALL BCO2(IGAU,MFR,IFOUR,NIFOUR,XEL,BPSS,SHPTOT,SHPWRK,
  897. . BGENE,DJAC,IRRT,IDIM)
  898. IF(IRRT.NE.0) THEN
  899. INTERR(1)=IB
  900. CALL ERREUR(764)
  901. GOTO 9985
  902. END IF
  903. C
  904. C ON CHERCHE LES CONTRAINTES -
  905. C
  906. MPTVAL=IVASTR
  907. DO 2723 ICOMP=1,NSTRS
  908. MELVAL=IVAL(ICOMP)
  909. IGMN=MIN(IGAU,VELCHE(/1))
  910. IBMN=MIN(IB ,VELCHE(/2))
  911. XSTRS(ICOMP)=VELCHE(IGMN,IBMN)
  912. 2723 CONTINUE
  913. C
  914. C
  915. C ON CHERCHE LES COORDONNEES DES NOEUDS DE L ELEMENT IB
  916. C
  917. CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  918. vx1= xe(1,2) - xe(1,1)
  919. vy1= xe(2,2) - xe(2,1)
  920. vz1=0.d0
  921. if(idim.eq.3) vz1= xe(3,2)-xe(3,1)
  922. xl= sqrt( vx1*vx1 + vy1*vy1 + vz1*vz1)
  923. vx1= vx1 / xl
  924. vy1= vy1/ xl
  925. if(idim.eq.3) THEN
  926. vz1=vz1 / xl
  927. ENDIF
  928. C
  929. C recuperation de la section et calcul du diamètre
  930. C
  931. MPTVAL=IVACAR
  932. DO 2729 ICOMP=1,NCARR
  933. MELVAL=IVAL(ICOMP)
  934. IGMN = VELCHE(/1)
  935. IBMN=MIN(IB,VELCHE(/2))
  936. SECA =VELCHE(IGMN,IBMN)
  937. 2729 CONTINUE
  938. diam = sqrt(4.d0*SECA/xpi)
  939. C
  940. C ON CALCULE B*EFFORTS
  941. C
  942. DJAC=DJAC*POIGAU(IGAU)
  943. DO i=1,LRE
  944. cc = 0.D0
  945. DO j=1,NSTRS
  946. cc = cc + (BGENE(j,i) * XSTRS(j))
  947. ENDDO
  948. WORK(i) = xl * cc * DJAC * diam
  949. ENDDO
  950.  
  951. IF (IDIM.EQ.2) THEN
  952. WORK(1) = WORK(1) * xpi
  953. WORK(3) = WORK(3) * xpi
  954. WORK(5) = WORK(5) * xpi
  955. WORK(7) = WORK(7) * xpi
  956. WORK(1+4*idim) = WORK(1)*vx1 - WORK(2)*vy1
  957. WORK(2+4*idim) = WORK(1)*vy1 + WORK(2)*vx1
  958. WORK(3+4*idim) = WORK(3)*vx1 - WORK(4)*vy1
  959. WORK(4+4*idim) = WORK(3)*vy1 + WORK(4)*vx1
  960. WORK(5+4*idim) = WORK(5)*vx1 - WORK(6)*vy1
  961. WORK(6+4*idim) = WORK(5)*vy1 + WORK(6)*vx1
  962. WORK(7+4*idim) = WORK(7)*vx1 - WORK(8)*vy1
  963. WORK(8+4*idim) = WORK(7)*vy1 + WORK(8)*vx1
  964. ELSE IF (IDIM.EQ.3) THEN
  965. WORK(1) = WORK(1) * xpi
  966. WORK(4) = WORK(4) * xpi
  967. WORK(7) = WORK(7) * xpi
  968. WORK(10) = WORK(10) * xpi
  969. IF (vy1.EQ.0.0D0.AND.vz1.EQ.0.0D0) THEN
  970. vx22 = 0.0D0
  971. vy22 = 1.0D0
  972. vz22 = 0.0D0
  973. ELSE IF ((vx1.EQ.0.0D0.AND.vy1.EQ.0.0D0).OR.
  974. . (vx1.EQ.0.0D0.AND.vz1.EQ.0.0D0)) THEN
  975. vx22 = 1.0D0
  976. vy22 = 0.0D0
  977. vz22 = 0.0D0
  978. ELSE IF (vy1.NE.0.0D0.AND.vz1.NE.0.0D0) THEN
  979. Vx22 = 0.0D0
  980. Vy22 = -vz1
  981. Vz22 = vy1
  982. LLL = 1
  983. ELSE IF (vx1.NE.0.0D0.AND.vz1.NE.0.0D0) THEN
  984. Vx22 = -vz1
  985. Vy22 = 0.0D0
  986. Vz22 = vx1
  987. LLL = 1
  988. ELSE IF (vy1.NE.0.0D0.AND.vx1.NE.0.0D0) THEN
  989. Vx22 = -vy1
  990. Vy22 = vx1
  991. Vz22 = 0.0D0
  992. LLL = 1
  993. END IF
  994. xl22 = sqrt((vx22*vx22) + (vy22*vy22)+ (vz22*vz22))
  995. vx22 = vx22/xl22
  996. vy22 = vy22/xl22
  997. vz22 = vz22/xl22
  998. vx3 = (vy1*Vz22) - (vz1*Vy22)
  999. vy3 = (vz1*Vx22) - (vx1*vz22)
  1000. vz3 = (vx1*vy22) - (vy1*Vx22)
  1001. xl3 = sqrt((vx3*vx3) + (vy3*vy3)+ (vz3*vz3))
  1002. vx3 = vx3/xl3
  1003. vy3 = vy3/xl3
  1004. vz3 = vz3/xl3
  1005. vx2 = (vy3*vz1) - (vz3*vy1)
  1006. vy2 = (vz3*vx1) - (vx3*vz1)
  1007. vz2 = (vx3*vy1) - (vy3*vx1)
  1008. xl2 = sqrt((vx2*vx2) + (vy2*vy2)+ (vz2*vz2))
  1009. vx2 = vx2/xl2
  1010. vy2 = vy2/xl2
  1011. vz2 = vz2/xl2
  1012. WORK(1+4*idim) = WORK(1)*vx1 + WORK(2)*vx2 + WORK(3)*vx3
  1013. WORK(2+4*idim) = WORK(1)*vy1 + WORK(2)*vy2 + WORK(3)*vy3
  1014. WORK(3+4*idim) = WORK(1)*vz1 + WORK(2)*vz2 + WORK(3)*vz3
  1015. WORK(4+4*idim) = WORK(4)*vx1 + WORK(5)*vx2 + WORK(6)*vx3
  1016. WORK(5+4*idim) = WORK(4)*vy1 + WORK(5)*vy2 + WORK(6)*vy3
  1017. WORK(6+4*idim) = WORK(4)*vz1 + WORK(5)*vz2 + WORK(6)*vz3
  1018. WORK(7+4*idim) = WORK(7)*vx1 + WORK(8)*vx2 + WORK(9)*vx3
  1019. WORK(8+4*idim) = WORK(7)*vy1 + WORK(8)*vy2 + WORK(9)*vy3
  1020. WORK(9+4*idim) = WORK(7)*vz1 + WORK(8)*vz2 + WORK(9)*vz3
  1021. WORK(10+4*idim) = WORK(10)*vx1 + WORK(11)*vx2 + WORK(12)*vx3
  1022. WORK(11+4*idim) = WORK(10)*vy1 + WORK(11)*vy2 + WORK(12)*vy3
  1023. WORK(12+4*idim) = WORK(10)*vz1 + WORK(11)*vz2 + WORK(12)*vz3
  1024. DO i=1,LRE
  1025. XFORC(i) = WORK(i+4*idim) + XFORC(i)
  1026. ENDDO
  1027. END IF
  1028. 2722 CONTINUE
  1029. c write(6,*) 'xforc coa2', (xforc(io),io = 1,lre)
  1030. C
  1031. C RANGEMENT DANS MELVAL
  1032. C
  1033. IE=0
  1034. MPTVAL=IVAFOR
  1035. C
  1036. C NODE=4= NOMBRE DE NOEUDS
  1037. C ICOMP=NSTRS= NOMBRE DE COMPOSANTES DE CONTRAINTES PAR NOEUD
  1038. C
  1039. DO 2724 NODE=1,4
  1040. DO 2724 ICOMP=1,NSTRS
  1041. IE=IE+1
  1042. MELVAL=IVAL(ICOMP)
  1043. IBMN=MIN(IB ,VELCHE(/2))
  1044. VELCHE(NODE,IBMN)=XFORC(IE)
  1045. 2724 CONTINUE
  1046. 2721 CONTINUE
  1047. 9985 CONTINUE
  1048. SEGSUP WRK1,WRK2,WRK3,WRK4
  1049. GO TO 510
  1050.  
  1051. *-----------------------------------------------
  1052. * element shb8
  1053. *------------------------------------------------
  1054. 260 CONTINUE
  1055. NBBB=NBNN
  1056. SEGINI WRK1,WRK7
  1057. DO 3260 IB=1,NBELEM
  1058. c
  1059. c on cherche les coordonnees des noeuds de l element ib
  1060. c
  1061. CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  1062. c
  1063. c on cherche les contraintes
  1064. c
  1065. IE=0
  1066. MPTVAL=IVASTR
  1067. DO 7260 IGAU=1,NBPGAU
  1068. DO 7260 ICOMP=1,NSTRS
  1069. IE=IE+1
  1070. MELVAL=IVAL(ICOMP)
  1071. IGMN=MIN(IGAU,VELCHE(/1))
  1072. IBMN=MIN(IB ,VELCHE(/2))
  1073. WORK1(IE)=VELCHE(IGMN,IBMN)
  1074. 7260 CONTINUE
  1075. c
  1076. c on cherche les materiaux
  1077. c
  1078. yyo=0.d0
  1079. ynu=0.d0
  1080. MPTVAL=IVAMAT
  1081. segact mptval
  1082. DO 6260 IGAU=1,NBPGAU
  1083. DO 6260 ICOMP=1,2
  1084. MELVAL=IVAL(ICOMP)
  1085. IF (MELVAL.NE.0) THEN
  1086. IGMN=MIN(IGAU,VELCHE(/1))
  1087. IBMN=MIN(IB ,VELCHE(/2))
  1088. if(icomp.eq.1) yyo=yyo+VELCHE(IGMN,IBMN)/nbpgau
  1089. if(icomp.eq.2)ynu=ynu+VELCHE(IGMN,IBMN)/nbpgau
  1090. ENDIF
  1091. 6260 CONTINUE
  1092. c
  1093. c on calcule b*sigma
  1094. c
  1095. d(1)=yyo
  1096. d(2)=ynu
  1097. d(3)=0
  1098. CALL shb8(8,XE,D,propel,work1,rel,xforc)
  1099. c
  1100. c rangement dans melval
  1101. c
  1102. IE=0
  1103. MPTVAL=IVAFOR
  1104. DO 9260 IGAU=1,Nbnn
  1105. DO 9260 ICOMP=1,3
  1106. IE=IE+1
  1107. MELVAL=IVAL(ICOMP)
  1108. IBMN=MIN(IB ,VELCHE(/2))
  1109. VELCHE(IGAU,IBMN)=XFORC(IE)
  1110. 9260 CONTINUE
  1111. 3260 CONTINUE
  1112. SEGSUP WRK1,WRk7
  1113. GO TO 510
  1114. C_______________________________________________________________________
  1115. C
  1116. C LIA2 : element de liaison a 2 noeuds (6 ddl par
  1117. C noeuds)
  1118. C_______________________________________________________________________
  1119. C
  1120. 125 CONTINUE
  1121. NBBB=NBNN
  1122. NSTN=3
  1123. LRN =3
  1124. SEGINI WRK1,WRK3,WRK5
  1125. C
  1126. KERRE=0
  1127. DO 3109 IB=1,NBELEM
  1128. C
  1129. C MISE A ZERO DES FORCES INTERNES
  1130. C
  1131. CALL ZERO(XFORC,1,LRE)
  1132. C
  1133. MPTVAL=IVACAR
  1134. DO IC=1,NCARR
  1135. IF(IVAL(IC).NE.0) THEN
  1136. MELVAL=IVAL(IC)
  1137. IBMN=MIN(IB,VELCHE(/2))
  1138. WORK(IC)=VELCHE(1,IBMN)
  1139. ELSE
  1140. WORK(IC)=0.D0
  1141. ENDIF
  1142. END DO
  1143. C
  1144. CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  1145. CALL MAPALI(XE,NBNN,WORK,XGENE,KERRE)
  1146. IF(KERRE.NE.0) THEN
  1147. INTERR(1)=ISOUS
  1148. INTERR(2)=IB
  1149. IF(KERRE.EQ.1) CALL ERREUR(128)
  1150. ENDIF
  1151. C
  1152. C ON CHERCHE LES CONTRAINTES -
  1153. C
  1154. MPTVAL=IVASTR
  1155. IE=0
  1156. DO 7109 IGAU=1,2
  1157. DO 7109 ICOMP=1,NSTRS
  1158. IE=IE+1
  1159. MELVAL=IVAL(ICOMP)
  1160. IBMN=MIN(IB ,VELCHE(/2))
  1161. IGMN=MIN(IGAU,VELCHE(/1))
  1162. XSTRS(IE)=VELCHE(IGMN,IBMN)
  1163. 7109 CONTINUE
  1164. C
  1165. C ON CALCULE LES FORCES INTERNES !!!! a ameliorer
  1166. C
  1167. CALL LIAINT(XSTRS,XGENE,XFORC,LRE,NSTRS)
  1168. C
  1169. C RANGEMENT DANS MELVAL
  1170. C
  1171. IE=0
  1172. MPTVAL=IVAFOR
  1173. DO 9109 IGAU=1,NBNN
  1174. DO 9109 ICOMP=1,NFOR
  1175. IE=IE+1
  1176. MELVAL=IVAL(ICOMP)
  1177. IBMN=MIN(IB ,VELCHE(/2))
  1178. VELCHE(IGAU,IBMN)=XFORC(IE)
  1179. 9109 CONTINUE
  1180. 3109 CONTINUE
  1181. SEGSUP WRK1,WRK3,WRK5
  1182. GO TO 510
  1183. C_______________________________________________________________________
  1184. C
  1185. C JOI1 : element de liaison a 2 noeuds (6 ddl par
  1186. C noeuds)
  1187. C_______________________________________________________________________
  1188. C
  1189. 265 CONTINUE
  1190. NBBB=NBNN
  1191. NSTN=3
  1192. LRN =3
  1193. SEGINI WRK1,WRK3,WRK4
  1194. C
  1195. KERRE=0
  1196. DO 3110 IB=1,NBELEM
  1197. C
  1198. C MISE A ZERO DES FORCES INTERNES
  1199. C
  1200. CALL ZERO(XFORC,1,LRE)
  1201. C
  1202. MPTVAL=IVAMAT
  1203. DO IC=1,NMATT
  1204. IF(IVAL(IC).NE.0) THEN
  1205. MELVAL=IVAL(IC)
  1206. IBMN=MIN(IB,VELCHE(/2))
  1207. WORK(IC)=VELCHE(1,IBMN)
  1208. ELSE
  1209. WORK(IC)=0.D0
  1210. ENDIF
  1211. END DO
  1212. C
  1213. CALL MAPALU(NMATT,WORK,BPSS,IDIM)
  1214. C
  1215. C ON CHERCHE LES CONTRAINTES -
  1216. C
  1217. MPTVAL=IVASTR
  1218.  
  1219. IE=0
  1220. DO 7110 ICOMP=1,NSTRS
  1221. IE=IE+1
  1222. MELVAL=IVAL(ICOMP)
  1223. IBMN=MIN(IB ,VELCHE(/2))
  1224. XSTRS(IE)=VELCHE(1,IBMN)
  1225. 7110 CONTINUE
  1226. C
  1227. C ON CALCULE LES FORCES INTERNES
  1228. C
  1229. CALL INTJOI(XSTRS,XFORC,NSTRS)
  1230. C
  1231. C ON PASSE LES CONTRAINTES DANS LE REPERE GLOBAL
  1232. C
  1233. IAW1 = 101
  1234. IAW2=IAW1+LRE
  1235. CALL JOIGLV(XFORC,BPSS,WORK(IAW1),WORK(IAW2),LRE,IDIM)
  1236. C
  1237. C RANGEMENT DANS MELVAL
  1238. C
  1239. IE=0
  1240. MPTVAL=IVAFOR
  1241. DO 9110 IGAU=1,NBNN
  1242. DO 9110 ICOMP=1,NFOR
  1243. IE=IE+1
  1244. MELVAL=IVAL(ICOMP)
  1245. IBMN=MIN(IB ,VELCHE(/2))
  1246. VELCHE(IGAU,IBMN)=XFORC(IE)
  1247. 9110 CONTINUE
  1248. 3110 CONTINUE
  1249. SEGSUP WRK1,WRK3,WRK4
  1250. GO TO 510
  1251. c_______________________________________________________________________
  1252. c
  1253. c element tuyo
  1254. c_______________________________________________________________________
  1255. c
  1256. 96 CONTINUE
  1257. NBBB=NBNN
  1258. SEGINI WRK1,WRK3
  1259. c
  1260. DO 3096 IB=1,NBELEM
  1261. c
  1262. c on cherche les coordonnees des noeuds de l elementib
  1263. c
  1264. CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  1265. c
  1266. c il faudrait aussi modifier le vecteur local de la poutre
  1267. c
  1268. c mise a zero des forces internes
  1269. c
  1270. CALL ZERO(XFORC,1,LRE)
  1271. c
  1272. c rangement des caracteristiques dans work
  1273. c
  1274. MPTVAL=IVACAR
  1275. DO 6096 IC=1,NCARR
  1276. IF (IVAL(IC).NE.0) THEN
  1277. MELVAL=IVAL(IC)
  1278. IBMN=MIN(IB,VELCHE(/2))
  1279. WORK(IC)=VELCHE(1,IBMN)
  1280. ELSE
  1281. WORK(IC)=0.D0
  1282. ENDIF
  1283. 6096 CONTINUE
  1284. c
  1285. c cas ou on a lu le mot vecteur
  1286. c
  1287. IF (IVECT.EQ.1) THEN
  1288. IF (IVAL(NCARR).NE.0) THEN
  1289. MELVAL=IVAL(NCARR)
  1290. IBMN=MIN(IB,IELCHE(/2))
  1291. IP=IELCHE(1,IBMN)
  1292. IREF=(IP-1)*(IDIM+1)
  1293. DO 6196 IC=1,IDIM
  1294. WORK(NCARR+IC-1)=XCOOR(IREF+IC)
  1295. 6196 CONTINUE
  1296. ELSE
  1297. DO 6296 IC=1,IDIM
  1298. WORK(NCARR+IC-1)=0.D0
  1299. 6296 CONTINUE
  1300. ENDIF
  1301. c
  1302. c cas du chamelem comverti
  1303. c
  1304. ELSE IF (IVECT.EQ.2) THEN
  1305. DO 6496 IC=1,IDIM
  1306. MELVAL=IVAL(NCARR+IC-3)
  1307. IF (MELVAL.NE.0) THEN
  1308. IBMN=MIN(IB,VELCHE(/2))
  1309. WORK(NCARR+IC-3)=VELCHE(1,IBMN)
  1310. ELSE
  1311. WORK(NCARR+IC-3)=0.D0
  1312. ENDIF
  1313. 6496 CONTINUE
  1314. ENDIF
  1315. c
  1316. c
  1317. c cas des tuyaux - on calcule les caracteristiques de la poutre
  1318. c equivalente
  1319. c
  1320. *Bizarre ici MELE = 96 ???
  1321. IF(MELE.EQ.42) THEN
  1322. CISA=WORK(4)
  1323. VX =WORK(5)
  1324. VY =WORK(6)
  1325. VZ =WORK(7)
  1326. CALL TUYCAR(WORK,CISA,VX,VY,VZ,KERRE,0)
  1327. IF (KERRE.EQ.77) THEN
  1328. CALL ERREUR(77)
  1329. GOTO 510
  1330. ENDIF
  1331. ENDIF
  1332. c
  1333. c on cherche les contraintes -
  1334. c
  1335. MPTVAL=IVASTR
  1336. IE=9
  1337. DO 7096 IGAU=1,2
  1338. DO 7096 ICOMP=1,NSTRS
  1339. IE=IE+1
  1340. MELVAL=IVAL(ICOMP)
  1341. IBMN=MIN(IB ,VELCHE(/2))
  1342. IGMN=MIN(IGAU,VELCHE(/1))
  1343. WORK(IE)=VELCHE(IGMN,IBMN)
  1344. 7096 CONTINUE
  1345. c
  1346. c on calcule les forces internes
  1347. c
  1348. CALL POUBSG(XFORC,WORK,XE,WORK(10),WORK(22),KERRE)
  1349. IF(KERRE.EQ.0) GO TO 5096
  1350. INTERR(1)=ISOUS
  1351. INTERR(2)=IB
  1352. SEGSUP WRK1,WRK3
  1353. CALL ERREUR(128)
  1354. GO TO 510
  1355. 5096 CONTINUE
  1356. c
  1357. c rangement dans melval
  1358. c
  1359. IE=0
  1360. MPTVAL=IVAFOR
  1361. DO 9096 IGAU=1,NBNN
  1362. DO 9096 ICOMP=1,NFOR
  1363. IE=IE+1
  1364. MELVAL=IVAL(ICOMP)
  1365. IBMN=MIN(IB ,VELCHE(/2))
  1366. VELCHE(IGAU,IBMN)=XFORC(IE)
  1367. 9096 CONTINUE
  1368. 3096 CONTINUE
  1369. SEGSUP WRK1,WRK3
  1370. GO TO 510
  1371.  
  1372. C_______________________________________________________________________
  1373. C
  1374. C ELEMENTS ZONE_COHESIVE ZOC2,ZOC3,ZOC4
  1375. C_______________________________________________________________________
  1376. C
  1377. 266 CONTINUE
  1378.  
  1379. NDIM = 2
  1380. IF(IFOUR.GT.0) NDIM = 3
  1381. NBBB=NBNN
  1382. SEGINI WRK1,WRK2,WRK4
  1383. C
  1384. DO 3266 IB=1,NBELEM
  1385. C
  1386. C ON CHERCHE LES COORDONNEES DES NOEUDS DE L ELEMENT IB
  1387. C
  1388. CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  1389. C
  1390. C MISE A ZERO DES FORCES INTERNES
  1391. C
  1392. CALL ZERO(XFORC,1,LRE)
  1393. C
  1394. C BOUCLE SUR LES POINTS DE GAUSS
  1395. C
  1396. DO 6266 IGAU=1,NBPGAU
  1397. C
  1398. CALL ZCOLOC(XE,SHPTOT,NBNN,MELE,IFOUR,IGAU,BPSS)
  1399. C
  1400. CALL BZCO(IGAU,MFR,IFOUR,NIFOUR,XE,BPSS,SHPTOT,
  1401. . NSTRS,NBNN,LRE,MELE,SHPWRK,BGENE,DJAC,IERT)
  1402. IF (IERT.NE.0) THEN
  1403. INTERR(1)=IB
  1404. CALL ERREUR(612)
  1405. GOTO 99266
  1406. ENDIF
  1407. C
  1408. C ON CHERCHE LES CONTRAINTES -
  1409. C
  1410. MPTVAL=IVASTR
  1411. DO 7266 ICOMP=1,NSTRS
  1412. MELVAL=IVAL(ICOMP)
  1413. IGMN=MIN(IGAU,VELCHE(/1))
  1414. IBMN=MIN(IB ,VELCHE(/2))
  1415. XSTRS(ICOMP)=VELCHE(IGMN,IBMN)
  1416. 7266 CONTINUE
  1417. C
  1418. C ON CALCULE B*EFFORTS
  1419. C
  1420. DJAC=DJAC*POIGAU(IGAU)
  1421. CALL BSIG(BGENE,XSTRS,NSTRS,LRE,DJAC,XFORC)
  1422. 6266 CONTINUE
  1423. C
  1424. C TRAITEMENT DE XFORC ET RANGEMENT DANS MELVAL
  1425. C
  1426. IE=0
  1427. MPTVAL=IVAFOR
  1428. C
  1429. DO 9266 NODE=1,NBNN
  1430. DO 9266 ICOMP=1,NDIM
  1431. IE=IE+1
  1432. MELVAL=IVAL(ICOMP)
  1433. IBMN=MIN(IB ,VELCHE(/2))
  1434. VELCHE(NODE,IBMN)=XFORC(IE)
  1435. 9266 CONTINUE
  1436. 3266 CONTINUE
  1437.  
  1438. 99266 CONTINUE
  1439. SEGSUP WRK1,WRK2,WRK4
  1440. GOTO 510
  1441. c_______________________________________________________________________
  1442. 99 CONTINUE
  1443. MOTERR(1:4)=NOMTP(MELE)
  1444. MOTERR(5:12)='BSIGMA'
  1445. CALL ERREUR(86)
  1446. *
  1447. 510 CONTINUE
  1448. RETURN
  1449. END
  1450.  
  1451.  
  1452.  
  1453.  
  1454.  
  1455.  

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