Sources Esope | Cast3M

manur4
C MANUR4    SOURCE    CB215821  22/12/14    12:37:37     11527          
      SUBROUTINE MANUR4 (IPELEM,IPDESC,MTEMP4,IPMATR,IANTI)
************************************************************************
*
*                             M A N U R 4
*                             -----------
*
* FONCTION:
* ---------
*
*     CONSTRUCTION DES MATRICES ELEMENTAIRES DE RIGIDITE POUR UN OBJET
*     'RIGIDITE' CREE MANUELLEMENT.
*     L'UTILISATION DE CE SOUS-PROGRAMME N'EST PAS UNIVERSELLE.
*
* MODE D'APPEL:
* -------------
*
*     CALL MANUR4 (IPELEM,IPDESC,MTEMP4,IPMATR,IANTI)
*
* PARAMETRES:   (E)=ENTREE   (S)=SORTIE
* -----------
*
*     IPELEM  ENTIER    (E)  POINTEUR DE L'OBJET 'MAILLAGE' SUR LEQUEL
*                            VA S'APPUYER LA 'RIGIDITE'.
*     IPDESC  ENTIER    (E)  POINTEUR SUR LE SEGMENT DESCRIPTEUR DE
*                            L'OBJET 'RIGIDITE'.
*     MTEMP4  SEGMENT   (E)  REGROUPEMENT DE POINTEURS SUR DES
*                            'LISTREEL'.
*                            SOIT IL N'Y A QU'1 'LISTREEL', QUI CONTIENT
*                            TOUS LES TERMES DE LA
*                            MATRICE ELEMENTAIRE DE RIGIDITE, ECRITS
*                            LIGNE PAR LIGNE,
*                            SOIT IL Y A AUTANT DE 'LISTREEL' QUE DE
*                            LIGNES DANS LA MATRICE ELEMENTAIRE DE
*                            RIGIDITE, LE N-IEME 'LISTREEL' DECRIVANT LA
*                            N-IEME LIGNE DE LA MATRICE .
*     IPMATR  ENTIER    (S)  POINTEUR SUR LE SEGMENT CONTENANT LA LISTE
*                            DES POINTEURS DES MATRICES ELEMENTAIRES DE
*                            RIGIDITE.
*
*     EXEMPLE DE PRESENTATION DE LA MATRICE ELEMENTAIRE :
*     | A B C |
*     | D E F |
*     | G H I |
*     Elle peut etre donnee par:  (PROG A B C D E F G H I )
*     ou bien par :  (PROG A B C ) (PROG D E F) (PROG G H I )
*     si la matrice est symetrique ou antisymetrique on peut aussi
*     la decrire par 1 LISTREEEL : (PROG A D E G H I )
*     ou bien par plusieurs LISTREEL : (PROG A ) (PROG D E) (PROG G H I)
*
* LEXIQUE: (ORDRE ALPHABETIQUE)
* --------
*
*     LONG    ENTIER    LONGUEUR DU 'LISTREEL' TRAITE.
*     NBLREE  ENTIER    NOMBRE DE 'LISTREEL' REFERENCES PAR "MTEMP4".
*
*     LES AUTRES VARIABLES IMPORTANTES SONT EXPLIQUEES DANS LES MODULES
*     INCLUS.
*
************************************************************************
*
      IMPLICIT INTEGER(I-N)
      LOGICAL ZTRI
-INC PPARAM
-INC CCOPTIO
-INC SMELEME
-INC SMLREEL
-INC SMRIGID
-INC CCREEL
*
      SEGMENT /MTEMP4/ (ILREEL(0))
 
************************************************************************
*     INITIALISATIONS ET OUVERTURE
************************************************************************
 
*
      MELEME = IPELEM
      SEGACT,MELEME
      NELRIG = NUM(/2)
      SEGDES,MELEME
*
      DESCR = IPDESC
      SEGACT,DESCR
