Télécharger la plaquette Cast3M


fr en

Télécharger devfb4.eso

Retour à la liste

Numérotation des lignes : 

devfb4

  1. C DEVFB4    SOURCE    BP208322  20/09/18    21:15:29     10718          
  2.       SUBROUTINE DEVFB4(ITYP,FTOTB,XPTB,IPALB,IPLIB,XPALB,XVALB,NLIAB,
  3.      &                  NPLB,IND,IND1,INDM1,NPA,NPAM1,IND2,PDT,PDTS2,
  4.      &                  FEXPSM,NPC1,I,IERRD,iannul)
  5.       IMPLICIT INTEGER(I-N)
  6.       IMPLICIT REAL*8(A-H,O-Z)
  7. *--------------------------------------------------------------------*
  8. *                                                                    *
  9. *     Op{rateur DYNE : algorithme de Fu - de Vogelaere               *
  10. *     ________________________________________________               *
  11. *                                                                    *
  12. *     Calcul des forces de choc sur base B pour les liaisons         *
  13. *     de type PROFIL_PROFIL_...                                      *
  14. *                                                                    *
  15. *     Param}tres:                                                    *
  16. *                                                                    *
  17. * e   ITYP    type de la liaison.                                    *
  18. * es  FTOTB   Forces ext{rieures totalis{es sur la base B.           *
  19. * e   XPTB    Tableau des d{placements des points                    *
  20. * e   IPALB   Renseigne sur la liaison.                              *
  21. * e   IPLIB   Tableau contenant les num{ros "DYNE" de la liaison.    *
  22. * e   XPALB   Tableau contenant les param}tres de la liaison.        *
  23. * es  XVALB   Tableau contenant les variables internes de liaisons.  *
  24. * e   NLIAB   Nombre de liaisons sur la base B.                      *
  25. * e   NPLB    Nombre total de points intervenant dans les liaisons.  *
  26. * e   IND     Indice du pas.                                         *
  27. * e   I       num{ro de la liaison.                                  *
  28. * es  IERRD   num{ro d'erreur.                                       *
  29. *                                                                    *
  30. *                                                                    *
  31. *     Auteur, date de cr{ation:                                      *
  32. *                                                                    *
  33. *     Lionel VIVAN      : le 1 f{vrier 1990                          *
  34. *                                                                    *
  35. *--------------------------------------------------------------------*
  36. *
  37.       INTEGER IPALB(NLIAB,*),IPLIB(NLIAB,*)
  38.       REAL*8  XPALB(NLIAB,*),XPTB(NPLB,2,*),FTOTB(NPLB,*)
  39.       REAL*8  XVALB(NLIAB,4,*),XFOR(3),FEXPSM(NPLB,NPC1,2,*)
  40.       LOGICAL LPETPP
  41.       PARAMETER ( ZERO = 0.D0 , PRECIS = 1.D-15 )
  42. *
  43. * --- choc {l{mentaire PROFIL_PROFIL_INTERIEUR
  44. * --- choc {l{mentaire PROFIL_PROFIL_EXTERIEUR
  45. *
  46.       IF (ITYP.EQ.31 .OR. ITYP.EQ.32) THEN
  47.          NPOI   = IPLIB(I,1)
  48.          IDIM   = IPALB(I,3)
  49.          NOMBN1 = IPALB(I,4)
  50.          NOMBN2 = IPALB(I,5)
  51.          XRAID  = XPALB(I,1)
  52.          XSECT  = XPALB(I,2)
  53.          XEXPO  = XPALB(I,3)
  54.          ID1 = 3
  55.          ID2 = ID1 + IDIM
  56.          ID3 = ID1 + 2*IDIM
  57.          ID4 = ID1 + 3*IDIM
  58.          ID5 = ID1 + 4*IDIM
  59.          ID6 = ID1 + 5*IDIM
  60.          ID7 = ID6 + IDIM*NOMBN1
  61.          ID8 = ID7 + IDIM*NOMBN2
  62.          IP1 = 5
  63.          IP2 = IP1 + NOMBN1
  64.          IP3 = IP2
  65.          IP4 = IP1 + NOMBN1 + NOMBN1*NOMBN2
  66.          SOMFN = ZERO
  67. *
  68. *   Recherche du nombre de point qui ont travers{ le maillage fixe,
  69. *   donc s'il y a eu choc.
  70.          IPOIN = 0
  71.          ITIN = ID7
  72.          DO 10 IN = 1,NOMBN2
  73.             IP5 = IP4 + NOMBN1*(IN - 1)
  74. *   D{placement du point IN du profil mobile dans le plan des profils
  75.             XN1 = ZERO
  76.             YN1 = ZERO
  77.             DO 12 ID = 1,IDIM
  78.                XXX = XPALB(I,ITIN+ID) + XPTB(NPOI,1,ID)
  79.      &                                - XPALB(I,ID2+ID)
  80.                XN1 = XN1 + ( XXX * XPALB(I,ID3+ID) )
  81.                YN1 = YN1 + ( XXX * XPALB(I,ID4+ID) )
  82.  12            CONTINUE
  83. *           end do
  84. *   On boucle sur le nombre d'{l{ment du profil fixe pour savoir si
  85. *   le point IN @ travers{ le profil fixe
  86.             NUME1 = ID6
  87.             NUME2 = ID6 + IDIM
  88.             ICOEF = ID8
  89.             DO 14 IE = 1,NOMBN1
  90.                IPALB(I,IP5+IE) = 0
  91.                IPOSI = IPALB(I,IP3+IE)
  92.                ICAS = IPALB(I,IP1+IE)
  93.                XAIJ = XPALB(I,ICOEF+1)
  94.                XBIJ = XPALB(I,ICOEF+2)
  95.                IF (ICAS.EQ.1) THEN
  96.                   IF (XN1.GT.XAIJ) THEN
  97.                      IPOS2 = 1
  98.                   ELSE IF (XN1.LT.XAIJ) THEN
  99.                      IPOS2 = -1
  100.                   ELSE
  101.                      IPOS2 = 0
  102.                   ENDIF
  103.                ELSE IF (ICAS.EQ.2) THEN
  104.                   IF (YN1.GT.XBIJ) THEN
  105.                      IPOS2 = 1
  106.                   ELSE IF (YN1.LT.XBIJ) THEN
  107.                      IPOS2 = -1
  108.                   ELSE
  109.                      IPOS2 = 0
  110.                   ENDIF
  111.                ELSE
  112.                   RES = ( XAIJ * XN1 ) + XBIJ
  113.                   IF (YN1.GT.RES) THEN
  114.                      IPOS2 = 1
  115.                   ELSE IF (YN1.LT.RES) THEN
  116.                      IPOS2 = -1
  117.                   ELSE
  118.                      IPOS2 = 0
  119.                   ENDIF
  120.                ENDIF
  121.                IF (IPOSI.NE.IPOS2) THEN
  122. *   Le point IN a travers{ le maillage fixe,
  123. *      quel {l{ment a {t{ travers{ par ce point ?
  124.                   IF (IE.EQ.NOMBN1) NUME2 = ID6
  125.                   XE1 = ZERO
  126.                   YE1 = ZERO
  127.                   XE2 = ZERO
  128.                   YE2 = ZERO
  129.                   DO 16 ID = 1,IDIM
  130.                      XXX = XPALB(I,NUME1+ID) - XPALB(I,ID2+ID)
  131.                      YYY = XPALB(I,NUME2+ID) - XPALB(I,ID2+ID)
  132.                      XE1 = XE1 + ( XXX * XPALB(I,ID3+ID) )
  133.                      YE1 = YE1 + ( XXX * XPALB(I,ID4+ID) )
  134.                      XE2 = XE2 + ( YYY * XPALB(I,ID3+ID) )
  135.                      YE2 = YE2 + ( YYY * XPALB(I,ID4+ID) )
  136.  16                  CONTINUE
  137. *                 end do
  138.                   IF (XE2.GT.XE1) THEN
  139.                      IF (XE1.LE.XN1 .AND. XN1.LT.XE2) THEN
  140.                         IPOIN = IPOIN + 1
  141.                         IPALB(I,IP5+IE) = 1
  142.                         GOTO 18
  143.                      ENDIF
  144.                   ELSE
  145.                      IF (XE2.LT.XN1 .AND. XN1.LE.XE1) THEN
  146.                         IPOIN = IPOIN + 1
  147.                         IPALB(I,IP5+IE) = 1
  148.                         GOTO 18
  149.                      ENDIF
  150.                   ENDIF
  151.                   IF (YE2.GT.YE1) THEN
  152.                      IF (YE1.LE.YN1 .AND. YN1.LT.YE2) THEN
  153.                         IPOIN = IPOIN + 1
  154.                         IPALB(I,IP5+IE) = 1
  155.                         GOTO 18
  156.                      ENDIF
  157.                   ELSE
  158.                      IF (YE2.LT.YN1 .AND. YN1.LE.YE1) THEN
  159.                         IPOIN = IPOIN + 1
  160.                         IPALB(I,IP5+IE) = 1
  161.                         GOTO 18
  162.                      ENDIF
  163.                   ENDIF
  164.                ENDIF
  165.                NUME1 = NUME2
  166.                NUME2 = NUME2 + IDIM
  167.                ICOEF = ICOEF + 2
  168.  14            CONTINUE
  169. *           end do
  170.  18         CONTINUE
  171.             IP3 = IP3 + NOMBN1
  172.             ITIN = ITIN + IDIM
  173.  10         CONTINUE
  174. *        end do
  175. *
  176. *     Si le nombre de points ayant travers{s le profil fixe est {gale
  177. *     au nombre de points du maillage mobile, alors on imprime une
  178. *     erreur et on arrete le calcul.
  179.          IF (IPOIN.EQ.NOMBN2) THEN
  180.             IERRD = 1
  181.             RETURN
  182.          ENDIF
  183. *
  184. *        s'il y a eu choc,
  185. *        il y a IPOIN qui ont travers{ le maillage fixe
  186. *
  187.          IF (IPOIN.NE.0) THEN
  188. *
  189. *   Il y a eu choc :
  190. *    - recherche du premier point qui a travers{ le maillage fixe, d'o:
  191. *      d{termination du point P, intersection maillage fixe-mobile
  192. *    - recherche du dernier point qui a travers{ le maillage fixe, d'o:
  193. *      d{termination du point PP, intersection maillage fixe-mobile
  194. *    - @ partir des points P et PP, d{termination de la normale de choc
  195. *    - calcul de l'aire de choc, @ partir du point P, des points du
  196. *      maillage fixe, du point PP et des points du maillage mobile qui
  197. *      n'ont pas travers{ le maillage fixe
  198. *    - calcul de la force de choc.
  199. *
  200.             ICOMP = 0
  201.             IN2 = 0
  202.  20         CONTINUE
  203. *   On boucle sur les points du profil mobile
  204.             IN2 = IN2 + 1
  205.             IP5 = IP4 + NOMBN1*(IN2 - 1)
  206.             IF (IN2.GT.NOMBN2) GOTO 22
  207.             LPETPP = .FALSE.
  208.             DO 24 IE1 = 1,NOMBN1
  209.                IETAT = IPALB(I,IP5+IE1)
  210. *   On regarde si le point IN2 a travers{ un {l{ment du profil fixe
  211.                IF (IETAT.EQ.1) THEN
  212.                   ICOMP = ICOMP + 1
  213.                   IF (IN2.EQ.1) THEN
  214. *   Si on traite le premier point, on cherche les points en partant du
  215. *   dernier
  216.                      DO 26 INN2 = NOMBN2,1,-1
  217.                         IP6 = IP4 + NOMBN1*(INN2 - 1)
  218.                         DO 28 IEE1 = 1,NOMBN1
  219.                            IETA2 = IPALB(I,IP6+IEE1)
  220.                            IF (IETA2.EQ.1) THEN
  221.                               ICOMP = ICOMP + 1
  222.                               IEEE1 = IEE1
  223.                               GOTO 26
  224.                            ENDIF
  225.  28                        CONTINUE
  226. *                       end do
  227. *   Le point INN2 n'a pas travers{ le maillage fixe
  228. *   Le point (INN2-1) et la point INN2 forment le segment de droite
  229. *   pour d{finir le point d'intersection P avec l'{l{ment IEEE1
  230.                         IF (INN2.EQ.NOMBN2) THEN
  231.                            IPOIN1 = NOMBN2
  232.                            IPOIN2 = 1
  233.                            IELEM1 = IE1
  234.                         ELSE
  235.                            IPOIN1 = INN2
  236.                            IPOIN2 = INN2 + 1
  237.                            IELEM1 = IEEE1
  238.                         ENDIF
  239.                         GOTO 30
  240.  26                     CONTINUE
  241. *                    end do
  242.  30                  CONTINUE
  243.                   ELSE
  244.                      IPOIN1 = IN2 - 1
  245.                      IPOIN2 = IN2
  246.                      IELEM1 = IE1
  247.                   ENDIF
  248. *   D{termination du point P
  249. *   Le point IPOIN1 n'a pas travers{ le maillage fixe
  250. *   Le point IPOIN2 a travers{ le maillage fixe par l'{l{ment IELEM1
  251.                   CALL DYNE31(IPOIN1,IPOIN2,IELEM1,XPALB,IPALB,XPTB,
  252.      &                        NLIAB,NPLB,I,NPOI,IND,ID1,IP1,XP,YP)
  253. *
  254. *   Recherche du point PP, on regarde les points suivant
  255.                   IF (IN2.EQ.NOMBN2) THEN
  256.                      IPOIN3 = IN2
  257.                      IPOIN4 = 1
  258.                      IELEM2 = IE1
  259.                      GOTO 32
  260.                   ENDIF
  261.                   IN22 = IN2 + 1
  262.                   IEEE1 = IE1
  263.                   DO 34 INN2 = IN22,NOMBN2
  264.                      IP7 = IP4 + NOMBN1*(INN2 - 1)
  265.                      DO 36 IEE1 = 1,NOMBN1
  266.                         IETA2 = IPALB(I,IP7+IEE1)
  267.                         IF (IETA2.EQ.1) THEN
  268.                            ICOMP = ICOMP + 1
  269.                            IEEE1 = IEE1
  270.                            IF (INN2.EQ.NOMBN2) THEN
  271.                               IPOIN3 = INN2
  272.                               IPOIN4 = 1
  273.                               IELEM2 = IEE1
  274.                            ENDIF
  275.                            GOTO 34
  276.                         ENDIF
  277.  36                     CONTINUE
  278. *                    end do
  279. *   Le point INN2 n'a pas travers{ le maillage fixe
  280. *   Le point (INN2-1) et la point INN2 forment le segment de droite
  281. *   pour d{finir le point d'intersection P avec l'{l{ment IEEE1
  282.                      IPOIN3 = INN2 - 1
  283.                      IPOIN4 = INN2
  284.                      IELEM2 = IEEE1
  285.                      GOTO 32
  286.  34                  CONTINUE
  287. *                 end do
  288.  32               CONTINUE
  289. *   D{termination du point PP
  290. *   Le point IPOIN3 a travers{ le maillage fixe par l'{l{ment IELEM2
  291. *   Le point IPOIN4 n'a pas travers{ le maillage fixe
  292.                   CALL DYNE31(IPOIN3,IPOIN4,IELEM2,XPALB,IPALB,XPTB,
  293.      &                 NLIAB,NPLB,I,NPOI,IND,ID1,IP1,XPP,YPP)
  294.                   LPETPP = .TRUE.
  295.                   GOTO 40
  296.                ENDIF
  297.  24            CONTINUE
  298. *           end do
  299.  40         CONTINUE
  300. *
  301.             IF ( LPETPP ) THEN
  302. *   On a d{termin{ les points P et PP, calcul de la normale
  303.                XNORM1 = ( XP + XPP ) * 0.5D0
  304.                YNORM1 = ( YP + YPP ) * 0.5D0
  305.                XXX = XP - XPP
  306.                YYY = YP - YPP
  307.                IPT3  = ID7 + IDIM*(IPOIN3-1)
  308.                XN4 = ZERO
  309.                YN4 = ZERO
  310.                DO 42 ID = 1,IDIM
  311.                   XXX = XPALB(I,IPT3+ID) + XPTB(NPOI,1,ID)
  312.      &                                   - XPALB(I,ID2+ID)
  313.                   XN4 = XN4 + ( XXX * XPALB(I,ID3+ID) )
  314.                   YN4 = YN4 + ( XXX * XPALB(I,ID4+ID) )
  315.  42               CONTINUE
  316. *              END DO
  317.                IF (ABS(YYY).LT.PRECIS) THEN
  318.                   YNN = YN4 - YNORM1
  319.                   XNORM = ZERO
  320.                   YNORM = YNN / ABS(YNN)
  321.                ELSE IF (ABS(XXX).LT.PRECIS) THEN
  322.                   XNN = XN4 - XNORM1
  323.                   XNORM = XNN / ABS(XNN)
  324.                   YNORM = ZERO
  325.                ELSE
  326. * le premier vecteur est P-PP
  327. * la normale de choc sera le produit vectoriel P_PP, N
  328.                   XPPP = XPP - XP
  329.                   YPPP = YPP - YP
  330.                   PS = SQRT( XPPP*XPPP + YPPP*YPPP )
  331.                   XNORM = YPPP * XPALB(I,ID5+3) / PS
  332.                   YNORM = -XPPP * XPALB(I,ID5+3) / PS
  333.                ENDIF
  334. *   Calcul de l'aire
  335.                CALL DYNE32(IPOIN1,IPOIN4,IELEM1,IELEM2,XP,YP,XPP,YPP,
  336.      &         XPALB,IPALB,XPTB,NLIAB,NPLB,I,NPOI,ID1,IP1,IND,SOMAIR)
  337.                AIRE = XSECT - SOMAIR
  338. *   Calcul de la force de choc
  339. *   Liaison conditionnelle
  340.                IF (iannul.EQ.0) THEN
  341.                   FN = -XRAID * ( AIRE ** XEXPO )
  342.                ELSE
  343.                   FN = 0
  344.                ENDIF
  345.                SOMFN = SOMFN + FN
  346. *   On revient dans le repere global
  347.                XFOR1 = FN * XNORM
  348.                XFOR2 = FN * YNORM
  349.                DO 38 ID = 1,IDIM
  350.                   XXX = XFOR1*XPALB(I,ID3+ID) + XFOR2*XPALB(I,ID4+ID)
  351.                   FTOTB(NPOI,ID) = FTOTB(NPOI,ID) + XXX
  352.  38               CONTINUE
  353. *              end do
  354.                IN2 = INN2
  355.             ENDIF
  356.             IF (ICOMP.EQ.IPOIN) GOTO 22
  357.             GOTO 20
  358.  22         CONTINUE
  359.          ENDIF
  360.          DO 100 ID = 1,IDIM
  361.             XVALB(I,IND,3+ID) = XPTB(NPOI,1,ID)
  362.  100        CONTINUE
  363. *        end do
  364.          XVALB(I,IND,1) = SOMFN
  365. *
  366. * --- choc ...........
  367. *
  368. *     ELSE IF (ITYP.EQ.  ) THEN
  369. *        .......
  370. *        .......
  371. *
  372.       ENDIF
  373. *
  374.       END
  375.  
  376.  
  377.  
  378.  
  379.  

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