Sources Esope | Cast3M

hbmrw
C HBMRW     SOURCE    OF166741  26/05/11    21:15:14     12538          
 
*=======================================================================
* Rw: calcul de la derivee: Rw = Z,w(w)*X - P,w
*     TODO? : terme P,w ?
* Zw = 2*w*kron(DL**2,M)+kron(DL,C)
* X : coefficients de Fourier
* w : frequence angulaire
*=======================================================================
 
      SUBROUTINE HBMRW(NT,NDDL,ISD,X,OMEG,XM,XC,XK,Rw,LOGAMO)
 
*----- Declarations ----------------------------------------------------
 
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8(A-H,O-Z)
 
      INTEGER NT,NDDL,ISD
      REAL*8 OMEG
c-      DIMENSION*8 X(NT),Rw(NT)
c-      DIMENSION*8 XM(NDDL,NDDL),XC(NDDL,NDDL),XK(NDDL,NDDL)
      REAL*8 X(*),Rw(*)
      REAL*8 XM(NDDL,*),XC(NDDL,*),XK(NDDL,*)
      LOGICAL LOGAMO
 
      INTEGER H,I,J,M,P,Q,NH
      REAL*8 CF,CF2
 
*       Nombre d'harmoniques
      NH = (NT/NDDL - 1)/2
 
*------ Matrices diagonales --------------------------------------------
      IF (ISD.EQ.1) THEN
 
*         Termes constants
        DO P=1,NDDL
          Rw(P) = 0.D0
        ENDDO
 
c        -cas non conservatif (i.e. avec amortissement)
        IF (LOGAMO) THEN
 
          DO 10 J=1,2*NH
*           Termes pairs (sin)
            IF (MOD(J,2).EQ.0) THEN
              cf = -0.5D0*J
              cf2 = -2.D0*OMEG*(cf**2)
              DO 11 I = 1,NDDL
                Rw(J*nddl+I) = cf2 * XM(I,1) * X(J*nddl+I)
     &                        + cf * XC(I,1) * X((J-1)*nddl+I)
 11           CONTINUE
*              Termes impairs (cos)
            ELSE
              cf = 0.5D0*(J+1)
              cf2 = -2.D0*OMEG*(cf**2)
              DO 12 I = 1,NDDL
                Rw(J*nddl+I) = cf2 * XM(I,1) * X(J*nddl+I)
     &                        + cf * XC(I,1) * X((J+1)*NDDL+I)
 12           CONTINUE
            ENDIF
10        CONTINUE
 
c        -cas conservatif (i.e. sans amortissement)
        ELSE
          DO 20 J=1,2*NH
*           Termes pairs (sin)
            IF (MOD(J,2).EQ.0) THEN
cbp           cf = -0.5D0*J/2
              cf2 = -0.5D0*OMEG*(J**2)
              DO 21 I = 1,NDDL
cbp             Rw(J*nddl+I)= -2.*OMEG*(cf**2)*XM(I,1)*X(J*nddl+I)
                Rw(J*nddl+I) = cf2 * XM(I,1) * X(J*nddl+I)
 21           CONTINUE
*           Termes impairs (cos)
            ELSE
cbp           cf = 0.5D0*(J+1)
              cf2 = -0.5D0*OMEG*((J+1)**2)
              DO 22 I = 1,NDDL
cbp             Rw(J*nddl+I)=-2.*OMEG*(cf**2)*XM(I,1)*X(J*nddl+I)
                Rw(J*nddl+I) = cf2 * XM(I,1) * X(J*nddl+I)
 22           CONTINUE
            ENDIF
 20       CONTINUE
 
        ENDIF
 
*------- Matrices generales --------------------------------------------
      ELSE
 
*           Termes constants
        DO 101 P=1,NDDL
          Rw(P) = 0.D0
 101    CONTINUE
 
c       -cas non conservatif (i.e. avec amortissement)
        IF (LOGAMO) THEN
 
          DO 111 J=2,2*NH,2
            cf2 = 2.D0*OMEG*(J**2)
            DO 112 M=1,NDDL
              Rw(J*NDDL+M) = 0.D0
              DO 113 H=1,NDDL
                Rw(J*NDDL+M) = Rw(J*NDDL+M)
     &                        - cf2 * XM(M,H) * X(J*NDDL+H)
     &                        + J   * XC(M,H) * X((J+1)*NDDL+H)
                Rw((J+1)*NDDL+M) = Rw((J+1)*NDDL+M)
     &                            - cf2 * XM(M,H) * X((J+1)*NDDL+H)
     &                            - J   * XC(M,H) * X(J*NDDL+H)
 113          CONTINUE
 112        CONTINUE
 111      CONTINUE
 
c       -cas conservatif (i.e. sans amortissement)
        ELSE
 
          DO 121 J=2,2*NH,2
            cf2 = 2.D0*OMEG*(J**2)
            DO 122 M=1,NDDL
              Rw(J*NDDL+M) = 0.D0
              DO 123 H=1,NDDL
                Rw(J*NDDL+M) = Rw(J*NDDL+M)
     &                        - cf2 * XM(M,H) * X(J*NDDL+H)
                Rw((J+1)*NDDL+M) = Rw((J+1)*NDDL+M)
     &                            - cf2 * XM(M,H) * X((J+1)*NDDL+H)
 123          CONTINUE
 122        CONTINUE
 121      CONTINUE
 
        ENDIF
 
      ENDIF
 
      RETURN
      END