C ... La distinction entre les deux nombres est un peu artificielle,
C     car manur3 vérifie si les listmots sont de longueurs égales,
C     on en a juste besoin pour initialiser XMATRI ...
      NLIGRP = NOELEP(/1)
      NLIGRD = NOELED(/1)
C ... LVAL = nombre de termes d'une matrice pleine ...
      LVAL = NLIGRP * NLIGRD
      SEGDES,DESCR
*
      SEGINI,xMATRI
      IPMATR = xMATRI
*      SEGINI,XMATRI
*      DO 100 IB100=1,NELRIG
*         IMATTT(IB100) = XMATRI
*  100    CONTINUE
*     END DO
*      SEGDES,IMATRI
*
*     RQ: "XMATRI" EST GARDE ACTIF.
*
      SEGACT,MTEMP4
      NBLREE = ILREEL(/1)
*
************************************************************************
*     Cas 1 seul LISTREEL
************************************************************************
      IF (NBLREE .EQ. 1) THEN
*
         MLREEL = ILREEL(1)
         SEGACT,MLREEL
C    ... LVA1 = nombre de termes d'une matrice carrée défini par
C        sa moitié ...
         LVA1=NLIGRP*(NLIGRP+1)/2
         IF (LVAL .NE. PROG(/1).AND.PROG(/1).NE.LVA1 ) THEN
            NUMERR = 199
            PRINT *,'On attend ',LVAL, ' termes',NLIGRP
            CALL ERREUR (NUMERR)
            RETURN
         END IF
*
C    ... ILA = N° de ligne de la matrice élémentaire ...
         ILA=1
C    ... ILC = N° de colonne de la matrice élémentaire ...
         ILC=1
 
C    ... ZTRI dit si toutes les composantes ont été données ou
C        juste le triangle inférieur ...
         ZTRI=.FALSE.
cbp         IF(PROG(/1).EQ.LVA1) ZTRI=.TRUE.
         IF(PROG(/1).EQ.LVA1.AND.PROG(/1).NE.LVAL) ZTRI=.TRUE.
         IF(ZTRI .AND. IANTI.EQ.2) THEN
            CALL ERREUR(731)
C       ... On laisse les cochonneries dans XMATRI ...
*            SEGDES,XMATRI
C       ... puis on s'en va ...
            RETURN
         ENDIF
 
         DO 200 IB200=1,PROG(/1)
 
            RE(ILA,ILC,1) = PROG(IB200)
 
            IF(ZTRI) THEN
               IF(IANTI.EQ.1) THEN
                  RE(ILC,ILA,1)=-PROG(IB200)
               ELSE
                  RE(ILC,ILA,1)=PROG(IB200)
               ENDIF
            ENDIF
 
            ILC=ILC+1
C       ... On passe à la ligne si on a traversé la diagonale (cas triangulaire) ...
            IF(ILC.GT.ILA.AND. ZTRI) THEN
               ILC=1
               ILA=ILA+1
            ENDIF
 
C       ... On passe à la ligne si on est au bout (cas plein) ...
            IF(ILC.GT.NLIGRP) THEN
               ILC=1
               ILA=ILA+1
            ENDIF
 
  200    CONTINUE
        do ib=2,nelrig
                 do io=1,nligrp
                   do iu=1,nligrd
                     re(iu,io,ib)=re(iu,io,1)
                   enddo
                 enddo
                enddo
*        END DO
*
         SEGDES,MLREEL
 
 
************************************************************************
*     Cas plusieurs LISTREEL
************************************************************************
      ELSE IF (NBLREE .GT. 1) THEN
*
         IF (NBLREE .EQ. NLIGRD) THEN
*
            MLREEL=ILREEL(1)
            SEGACT MLREEL
 
C       ... Cas triangulaire ? ...
            ZTRI=.FALSE.
            IF(PROG(/1).EQ.1) ZTRI=.TRUE.
            IF(ZTRI .AND. IANTI.EQ.2) THEN
               CALL ERREUR(731)
