Sources Esope | Cast3M

ktange
C KTANGE    SOURCE    PV        20/04/02    21:15:26     10567          
 
      SUBROUTINE KTANGE
 
*_______________________________________________________________________
*
*           Operateur matrice de rigidite tangente KTAN
*
*    RI1 = 'KTAN' MOD1 SI1 VA1 CA1
*                 ('PREC' FLO1) ('DT  ' FLO2 ) ('SYME') ;
*
*     MOD1  MMODEL   Modele de calcul
*     SI1   MCHAML   Champ par element de contraintes
*     VA1   MCHAML   Champ par element de variables internes
*     CA1   MCHAML   Champ par element de caracteristiques
*                    (materielles et/ou geometriques)
*     FLO1  FLOTTANT (par defaut 1.D-3)
*                    précision avec laquelle on cherche si un état de
*                    contraintes est plastique ou non
*     FLO2  FLOTTANT nécessaire que pour les modèles visqueux
*                    pas de temps servant à calculer la matrice tangente
*
*_______________________________________________________________________
*
*     Calcul de la matrice de rigidite tangente par perturbation
*
*  RI1 = 'KTAN' 'PERT' MOD1 CHE1 CHE2
*                      ('C1' FLO1) ('C2' FLO2) ('SYME') ;
*
*        MOD1  MMODEL   Modele de calcul
*        CHE1  MCHAML   Champ par element a t pour COMP
*        CHE2  MCHAML   Champ par element a t+dt
*                       (champs pour COMP et les contraintes)
*        FLO1  FLOTTANT positif (par defaut 1.D-3)
*                       coefficient de perturbation de l increment de
*                       deformation
*        FLO2  FLOTTANT positif (par defaut FLO1/100.)
*                       perturbation minimale
*
*_______________________________________________________________________
*
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8(A-H,O-Z)
*
-INC PPARAM
-INC CCOPTIO
-INC SMCOORD
*
      LOGICAL BC2
      CHARACTER*4 DEUXMOT(3),PERTURB(4)
*
      DATA DEUXMOT / 'PREC','DT  ','SYME' /
      DATA PERTURB / 'PERT','C1  ','C2  ','SYME' /
*
      CALL LIRMOT(PERTURB,1,IPLA,0)
      IF (IERR.NE.0) RETURN
*
      IPRIGI=0
      segact mcoord
*
************************************************************************
*     MATRICE TANGENTE ANALYTIQUE SI POSSIBLE, SINON ELASTIQUE
************************************************************************
      IF (IPLA.EQ.0) THEN
*
*       Lecture d'un modele
*
        CALL LIROBJ('MMODEL',IPMODL,1,IRET)
        CALL ACTOBJ('MMODEL',IPMODL,1)
        IF (IERR.NE.0) RETURN
*
*       Lecture du MCHAML de CONTRAINTES
*
        CALL LIROBJ('MCHAML',IPIN,1,IRET)
        CALL ACTOBJ('MCHAML',IPIN,1)
        IF (IERR.NE.0) RETURN
          CALL REDUAF(IPIN,IPMODL,IPCHE1,0,IR,KER)
          IF(IR   .NE. 1) CALL ERREUR(KER)
          IF(IERR .NE. 0) RETURN
*
*       Lecture du MCHAML de VARIABLES INTERNES
*
        CALL LIROBJ('MCHAML',IPIN,1,IRET)
        CALL ACTOBJ('MCHAML',IPIN,1)
        IF (IERR.NE.0) RETURN
          CALL REDUAF(IPIN,IPMODL,IPCHE2,0,IR,KER)
          IF(IR   .NE. 1) CALL ERREUR(KER)
          IF(IERR .NE. 0) RETURN
*
*       Lecture du MCHAML de CARACTERISTIQUES
*
        CALL LIROBJ('MCHAML',IPIN,1,IRET)
        CALL ACTOBJ('MCHAML',IPIN,1)
        IF (IERR.NE.0) RETURN
          CALL REDUAF(IPIN,IPMODL,IPCHE3,0,IR,KER)
          IF(IR   .NE. 1) CALL ERREUR(KER)
          IF(IERR .NE. 0) RETURN
