Télécharger la plaquette Cast3M


fr en

Télécharger pron.eso

Retour à la liste

Numérotation des lignes : 

pron

  1. C PRON      SOURCE    MB234859  25/09/08    21:16:01     12358          
  2.  
  3.       SUBROUTINE PRON(IPMODE,IPCHE,IPCHT,ISUP,IPOUT)
  4.  
  5.       IMPLICIT INTEGER(I-N)
  6.       IMPLICIT REAL*8 (A-H,O-Z)
  7.  
  8. *-------------------------------------------------------------------- *
  9. *                                                                     *
  10. *   Sous-programme appelé par  l'opérateur   PROI                     *
  11. *                                                                     *
  12. *   Projection d'un chamelem aux noeuds ou aux points d'integration   *
  13. *   des elements support du modele                                    *
  14. *   Evaluation  faite au point courant par utilisation des fonctions  *
  15. *   de formes du modele associe au champ d origine (lu en option      *
  16. *   dans PRO2) ou par defaut a partir des fonctions geometriques      *
  17. *   associées à l'élément d origine                                   *
  18. *                                                                     *
  19. *   __________________________________________________________________*
  20. *                                                                     *
  21. *   ENTREE :                                                          *
  22. *          IPMODE  modele sur lequel on va projeter                   *
  23. *          IPCHE   MCHAML de caracteristiques (si le modele comprend  *
  24. *                  des coques)                                        *
  25. *          IPCHT  MCHAML a projeter (obligatoirement defini aux noeuds)
  26. *          ISUP   support du resultat  1 noeuds                       *
  27. *                                      2 pts d integration contraintes*
  28. *   SORTIE :                                                          *
  29. *          IPOUT   MCHAML resultat                                    *
  30. *                                                                     *
  31. *   Projection d'un chamelem sur le support d un modele               *
  32. *                                                                     *
  33. *---------------------------------------------------------------------*
  34. *
  35. -INC PPARAM
  36. -INC CCOPTIO
  37. -INC CCREEL
  38. *
  39. -INC SMCOORD
  40. -INC SMMODEL
  41. -INC SMCHAML
  42. -INC SMELEME
  43. -INC SMCHPOI
  44. -INC SMLCHPO
  45. -INC SMLMOTS
  46. -INC TMTRAV
  47. -INC SMINTE
  48.         character*(LOCOMP) nomm
  49.         segment snomip
  50.           character*(LOCOMP) nomip(0)
  51.         endsegment
  52.         segment ipcb(nsschp)
  53.         segment  icpr(2,nbpts)
  54. C
  55.         segment sicoq
  56.         integer icocoq(nbsous),imfr(nbsous),iiele(nbsous)
  57.         integer iminte(nbsous),imint1(nbsous),ibnap(nbsous)
  58.         integer inbg(nbsous),ibbas(nbsous)
  59.         endsegment
  60. C
  61.       SEGMENT WRK4
  62.        REAL*8 BPSS(3,3),XEL(3,NBN1) ,XU(3),SHP(6,NBN1),S(6)
  63.        REAL*8 TXR(3,3,NBN1),XE(3,NBG),XEG(3,NBG3)
  64.       ENDSEGMENT
  65.  
  66.       SEGMENT VECT
  67.         REAL*8  VEC1(IDIM)
  68.         REAL*8  VEC2(IDIM)
  69.         REAL*8  VECN(IDIM,NBN1*2)
  70.       ENDSEGMENT
  71.  
  72.       segment phachm
  73.         integer iphach(nsouz)
  74.         integer mphach(nsouz)
  75.       endsegment
  76.       LOGICAL logach
  77.  
  78. c     Compteur du nombre de sous modeles de sure
  79.       ISURE = 0
  80.  
  81.       IPOUT = 0
  82.  
  83.       MMODEL=IPMODE
  84.       NSOUS=KMODEL(/1)
  85. *
  86.       if(IPCHE.NE.0) THEN
  87.         MCHELM=IPCHE
  88.         NCHAM=ICHAML(/1)
  89.       endif
  90. *
  91.       segact mcoord*mod
  92.       NBNOE = nbpts       
  93. *
  94. *     Création du maillage provisoire support pour interpolation
  95. *
  96. C ipt2  contiendra les super couche au point de gauss
  97.       NBSOUS = NSOUS
  98.       NBREF = 0
  99.       NBELEM = 0
  100.       NBNN =  0
  101.       SEGINI IPT2
  102.       SEGINI SICOQ
  103. c listmots des phases
  104.       ilphmo = -1
  105.       jgn = 8
  106.       jgm = nsous
  107.       segini mlmots
  108.       ilphmo = mlmots
  109.       jgm = 1
  110. *
  111. * Boucle sur l'ensemble des sous zones du modeles sur lequel on va
  112. *   projetter   pour   initialisations
  113.       nponou=0
  114.       npcoq =0
  115. *
  116.       DO 30 ISOUS = 1, NSOUS
  117. *
  118.         IMODEL = KMODEL(ISOUS)
  119. *
  120.         MELE = NEFMOD
  121.         MELEME = IMAMOD
  122.         NBELE1= NUM(/2)
  123.         NBN1 = NUM(/1)
  124.       IF(INFMOD(/1).NE.0)THEN
  125.       IF(INFMOD(1).NE.0)THEN
  126.           NPINT=INFMOD(1)
  127.           IBNAP(isous) =NPINT
  128.       ELSE
  129.           NPINT=0
  130.           ibnap(isous)=1
  131.       ENDIF
  132.       ELSE
  133.           NPINT=0
  134.           ibnap(isous)=1
  135.       ENDIF
  136. C
  137.       mfr  =infele(13)
  138.       iele =infele(14)
  139.       minte=infmod(2+isup)
  140.       minte1=infmod(3)
  141.       nbg=nbn1
  142.       if(isup.eq.2) then
  143.       nbg=shptot(/3)
  144.       elseif(isup.eq.3) then
  145.       nbg=infele(6)
  146.       elseif(isup.eq.4) then
  147.       nbg=infele(3)
  148.       elseif(isup.eq.5) then
  149.       nbg=infele(4)
  150.       endif
  151. C
  152.       imfr(isous)=mfr
  153.       iiele(isous)=iele
  154.       iminte(isous)=minte
  155.       imint1(isous)=minte1
  156.       inbg(isous)=nbg
  157.       NBSOUS = 0
  158.       NBREF = 0
  159.       NBELEM = NBELE1
  160.       IF (mfr.EQ.3.OR.mfr.EQ.5.OR.mfr.EQ.9) THEN
  161. C
  162.       if(IPCHE.eq.0) then
  163.         call erreur(404)
  164.         return
  165.       endif
  166. C
  167. C----------------------------coque ---------------------
  168.         icocoq(isous) = 1
  169.         npcoq = npcoq + 1
  170. C  nombre de points par nappes
  171.         if(isup.GT.2) then
  172.                 if(mfr.eq.5) then
  173. C -----------------      cas des  coques epaisses
  174.                         if(nbg.eq.8.or.nbg.eq.6) then
  175.                         IBNAP(isous)=2
  176.                         ibbas(isous)=nbg/2
  177.                         else
  178.                         IBNAP(isous)=nbg/nbn1*2
  179.                         endif
  180.                 else
  181.                 ibbas(isous) = nbg/ibnap(isous)
  182.                 endif
  183.         elseif(isup.eq.1) then
  184.         ibnap(isous)=1
  185.         ibbas(isous) = nbn1
  186.         elseif(isup.eq.2) then
  187.         ibnap(isous)= 1
  188.         ibbas(isous) = 1
  189.         endif
  190. *
  191. *      write(6,*) 'nombre de noeuds par nappe ',
  192. *     &       ibbas(isous),' nappes ',ibnap(isous)
  193. *
  194. * création du maillage  des points supports pour cette sz
  195. *
  196.         NBNN = ibbas(isous)*ibnap(isous)
  197.         nponou = nponou + nbnn*nbele1
  198.         SEGINI IPT1
  199.         IPT1.ITYPEL = 28
  200.         ipt2.lisous(isous) = ipt1
  201. *write(6,*) 'sz ',isous ,'nbnn nbelem  ipt1' ,nbnn,nbelem,ipt1
  202.       ELSE
  203. C-----------------------------massif ------------------------
  204.         if(isup.GT.2) then
  205.         NBNN = NBG
  206.         elseif(isup.eq.2) then
  207.         NBNN=1
  208.         else
  209.         NBNN= nbn1
  210.         endif
  211. C
  212.         SEGINI IPT1
  213.         if(isup.eq.2) then
  214.         else
  215.         IPT1.ITYPEL = 28
  216.         endif
  217.         ipt2.lisous(isous) = ipt1
  218.         nponou=nponou+nbnn*nbele1
  219. C
  220.       ENDIF
  221.  
  222.         if (isous.eq.1) then
  223.          mots(1) = conmod(17:24)
  224.         else
  225.          do ipl = 1,jgm
  226.           if (mots(ipl).eq.conmod(17:24)) goto 27
  227.          enddo
  228.          jgm = jgm + 1
  229.          mots(jgm) = conmod(17:24)
  230.  27      continue
  231.         endif
  232. C
  233. 30    CONTINUE
  234.  
  235.       segadj mlmots
  236. C-----------------------------------------------------------
  237. C
  238. C     si il y a des coques ou si l'on veut le champ aux points
  239. C  d integration, on va fabriquer des points support provisoires
  240. C
  241. C----------------------------------------------------------
  242. C         write(6,*) ' nombre de points crees ' ,nponou
  243. C (fdp) ajustement de MCOORD pour les noeuds provisoires
  244. C       il sera remis a sa taille initiale a la fin du programme
  245.          NBPTS = NBNOE + nponou
  246.          SEGADJ MCOORD
  247. *-----------------------------------------------------------
  248. *     Boucle sur l'ensemble des sous zones du modeles
  249. *------------------------------------------------------------
  250.         inupo =NBNOE
  251.       DO 100 ISOUS = 1,NSOUS
  252.         IMODEL = KMODEL(ISOUS)
  253.         MELEME = IMAMOD
  254.         IF (itypel.eq.48) then
  255.             ISURE = ISURE+1
  256.             goto 100
  257.         ENDIF
  258.         nbele1= num(/2)
  259.         nbn1 = num(/1)
  260.         ipt1=ipt2.lisous(isous)
  261.         segact ipt1*mod
  262.         NNB = ibbas(isous)
  263.         NBG =inbg(isous)
  264.         NBG3=NBG*3
  265.         SEGINI WRK4
  266.         SEGINI VECT
  267.       IF(ICOCOQ(ISOUS).NE.0) THEN
  268. C-----------------------------------------------------------
  269. C  -------------------------------LA SZ EST UNE COQUE ------
  270. C-----------------------------------------------------------
  271.            NUCHA = 0
  272.        DO 15, NUCH = 1, NCHAM
  273. *
  274.             IF ( CONCHE(NUCH).EQ.CONMOD.AND.
  275.      &           IMACHE(NUCH).EQ.IMAMOD) NUCHA = NUCH
  276. *
  277.   15      CONTINUE
  278. *
  279.           IF (NUCHA.NE.0) THEN
  280.             MCHAML=ICHAML(NUCHA)
  281. *
  282.             MELVA1 = 0
  283.             MELVA2 = 0
  284.             NCOMP = IELVAL(/1)
  285.             DO 20, I = 1, NCOMP
  286.               IF (NOMCHE(I).EQ.'EPAI') THEN
  287.                 MELVA1 = IELVAL(I)
  288.               ELSEIF (NOMCHE(I).EQ.'EXCE') THEN
  289.                 MELVA2 = IELVAL(I)
  290.               ENDIF
  291.   20        CONTINUE
  292.           ENDIF
  293. C
  294.       MELE = NEFMOD
  295. c   coque integree ou pas ?
  296.         npint=infmod(1)
  297.         nbnap=ibnap(isous)
  298.         mfr  =imfr(isous)
  299.         iele =iiele(isous)
  300.         minte=iminte(isous)
  301.         minte1=imint1(isous)
  302.         nbg = inbg(isous)
  303.         nnb = ibbas(isous)
  304. *-------------------------------------------------------------------
  305. *  Boucle sur les elements de la sous zone du modele  (c est une coque)
  306. *-------------------------------------------------------------------
  307.         DO 95 iel=1, NBELE1
  308.         ipo=0
  309.         call doxe(xcoor,idim,nbn1,num,iel,xel)
  310. *    write(6,*) ' coordonnees des noeuds '
  311. *    write(6,2000) ((xel(ik,jk),ik=1,idim),jk=1,nbn1)
  312. *    write(6,*) '-------------------- '
  313. 2000     format(3E12.5)
  314. * on veut les coordonnees aux point d integration et normales
  315. *
  316.       IF(ISUP.GT.1) THEN
  317. C  --- projection des points d'integration sur la surface moyenne--
  318.       DO 601 ip=1,nnb
  319.       do 600 ik=1,idim
  320.       xe(ik,ip)=0.D0
  321.       DO 602 jk = 1,NBN1
  322.       XE(Ik,ip) = XE(Ik,ip) + SHPTOT(1,JK,ip)*(XEL(IK,jk))
  323. 602   CONTINUE
  324. 600   CONTINUE
  325. 601   CONTINUE
  326.       do i=1,nnb
  327.       do j=1,idim
  328.       xeg(j,i)=xe(j,i)
  329.       enddo
  330.       enddo
  331. C      ELSEIF(isup.eq.2) THEN
  332. C----------------------------          on travaille aux  CDG
  333. C      do  ik=1,idim
  334. C      xeg(ik,1)=0.D0
  335. C      do  ip=1,nbg
  336. C      DO  jk = 1,NBN1
  337. C      XEG(Ik,1) = XEG(Ik,jk) + SHPTOT(1,JK,ip)*(XEL(IK,jk))
  338. C      enddo
  339. C      enddo
  340. C      enddo
  341.       ELSE
  342. C  -------   on travaille aux noeuds
  343.       do i=1,nbn1
  344.       do j=1,idim
  345.       xeg(j,i)=xel(j,i)
  346.       enddo
  347.       enddo
  348.       ENDIF
  349. C      write(6,*) '  coordonnes des points de  base '
  350. C     write(6,2576) (ip,(xeg(il,ip),il=1,3),ip=1,nnb)
  351. C      write(6,*) ' ----------------- '
  352. 2576  format(I4,2X,3E12.5)
  353.  
  354. C----------- les normales aux points ad hoc ------------------
  355. C
  356. *---------------    coque minces   3 noeuds ---------------
  357.        if(iele.eq.4)  then
  358. C    --------coq3   dkt  dst   ----------------
  359.       do 98 i=1,idim
  360.       VEC1(i)=xel(i,2)-xel(i,1)
  361.       VEC2(i)=xel(i,3)-xel(i,1)
  362. 98    continue
  363.       vecn(1,1) = vec1(2)*vec2(3)-vec1(3)*vec2(2)
  364.       vecn(2,1) = vec1(3)*vec2(1)-vec1(1)*vec2(3)
  365.       vecn(3,1) = vec1(1)*vec2(2)-vec1(2)*vec2(1)
  366.       vnor =sqrt(vecn(1,1)*vecn(1,1)+vecn(2,1)*vecn(2,1)
  367.      & +vecn(3,1)*vecn(3,1))
  368.       do 91 i=1,idim
  369.       vecn(i,1)=vecn(i,1)/vnor
  370. 91    continue
  371.       do 999 j=2,nbn1
  372.       do 99 i=1,idim
  373.       vecn(i,j)=vecn(i,1)
  374.  99   continue
  375.  999  continue
  376. C  petite rectif necessaire pour le tri3
  377.       if(isup.eq.2) then
  378.        do j=1,idim
  379.        xeg(j,1) = XZERO
  380.                do ip=1,3
  381.                xeg(j,1) = xeg(j,1) + xel(j,ip)
  382.                enddo
  383.        xeg(j,1)=xeg(j,1)/3.D0
  384.        enddo
  385.        endif
  386.        elseif (iele.eq.8) then
  387. C    --------------- coq4 ----------------
  388.       call cq4loc(xel,xe,BPSS,irrt,1)
  389. C  ici c est une estimation du plan moyen
  390. *     write(6,fmt='(3E12.5)') ((bpss(i,j),i=1,3),j=1,3)
  391.       do ip=1,ibbas(isous)
  392.       do i=1,3
  393.       vecn(i,ip)=bpss(3,i)
  394.       enddo
  395.       enddo
  396.         elseif (iele.eq.2) then
  397. C   ---------------- coq2 --------------
  398.         vnor=0.D0
  399.         do 92 i=1,idim
  400.         VEC1(i)=xel(i,2)-xel(i,1)
  401.         vnor=vnor + vec1(i)*vec1(i)
  402. 92    continue
  403.         vnor = sqrt(vnor)
  404.         vecn(1,1)=-vec1(2)/vnor
  405.         vecn(2,1)=vec1(1)/vnor
  406.         vecn(1,2)=vecn(1,1)
  407.         vecn(2,2)=vecn(2,1)
  408. C
  409.         elseif (iele.eq.5.or.iele.eq.10) then
  410. *-----------    coques epaisses  coq8 coq6
  411.        if(isup.eq.1) then
  412.       call cq8loc(xel,nbn1,MINTE1.SHPTOT,TXR,irr)
  413.        else
  414.       call cq8loc(xel,nbn1,SHPTOT,TXR,irr)
  415.        endif
  416.  
  417. C      write(6,*) ' normales aux points ad hoc  coques epaisses'
  418.         do i=1,nbn1
  419.         do k=1,3
  420.         vecn(k,i)=txr(k,3,i)
  421.         enddo
  422.         enddo
  423. 2222  format(I4,3E12.5)
  424. 2001  format(3E10.3)
  425.         endif
  426. C--------------------------------------------------
  427. *     write(6,*)  ' normales   a l element  '
  428. *     write(6,fmt='(3E12.5)') ((vecn(l,ip),l=1,3),ip=1,nbn1)
  429. C--------------------------------------------------
  430. C
  431.         do 604 ipoi=1,nnb
  432.         if(melva1.ne.0) then
  433.         ibmn =MIN(iel,melva1.velche(/2))
  434.         igmn =MIN(ipoi,melva1.velche(/1))
  435.         epai = melva1.velche(igmn,ibmn)
  436.         else
  437.         epai=0.D0
  438.         endif
  439. C
  440.         if(melva2.ne.0) then
  441.         ibmn =MIN(iel,melva2.velche(/2))
  442.         igmn =MIN(ipoi,melva2.velche(/1))
  443.         exce = melva2.velche(igmn,ibmn)
  444.         else
  445.         exce =0.D0
  446.         endif
  447. C  si on projette aux points d integration  calcul des coordonnes
  448. C                sur l element reel
  449. *  on range les points sur les normales et non   pas par nappes
  450. *  pour ne pas chercher n fois exce et epai attention au MCHAML
  451. *
  452.         do 603 inap =1,ibnap(isous)
  453.         ipo=ipo+1
  454.         inupo = inupo+1
  455.                 if(ibnap(isous).gt.1)  then
  456. C                if(npint.ne.0.or.mfr.eq.5) then
  457.                 igau=(inap-1)*inbg(isous)/ibnap(isous)+1
  458.                 zzz = dzegau(igau)
  459.                 else
  460. C               zzz = (1.D0*inap-2.D0)
  461.                 zzz = XZERO
  462.                 endif
  463.         dnor=exce+epai*zzz/2.D0
  464.         ipt1.num(ipo,iel) =inupo
  465.         do 97 i=1,idim
  466.         XCOOR((inupo-1)*(IDIM+1)+i) = xeg(i,ipoi)+vecn(i,ipoi)*dnor
  467. 97    continue
  468. C       write(6,2003) iel,ipoi,inupo-nbnoe,
  469. C     & ( XCOOR((inupo-1)*(IDIM+1)+i),i=1,3),dnor,igau
  470. 603     continue
  471. 604     continue
  472. C        write(6,*) 'ipt1 iel ipo nbp crees ' ,
  473. C     & ipt1,iel,inap,(inupo-nbnoe)
  474.  
  475. 95     CONTINUE
  476. 2003   format(3I4,2X,4e12.5,I4)
  477. C
  478. C
  479.       ELSEIF(ICOCOQ(ISOUS).EQ.0) THEN
  480. C---------------------------------------------------------
  481. C ---------------------------- LA SZ  EST UN MASSIF ------
  482. C---------------------------------------------------------
  483.         IF(ISUP.GT.1) THEN
  484. C ----------------     on travaille aux points d integration
  485.         minte=iminte(isous)
  486. 2234    format(6e12.5)
  487. C
  488.       DO 195 iel=1,nbele1
  489.        call doxe(xcoor,idim,nbn1,num,iel,xel)
  490. C  ------ PROJECTION
  491.       DO 700 ip=1,shptot(/3)
  492.                 DO ik=1,idim
  493.                 XU(ik)=0.D0
  494.                 DO  jk = 1,NBN1
  495.                 XU(Ik) = XU(Ik) + SHPTOT(1,JK,ip)*(XEL(IK,jk))
  496.                 ENDDO
  497.                 ENDDO
  498.        inupo=inupo+1
  499.         ipt1.num(ip,iel) = inupo
  500.                 do  i=1,idim
  501.                 XCOOR((inupo-1)*(IDIM+1)+i) = xu(i)
  502.                 enddo
  503. C       write(6,2223) inupo,(xu(i),i=1,idim)
  504. 700   CONTINUE
  505. 195   CONTINUE
  506. 2223  format(I6,2X,3E12.5)
  507.        ELSE
  508. C     on travaille aux noeuds
  509.        do 197 iel=1,nbele1
  510.        call doxe(xcoor,idim,nbn1,num,iel,xel)
  511.                 do j=1,nbn1
  512.                 inupo=inupo+1
  513.                 ipt1.num(j,iel) = inupo
  514.                 do i=1,idim
  515.                 XCOOR((inupo-1)*(IDIM+1)+i) = xel(i,j)
  516.                 enddo
  517.                 enddo
  518. 197    continue
  519.        ENDIF
  520. C    fin de  l'alternative   coques ou massifs
  521.       ENDIF
  522. C--------------------------------
  523.       SEGSUP WRK4,VECT
  524. *
  525.  100  CONTINUE
  526. C   fin de creation des noeuds supports provisoires
  527. C        ENDIF
  528. C--------------//////////////////////////////////------------------
  529. C  maintenant il faut un meleme poi1 de tous les points sur lesquels
  530. C          on doit interpoler une valeur
  531.          segini,ipt4=ipt2
  532.         ipgeom=ipt4
  533.         if( isup.ne.2) then
  534.         call change(ipgeom,1)
  535.         segsup ipt4
  536.         endif
  537. C       call ecrobj('MAILLAGE',ipgeom)
  538. C----------------------------------------------------------------
  539. C----------------------------------------------------------------
  540. *-------------on est pret a faire l interpolation sur ipgeom ----
  541. C
  542. *write(6,*) ' maillage apres changer poi1 ' ,ipgeom
  543. C       call ecmail(ipgeom)
  544.         isort=0
  545.  
  546. *        write(6,*) ' dans pron appel a pro2'
  547.           CALL PRO2(IPGEOM,IPCHT,isort,IPOUT,ilphmo)
  548. *      write(6,*) ' dans pron sortie de pro2',nsous,ilphmo
  549.          if (ierr.ne.0) return
  550. C----------------------------------------------------------------
  551. C----------------------------------------------------------------
  552. C----------------------------------------------------------------
  553. C  ---      il faut maintenant reconstituer les MCHAML
  554. C  ---       a partir du chpo  construit  sur ipgeom ----------------
  555. C------------ INITIALISATION   du MCHAML RESULTAT ------------------
  556.        mlchpo = ipout
  557.        nsouz  = nsous
  558.        segini phachm
  559.         N1=NSOUS - ISURE
  560.         L1=12
  561.         N3=6
  562.         SEGINI MCHEL1
  563.         MCHEL1.TITCHE='SCALAIRE'
  564.         MCHEL1.IFOCHE=IFOUR
  565.  
  566. * boucle phases modele
  567.        DO 7000 IPHAS = 1,JGM
  568. *         write(6,*) 'iphas ',jgm,iphas, ' mots ',mots(iphas),ichpoi(/1)
  569.         segini icpr,snomip
  570. *        mdata=isort
  571.         mchpoi=ichpoi(iphas)
  572.          if (mchpoi.le.0) goto 7000
  573.         nsschp=ipchp(/1)
  574.         segini ipcb
  575. C
  576.         do i=1,nsschp
  577. *        mchpoi=ipca(i)
  578. C       call ecchpo( mchpoi)
  579.         msoupo=ipchp(i)
  580.         ipcb(i)=msoupo
  581. C       write(6,*) i,ipca(i),msoupo
  582.         mpoval=ipoval
  583.         nc = nocomp(/2)
  584. C
  585. C tableau des general des composantes
  586.         do ic = 1,nc
  587.         nomm= nocomp(ic)
  588.         if(nomip(/2).eq.0) then
  589.           nomip(**)=nomm
  590.         else
  591.           do  k=1,nomip(/2)
  592.           if(nomm.eq.nomip(k)) goto 4322
  593.         enddo
  594.         nomip(**)=nomm
  595. 4322    continue
  596.         endif
  597.         enddo
  598. *       write(6,*) (nomip(l),l=1,nomip(/2))
  599. *    reperage de la position  ses points dans le chpo
  600.         ipt5=igeoc
  601.         do j=1,ipt5.num(/2)
  602.         icpr(1,ipt5.num(1,j))=j
  603.         icpr(2,ipt5.num(1,j))= i
  604.         enddo
  605.         enddo
  606. C  i=icpr(2,k)   point k venant du msoupo du msoupo i  dans ipcb
  607.  
  608.         kphach = 0
  609.         inupo=0
  610.         DO 2100 isous = 1,nsous
  611.         imodel= KMODEL(isous)
  612.          if (conmod(17:24).ne.mots(iphas)) goto 2100
  613.         meleme = imamod
  614.          if (itypel.eq.48) goto 2100; 
  615. C
  616.         ipt5=ipt2.lisous(isous)
  617.         N2 = nomip(/2)
  618.         mfr=imfr(isous)
  619. * pour une phase donnée ne cree le mchaml qu une fois par maillage
  620. * distinct
  621.  
  622.         if (kphach.gt.0) then
  623.           do kpha = 1,kphach
  624.             if (imamod.eq.iphach(kpha)) then
  625.               mchaml = mphach(kpha)
  626.               logach = .true.
  627.               goto 2105
  628.             endif
  629.           enddo
  630. *             création du nouveau segment MCHAML
  631.               SEGINI MCHAML
  632.               kphach = kphach + 1
  633.               iphach(kphach) = imamod
  634.               mphach(kphach) = mchaml
  635.               logach = .false.
  636.         else
  637. *             création du nouveau segment MCHAML
  638.               SEGINI MCHAML
  639.               kphach = kphach + 1
  640.               iphach(kphach) = imamod
  641.               mphach(kphach) = mchaml
  642.               logach = .false.
  643.         endif
  644.  
  645.  2105         MCHEL1.IMACHE(isous)=MELEME
  646.               MCHEL1.ICHAML(isous)=MCHAML
  647.               MCHEL1.CONCHE(isous)=CONMOD
  648.               MCHEL1.INFCHE(isous,1)=0
  649.               MCHEL1.INFCHE(isous,2)=0
  650.               MCHEL1.INFCHE(isous,3)=NIFOUR
  651.           IF(isup.eq.1) then
  652.               MCHEL1.INFCHE(isous,4)=0
  653.               MCHEL1.INFCHE(isous,5)=0
  654.               MCHEL1.INFCHE(isous,6)=1
  655.           ELSE
  656.               MCHEL1.INFCHE(isous,4)=iminte(isous)
  657.               MCHEL1.INFCHE(isous,5)=0
  658.               MCHEL1.INFCHE(isous,6)=isup
  659.           ENDIF
  660. *
  661.         if (logach) goto 2100
  662.  
  663.         N1EL=NUM(/2)
  664.         N2PTEL=0
  665.         N2EL=0
  666. C
  667.         IF(ICOCOQ(ISOUS).EQ.1) THEN
  668. C--------------------------------------------------------
  669. C       boucle sur les composantes   et c est une coque
  670. C--------------------------------------------------------
  671.         inocom=0
  672.         DO  2200 icomp=1,nomip(/2)
  673.         inocom=inocom+1
  674. *write(6,*) 'composante ', nomip(icomp)
  675.         idebco = inupo
  676.  
  677. C         N1PTEL= NUM(/1)*nbnap
  678.           N1PTEL= nnb*nbnap
  679.           SEGINI MELVAL
  680.           IELVAL(inocom)= MELVAL
  681.           NOMCHE(inocom)= NOCOMP(icomp)
  682.           TYPCHE(inocom)= 'REAL*8'
  683. C
  684.           DO 161 NUEL=1,num(/2)
  685.           DO 162 NUPT=1,nnb
  686.              DO 172 IPOS = 1,nbnap
  687.              ipop = (ipos-1)*nnb+nupt
  688.              inupo = ipt5.num(ipop,nuel)
  689.        jh=icpr (1,inupo)
  690.        ipa=icpr(2,inupo)
  691.        if(jh.eq.0) go to 172
  692. C  il faut verifier si vpocha existe pour ce point
  693.        msoupo = ipcb(ipa)
  694.        mpoval=ipoval
  695.        do l=1,nocomp(/2)
  696.        if(nocomp(l).eq.nomip(icomp)) then
  697.        vvv=vpocha(jh,icomp)
  698.        goto  4557
  699.        endif
  700.        enddo
  701.        vvv=XZERO
  702. 4557   continue
  703. C       write(6,2004) NUEL,inupo,jh,ipop,
  704. C     & ( XCOOR((inupo+nbnoe-1)*(IDIM+1)+i),i=1,3),vvv
  705.               VELCHE(ipop,NUEL) = vvv
  706.  172      CONTINUE
  707.  162      CONTINUE
  708.  161      CONTINUE
  709.  
  710.         if(icomp.lt.nocomp(/2)) inupo =idebco
  711. 2200    continue
  712.  
  713.         ELSE
  714. C--------------------------------------------------------
  715. C       boucle sur les composantes   et c est un massif
  716. C--------------------------------------------------------
  717.  
  718.         if(isup.eq.1) then
  719.         N1PTEL=NUM(/1)
  720.         elseif(isup.eq.2) then
  721.         N1PTEL=1
  722.         else
  723.         N1PTEL= inbg(isous)
  724.         endif
  725. C
  726.         DO  2220 icomp=1,NOMIP(/2)
  727. *        write(6,*) 'composante ', nomip(icomp)
  728.         idebco = inupo
  729. C
  730.  
  731.         SEGINI MELVAL
  732.         NOMCHE(icomp)=nomip(icomp)
  733.         TYPCHE(icomp)='REAL*8'
  734.         IELVAL(icomp)=MELVAL
  735. C
  736.           DO 164 NUEL=1,N1EL
  737.           DO 163 NUPT=1,N1PTEL
  738.               inupo=ipt5.num(nupt,nuel)
  739.        jh=icpr(1,inupo)
  740.        ipa=icpr(2,inupo)
  741.        if(jh.eq.0) go to 163
  742. C  il faut verifier si vpocha existe pour ce point
  743.        msoupo = ipcb(ipa)
  744.        mpoval=ipoval
  745.        do l=1,nocomp(/2)
  746.        if(nocomp(l).eq.nomip(icomp)) then
  747.        vvv=vpocha(jh,icomp)
  748.        goto  4558
  749.        endif
  750.        enddo
  751.        vvv=XZERO
  752. 4558   continue
  753.  
  754. C        write(6,2003) NUEL,inupo,jh,
  755. C     & ( XCOOR((inupo+nbnoe-1)*(IDIM+1)+i),i=1,3),vvv
  756.               VELCHE(NUPT,NUEL) = vvv
  757.  163      CONTINUE
  758.  164      CONTINUE
  759.         ipcham= mchaml
  760.         ipc1 = melval
  761.         call comred(ipc1)
  762.         IELVAL(icomp)=ipc1
  763.         melval = ipc1
  764.         if(icomp.lt.nocomp(/2)) inupo =idebco
  765. C  ----  fin de la boucle sur les composantes
  766. 2220    continue
  767.  
  768. C----------------  fin du traitement des massifs
  769.         ENDIF
  770.         segsup ipt5
  771. C    fin de la boucle sur les sous zones du modele
  772. 2100    continue
  773. C   destruction des chpo intermediaires
  774.         do i=1,ipcb(/1)
  775.           msoupo=ipcb(i)
  776. cgf          ipt5= igeoc (correction 7284)
  777.           mpoval= ipoval
  778. *          mchpoi=ipca(i)
  779.           segsup mpoval,msoupo
  780. cgf          segsup ipt5 (correction 7284)
  781.         enddo
  782.         segsup mchpoi
  783.         segsup ipcb,snomip
  784.         segsup icpr
  785.  7000   CONTINUE
  786. C (fdp) re-ajustement de MCOORD a sa taille initiale
  787.         NBPTS = NBNOE
  788.         SEGADJ MCOORD
  789. C       retrait des maillages temporaires du pre-conditionnement
  790. c       (leurs numero de noeuds depasse la taille de MCOOR)
  791.         call crech1b
  792. C (fdp) suppression du maillage temporaire IPT1
  793.         ipt1=ipgeom
  794.         segsup ipt1
  795.         segsup ipt2
  796.         segsup phachm
  797.         IPOUT=MCHEL1
  798.         segsup sicoq
  799.  
  800.         return
  801.       END
  802.  
  803.  
  804.  
  805.  
  806.  
  807.  
  808.  

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