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

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