Télécharger la plaquette Cast3M


fr en

Télécharger prot.eso

Retour à la liste

Numérotation des lignes : 

prot

  1. C PROT      SOURCE    CB215821  24/04/12    21:16:58     11897          
  2.       SUBROUTINE PROT(IPMODE,IPCHT,IPCHE,ITPR)
  3.       IMPLICIT INTEGER(I-N)
  4.       IMPLICIT REAL*8 (A-H,O-Z)
  5. *--------------------------------------------------------------------*
  6. *                                                                    *
  7. *   Sous-programme associé à l'opérateur   CALP                      *
  8. *   ____________________________________________                     *
  9. *                                                                    *
  10. *   Projection d'un chamelem de temperature sur une géometrie
  11. *                                              constituée de coques  *
  12. *                                                                    *
  13. *                                                                    *
  14. *   Auteur, date de création:                                        *
  15. *   -------------------------                                        *
  16. *                                                                    *
  17. *   Bruno VIGAN, le 26 février 1997.                                 *
  18. *                                                                    *
  19. *--------------------------------------------------------------------*
  20. *
  21.  
  22. -INC PPARAM
  23. -INC CCOPTIO
  24. -INC SMCOORD
  25. -INC SMMODEL
  26. -INC SMCHAML
  27. -INC SMELEME
  28. -INC SMCHPOI
  29. -INC SMLCHPO
  30. -INC SMLMOTS
  31. -INC TMTRAV
  32.  
  33.       SEGMENT VECT
  34.         REAL*8  VEC1(IDIM)
  35.         REAL*8  VEC2(IDIM)
  36.         REAL*8  VECN(IDIM)
  37.       ENDSEGMENT
  38.       SEGINI VECT
  39. *
  40.       SEGMENT ICPR(NBNOE,NCHAM)
  41.       SEGMENT NKON(IKOUR)
  42.       SEGMENT NUIN(IKOUR)
  43. *
  44.       SEGMENT ICARAC
  45.         REAL*8  XEPAI(NCHAM)
  46.         REAL*8  XEXCE(NCHAM)
  47.       ENDSEGMENT
  48.       SEGMENT NCARAC(NCHAM)
  49.  
  50.  
  51. *
  52.       MMODEL=IPMODE
  53.       SEGACT,MMODEL
  54.       NMOD=KMODEL(/1)
  55. *
  56.       MCHELM=IPCHE
  57.       SEGACT,MCHELM
  58.       NCHAM=ICHAML(/1)
  59. *
  60.       segact mcoord*mod
  61.       NBNOE = nbpts 
  62.       SEGINI ICPR
  63.       SEGINI ICARAC
  64.       SEGINI NCARAC
  65. *
  66.       DO 10, I = 1, NCHAM
  67.         ICARAC.XEPAI(I)= 0.
  68.         ICARAC.XEXCE(I)= 0.
  69.         DO  10, J=1, NBNOE
  70.           ICPR(J,I)=0
  71.         CONTINUE
  72.   10  CONTINUE
  73.       NBCAR = 0
  74. *
  75. *     Création du maillage principal
  76. *
  77.       NBSOUS = 0
  78.       NBREF = 0
  79.       NBELEM = 0
  80.       NBNN =  0
  81.       SEGINI IPT2
  82.       IKOUR=0
  83. c listmots des phases
  84.       ilphmo = -1
  85.       jgn = 8
  86.       jgm = nmod
  87.       segini mlmots
  88.       ilphmo = mlmots
  89.       jgm = 1
  90. *
  91. *     Boucle sur l'ensemble des sous zones du modeles
  92. *
  93.       DO 30, NUMO = 1, NMOD
  94. *
  95.         IMODEL = KMODEL(NUMO)
  96.         SEGACT,IMODEL
  97. *
  98. *       Test si le modele est une coque
  99. *
  100.         NUFOR = NUMMFR(NEFMOD)
  101.         IF (NUFOR.EQ.3 .OR. NUFOR.EQ.5 .OR. NUFOR.EQ.9)THEN
  102. *
  103. *         Recherche du chamemlem de caracteristique assossiée
  104. *
  105.           NUCHA = 0
  106.           DO 15, NUCH = 1, NCHAM
  107. *
  108.             IF ( CONCHE(NUCH).EQ.CONMOD.AND.
  109.      C           IMACHE(NUCH).EQ.IMAMOD) NUCHA = NUCH
  110. *
  111.   15      CONTINUE
  112. *
  113.           IF (NUCHA.NE.0) THEN
  114.             MCHAML=ICHAML(NUCHA)
  115.             SEGACT,MCHAML
  116. *
  117.             XEXCE1 = 0.
  118.             XEPAI1 = 0.
  119.             NCOMP = IELVAL(/1)
  120.             DO 20, I = 1, NCOMP
  121.               IF (NOMCHE(I).EQ.'EPAI') THEN
  122.                 MELVAL = IELVAL(I)
  123.                 SEGACT, MELVAL
  124.                 XEPAI1 = VELCHE(1,1)
  125.               ELSEIF (NOMCHE(I).EQ.'EXCE') THEN
  126.                 MELVAL = IELVAL(I)
  127.                 SEGACT, MELVAL
  128.                 XEXCE1 = VELCHE(1,1)
  129.               ENDIF
  130.   20        CONTINUE
  131. *
  132. *           recherche du numero de caracteristique associe
  133. *           a l'epaisseur et l'excentricitee
  134. *
  135.             NUCAR = 0
  136.             DO 22, I = 1, NBCAR
  137.               IF (ICARAC.XEPAI(I).EQ.XEPAI1.AND.
  138.      C            ICARAC.XEXCE(I).EQ.XEXCE1) NUCAR = I
  139.   22        CONTINUE
  140. *
  141.             IF (NUCAR.EQ.0) THEN
  142.               NUCAR = NBCAR+1
  143.               ICARAC.XEPAI(NUCAR)=XEPAI1
  144.               ICARAC.XEXCE(NUCAR)=XEXCE1
  145.               NBCAR = NUCAR
  146.             ENDIF
  147.             NCARAC(NUCHA)=NUCAR
  148. *
  149.             MELEME = IMAMOD
  150.             SEGACT MELEME
  151. *
  152. *           recherche du nombre de noeuds
  153. *
  154.             DO 25 I=1, NUM(/1)
  155.               DO 25 J=1, NUM(/2)
  156.                 ITH= NUM(I,J)
  157.                 IF (ICPR(ITH,NUCAR).EQ.0) THEN
  158.                   IKOUR=IKOUR+1
  159.                   ICPR(ITH,NUCAR)=IKOUR
  160.                 ENDIF
  161.   25        CONTINUE
  162.  
  163.           ENDIF
  164.         ENDIF
  165. *
  166.         if (numo.eq.1) then
  167.          mots(1) = conmod(17:24)
  168.         else
  169.          do ipl = 1,jgm
  170.           if (mots(ipl).eq.conmod(17:24)) goto 27
  171.          enddo
  172.          jgm = jgm + 1
  173.          mots(jgm) = conmod(17:24)
  174.  27      continue
  175.         endif
  176. C
  177.   30  CONTINUE
  178. *
  179.       segadj mlmots
  180. *     Augmentation du tableau de coordonnées
  181. *
  182.       NBPTS = NBNOE+3*IKOUR
  183.       SEGADJ MCOORD
  184. *
  185.       NNNO = IKOUR
  186.       SEGINI NKON
  187.       SEGINI NUIN
  188. *
  189.       DO 40, I = 1, NNNO
  190.         NKON(I)=0
  191.         DO 40, K = 1, IDIM
  192.           XCOOR((NBNOE+I-1)*(IDIM+1)+K) = 0.
  193.           XCOOR((NBNOE+I-1+IKOUR)*(IDIM+1)+K) = 0.
  194.           XCOOR((NBNOE+I-1+2*IKOUR)*(IDIM+1)+K) = 0.
  195.   40  CONTINUE
  196. *
  197. *     Boucle sur l'ensemble des sous zones du modeles
  198. *
  199.       DO 100, NUMO = 1, NMOD
  200.         IMODEL = KMODEL(NUMO)
  201. *
  202. *       Test si le modele est une coque
  203. *
  204.         NUFOR = NUMMFR(NEFMOD)
  205.         IF (NUFOR.EQ.3 .OR. NUFOR.EQ.5 .OR. NUFOR.EQ.9) THEN
  206. *
  207. *         Recherche du chamemlem de caracteristique assossiée
  208. *
  209.           NUCHA = 0
  210.           DO 50, NUCH = 1, NCHAM
  211. *
  212.             IF ( CONCHE(NUCH).EQ.CONMOD.AND.
  213.      C           IMACHE(NUCH).EQ.IMAMOD) NUCHA = NUCH
  214. *
  215.   50      CONTINUE
  216. *
  217.           IF (NUCHA.NE.0) THEN
  218. *
  219.             NUCAR = NCARAC(NUCHA)
  220.             MELEME = IMAMOD
  221. *
  222. *           création du nouveau maillage
  223. *
  224.             NBSOUS = 0
  225.             NBREF = 0
  226.             NBELE1 = NUM(/2)
  227.             NBELEM = 3* NBELE1
  228.             NBNN =  NUM(/1)
  229.  
  230.             SEGINI IPT1
  231.             IPT1.ITYPEL = ITYPEL
  232. *
  233.             DO 95 J=1, NBELE1
  234.               IPT1.ICOLOR(J) = ICOLOR(J)
  235.               IPT1.ICOLOR(J+NBELE1) = ICOLOR(J)
  236.               IPT1.ICOLOR(J+2*NBELE1) = ICOLOR(J)
  237. *
  238. *             Recherche d'une normale a l'element courant
  239. *
  240.               XNORM = 0.
  241.               DO 55, K = 1, IDIM
  242.                 VECN(K) = 0.
  243.   55          CONTINUE
  244.               IF (IDIM.EQ.2) THEN
  245.                 ICO1 = NUM(NBNN,J)
  246.                 ICO2 = NUM(1,J)
  247.                 DO 57, K = 1, IDIM
  248.                   VEC1(K) = XCOOR((ICO1-1)*(IDIM+1)+K)-
  249.      C                      XCOOR((ICO2-1)*(IDIM+1)+K)
  250.                   K1 = K+1
  251.                   IF (K1.GT.IDIM) K1 = 1
  252.                   VECN(K) =  VEC1(K1)*(-1)**K
  253.                   XNORM = XNORM +VECN(K)*VECN(K)
  254.  57             CONTINUE
  255.               ENDIF
  256.               IF (IDIM.EQ.3) THEN
  257.                 ICO1 = NUM(NBNN-1,J)
  258.                 ICO2 = NUM(NBNN,J)
  259. *
  260.                 DO 65 I=1, NBNN
  261.                   ICO3 = NUM(I,J)
  262.                   DO 60, K = 1, IDIM
  263.                     VEC1(K) = XCOOR((ICO1-1)*(IDIM+1)+K)-
  264.      C                        XCOOR((ICO2-1)*(IDIM+1)+K)
  265.                     VEC2(K) = XCOOR((ICO2-1)*(IDIM+1)+K)-
  266.      C                        XCOOR((ICO3-1)*(IDIM+1)+K)
  267.   60              CONTINUE
  268. *
  269.                   ICO1 = ICO2
  270.                   ICO2 = ICO3
  271.                   DO 65, K = 1, IDIM
  272.                     K1 = K+1
  273.                     K2 = K+2
  274.                     IF (K1.GT.IDIM) K1 = K1 - IDIM
  275.                     IF (K2.GT.IDIM) K2 = K2 - IDIM
  276.                     VECN(K) = VEC1(K1)*VEC2(K2) -VEC2(K1)*VEC1(K2)
  277.      C                      + VECN(K)
  278.                     IF (I.EQ.NBNN) XNORM = XNORM + VECN(K)*VECN(K)
  279.   65            CONTINUE
  280.               ENDIF
  281.               XNORM = SQRT(XNORM)
  282. *
  283.               DO 70, K = 1, IDIM
  284.                 VECN(K) = VECN(K)/XNORM
  285.   70          CONTINUE
  286.  
  287.               DO 95 I=1, NBNN
  288. *
  289.                 ICOU = NUM(I,J)
  290.                 IKOUR = ICPR(ICOU,NUCAR)
  291.                 NKON(IKOUR) = NKON(IKOUR)+1
  292.                 NUIN(IKOUR) = ICOU
  293.                 IPT1.NUM(I,J)= NBNOE+IKOUR
  294.                 IPT1.NUM(I,J+NBELE1)= NBNOE+IKOUR+NNNO
  295.                 IPT1.NUM(I,J+2*NBELE1)= NBNOE+IKOUR+2*NNNO
  296. *
  297. *               Calcul des coordonées des nouveaux points
  298. *
  299.                 DO 90, K = 1, IDIM
  300.                   XCOOR((IPT1.NUM(I,J)-1)*(IDIM+1)+K) =
  301.      C            XCOOR((IPT1.NUM(I,J)-1)*(IDIM+1)+K) +
  302.      C            VECN(K)*ICARAC.XEXCE(NUCAR)
  303.                   XCOOR((IPT1.NUM(I,J)+NNNO-1)*(IDIM+1)+K) =
  304.      C            XCOOR((IPT1.NUM(I,J)+NNNO-1)*(IDIM+1)+K) +
  305.      C            VECN(K)*(ICARAC.XEXCE(NUCAR)+ICARAC.XEPAI(NUCAR)/2)
  306.                   XCOOR((IPT1.NUM(I,J)+2*NNNO-1)*(IDIM+1)+K) =
  307.      C            XCOOR((IPT1.NUM(I,J)+2*NNNO-1)*(IDIM+1)+K) +
  308.      C            VECN(K)*(ICARAC.XEXCE(NUCAR)-ICARAC.XEPAI(NUCAR)/2)
  309.  
  310.   90            CONTINUE
  311.  
  312.   95        CONTINUE
  313. *
  314. *           Ajustement du pointeur maillage principal
  315. *
  316.             NBSOUS = IPT2.LISOUS(/1)+1
  317.             NBNN = 0
  318.             NBREF = 0
  319.             NBELEM = 0
  320.             SEGADJ IPT2
  321.             IPT2.LISOUS(NBSOUS) = IPT1
  322.           ENDIF
  323.         ENDIF
  324. *
  325.  100  CONTINUE
  326. *
  327.       DO 110 I=1, NNNO
  328.         DO 110, K=1, IDIM
  329.           XCOOR((NBNOE+I-1)*(IDIM+1)+K) =
  330.      C              XCOOR((NBNOE+I-1)*(IDIM+1)+K)/NKON(I) +
  331.      C              XCOOR((NUIN(I)-1)*(IDIM+1)+K)
  332.           XCOOR((NBNOE+I+NNNO-1)*(IDIM+1)+K) =
  333.      C              XCOOR((NBNOE+I+NNNO-1)*(IDIM+1)+K)/NKON(I) +
  334.      C              XCOOR((NUIN(I)-1)*(IDIM+1)+K)
  335.           XCOOR((NBNOE+I+2*NNNO-1)*(IDIM+1)+K) =
  336.      C              XCOOR((NBNOE+I+2*NNNO-1)*(IDIM+1)+K)/NKON(I) +
  337.      C              XCOOR((NUIN(I)-1)*(IDIM+1)+K)
  338.  110  CONTINUE
  339. *
  340.       SEGSUP ICARAC
  341.       SEGSUP NKON
  342.       SEGSUP NUIN
  343.       SEGSUP VECT
  344.  
  345.       NMAILL = IPT2.LISOUS(/1)
  346.       IF (NMAILL.GE.1) THEN
  347.         IF (NMAILL.EQ.1) THEN
  348.           IPT3 = IPT2.LISOUS(1)
  349.           SEGSUP IPT2
  350.           IPT2 = IPT3
  351.         ENDIF
  352. *
  353. *       appel a PRO2 pour projeter les temperature sur le maillage
  354. *       cree.
  355.          isort= 1
  356. *
  357.        CALL PRO2(IPT2,IPCHT,isort,IPOUT,ilphmo)
  358.          if (ierr.ne.0) return
  359. *
  360. *       Recopie des valeurs du champoint dans un Chamelem image
  361. *       de la geometrie initiale de la coque
  362. *
  363.        mlchpo = ipout
  364.        segact mlchpo
  365. * kich : pour la projection du champ de temperature on n attend qu une phase
  366.        if (ichpoi(/1).ne.1) call erreur(5)
  367.         MCHPOI = ICHPOI(1)
  368.         SEGACT MCHPOI
  369. *
  370. *       Creation du Chamelem
  371. *
  372.         N1 = NMAILL
  373.         N3 = 6
  374.         L1 = 12
  375.         SEGINI MCHEL1
  376.         MCHEL1.TITCHE='SCALAIRE'
  377.         MCHEL1.IFOCHE=IFOUR
  378.         NUCHAM = 0
  379. *
  380. *       Boucle sur l'ensemble des sous zones du modeles
  381. *
  382.         DO 200, NUMO = 1, NMOD
  383. *
  384.           IMODEL = KMODEL(NUMO)
  385.           SEGACT IMODEL
  386. *
  387. *         Test si le modele est une coque
  388. *
  389.           NUFOR = NUMMFR(NEFMOD)
  390.           IF (NUFOR.EQ.3 .OR. NUFOR.EQ.5 .OR. NUFOR.EQ.9) THEN
  391. *
  392. *           Recherche du chamemlem de caracteristique assossiée
  393. *
  394.             NUCHA = 0
  395.             DO 120, NUCH = 1, NCHAM
  396. *
  397.               IF ( CONCHE(NUCH).EQ.CONMOD.AND.
  398.      C             IMACHE(NUCH).EQ.IMAMOD) NUCHA = NUCH
  399. *
  400.  120        CONTINUE
  401. *
  402.             IF (NUCHA.NE.0) THEN
  403. *
  404.               NUCAR = NCARAC(NUCHA)
  405.               MELEME = IMAMOD
  406.               SEGACT MELEME
  407. *
  408. *             création du nouveau segment MCHAML
  409. *
  410.               N2 = 3
  411.               SEGINI MCHAML
  412.               NUCHAM = NUCHAM+1
  413.               MCHEL1.IMACHE(NUCHAM)=MELEME
  414.               MCHEL1.ICHAML(NUCHAM)=MCHAML
  415.               MCHEL1.CONCHE(NUCHAM)=CONMOD
  416.               MCHEL1.INFCHE(NUCHAM,1)=0
  417.               MCHEL1.INFCHE(NUCHAM,2)=0
  418.               MCHEL1.INFCHE(NUCHAM,3)=0
  419.               MCHEL1.INFCHE(NUCHAM,4)=0
  420.               MCHEL1.INFCHE(NUCHAM,5)=0
  421.               MCHEL1.INFCHE(NUCHAM,6)=1
  422. *
  423.               N1PTEL = NUM(/1)
  424.               N1EL = NUM(/2)
  425.               N2PTEL = 0
  426.               N2EL = 0
  427. *
  428.               DO 170, IPOS = 1, N2
  429. *
  430.                 SEGINI MELVAL
  431.                 IF (IPOS.EQ.1) THEN
  432.                   NOMCHE(IPOS) = 'T'
  433.                   IMUL = 0
  434.                 ELSEIF (IPOS.EQ.2) THEN
  435.                   NOMCHE(IPOS) = 'TSUP'
  436.                   IMUL = 1
  437.                 ELSEIF (IPOS.EQ.3) THEN
  438.                   NOMCHE(IPOS) = 'TINF'
  439.                   IMUL = 2
  440.                 ENDIF
  441.                 IELVAL(IPOS) =  MELVAL
  442.                 TYPCHE(IPOS) = 'REAL*8'
  443. *
  444.                 DO 160 NUEL=1, N1EL
  445. *
  446.                   DO 160 NUPT=1, N1PTEL
  447. *
  448.                     ICO3 = NUM(NUPT,NUEL)
  449.                     IKOUR = ICPR(ICO3,NUCAR)
  450. *
  451. *
  452. *                   Boucle sur les sous-zones du champoint
  453. *
  454.                     DO 150, I = 1, IPCHP(/1)
  455. *
  456.                       MSOUPO = IPCHP(I)
  457.                       SEGACT MSOUPO
  458.                       MPOVAL = IPOVAL
  459.                       SEGACT MPOVAL
  460.                       IPT1 = IGEOC
  461.                       SEGACT IPT1
  462. *
  463. *                     Boucle sur les composantes du champoint
  464. *
  465.                       DO 140, J = 1, NOCOMP(/2)
  466. *
  467.                         IF (NOCOMP(J).EQ.'T   ') THEN
  468. *
  469. *                         Boucle sur les points
  470. *
  471.                           DO  130, K = 1, IPT1.NUM(/2)
  472. *
  473. *                           Comparaison des numeros de points
  474. *                           entre le champoint et la geometrie creee
  475. *
  476.                             IF (IPT1.NUM(1,K).EQ.IKOUR+NBNOE+IMUL*NNNO)
  477.      C                      THEN
  478.                               VELCHE(NUPT,NUEL) = VPOCHA(K,J)
  479.                               GOTO 160
  480.                             ENDIF
  481. *
  482.  130                      CONTINUE
  483.                         ENDIF
  484.  140                  CONTINUE
  485.  150                CONTINUE
  486.  160            CONTINUE
  487.  170          CONTINUE
  488.             ENDIF
  489.           ENDIF
  490. *
  491.  200    CONTINUE
  492. *
  493. *       Suppression du champoint
  494. *
  495.         DO 220, I = 1, IPCHP(/1)
  496. *
  497.           MSOUPO = IPCHP(I)
  498.           MPOVAL = IPOVAL
  499.           IPT1 = IGEOC
  500. *****     SEGSUP IPT1
  501.           SEGSUP MPOVAL
  502.           SEGSUP MSOUPO
  503. *
  504.  220    CONTINUE
  505.         SEGSUP MCHPOI
  506. *
  507. *       Suppression du maillage intermediaire
  508. *
  509.         SEGACT IPT2
  510. *
  511.         DO 240, IOB =1, IPT2.LISOUS(/1)
  512. *
  513.           IPT1 = IPT2.LISOUS(IOB)
  514.           SEGSUP IPT1
  515. *
  516.  240    CONTINUE
  517. *****   SEGSUP IPT2
  518.  
  519. *
  520. *       Reajustement du tableau de coordonées
  521. *
  522.         NBPTS = NBNOE
  523.         SEGADJ MCOORD
  524. *
  525. *       RESTITUTION DU CHAMP DE SORTIE
  526. *
  527.         ITPR= MCHEL1
  528. *
  529.       ELSE
  530.         CALL ERREUR(704)
  531.       ENDIF
  532. *
  533.       SEGSUP ICPR
  534.       SEGSUP NCARAC
  535.       END
  536.  
  537.  
  538.  
  539.  
  540.  
  541.  
  542.  

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