Sources Esope | Cast3M

cfl2
C CFL2      SOURCE    OF166741  25/02/21    21:15:29     12166          
      SUBROUTINE CFL2(ICAS,IPMAIL,MELE,IVAM1,IVAM2,MELV1,MELV2,N2)
*---------------------------------------------------------------------*
*
*          calcul du pas de temps CFL
*
* elements massifs formulation mécanique
*
*
*  entree
*       icas : cas à traiter
*        = 1 calcul du pas de temps complet ivam1 avec cara
*        = 2 calcul du pas de temps lorsque cson est donne ivam2
*        = 3 calcul du pas de temps lorsque la taille est donnée ivam1 si cara
*        = 4 calcul de la vitesse du son ivam1 donné
*        = 5 calcul du parametre de taille ivam1 si cara
*
*       ipmail : pointeur vers le maillage a traiter
*       mele   : numero de l'élément finis dans nomtp
*       ivam1  : pointeur vers mptval du cham1 actif
*       ivam2  : pointeur vers mptval du cham2 actif
*        n2    : nombre de comosante en sortie
*
*  sortie
*        melv1 : melval de la première composante du chamelem resultat
*                inactif en sortie
*        melv2 : melval de la deuxième composante du chamelem resultat
*                inactif en sortie
*
*
*---------------------------------------------------------------------*
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8(A-H,O-Z)
 
-INC PPARAM
-INC CCOPTIO
-INC CCHAMP
 
-INC SMCHAML
-INC SMELEME
 
-INC TMPTVAL
 
      MPTVA1 = IVAM1
      MPTVA2 = IVAM2
*
*    branchement en fonction de l'élément fini
*
*           0            5              0              5             0
      GOTO (99,99,99, 4,99,99,99, 4,99),MELE
      GOTO 99
*
* éléments massifs
*
 4    CONTINUE
*       ================  calcul de la vitesse du son
*       le resultat est stocké dans  melval  avec n1ptel =1
        IF (ICAS.EQ.1.OR.ICAS.EQ.3.OR.ICAS.EQ.4) THEN
*          recherche des paramètre matériau
*            module d'young
           MELVA3 = MPTVA1.IVAL(1)
*            densite
           MELVA4 = MPTVA1.IVAL(3)
           SEGACT MELVA3,MELVA4
*
           N1EL = MIN(MELVA3.VELCHE(/2),MELVA4.VELCHE(/2))
           N1PTEL = 1
           N2PTEL = 0
           N2EL = 0
           SEGINI MELVAL
           MCSON = MELVAL
*    boucle sur les éléments pour calculer la vitesse du son
           DO 403 I=1,N1EL
*            on prend les valeurs moyennes sur les éléments
               YOU1 = 0.D0
               I3 = MIN(I,MELVA3.VELCHE(/2))
               DO 401 J=1,MELVA3.VELCHE(/1)
                  YOU1 = YOU1 + MELVA3.VELCHE(J,I3)
 401          CONTINUE
               YOU1 = YOU1 / MELVA3.VELCHE(/1)
*
               RO1 = 0.D0
               I4 = MIN(I,MELVA4.VELCHE(/2))
               DO 402 J=1,MELVA4.VELCHE(/1)
                  RO1 = RO1 + MELVA4.VELCHE(J,I4)
 402          CONTINUE
               RO1 = RO1 / MELVA4.VELCHE(/1)
               IF (RO1.EQ.0.D0) THEN
                 SEGDES MELVA4,MELVA3
                 SEGSUP MELVAL
                 CALL ERREUR(855)
                 RETURN
               ENDIF
*
               IF (YOU1.EQ.0.D0) THEN
                 SEGDES MELVA4,MELVA3
                 SEGSUP MELVAL
                 CALL ERREUR(856)
                 RETURN
               ENDIF
               VELCHE(1,I) = SQRT(YOU1/RO1)
*               write(6,*) 'Element', i , 'Cson' , VELCHE(1,i)
 403      CONTINUE
           SEGDES MELVA4,MELVA3
           IF (ICAS.EQ.4) THEN
*          cas ou seule la vitesse du son est demandée
             MELVA2 = 0
             MELV1 = MELVAL
             SEGDES MELVAL
             RETURN
           ENDIF
        ELSE IF (ICAS.EQ.2) THEN
*            recuperation du champ
             SEGACT MPTVA2
             MELVA1 = MPTVA2.IVAL(1)
             SEGACT MELVA1
             MCSON = MELVA1
             SEGDES MPTVA2
        ENDIF
*       ================ paramètre geometrique
*       stocké dans un melval mtaille
        IF (ICAS.EQ.1.OR.ICAS.EQ.2.OR.ICAS.EQ.5) THEN
           MELEME = IPMAIL
           SEGACT MELEME
           N1EL = NUM(/2)
           N1PTEL = 1
           N2PTEL = 0
           N2EL = 0
           SEGINI MELVAL
           MTAIL1 = MELVAL
           MTAIL2 = 0
           CALL CFLTAI(MTAIL1,IPMAIL,MELE)
           SEGDES MELEME
           IF (ICAS.EQ.5) THEN
             MELV1 = MTAIL1
             MELV2 = 0
             SEGDES MELVAL
             RETURN
           ENDIF
        ELSE IF (ICAS.EQ.3) THEN
*            recuperation du champ
             SEGACT MPTVA2
             MELVA1 = MPTVA2.IVAL(1)
             SEGACT MELVA1
             MTAIL1 = MELVA1
             SEGDES MPTVA2
         ENDIF
*       ================ pas de temps cfl
        IF (ICAS.EQ.1.OR.ICAS.EQ.2.OR.ICAS.EQ.3) THEN
*          recuperation de la vitesse du son
*          et du paramètre de taille
           MELVA1 = MCSON
           MELVA2 = MTAIL1
*          creation du melval résultat
           N1EL = MAX(MELVA1.VELCHE(/2),MELVA2.VELCHE(/2))
           N1PTEL = 1
           N2EL = 0
           N2PTEL = 0
           SEGINI MELVAL
*
           DO 404 I=1,N1EL
              I1 = MIN(I,MELVA1.VELCHE(/2))
              I2 = MIN(I,MELVA2.VELCHE(/2))
              VELCHE(1,I)=MELVA2.VELCHE(1,I2)/MELVA1.VELCHE(1,I1)
*              write(6,*) 'Element', i , 'Dtcfl' , VELCHE(1,i)
 404       CONTINUE
           MELV1 = MELVAL
           MELV2 = 0
           SEGDES MELVAL
           IF (ICAS.EQ.1) THEN
                SEGSUP MELVA1,MELVA2
           ELSE IF (ICAS.EQ.2) THEN
                SEGSUP MELVA2
                SEGDES MELVA1
           ELSE
                SEGSUP MELVA1
                SEGDES MELVA2
           ENDIF
           RETURN
        ENDIF
 
   99 CONTINUE
      MOTERR(1:4)=NOMTP(MELE)
      MOTERR(9:12)='CFL5'
      CALL ERREUR(86)
 
      RETURN
      END