Télécharger pjmode.eso

Retour à la liste

Numérotation des lignes :

  1. C PJMODE SOURCE CB215821 19/07/30 21:17:31 10273
  2. SUBROUTINE PJMODE(ipmode)
  3. IMPLICIT INTEGER(I-N)
  4. IMPLICIT REAL*8 (A-H,O-Z)
  5. C=======================================================================
  6. C OPERATEUR PJBA :
  7. C PROJECTION D'UN CHPOINT, D'UN CHARGEMENT OU D'UNE RIGIDITE
  8. C SUR LES ELEMENTS D'UNE BASE MODALE B.
  9. C LE RESULTAT EST DU MEME TYPE.
  10. C
  11. C SYNTAXE :
  12. C * FN = PJBA B OBJET ; SI BASE ELEMENTAIRE
  13. C * FN = PJBA B STR1 (N) OBJET ; SI BASE COMPLEXE
  14. C
  15. C OBJET POUVANT ETRE UNE FORCE OU UN CHARGEMENT,
  16. C OU UNE RIGIDITE DANS LE PREMIER CAS.
  17. C=======================================================================
  18. ***********************************************************
  19. * PROJECTION D'UNE MATRICE SUR UNE BASE DE MODES *
  20. * _______________________________________________________ *
  21. * *
  22. * DATE : le 11 Avril 1995 *
  23. * AUTEUR : Nicolas BENECH *
  24. * _______________________________________________________ *
  25. * *
  26. * MODULE(S) APPELANT(S) : PJBA *
  27. * *
  28. * MODULE(S) APPELE(S) : ACCTAB, YTMX *
  29. * _______________________________________________________ *
  30. * *
  31. * EN ENTREE : *
  32. * MRIGID : Matrice a projeter *
  33. * MTAB1 : Base de modes, reels ou complexes *
  34. * 'REEL' : indique que l'on utilise le produit *
  35. * scalaire reel (pas de conjugaison) *
  36. * *
  37. * EN SORTIE : *
  38. * RI1 : Matrice projetee (partie reelle) *
  39. * RI2 : Matrice projetee (partie imaginaire) *
  40. * _______________________________________________________ *
  41. * *
  42. * REMARQUE : *
  43. * L'operation realisee est : *
  44. * (MTAB1)t . MRIGID. MTAB1 *
  45. * Dans le cas complexe, la transposition est accompagnee *
  46. * de la conjugaison (si REEL n'est pas mentionne). *
  47. *
  48. * voir aussi PROJTA
  49. ***********************************************************
  50. *
  51. -INC SMCHPOI
  52. -INC SMCHARG
  53. -INC SMLCHPO
  54. -INC CCOPTIO
  55. -INC CCGEOME
  56. -INC CCREEL
  57. -INC SMELEME
  58. -INC CCHAMP
  59. -INC SMCHAML
  60. -INC SMMODEL
  61. -INC SMRIGID
  62. -INC SMLMOTS
  63. -INC SMLENTI
  64.  
  65. C
  66. * Declarations
  67. *
  68. PARAMETER(ZERO=0.D0)
  69. REAL*8 XVAL, RMAX
  70. CHARACTER*8 LETYPE
  71. CHARACTER*8 TYPMOD,cmate
  72. LOGICAL MODCOM,dedans,dchpo,l3,lr2,lirl
  73. INTEGER I, J, NBMOD, POS, IREEL, IVALRE, IOBRE
  74. REAL*8 XVALRE
  75. LOGICAL LOGRE
  76. segment plcf
  77. integer lpref(ldepl),ldefo(ldepl),lmade(ldepl)
  78. real*8 prmas(ldepl)
  79. endsegment
  80. segment pmod
  81. integer kdefo(nbmod)
  82. endsegment
  83. segment prigmat
  84. integer lrigmat(nrigmat,2+9)
  85. endsegment
  86. segment pmapmo
  87. integer lmapmo(nmapmo),defpmo(nmapmo),dimpmo(nmapmo)
  88. character*4 compmo(nmapmo)
  89. real*8 coepmo(nmapmo)
  90. endsegment
  91. segment pcompo
  92. character*4 mcol
  93. real*8 valmod(nipmod)
  94. endsegment
  95. LOGICAL L0,L1,lcf
  96. PARAMETER (ncod=8)
  97. CHARACTER*4 IDDL,lcod(ncod),lcof(ncod),motinc
  98. CHARACTER*8 TYPRET,CHARRE
  99. data xlopre/1.d-11/
  100. DATA KZERO/0/
  101. data lcod/'UX','UY','UZ','RX','RY','RZ','UR','UT'/
  102. data lcof/'FX','FY','FZ','MX','MY','MZ','FR','FT'/
  103.  
  104. plcf = 0
  105. jgn = 4
  106. jgm = ncod
  107. segini mlmot5
  108. segini mlmot6
  109. do io = 1,ncod
  110. mlmot5.mots(io) = lcod(io)
  111. mlmot6.mots(io) = lcof(io)
  112. enddo
  113.  
  114. modcom = .false.
  115. dchpo = .false.
  116. iriout = 0
  117. iriout1 = 0
  118. iriout2 = 0
  119. mmodel = ipmode
  120. n1 = kmodel(/1)
  121. segini mmode1
  122. jn = 0
  123. do im = 1, n1
  124. imodel = kmodel(im)
  125. if (formod(1).eq.'MECANIQUE'.and.MATMOD(1).eq.'ELASTIQUE'
  126. &.and.(MATMOD(2).eq.'MODAL'.or.MATMOD(2).eq.'STATIQUE')) then
  127. jn = jn + 1
  128. mmode1.kmodel(jn) = imodel
  129. endif
  130. enddo
  131. if (jn.ne.0) then
  132. n1 = jn
  133. segadj mmode1
  134. ipmode = mmode1
  135. else
  136. segsup mmode1
  137. * cas de projection non pr�vue
  138. call erreur(5)
  139. return
  140. endif
  141.  
  142. call lirobj('MCHAML ',IPCAR1,1,iretou)
  143. call actobj('MCHAML ',IPCAR1,1)
  144. if (ierr.ne.0) return
  145.  
  146. ipchpo = 0
  147. iprigi = 0
  148. call lirobj('CHARGEME',IPCHAR,0,iretou)
  149. if (iretou.eq.0) then
  150. call lirobj('CHPOINT ',IPCHPO,0,iretou)
  151. if(iretou .EQ. 1)call actobj('CHPOINT ',IPCHPO,1)
  152. endif
  153. if (iretou.eq.0) call lirobj('RIGIDITE',IPRIGI,0,iretou)
  154.  
  155. if (iretou.eq.0) then
  156. * manque un op�rande
  157. call erreur(5)
  158. return
  159. endif
  160.  
  161. call reduaf (ipcar1,ipmode,IPCARA,1,iretr,kerre)
  162. if (ierr.ne.0) return
  163. if( iretr.ne.1) then
  164. call erreur (kerre)
  165. return
  166. endif
  167.  
  168. lcf = .false.
  169. mmodel = ipmode
  170. mchelm = ipcara
  171. if (ipchar.ne.0) goto 100
  172. if (iprigi.ne.0) goto 200
  173. if (ipchpo.ne.0) then
  174. n = 1
  175. segini mcharg
  176. ipchar = mcharg
  177. segini icharg
  178. kcharg(1) = icharg
  179. ichpo1 = ipchpo
  180. goto 100
  181. endif
  182.  
  183.  
  184. 100 continue
  185. MCHAR1=IPCHAR
  186. SEGINI,MCHARG=MCHAR1
  187. NBCHG=KCHARG(/1)
  188. DO 10 INCHA=1,NBCHG
  189. ICHAR1=KCHARG(INCHA)
  190. SEGINI,ICHARG=ICHAR1
  191. KCHARG(INCHA)=ICHARG
  192. IP1=ICHPO1
  193. c
  194. IRET = 0
  195. c
  196. c deplacement impose => idepi=1
  197. c force imposee => idepi=0
  198. c
  199. IDEPI = 0
  200. c idepi = -1
  201. KDEPI = 0
  202. MCHPOI = IP1
  203. IF (MTYPOI.EQ.'FLX ') IDEPI = 1
  204. c if (idepi.lt.0) then
  205. c moterr(1:8) = 'chpoint'
  206. c call erreur(302)
  207. c return
  208. c endif
  209. c
  210. NBNN = 1
  211. NBREF = 0
  212. NBSOUS = 0
  213. *
  214. LDEPL = kmodel(/1)
  215. if (.not.lcf) segini plcf
  216. c
  217. c
  218. c **** on initialise le chpoint
  219. c
  220. NSOUPO = 1
  221. NAT=1
  222. SEGINI,MCHPOI
  223. IRET = MCHPOI
  224. MTYPOI = ' '
  225. MOCHDE=' J''AI ETE FABRIQUE PAR L''OPERATEUR PJBA'
  226. IFOPOI = IFOUR
  227. * champ de force nodal: nature discrete
  228. JATTRI(1)=2
  229. NC = 1
  230. SEGINI,MSOUPO
  231. IPCHP(1) = MSOUPO
  232. NOHARM(1) = NIFOUR
  233. NOCOMP(1) = 'FALF'
  234.  
  235. do 101 inocomp=1,2
  236.  
  237. N = LDEPL
  238. SEGINI MPOVAL
  239. IPOVAL = MPOVAL
  240. *
  241. NBNN = 1
  242. NBELEM = LDEPL
  243. NBSOUS = 0
  244. NBREF = 0
  245. SEGINI MELEME
  246. IGEOC = MELEME
  247. ITYPEL = 1
  248.  
  249. knum = 0
  250. c
  251. c ****boucle sur les chpoints de depl
  252. c
  253. DO 11 IM = 1,kmodel(/1)
  254. imodel = kmodel(im)
  255. nomid = lnomid(2)
  256. if (.not.lcf) then
  257. ipt1 = imamod
  258. iptr = ipt1.num(1,1)
  259. lpref(im) = iptr
  260.  
  261. indc = 1
  262. 34 if (imache(indc).ne.imamod.or.conche(indc).ne.conmod) then
  263. indc = indc + 1
  264. if (indc.gt.imache(/1)) then
  265. * champ de caracteristiques incomplet
  266. goto 99
  267. endif
  268. goto 34
  269. endif
  270.  
  271. mchaml = ichaml(indc)
  272. do iij = 1, nomche(/2)
  273. if (nomche(iij).eq.'DEFO') then
  274. melval = ielval(iij)
  275. ipp1 = ielche(1,1)
  276. ldefo(im) = ipp1
  277. endif
  278. if (nomche(iij).eq.'MADE') then
  279. melval = ielval(iij)
  280. ipp2 = ielche(1,1)
  281. lmade(im) = ipp2
  282. endif
  283. if (nomche(iij).eq.'MASS') then
  284. melval = ielval(iij)
  285. ymass = velche(1,1)
  286. prmas(im) = ymass
  287. endif
  288. if(ldefo(im).gt.0.and.lmade(im).gt.0.and.
  289. &prmas(im).gt.0) goto 35
  290. enddo
  291. 35 continue
  292. if (ldefo(im).eq.0) goto 99
  293. if (prmas(im).le.0.and.cmatee(1:5).eq.'MODAL') goto 99
  294. if (lmade(im).eq.0.and.cmatee(1:8).eq.'STATIQUE') goto 99
  295.  
  296. endif
  297.  
  298. if (NOCOMP(1).ne.lesobl(1)(1:4)) goto 11
  299. knum = knum + 1
  300.  
  301. iptr = lpref(im)
  302. ipp1 = ldefo(im)
  303. NUM(1,knum) = IPTR
  304. ICOLOR(knum) = IDCOUL
  305. XRET = 0.D0
  306. call xty1(ipp1,ip1,mlmot5,mlmot6,xret)
  307. IF (IDEPI.NE.1) THEN
  308. ELSE
  309. * ??
  310. indn = 1
  311. 45 if (nomche(indn).ne.'FREQ') then
  312. indn = indn + 1
  313. if (indn.gt.nomche(/2)) then
  314. * pas la composante FREQ
  315. goto 99
  316. endif
  317. goto 45
  318. endif
  319.  
  320. melval = ielval(indn)
  321. x1 = velche(1,1)
  322. OM = X1
  323. OM = 2.D0 * XPI * OM
  324. OM = OM * OM
  325. XRET = -XRET / OM
  326. ENDIF
  327. IF (IFOUR .EQ. 1) THEN
  328. IF (NIFOUR .NE. 0) THEN
  329. XRET = XRET*XPI
  330. ELSE
  331. XRET = XRET*2.D0*XPI
  332. ENDIF
  333. ENDIF
  334. VPOCHA(knum,1) = XRET
  335. *
  336. if (cmatee(1:5).eq.'MODAL') then
  337. ymass = prmas(im)
  338. elseif (cmatee(1:8).eq.'STATIQUE') then
  339. ipp2 = lmade(im)
  340. call xty1(ipp1,ipp2,mlmot5,mlmot6,ymass)
  341. else
  342. endif
  343. if (lmade(im).gt.0.and.ABS(XRET).gt.(1.d-10*ymass).and.
  344. & ymass.gt.0.and.cmatee(1:5).eq.'MODAL') then
  345. * kich : on enleve la projection sur base modale - a creuser pour statique !
  346. CALL ADCHPO(IP1,IPP2,IP2,1.d0,(XRET/ymass*(-1.d0)))
  347. IP1 = IP2
  348. endif
  349. *
  350. 11 CONTINUE
  351. *
  352. lcf = .true.
  353. *
  354. *
  355. if (knum.eq.kmodel(/1)) goto 102
  356. if (inocomp.eq.1) then
  357. if (knum.eq.0) then
  358. NOCOMP(1) = 'FBET'
  359. else
  360. N = knum
  361. NBELEM = knum
  362. segadj MPOVAL,MELEME
  363. NSOUPO = 2
  364. segadj MCHPOI
  365. SEGINI,MSOUPO
  366. IPCHP(2) = MSOUPO
  367. NOCOMP(1) = 'FBET'
  368. endif
  369. endif
  370. 101 continue
  371.  
  372. 102 continue
  373. N = knum
  374. NBELEM = knum
  375. segadj MPOVAL,MELEME
  376.  
  377. IF(IERR.NE.0) RETURN
  378. ICHPO1=IRET
  379. SEGDES,ICHARG
  380. 10 CONTINUE
  381. segsup mlmot5,mlmot6,plcf
  382. if (ipchpo.gt.0) then
  383. segsup icharg,mcharg
  384. call actobj('CHPOINT ',iret,1)
  385. call ecrobj('CHPOINT ',iret)
  386. goto 999
  387. endif
  388. SEGDES,MCHARG
  389. CALL ECROBJ('CHARGEME',MCHARG)
  390.  
  391. goto 999
  392. 99 segsup mpoval,msoupo,mchpoi
  393. call erreur(26)
  394. return
  395.  
  396.  
  397. 200 continue
  398. ipri1 = iprigi
  399. call SEPA(ipri1,1)
  400. ipri2 = iprigi
  401. call SEPA(ipri2,2)
  402. *
  403. *
  404. *
  405. *
  406. nmapmo = 100
  407. kpmo = 0
  408. segini pmapmo
  409. do isous = 1,kmodel(/1)
  410. imodel = kmodel(isous)
  411. cmate = cmatee
  412. meleme = imamod
  413. if (itypel.ne.1) call erreur(5)
  414. if (num(/1).ne.1) call erreur(5)
  415. if (cmate.eq.'STATIQUE'.or.cmate.EQ.'MODAL') then
  416. do ilp = 1,num(/2)
  417. kpmo = kpmo + 1
  418. if (kpmo.gt.nmapmo) then
  419. nmapmo = nmapmo + 100
  420. segadj pmapmo
  421. endif
  422. lmapmo(kpmo) = num(1,ilp)
  423. if (cmate.eq.'STATIQUE') then
  424. compmo(kpmo) = 'BETA'
  425. elseif (cmate.eq.'MODAL') then
  426. compmo(kpmo) = 'ALFA'
  427. endif
  428. do im = 1 , imache(/1)
  429. if (imache(im).eq.imamod) then
  430. if (conche(im).eq.conmod) then
  431. mchaml = ichaml(im)
  432. do iv = 1,ielval(/1)
  433. if (nomche(iv).eq.'DEFO') then
  434. melval = ielval(iv)
  435. ibmn = min(ilp,ielche(/2))
  436. defpmo(kpmo) = ielche(1,ibmn)
  437. endif
  438. if (nomche(iv).eq.'IDEF') then
  439. melval = ielval(iv)
  440. ibmn = min(ilp,ielche(/2))
  441. dimpmo(kpmo) = ielche(1,ibmn)
  442. endif
  443. enddo
  444. endif
  445. endif
  446. enddo
  447.  
  448. enddo
  449. endif
  450. enddo
  451.  
  452. nmapmo = kpmo
  453. segadj pmapmo
  454. nbmod = nmapmo
  455. *
  456. N1 = NBMOD
  457. nbcod = 8
  458. segini pmod
  459. SEGINI, MLCHP1
  460. SEGINI, MLCHP2
  461. jgm = 1
  462. jgn = 4
  463. segini mlmot4
  464. *
  465. * Constitution du maillage support et du segment descriptif
  466. *
  467. NBNN = NBMOD
  468. NBELEM = 1
  469. NBSOUS = 0
  470. NBREF = 0
  471. SEGINI, MELEME
  472. ITYPEL = 1
  473. *
  474. NLIGRD=NBMOD
  475. NLIGRP=NBMOD
  476. SEGINI, DESCR
  477. *
  478. mrigid = ipri1
  479. segact mrigid
  480. nrigel = coerig(/1)
  481. if (nrigel.lt.1) goto 250
  482. typmod = ' '
  483. IREEL = -1
  484. C* POS ? IF (POS.EQ.1) IREEL = 1
  485. *
  486. LETYPE = ' '
  487. DO 210 IM=1,NBMOD
  488. IPT1 = 0
  489. *
  490. imodel = kmodel(im)
  491. ipt1 = imamod
  492. iptr = ipt1.num(1,1)
  493. * Cas reel ou cas complexe ?
  494. *
  495. if (dimpmo(im).gt.0) TYPMOD = 'MODE_COM'
  496.  
  497. IF (TYPMOD .EQ. 'MODE_COM') THEN
  498. MODCOM=.TRUE.
  499. mchpoi = defpmo(im)
  500. MLCHP1.ICHPOI(IM) = MCHPOI
  501. mchpoi = dimpmo(im)
  502. MLCHP2.ICHPOI(IM) = MCHPOI
  503. ELSE
  504. MODCOM = .FALSE.
  505. mchpoi = defpmo(im)
  506. MLCHP1.ICHPOI(IM) = MCHPOI
  507. ENDIF
  508. *
  509. ipt1 = iptr
  510. *
  511. MELEME.NUM(IM,1)=IPT1
  512. *
  513. if (cmatee.eq.'MODAL') then
  514. DESCR.LISINC(IM) = 'ALFA'
  515. DESCR.LISDUA(IM) = 'FALF'
  516. else if (cmatee.eq.'STATIQUE') then
  517. DESCR.LISINC(IM) = 'BETA'
  518. DESCR.LISDUA(IM) = 'FBET'
  519. endif
  520. DESCR.NOELEP(IM) = IM
  521. DESCR.NOELED(IM) = IM
  522. *
  523. 210 CONTINUE
  524. *
  525. * Constitution des segments XMATRI
  526. *
  527. NLIGRD=NBMOD
  528. NLIGRP=NBMOD
  529. nelrig=1
  530. *
  531. IF (LETYPE .EQ. 'ANNULE') THEN
  532. SEGINI, XMATR1
  533. IF (MODCOM) THEN
  534. SEGINI, XMATR2
  535. SEGDES, XMATR2
  536. ENDIF
  537. SEGDES, XMATR1
  538. GOTO 55
  539. ENDIF
  540. *
  541. * Cas reel
  542. *
  543. SEGINI, XMATR1
  544. DO 20, I=1, NBMOD
  545. MCHPO1 = MLCHP1.ICHPOI(I)
  546. DO 20, J=1, NBMOD
  547. MCHPO2 = MLCHP1.ICHPOI(J)
  548. CALL YTMX (MCHPO2, MCHPO1, MRIGID, XVAL)
  549. XMATR1.RE(I,J,1)=XVAL
  550. 20 CONTINUE
  551. SEGDES, XMATR1
  552. *
  553. * Cas complexe : calcul de termes complementaires
  554. *
  555. IF (MODCOM) THEN
  556. SEGACT, XMATR1*mod
  557. DO 30, I=1, NBMOD
  558. MCHPO1 = MLCHP2.ICHPOI(I)
  559. DO 30, J=1, NBMOD
  560. MCHPO2 = MLCHP2.ICHPOI(J)
  561. CALL YTMX (MCHPO1, MCHPO2, MRIGID, XVAL)
  562. XMATR1.RE(I,J,1)=XMATR1.RE(I,J,1)-IREEL*XVAL
  563. 30 CONTINUE
  564. SEGDES, XMATR1
  565. *
  566. SEGINI, XMATR2
  567. DO 40, I=1, NBMOD
  568. MCHPO1 = MLCHP1.ICHPOI(I)
  569. DO 40, J=1, NBMOD
  570. MCHPO2 = MLCHP2.ICHPOI(J)
  571. CALL YTMX (MCHPO1, MCHPO2, MRIGID, XVAL)
  572. XMATR2.RE(I,J,1)=XVAL
  573. 40 CONTINUE
  574. DO 50, I=1, NBMOD
  575. MCHPO1 = MLCHP2.ICHPOI(I)
  576. DO 50, J=1, NBMOD
  577. MCHPO2 = MLCHP1.ICHPOI(J)
  578. CALL YTMX (MCHPO1, MCHPO2, MRIGID, XVAL)
  579. XMATR2.RE(I,J,1)=XMATR2.RE(I,J,1)+IREEL*XVAL
  580. 50 CONTINUE
  581. SEGDES, XMATR2
  582. ENDIF
  583. *
  584. SEGACT, MRIGID
  585. LETYPE = MRIGID.MTYMAT
  586. SEGDES, MRIGID
  587. *
  588. * Creation des segments IMATRI
  589. *
  590. 55 NELRIG = 1
  591. * SEGINI, IMATR1
  592. * IMATR1.IMATTT(1) = XMATR1
  593. SEGDES, xMATR1
  594. IF (MODCOM) THEN
  595. * SEGINI, IMATR2
  596. * IMATR2.IMATTT(1) = XMATR2
  597. SEGDES, xMATR2
  598. ENDIF
  599. *
  600. * Creation des rigidites calculees
  601. *
  602. NRIGE=7
  603. NRIGEL=1
  604. SEGINI, RI1
  605. RI1.MTYMAT = LETYPE
  606. RI1.IFORIG = IFOMOD
  607. RI1.IMGEO1 = 0
  608. RI1.IMGEO2 = 0
  609. RI1.COERIG = 1.D0
  610. RI1.IRIGEL(1,1) = MELEME
  611. RI1.IRIGEL(2,1) = 0
  612. RI1.IRIGEL(3,1) = DESCR
  613. RI1.IRIGEL(4,1) = xMATR1
  614. RI1.IRIGEL(5,1) = NIFOUR
  615. RI1.IRIGEL(6,1) = 0
  616. RI1.IRIGEL(7,1) = 2
  617. segact xmatr1*mod
  618. xmatr1.symre=2
  619. segdes xmatr1
  620. SEGDES, RI1
  621. IF (MODCOM) THEN
  622. SEGINI, RI2 = RI1
  623. RI2.IRIGEL(4,1) = xMATR2
  624. SEGDES, RI2
  625. ELSE
  626. RI2 = 0
  627. SEGSUP, MLCHP2
  628. ENDIF
  629. *
  630. iriout1 = ri1
  631. iriout2 = ri2
  632.  
  633. 250 continue
  634. mrigid = ipri2
  635. segact mrigid
  636. nrigel = coerig(/1)
  637. if (nrigel.lt.1) goto 290
  638. typmod = ' '
  639.  
  640. nrigmat =100
  641. kgmat = 0
  642. segini prigmat
  643.  
  644. KRIGEL = 0
  645. nrigel = irigel(/2)
  646. nrige = irigel(/1)
  647. segini ri1
  648. ri1.mtymat = mtymat
  649. nrige0 = nrigel
  650.  
  651. kige = 0
  652. kige1 = 100
  653. nrigel = kige1
  654. segini ri2
  655. ri2.mtymat = mtymat
  656.  
  657. DO ire = 1,nrige0
  658. meleme = irigel (1,ire)
  659. segact meleme
  660. if (itypel.ne.22) then
  661. call erreur(977)
  662. return
  663. endif
  664. nbelem = num(/2)
  665. nbele0 = nbelem
  666. descr = irigel(3,ire)
  667. segact descr
  668. nligrp0 = noelep(/1)
  669. nligrd0 = noeled(/1)
  670. nligrp = nligrp0 + nmapmo
  671. nligrd = nligrd0 + nmapmo
  672.  
  673. nbnn = num(/1)
  674. nbsous = 0
  675. nbref = 0
  676. segini ipt2
  677. ipt2.itypel = itypel
  678. nbelem = 1
  679. nbnn = nligrd
  680. segini ipt1
  681. ipt1.itypel = itypel
  682. ri1.coerig(ire) = coerig(ire)
  683. kele = 0
  684.  
  685. xmatr1 = irigel(4,ire)
  686. segact xmatr1
  687. nelrig0 = xmatr1.re(/3)
  688. nelrig = nelrig0 + nmapmo
  689. segini xmatr2
  690. DO iele = 1,nbele0
  691. ie2 = min(iele,nelrig0)
  692. * xmatr1 = imatr1.imattt(ie2)
  693. * segact xmatr1
  694. nligrp = nligrp0 + nmapmo
  695. nligrd = nligrd0 + nmapmo
  696. nelrig=1
  697. segini des2,xmatri
  698. des2.lisinc(1) = 'LX'
  699. des2.lisdua(1) = 'FLX'
  700. des2.noelep(1) = 1
  701. des2.noeled(1) = 1
  702. * le premier point correspond aux multiplicateurs
  703. CALL CREPO1 (ZERO, ZERO, ZERO, IPTS)
  704. ipt1.num(1,1) = ipts
  705. kgrp = 1
  706. kirp = 1
  707. do ipmo = 1,nmapmo
  708. coepmo(ipmo) = 0.d0
  709. enddo
  710. do igrp = 2,nligrp0
  711. jno = noelep(igrp)
  712. motinc = lisinc(igrp)
  713. IP1 = num(jno,iele)
  714. * recherche association noeud physique - points support déformée
  715. do ilmat = 1,kgmat
  716. if(lrigmat(ilmat,1).eq.ip1) goto 315
  717. enddo
  718.  
  719. kgmat = kgmat+1
  720. ilmat = kgmat
  721. if (kgmat.gt.nrigmat) then
  722. nrigmat = nrigmat + 100
  723. segadj prigmat
  724. endif
  725. kpb = 0
  726. jg = 100
  727. segini mlent3
  728. lrigmat(kgmat,1) = ip1
  729. do ikmo = 1, nmapmo
  730. ichp1 = defpmo(ikmo)
  731. call ecrcha('NOMU')
  732. call ecrcha('MAIL')
  733. call ecrobj('CHPOINT ',ichp1)
  734. call extrai
  735. call ecrobj('POINT ',IP1)
  736. call DANS
  737. call lirlog(l3,1,iretou)
  738. if(l3) then
  739. kpb = kpb + 1
  740. if (kpb.gt.jg) then
  741. jg = jg + 100
  742. segadj mlent3
  743. endif
  744. mlent3.lect(kpb) = ikmo
  745. endif
  746. enddo
  747. jg = kpb
  748. segadj mlent3
  749. if (kpb.gt.0) then
  750. lrigmat(ilmat,2) = mlent3
  751. else
  752. lrigmat(ilmat,2) = 0
  753. segsup mlent3
  754. endif
  755.  
  756. 315 continue
  757. ilr3 = lrigmat(ilmat,2)
  758. if (ilr3.eq.0) goto 253
  759. mlent3 = ilr3
  760. segact mlent3
  761. * selection selon nom composante
  762. mlmat = 0
  763. do lmo = 1,9
  764. if (motinc.eq.lcod(lmo)) mlmat = lmo+2
  765. enddo
  766. if (mlmat.eq.0) then
  767. * WRITE(6,*) 'coefs pour cette composante non trouves'
  768. goto 253
  769. endif
  770. if (lrigmat(ilmat,mlmat).ne.0) then
  771. pcompo = lrigmat(ilmat,mlmat)
  772. segact pcompo
  773. nipmod = valmod(/1)
  774. do ilg = 1,nipmod
  775. lkmo = mlent3.lect(ilg)
  776. coepmo(lkmo) = (valmod(ilg)* xmatr1.re(1,igrp,ie2))
  777. & + coepmo(lkmo)
  778. enddo
  779. else
  780. jg = mlent3.lect(/1)
  781. nipmod = jg
  782. segini pcompo
  783. mcol = motinc
  784. do ilg = 1,nipmod
  785. lkmo = mlent3.lect(ilg)
  786. ichp1 = defpmo(lkmo)
  787. CALL EXTRA9(ICHP1,ip1,motinc,0,.false.,XFLOT,IRET)
  788. coepmo(lkmo) = (xflot * xmatr1.re(1,igrp,ie2))
  789. & + coepmo(lkmo)
  790. valmod(ilg) = xflot
  791. enddo
  792. lrigmat(ilmat,mlmat) = pcompo
  793. endif
  794.  
  795. 253 continue
  796. enddo
  797.  
  798. xmaut1 = 0.d0
  799. do kpmo = 1,nmapmo
  800. xmaut1 = max(xmaut1,ABS(coepmo(kpmo)))
  801. enddo
  802.  
  803. * synthèse
  804. do igrp = 2,nligrp0
  805. jno = noelep(igrp)
  806. motinc = lisinc(igrp)
  807. IP1 = num(jno,iele)
  808. lr2 = .false.
  809. do jgmat = 1,kgmat
  810. if(lrigmat(jgmat,1).eq.ip1) goto 325
  811. enddo
  812. c WRITE(6,*) 'bizarre, point dans l element non repertorie'
  813. call erreur(5)
  814. return
  815. 325 continue
  816. mlmat = 0
  817. do lmo = 1,9
  818. if (motinc.eq.lcod(lmo)) mlmat = lmo+2
  819. enddo
  820. if (mlmat.eq.0) lr2 = .true.
  821. if (lrigmat(jgmat,mlmat).eq.0) lr2 = .true.
  822. if(lr2) then
  823. jirp = 0
  824. do iirp = 1,kgrp
  825. if (ipt1.num(iirp,1).eq.ip1) jirp = iirp
  826. enddo
  827. c recopie
  828. kgrp = kgrp + 1
  829. if (jirp.ne.0) then
  830. des2.noelep(kgrp) = des2.noelep(jirp)
  831. des2.noeled(kgrp) = des2.noeled(jirp)
  832. else
  833. kirp = kirp + 1
  834. ipt1.num(kirp,1) = ip1
  835. des2.noelep(kgrp) = kirp
  836. des2.noeled(kgrp) = kirp
  837. endif
  838. des2.lisinc(kgrp) = lisinc(igrp)
  839. des2.lisdua(kgrp) = lisdua(igrp)
  840. re(1,kgrp,1) = xmatr1.re(1,igrp,ie2)
  841. re(kgrp,1,1) = re(1,kgrp,1)
  842. endif
  843. *
  844. enddo
  845.  
  846. do kpmo = 1,nmapmo
  847. if (ABS(coepmo(kpmo)).gt.xlopre*xmaut1) then
  848. kirp = kirp + 1
  849. kgrp = kgrp + 1
  850. ipt1.num(kirp,1) = lmapmo(kpmo)
  851. des2.noelep(kgrp) = kirp
  852. des2.noeled(kgrp) = kirp
  853. des2.lisinc(kgrp) = compmo(kpmo)
  854. if (compmo(kpmo).eq.'ALFA') des2.lisdua(kgrp) = 'FALF'
  855. if (compmo(kpmo).eq.'BETA') des2.lisdua(kgrp) = 'FBET'
  856. re(1,kgrp,1) = coepmo(kpmo)
  857. re(kgrp,1,1) = re(1,kgrp,1)
  858. endif
  859. enddo
  860. *
  861. lirl = .false.
  862. if (kirp.ne.num(/1)) then
  863. lirl = .true.
  864. else
  865. do io = 1,kirp
  866. if (num(io,iele).ne.ipt1.num(io,1)) lirl=.true.
  867. enddo
  868. endif
  869. c creation d'un irigel
  870. if (lirl) then
  871. kige = kige + 1
  872. if (kige.gt.kige1) then
  873. nrigel = kige1 + 100
  874. segadj ri2
  875. kige1 = nrigel
  876. endif
  877. nbelem = 1
  878. nbnn = kirp
  879. segini ipt3
  880. ipt3.itypel = itypel
  881. do io =1,nbnn
  882. ipt3.num(io,1) = ipt1.num(io,1)
  883. enddo
  884. nligrp = kgrp
  885. nligrd = kgrp
  886. nelrig=1
  887. segadj xmatri,des2
  888. * segini imatr3
  889. * imatr3.imattt(1) = xmatri
  890. segdes ipt3,des2,xmatri
  891. RI2.IRIGEL(1,kige) = IPT3
  892. RI2.IRIGEL(3,kige) = DES2
  893. RI2.IRIGEL(4,kige) = xmatri
  894. RI2.IRIGEL(2,kige) = 0
  895. RI2.IRIGEL(5,kige) = irigel(5,ire)
  896. RI2.IRIGEL(6,kige) = irigel(6,ire)
  897. ri2.coerig(kige) = coerig(ire)
  898. else
  899. * relation non modifiee pour cet element
  900. kele = kele + 1
  901. do ig = 1,nligrp0
  902. ipt2.num(ig,kele) = ipt1.num(ig,1)
  903. enddo
  904. * imatr2.imattt(kele) = xmatr1
  905. * kich : a tester
  906. do ju = 1,kgrp
  907. xmatr2.re(1,ju,kele) = re(1,ju,1)
  908. xmatr2.re(ju,1,kele) = re(ju,1,1)
  909. enddo
  910. segsup xmatri,des2
  911. endif
  912. ENDDO
  913.  
  914. nbelem = kele
  915. nelrig = kele
  916. nligrd=xmatr2.re(/1)
  917. nligrp=xmatr2.re(/2)
  918. if (nbelem.gt.0) then
  919. segadj ipt2
  920. segadj xmatr2
  921. krigel = krigel + 1
  922. RI1.IRIGEL(1,krigel) = IPT2
  923. RI1.IRIGEL(3,krigel) = irigel(3,ire)
  924. RI1.IRIGEL(4,krigel) = xmatr2
  925. RI1.IRIGEL(2,krigel) = 0
  926. RI1.IRIGEL(5,krigel) = irigel(5,ire)
  927. RI1.IRIGEL(6,krigel) = irigel(6,ire)
  928. segdes ipt2,xmatr2
  929. else
  930. segsup ipt2
  931. endif
  932. segsup ipt1
  933. ENDDO
  934.  
  935. iriout = 0
  936. nrigel = krigel
  937. segadj ri1
  938. nrigel = kige
  939. segadj ri2
  940. segdes mrigid,ri1,ri2
  941. if (kige.eq.0) segsup ri2
  942. if (krigel.eq.0) segsup ri1
  943. if (kige.gt.0.and.krigel.gt.0) then
  944. c WRITE(6,*) 'fus', ri1,ri2,kige,krigel
  945. call fusrig(ri1,ri2,iriout)
  946. segsup ri1, ri2
  947. return
  948. endif
  949. if (kige.gt.0) iriout = ri2
  950. if (krigel.gt.0) iriout = ri1
  951. if (iriout.eq.0) call erreur(-5)
  952. c WRITE(6,*) 'iriout', iriout
  953.  
  954. 290 continue
  955. if (iriout.ne.0) iriout3 = iriout
  956. if (iriout1.ne.0) iriout3 = iriout1
  957. if (iriout.ne.0.and.iriout1.ne.0) then
  958. call fusrig(iriout, iriout1,iriout3)
  959. ri1 = iriout
  960. ri2 = iriout1
  961. segsup ri1,ri2
  962. endif
  963.  
  964. call ecrobj('RIGIDITE',iriout3)
  965. if (modcom) call ecrobj('RIGIDITE',iriout2)
  966.  
  967. goto 999
  968.  
  969. 199 continue
  970. segsup descr,meleme,mlchp1,mlchp2
  971. call erreur(5)
  972. return
  973.  
  974. 999 continue
  975.  
  976. if (plcf.ne.0) segsup plcf
  977.  
  978. END
  979.  
  980.  
  981.  

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