C          ... On laisse les cochonneries dans XMATRI ...
*               SEGDES,XMATRI
C          ... puis on s'en va ...
               RETURN
            ENDIF
 
C       ... Boucle sur les lignes (IB300 = N° de la ligne) ...
            DO 300 IB300=1,NBLREE
*
               MLREEL = ILREEL(IB300)
               SEGACT,MLREEL
               LONG = PROG(/1)
C          ... Cas lignes pleines : longueur doit être NLIGRP ...
               IF (.not.ZTRI .AND. LONG.NE.NLIGRP) THEN
                  NUMERR = 200
                  CALL ERREUR (NUMERR)
                  RETURN
               END IF
C          ... Cas triangulaire : longueur doit être N° de la ligne ...
               IF (ZTRI .AND. LONG.NE.IB300) THEN
                  NUMERR = 200
                  CALL ERREUR (NUMERR)
                  RETURN
               END IF
 
C          ... Boucle sur les colonnes (IB310 = N° de la colonne) ...
               DO 310 IB310=1,LONG
 
                  RE(IB300,IB310,1) = PROG(IB310)
       IF(IANTI.EQ.1.AND.ZTRI) RE(IB310,IB300,1)=-RE(IB300,IB310,1)
       IF(IANTI.EQ.0.AND.ZTRI) RE(IB310,IB300,1)=RE(IB300,IB310,1)
 
  310          CONTINUE
*              END DO
*
               SEGDES,MLREEL
*
  300       CONTINUE
            do ib=2,nelrig
              do io=1,nligrp
                do iu=1,nligrd
                   re(iu,io,ib)=re(iu,io,1)
                enddo
              enddo
            enddo
*           END DO
 
C    ... c.à.d. le nombre de LISTREEL est différent du nombre de variables duales ...
         ELSE
*
            NUMERR = 201
            CALL ERREUR (NUMERR)
            RETURN
*
         END IF
 
 
************************************************************************
*     Cas aucun LISTREEL !
************************************************************************
      ELSE
*
*        AUCUN 'LISTREEL' N'A ETE FOURNI EN DONNEE.
         MOTERR(1:8) = 'LISTREEL'
         NUMERR = 37
         CALL ERREUR (NUMERR)
         RETURN
*
      END IF
 
 
************************************************************************
*     VERIFICATION EN FONCTION DES CAS (ajout, bp 2020)
************************************************************************
 
      IF(IANTI.EQ.0) THEN
        XMATRI.SYMRE=0
        IF(.NOT.ZTRI) THEN
*         SI SYM      ET COMPLET : VERIF DES TERMES EXTRA-DIAGONAUX Aij=Aji
          call versym(re,re(/1),re(/2),re(/3),IANTI)
          if (ierr.ne.0) return
 
*       ELSE
*         SI SYM      ET TRIANGULAIRE : PAS DE VERIF
        ENDIF
 
      ELSEIF(IANTI.EQ.1) THEN
        XMATRI.SYMRE=1
        IF(ZTRI) THEN
*         SI ANTI-SYM ET TRIANGULAIRE : VERIF DE LA DIAGONALE
          xzref=(xpetit/xzprec)
          do iel=1,re(/3)
            do ir=1,re(/1)
              re1=re(ir,ir,iel)
              if (abs(re1).gt.xzref) then
                MOTERR(1:15)='ANTI-SYMETRIQUE'
                reaerr(1)=re1
                reaerr(2)=0.D0
                reaerr(3)=abs(re1)
                call erreur(1044)
                return
              endif
            enddo
          enddo
 
        ELSE
*         SI ANTI-SYM ET COMPLET : VERIF DE TOUS LES TERMES Aij=-Aji
          call versym(re,re(/1),re(/2),re(/3),IANTI)
          if (ierr.ne.0) return
        ENDIF
 
      ELSE
*       SI QUELCONQUE          : On place SYMRE
        XMATRI.SYMRE=2
      ENDIF
 
      SEGDES,XMATRI
      SEGDES,MTEMP4
 
      END