majzz.f

Go to the documentation of this file.
00001 C:\opentelemac\v7p0\sources\telemac2d\majzz.f
00002 !
00058                         SUBROUTINE MAJZZ
00059 !                       ****************
00060 !
00061      &(W,FLUX,FLUX_OLD,AIRS,DT,NPOIN,CF,KFROT,SMH,
00062      & HN,QU,QV,LT,GAMMA,
00063      & NPTFR,NBOR,LIMPRO,XNEBOR,YNEBOR,KNEU,G)
00064 !
00065 !***********************************************************************
00066 ! TELEMAC2D   V6P3                                         01/07/2013
00067 !***********************************************************************
00068 !
00069 !
00070 !
00071 !
00072 !
00073 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
00074 !|  AIRS          |-->|  CELL AREAS
00075 !|  CF            |-->|  FRICTION COEFFICIENTS
00076 !|  DT            |-->|  TIME STEP
00077 !|  FLUX          |-->|  FLUX AT TN+1
00078 !|  FLUX_OLD      |-->|  FLUX AT TN
00079 !|  G             |-->|  GRAVITY
00080 !|  GAMMA         |-->|  NEWMARK PARAMETER (SEE BELOW)
00081 !|  HN,QU;QV      |-->|  H, HU AND HV AT TN
00082 !|  KFROT         |-->|  LOGICAL! FRICTION OR NO FRICTION
00083 !|  KNEU          |-->|  CONVENTION FOR NEUMANN POINT
00084 !|  LIMPRO        |-->|  BC TYPE
00085 !|  LT            |-->|  CURRENT TIME ITERATION
00086 !|  NBOR          |-->|  GLOBAL INDEX OF BOUNDARY NODES
00087 !|  NPOIN         |-->|  TOTAL NUMBER OF NODES
00088 !|  NPTFR         |-->|  TOTAL NUMNER OF BOUNDARY NODES
00089 !|  SMH           |-->|  MASS SOURCE
00090 !|  XNEBOR,YNEBOR |-->|  X AND Y COMPONENT OF THE OUTWARD UNIT NORMAL
00091 !|  W             |<--|  (H,HU,HV)
00092 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
00093 !
00094       USE INTERFACE_TELEMAC2D, EX_MAJZZ => MAJZZ
00095       IMPLICIT NONE
00096 !
00097       INTEGER LNG,LU
00098       COMMON/INFO/LNG,LU
00099 !
00100 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
00101 !
00102       INTEGER, INTENT(IN)             :: NPOIN,KFROT,LT,NPTFR,KNEU
00103       INTEGER, INTENT(IN)             :: NBOR(NPTFR),LIMPRO(NPTFR,6)
00104       DOUBLE PRECISION, INTENT(IN)    :: XNEBOR(2*NPTFR),YNEBOR(2*NPTFR)
00105       DOUBLE PRECISION, INTENT(INOUT) :: W(3,NPOIN)
00106       DOUBLE PRECISION, INTENT(IN)    :: FLUX(NPOIN,3),DT
00107       DOUBLE PRECISION, INTENT(IN)    :: FLUX_OLD(NPOIN,3),GAMMA
00108       DOUBLE PRECISION, INTENT(IN)    :: AIRS(NPOIN)
00109       DOUBLE PRECISION, INTENT(IN)    :: CF(NPOIN),SMH(NPOIN)
00110       DOUBLE PRECISION, INTENT(IN)    :: HN(NPOIN),QU(NPOIN),QV(NPOIN)
00111       DOUBLE PRECISION, INTENT(IN)    :: G
00112 !
00113 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
00114 !
00115       INTEGER I
00116       DOUBLE PRECISION FACT,UNMGAMMA
00117 !
00118 !
00119 !=======================
00120 !---- TEST FOR DT
00121 !=======================
00122 !
00123       IF(DT.LE.0.D0) THEN
00124         WRITE(LU,*)'*********************************************'
00125         WRITE(LU,*)'          WARNING: TIME STEP =0'
00126         WRITE(LU,*)'          IN MAJZZ SUBROUTINE...'
00127         WRITE(LU,*)'*********************************************'
00128         CALL PLANTE(1)
00129         STOP
00130       ENDIF
00131 !
00132 !++++++++++++++++++++++++++++++++++++
00133 ! TIME INTEGRATION
00134 !++++++++++++++++++++++++++++++++++++
00135 !
00136       IF(GAMMA.EQ.1.D0) THEN
00137 !
00138 !==========================
00139 !---- EULER EXPLICIT SCHEME
00140 !==========================
00141 !
00142         DO I=1,NPOIN
00143           FACT=DT/AIRS(I)
00144           W(1,I) = HN(I) + FACT*(FLUX(I,1)+SMH(I))
00145           W(2,I) = QU(I) + FACT* FLUX(I,2)
00146           W(3,I) = QV(I) + FACT* FLUX(I,3)
00147         ENDDO
00148 !
00149       ELSEIF(GAMMA.GE.0.D0.AND.GAMMA.LT.1.D0) THEN
00150 !
00151 !==========================
00152 !---- NEWMARK SCHEME
00153 !==========================
00154 !
00155         UNMGAMMA = 1.D0-GAMMA
00156 !
00157 !       - FOR GAMMA=0.5, THIS CHOICE GIVES ORDER 2 ACCURACY AND
00158 !         THE SCHEME IS UNCONDITIALLY STABLE
00159 !       - FOR USER WHO PREFERS (EULER) EXPLICIT SCHEME,
00160 !         YOU HAVE TO PUT GAMMA=1
00161         DO I=1,NPOIN
00162           FACT=DT/AIRS(I)
00163           !--- FIRST TIME STEP
00164           IF(LT.EQ.1)THEN
00165             W(1,I) = HN(I) + FACT*(FLUX(I,1)+SMH(I) )
00166             W(2,I) = QU(I) + FACT* FLUX(I,2)
00167             W(3,I) = QV(I) + FACT* FLUX(I,3)
00168           ELSE
00169             W(1,I) = HN(I) + FACT*(UNMGAMMA*FLUX_OLD(I,1) +
00170      &                             GAMMA*FLUX(I,1)+SMH(I))
00171             W(2,I) = QU(I) + FACT*(UNMGAMMA*FLUX_OLD(I,2) +
00172      &                             GAMMA*FLUX(I,2))
00173             W(3,I) = QV(I) + FACT*(UNMGAMMA*FLUX_OLD(I,3) +
00174      &                             GAMMA*FLUX(I,3))
00175           ENDIF
00176         ENDDO
00177 !
00178       ELSE
00179         IF(LNG.EQ.1) THEN
00180            WRITE(LU,*) 'MAJZZ: ERREUR: COEFFICIENT DE NEWMARK DOIT ...'
00181            WRITE(LU,*) '... ETRE ENTRE 0 ET 1: ',GAMMA
00182         ELSEIF(LNG.EQ.2) THEN
00183            WRITE(LU,*) 'MAJZZ: ERROR:NEWMARK COEFFICIENT MUST...'
00184            WRITE(LU,*) '... BE BETWEEN 0 AND 1: ',GAMMA
00185         ENDIF
00186         CALL PLANTE(1)
00187         STOP
00188       ENDIF
00189 !
00190 !     PROJECTION ON THE SLIPPING BOUNDARY CONDITIONS
00191 !     **********************************************
00192 !
00193       CALL CDLPROJ(NPOIN,NPTFR,NBOR,LIMPRO,XNEBOR,YNEBOR,KNEU,W)
00194 !
00195       DO I =1,NPOIN
00196         IF(W(1,I).LE.1.D-12) W(1,I)=0.D0
00197         IF(ABS(W(2,I)).LE.1.D-12) W(2,I)=0.D0
00198         IF(ABS(W(3,I)).LE.1.D-12) W(3,I)=0.D0
00199       ENDDO
00200 !
00201 !     SEMI IMPLICIT FRICTION INTRODUCTION
00202 !     ***********************************
00203 !
00204       IF(KFROT.NE.0) CALL FRICTION(NPOIN,G,DT,W,HN,QU,QV,CF)
00205 !
00206 !-----------------------------------------------------------------------
00207 !
00208       RETURN
00209       END

Generated on Fri Aug 31 2013 18:12:58 by S.E.Bourban (HRW) using doxygen 1.7.0