C TENS2     SOURCE    GOUNAND   26/01/16    21:15:18     12450          
      SUBROUTINE TENS2(IOTENS,IKAS,A,D,R)
      IMPLICIT REAL*8 (A-H,O-Z)
      IMPLICIT INTEGER (I-N)
C***********************************************************************
C NOM         : TENS2
C DESCRIPTION : Opérations sur des tenseurs (unaires pour l'instant)
C               Ici on gère les opérations à l'échelle du noeud
C
C
C LANGAGE     : ESOPE
C AUTEUR      : Stephane GOUNAND (CEA/DES/ISAS/DM2S/SEMT/LTA)
C               mel : gounand@semt2.smts.cea.fr
C***********************************************************************
C APPELES          :
C APPELES (E/S)    :
C APPELES (BLAS)   :
C APPELES (CALCUL) :
C APPELE PAR       :
C***********************************************************************
C SYNTAXE GIBIANE    :
C ENTREES            :  A si IOTENS.NE.12
C                       D et R si IOTENS.EQ.12
C ENTREES/SORTIES    :
C SORTIES            :  R si IOTENS.NE.12
C                       D si IOTENS.EQ.11
C                       A si IOTENS.EQ.12
C***********************************************************************
C VERSION    : v1, 28/08/2024, version initiale
C HISTORIQUE : v1, 28/08/2024, creation
C HISTORIQUE :
C HISTORIQUE :
C***********************************************************************
-INC PPARAM
-INC CCOPTIO
-INC CCREEL
*
C   ENTREEES
C        A(3,3) = MATRICE A DIAGONALISER
C   SORTIES
C        D(3)   =  VALEURS PROPRES D(1) > D(2) ET D(3)=A(3,3)
C        S(3,3) =  VECTEURS PROPRES
C        R(3,3) =  RESULTAT DES OPERATIONS TENSORIELLES
      DIMENSION A(3,3),D(3),S(3,3),R(3,3)
      LOGICAL LVERIF,LDBG
      DATA ((S(I,J),I=1,3),J=1,3) /9*0.D0/
*                     1       2        3     4        5        6      7
*     DATA MOTENS/'NORM2','NORMINF','DET','TRACE','INVERSE','IDEN' ,'LOG'
*              8     9      10        11         12         13
*     $     ,'EXP','INVS','ABSOLU','PRINCIPA','RECOMPOS','TRANSPOS'
*              14
*     $     ,'VALP'/
*
*
* Executable statements
*
      lverif=.false.
      ldbg=.false.
      if (ldbg) then
         WRITE(IOIMP,*) 'TENS2 : IOTENS,IKAS=',IOTENS,IKAS
         IF (IOTENS.NE.12) THEN
            write(ioimp,*) 'Input  A,I,J=',IDIM,IDIM
            write(ioimp,*) ((A(I,J),I=1,IDIM),J=1,IDIM)
         ELSE
            write(ioimp,*) 'Input  D,J=',IDIM
            write(ioimp,*) (D(J),J=1,IDIM)
            write(ioimp,*) 'Input  R,I,J=',IDIM,IDIM
            write(ioimp,*) ((R(I,J),I=1,IDIM),J=1,IDIM)
         ENDIF
      endif
*
      IF (IOTENS.GT.6.AND.IOTENS.LT.12) THEN
*
*     Décomposition spectrale (positive pour les cas scalaires et
*     vecteurs)
*
         IF (IKAS.EQ.1.OR.(IKAS.EQ.3.AND.IDIM.EQ.1)) THEN
* Difficile de donner un sens...
*            D(1)=ABS(A(1,1))
*            S(1,1)=SIGN(1.D0,A(1,1))
            D(1)=A(1,1)
            S(1,1)=1.D0
            JD=1
         ELSEIF (IKAS.EQ.2) THEN
* Pas de sens pour les vecteurs
c$$$            XNO=0.D0
c$$$            DO I=1,IDIM
c$$$               XNO=XNO+A(I,1)**2
c$$$            ENDDO
c$$$            XNA=SQRT(XNO)
c$$$            D(1)=XNA
c$$$            DO I=1,IDIM
c$$$               S(I,1)=A(I,1)/XNA
c$$$            ENDDO
            CALL ERREUR(980)
            RETURN
c$$$         write(ioimp,*) 'A,I,J=',IDIM,1
c$$$         write(ioimp,*) ((A(I,J),I=1,IDIM),J=1,1)
c$$$         write(ioimp,*) 'D,J=',JDIM
c$$$         write(ioimp,*) (D(J),J=1,JDIM)
c$$$         write(ioimp,*) 'S,I,J=',IDIM,1
c$$$         write(ioimp,*) ((S(I,J),I=1,IDIM),J=1,1)

         ELSEIF (IKAS.EQ.3) THEN
            CALL JACOB3(A,IDIM,D,S)
            IF (IERR.NE.0) RETURN
            JD=IDIM
            if (lverif) then
               write(ioimp,*) 'A,I,J=',IDIM,IDIM
               write(ioimp,*) ((A(I,J),I=1,IDIM),J=1,IDIM)
               write(ioimp,*) 'D,J=',IDIM
               write(ioimp,*) (D(J),J=1,IDIM)
               write(ioimp,*) 'S,I,J=',IDIM,IDIM
               write(ioimp,*) ((S(I,J),I=1,IDIM),J=1,IDIM)
*     Verif
               DO K=1,IDIM
*                  DO I=K,IDIM
                  DO I=1,IDIM
                     R(I,K)=0.D0
                     DO J=1,IDIM
*                        R(I,K)=R(I,K)+D(J)*S(J,I)*S(J,K)
                        R(I,K)=R(I,K)+D(J)*S(I,J)*S(K,J)
                     ENDDO
                  ENDDO
               ENDDO
               write(ioimp,*) 'A2,I,J=',IDIM,IDIM
               write(ioimp,*) ((R(I,J),I=1,IDIM),J=1,IDIM)
c$$$               DO J=1,IDIM
c$$$                  DO I=1,IDIM
c$$$                     V(I)=0.D0
c$$$                     DO K=1,IDIM
c$$$                        V(I)=V(I)+A(I,K)*S(K,J)
c$$$                     ENDDO
c$$$                  ENDDO
c$$$                  write(ioimp,*) 'VJ(I)=AIK.SKJ,J',J
c$$$                  write(ioimp,*) (V(I),I=1,IDIM)
c$$$                  write(ioimp,*) 'VJ(I)=D(J).SIJ,J',J
c$$$                  write(ioimp,*) (S(I,J)*D(J),I=1,IDIM)
c$$$               ENDDO
            endif
         ELSE
*     D'abord faire la décomposition polaire avant la décomposition
*     spectrale dans le cas non-symetrique (mais inversible).
*     Pas encore fait
            CALL ERREUR(969)
            RETURN
         ENDIF
      ENDIF
*
* Les opérations tensorielles proprement dites
*
      IF (IOTENS.EQ.1) THEN
         IF (IKAS.EQ.1) THEN
            XNORM=ABS(A(1,1))
         ELSE
            IF (IKAS.EQ.2) THEN
               ID=IDIM
               JD=1
            ELSE
               ID=IDIM
               JD=IDIM
            ENDIF
            XNORM=0.D0
            DO J=1,JD
               DO I=1,ID
                  XNORM=XNORM+A(I,J)**2
               ENDDO
            ENDDO
            XNORM=SQRT(XNORM)
         ENDIF
         R(1,1)=XNORM
      ELSEIF (IOTENS.EQ.2) THEN
         IF (IKAS.EQ.1) THEN
            XNORM=ABS(A(1,1))
         ELSE
            IF (IKAS.EQ.2) THEN
               ID=IDIM
               JD=1
            ELSE
               ID=IDIM
               JD=IDIM
            ENDIF
            XNORM=0.D0
            DO J=1,JD
               DO I=1,ID
                  XNORM=max(xnorm,abs(A(I,J)))
               ENDDO
            ENDDO
         ENDIF
         R(1,1)=XNORM
      ELSEIF (IOTENS.EQ.3.OR.IOTENS.EQ.5) THEN
*
         IF (IKAS.EQ.1) THEN
            JD=1
         ELSE
            JD=IDIM
         ENDIF
         DO J=1,JD
            DO I=1,JD
               R(I,J)=A(I,J)
            ENDDO
         ENDDO
* D sert de vecteur de travail ici
         CALL INVERE(R,JD,3,D,DET)
         IF(IOTENS.EQ.3) THEN
            R(1,1)=DET
         ELSE
            IF (DET.EQ.XZERO) THEN
* Matrice singuliere. Numero de ligne =%i1
               INTERR(1)=0
               CALL ERREUR(49)
*            write(ioimp,*) ((R(I,J),I=1,IDIM),J=1,IDIM)
               RETURN
            ENDIF
         ENDIF
      ELSEIF (IOTENS.EQ.4) THEN
         IF (IKAS.EQ.1) THEN
            XTRACE=A(1,1)
         ELSE
            XTRACE=0.D0
            DO I=1,IDIM
               XTRACE=XTRACE+A(I,I)
            ENDDO
         ENDIF
         R(1,1)=XTRACE
      ELSEIF (IOTENS.EQ.6) THEN
         DO i=1,idim
            DO j=1,idim
               R(J,I)=A(J,I)
           ENDDO
         ENDDO
      ELSEIF (IOTENS.EQ.7) THEN
         DO J=1,JD
            IF (D(J).GE.XPETIT) THEN
               D(J)=LOG(D(J))
            ELSE
               REAERR(1)=D(J)
               REAERR(2)=XPETIT
               CALL ERREUR(1063)
               RETURN
            ENDIF
         ENDDO
      ELSEIF (IOTENS.EQ.8) THEN
*         write(ioimp,*) 'TENS2 D(J) avant:',(D(idx),idx=1,jd)
         DO J=1,JD
            D(J)=EXP(D(J))
         ENDDO
*         write(ioimp,*) 'TENS2 D(J) apres:',(D(idx),idx=1,jd)
      ELSEIF (IOTENS.EQ.9) THEN
         DO J=1,JD
            IF (ABS(D(J)).GE.XPETIT) THEN
               D(J)=1.D0/D(J)
            ELSE
               CALL ERREUR(835)
               RETURN
            ENDIF
         ENDDO
      ELSEIF (IOTENS.EQ.10) THEN
         DO J=1,JD
            D(J)=ABS(D(J))
         ENDDO
      ELSEIF (IOTENS.EQ.11) THEN
            DO J=1,JD
               DO I=1,JD
               R(I,J)=S(I,J)
            ENDDO
         ENDDO
      ELSEIF (IOTENS.EQ.12) THEN
         DO i=1,idim
            DO j=1,idim
               S(J,I)=R(J,I)
           ENDDO
         ENDDO
      ELSEIF (IOTENS.EQ.13) THEN
         DO i=1,idim
            DO j=1,idim
               R(I,J)=A(J,I)
           ENDDO
         ENDDO
      ELSEIF (IOTENS.EQ.14) THEN
         IF (IDIM.EQ.1) THEN
            D(1)=A(1,1)
         ELSEIF (IDIM.EQ.2) THEN
            CALL JACOD2(A,D)
         ELSEIF (IDIM.EQ.3) THEN
            CALL JACOD3(A,3,D)
         ELSE
            write(ioimp,*) 'idim=',idim,' not implemented'
            MOTERR(1:8)='TENS2'
            call erreur(1039)
            return
         ENDIF
      ELSE
         write(ioimp,*) 'iotens=',iotens,' not implemented'
         MOTERR(1:8)='TENS2'
         call erreur(1039)
         return
      ENDIF
* Recomposition spectrale
      IF (IOTENS.GT.6.AND.IOTENS.NE.11.AND.IOTENS.LT.13) THEN
c$$$            write(ioimp,*) 'D,J=',JDIM
c$$$            write(ioimp,*) (D(J),J=1,JDIM)
c$$$            write(ioimp,*) 'IDIM,IKAS=',IDIM,IKAS
         IF (IKAS.EQ.1) THEN
            R(1,1)=D(1)*S(1,1)
         ELSEIF (IKAS.EQ.2) THEN
c$$$            DO I=1,IDIM
c$$$               R(I,1)=D(1)*S(I,1)
c$$$            ENDDO
            CALL ERREUR(980)
            RETURN
         ELSEIF (IKAS.EQ.3) THEN
            DO K=1,IDIM
* Tenir compte de la symétrie               DO I=1,IDIM
               DO I=K,IDIM
                  R(I,K)=0.D0
                  DO J=1,IDIM
*     faux!                        R(I,K)=R(I,K)+D(J)*S(J,I)*S(J,K)
                     R(I,K)=R(I,K)+D(J)*S(I,J)*S(K,J)
                  ENDDO
                  IF (I.NE.K) R(K,I)=R(I,K)
               ENDDO
            ENDDO
         ELSEIF (IKAS.EQ.5) THEN
            DO K=1,IDIM
* Tenir compte de la symétrie               DO I=1,IDIM
               DO I=K,IDIM
                  A(I,K)=0.D0
                  DO J=1,IDIM
                     A(I,K)=A(I,K)+D(J)*S(I,J)*S(K,J)
                  ENDDO
                  IF (I.NE.K) A(K,I)=A(I,K)
               ENDDO
            ENDDO
         ELSE
*     Faire la recomposition polaire après la recomposition
*     spectrale dans le cas non-symetrique (mais inversible).
*     Pas encore fait
            CALL ERREUR(969)
            RETURN
         ENDIF
      ENDIF
      if (ldbg) then
         IF (IOTENS.NE.12) THEN
            write(ioimp,*) 'Output R,I,J=',IDIM,IDIM
            write(ioimp,*) ((R(I,J),I=1,IDIM),J=1,IDIM)
         ELSE
            write(ioimp,*) 'Output A,I,J=',IDIM,IDIM
            write(ioimp,*) ((A(I,J),I=1,IDIM),J=1,IDIM)
         ENDIF
         IF (IOTENS.EQ.11) THEN
            write(ioimp,*) 'Output D,J=',IDIM
            write(ioimp,*) (D(J),J=1,IDIM)
         ENDIF
      endif
*
* Normal termination
*
      RETURN
*
* Format handling
*
*
* Error handling
*
*
* End of subroutine TENS2
*
      END