Télécharger pjmode.eso

Retour à la liste

Numérotation des lignes :

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

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