*
*       On lit un ou deux flottants xprec et dt
*
        XPREC  = 1.D-3
        DTPS   = 0.D0
        IKTSYM = 0
 1      CONTINUE
        CALL LIRMOT(DEUXMOT,3,IVAL,0)
        IF (IERR.NE.0) RETURN
        IF (IVAL .EQ. 1) THEN
          CALL LIRREE(XPREC,1,IRET)
          IF (IERR.NE.0) RETURN
          GOTO 1
        ELSE IF (IVAL .EQ. 2) THEN
          CALL LIRREE(DTPS,1,IRET)
          IF (IERR.NE.0) RETURN
          GOTO 1
        ELSE IF (IVAL .EQ. 3) THEN
          IKTSYM = 1
          GOTO 1
        ENDIF
*
*       Calcul de la matrice tangente
*
        CALL KTANGA(IPMODL,IPCHE1,IPCHE2,IPCHE3,XPREC,DTPS,IKTSYM,
     &              IPRIGI)
*
************************************************************************
*         MATRICE DE RAIDEUR TANGENTE PAR PERTURBATION
************************************************************************
* Appel a 'KTAN' 'PERT'
      ELSE IF (IPLA.EQ.1) THEN
*
*       Lecture d'un modele
*
        CALL LIROBJ('MMODEL',IPMODL,1,IRET)
        CALL ACTOBJ('MMODEL',IPMODL,1)
        IF (IERR.NE.0) RETURN
*
*       Lecture du MCHAML initial
*
        CALL LIROBJ('MCHAML',IPIN,1,IRET)
        IF (IERR.NE.0) RETURN
        CALL REDUAF(IPIN,IPMODL,IPCHE1,0,IR,KER)
        IF(IR   .NE. 1) CALL ERREUR(KER)
        IF(IERR .NE. 0) RETURN
*
*       Lecture du MCHAML final
*
        CALL LIROBJ('MCHAML',IPIN,1,IRET)
        CALL ACTOBJ('MCHAML',IPIN,1)
        IF (IERR.NE.0) RETURN
        CALL REDUAF(IPIN,IPMODL,IPCHE2,0,IR,KER)
        IF(IR   .NE. 1) CALL ERREUR(KER)
        IF(IERR .NE. 0) RETURN
*
*       On lit un ou deux flottants CO1 et CO2
*       par defaut CO1 = 1.D-3 et CO2 = CO1/100.
*
        CO1 = 1.D-3
        BC2 = .TRUE.
        IKTSYM = 0
 100    CONTINUE
        CALL LIRMOT(PERTURB(2),3,IVAL,0)
        IF (IERR.NE.0) RETURN
        IF (IVAL .EQ. 1) THEN
          CALL LIRREE(CO1,1,IRET)
          IF (IERR.NE.0) RETURN
          GOTO 100
        ELSE IF (IVAL .EQ. 2) THEN
          CALL LIRREE(CO2,1,IRET)
          IF (IERR.NE.0) RETURN
          BC2 = .FALSE.
          GOTO 100
        ELSE IF (IVAL .EQ. 3) THEN
          IKTSYM = 1
          GOTO 100
        ENDIF
*
*   Verification des coefficients de perturbation C1 et C2
*
        IF (CO1.LE.0.) THEN
          CO1 = 1.D-3
          WRITE (IOIMP,*) 'Valeur de C1 negative ==> ',
     &                    'valeur par defaut (1.D-3) utilisee'
        ENDIF
        IF (CO2.LE.0.) THEN
          BC2 = .TRUE.
          WRITE (IOIMP,*) 'Valeur de C2 negative ==> ',
     &                    'valeur par defaut (C1/100.) utilisee'
        ENDIF
        IF (BC2) CO2 = CO1/1.D2
*
*   Calcul de la matrice par perturbation
*
        CALL KTAPER(IPMODL,IPCHE1,IPCHE2,CO1,CO2,IKTSYM, IPRIGI)
*
      ENDIF
*
************************************************************************
*       Ecriture de la RIGIDITE
************************************************************************
      IF (IPRIGI.GT.0) CALL ECROBJ('RIGIDITE',IPRIGI)
 
      END