sisyphe.f

Go to the documentation of this file.
00001 C:\opentelemac\v7p0\sources\sisyphe\sisyphe.f
00002 !
00447                      SUBROUTINE SISYPHE
00448 !                    ******************
00449 !
00450      &(PART,LOOPCOUNT,GRAFCOUNT,LISTCOUNT,TELNIT,
00451      & U_TEL,V_TEL,H_TEL,HN_TEL,ZF_TEL,UETCAR,CF_TEL,KS_TEL,
00452      & CONSTFLOW,NSIS_CFD,SISYPHE_CFD,CODE,PERICOU,
00453      & U3D,V3D,T_TEL,VISC_TEL,DT_TEL,CHARR_TEL,SUSP_TEL,
00454      & FLBOR_TEL,SOLSYS,DM1,UCONV_TEL,VCONV_TEL,ZCONV,
00455      & THETAW_TEL,HW_TEL,TW_TEL,UW_TEL)
00456 !
00457 !***********************************************************************
00458 ! SISYPHE   V7P0                                   02/01/2014
00459 !***********************************************************************
00460 !
00461 !
00462 !
00463 !
00464 !
00465 !
00466 !
00467 !
00468 !
00469 !
00470 !
00471 !
00472 !
00473 !
00474 !
00475 !
00476 !
00477 !
00478 !
00479 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
00480 !| CF_TEL         |<->| QUADRATIC FRICTION COEFFICIENT FROM TELEMAC
00481 !| CHARR_TEL      |<->| LOGICAL, BED LOAD OR NOT: Sent to TELEMAC-2D
00482 !| CODE           |-->| NAME OF CALLING PROGRAMME (TELEMAC2D OR 3D)
00483 !| CONSTFLOW      |<->| LOGICAL, CONSTANT FLOW DISCHARGE OR NOT (A SUPPRIMER)
00484 !| CSRATIO        |<->| EQUILIBRIUM CONCENTRATION FOR SOULSBY-VAN RIJN EQ.
00485 !| DM1            |-->| THE PIECE-WISE CONSTANT PART OF ADVECTION FIELD
00486 !|                |   | IS DM1*GRAD(ZCONV)
00487 !| DT_TEL         |-->| TIME STEP FROM TELEMAC
00488 !| FLBOR_TEL      |-->| FLOW FLUXES AT BOUNDARIES
00489 !| GRAFCOUNT      |-->| PERIOD OF GRAPHICAL OUTPUTS
00490 !| HN_TEL         |-->| WATER DEPTH FROM TEL HN
00491 !| H_TEL          |-->| WATER DEPTH FROM TEL H (N+1)
00492 !| KS_TEL         |-->| BED ROUGHNESS SENT TO TELEMAC
00493 !| LISTCOUNT      |-->| PERIODE DE SORTIE LISTING
00494 !| LOOPCOUNT      |-->| NUMERO DE L'ITERATION
00495 !| UW_TEL         |-->| ORBITAL VELOCITY
00496 !| NSIS_CFD       |---| NUMBER OF ITERATIONS FOR TELEMAC (CONSTANT FLOW DISCHARGE OPTION-TO BE SUPRESSED)
00497 !| PART           |-->| IF -1, NOT COUPLING : ON PASSE TOUTE LA
00498 !|                |   | SUBROUTINE. SINON, INDIQUE LA PARTIE DE LA
00499 !|                |   | SUBROUTINE DANS LAQUELLE ON PASSE
00500 !| PERICOU        |-->| COUPLING PERIOD FOR BEDLOAD
00501 !| SISYPHE_CFD    |<->| LOGICAL, CONSTANT FLOW DISCHARGE OR NOT
00502 !| SOLSYS         |-->|1 OR 2. IF 2 ADVECTION FIELD IS UCONV + DM1*GRAD(ZCONV)
00503 !| SUSP_TEL       |<->|LOGICAL, SUSPENDED LOAD OR NOT: Sent to TELEMAC-2D
00504 !| TELNIT         |-->| NUMBER OF TELEMAC ITERATIONS
00505 !| T_TEL          |-->| CURRENT TIME IN CALLING PROGRAMME
00506 !| U3D,V3D        |-->| 3D VELOCITY SENT BY TELEMAC 3D
00507 !| UCONV_TEL      |-->| ADVECTION VELOCITY FROM TELEMAC2D (X-DIRECTION)
00508 !| UETCAR         |<->| SQUARE OF THE FRICTION VELOCITY (COUPLING WITH TEL 3D)
00509 !| U_TEL          |-->| U VELOCITY FROM TELEMAC
00510 !| VCONV_TEL      |-->| ADVECTION VELOCITY FROM TELEMAC2D (Y-DIRECTION)
00511 !| VISC_TEL       |-->| VELOCITY DIFFUSIVITY (TELEMAC-2D)
00512 !| V_TEL          |-->| V VELOCITY FROM TELEMAC
00513 !| ZCONV          |-->| THE PIECE-WISE CONSTANT PART OF ADVECTION FIELD
00514 !|                |   | IS DM1*GRAD(ZCONV), SEE SOLSYS.
00515 !| ZF_TEL         |<->| BOTTOM ELEVATION OF THE CALLING TELEMAC
00516 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
00517 !
00518       USE INTERFACE_SISYPHE, EX_SISYPHE => SISYPHE
00519       USE BIEF
00520       USE DECLARATIONS_TELEMAC
00521       USE DECLARATIONS_SISYPHE
00522 !
00523       IMPLICIT NONE
00524       INTEGER LNG,LU
00525       COMMON/INFO/LNG,LU
00526 !
00527 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
00528 !
00529       INTEGER,           INTENT(IN)    :: PART,LOOPCOUNT,GRAFCOUNT
00530       INTEGER,           INTENT(IN)    :: LISTCOUNT,TELNIT,PERICOU
00531       CHARACTER(LEN=24), INTENT(IN)    :: CODE
00532       TYPE(BIEF_OBJ),    INTENT(IN)    :: U_TEL,V_TEL,H_TEL,HN_TEL
00533       TYPE(BIEF_OBJ),    INTENT(INOUT) :: ZF_TEL,UETCAR,KS_TEL
00534       INTEGER,           INTENT(INOUT) :: NSIS_CFD
00535       LOGICAL,           INTENT(INOUT) :: CONSTFLOW,SISYPHE_CFD
00536       TYPE(BIEF_OBJ),    INTENT(IN)    :: U3D,V3D,VISC_TEL
00537       TYPE(BIEF_OBJ),    INTENT(INOUT) :: CF_TEL
00538       DOUBLE PRECISION,  INTENT(IN)    :: T_TEL
00539       LOGICAL,           INTENT(INOUT) :: CHARR_TEL,SUSP_TEL
00540       DOUBLE PRECISION,  INTENT(IN)    :: DT_TEL
00541       INTEGER,           INTENT(IN)    :: SOLSYS
00542       TYPE(BIEF_OBJ),    INTENT(IN)    :: FLBOR_TEL,DM1,ZCONV
00543       TYPE(BIEF_OBJ),    INTENT(IN)    :: UCONV_TEL,VCONV_TEL
00544       TYPE(BIEF_OBJ),    INTENT(IN)    :: THETAW_TEL,HW_TEL,TW_TEL
00545       TYPE(BIEF_OBJ),    INTENT(IN)    :: UW_TEL
00546 !
00547 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
00548 !
00549       INTEGER                        P_IMAX
00550       DOUBLE PRECISION P_DMAX,P_DMIN
00551       EXTERNAL         P_DMAX,P_DMIN,P_IMAX
00552 !
00553       INTEGER, PARAMETER :: NHIST = 0
00554       INTEGER, PARAMETER :: NSOR = 100
00555       INTEGER            :: VALNIT,NLISS
00556       INTEGER            :: I,J,K,MN,MT,ISOUS,NIDT,NIT,IMA,IMI,IELEB,KP1
00557       INTEGER            :: IMIN,IMAX,NCALCU,NUMEN,NUMENX,NUMDEB
00558       INTEGER            :: ALIRE(MAXVAR),ALIRV(MAXVAR),ALIRH(MAXVAR)
00559       INTEGER            :: ALIR0(MAXVAR)
00560       INTEGER            :: TROUVE(MAXVAR+10)
00561       DOUBLE PRECISION   :: AT0,DTS,BID,XMA,XMI
00562       DOUBLE PRECISION   :: XMIN,XMAX
00563       DOUBLE PRECISION   :: AT,VCUMU,MASS_GF
00564       DOUBLE PRECISION   :: HIST(1)
00565       LOGICAL            :: PASS,PASS_SUSP
00566       LOGICAL            :: ENTETS,YAZR
00567 !
00568       DOUBLE PRECISION, POINTER, DIMENSION(:) :: SAVEZF,SAVEQU,SAVEQV
00569       DOUBLE PRECISION, POINTER, DIMENSION(:) :: SAVEZ
00570       DOUBLE PRECISION, POINTER, DIMENSION(:) :: SAVEUW
00571 !
00572 !     SAVES LOCAL VARIABLES
00573 !
00574       SAVE VCUMU       ! FOR THE BALANCE
00575       SAVE MASS_GF             ! FOR GRAIN-FEEDING
00576       SAVE PASS, PASS_SUSP     ! IDENTIFIES 1ST TIMESTEP
00577       SAVE NIDT, NCALCU, NUMEN, NIT, VALNIT !
00578       SAVE AT0                 ! TIME
00579 !     NUMEN0 : 1ST RECORD TO READ
00580       INTEGER NUMEN0
00581 !
00582 !     VARIABLES TO READ IF COMPUTATION IS CONTINUED
00583 !     --------------------------------
00584 !     0 : DISCARD
00585 !     1 : READ  (SEE SUBROUTINE NOMVAR)
00586 !
00587 !     HYDRO + EVOLUTION
00588       DATA ALIRE /1,1,1,1,1,1,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,
00589      &            0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,
00590      &            0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,
00591      &            0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,400*0/
00592 !   WAVES ONLY
00593       DATA ALIRH /0,0,0,0,0,0,0,0,0,0,0,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,
00594      &            0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,
00595      &            0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,
00596      &            0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,400*0/
00597 !   NOTHING TO READ
00598       DATA ALIR0 /0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,
00599      &            0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,
00600      &            0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,
00601      &            0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,400*0/
00602 !
00603 !    FOR VALIDATION, EACH VARIABLE IN THE FILE IS COMPARED
00604 !
00605       DATA ALIRV /500*1/
00606 !
00607 !======================================================================!
00608 !======================================================================!
00609 !                               PROGRAM                                !
00610 !======================================================================!
00611 !======================================================================!
00612 !
00613 !------------------------------------------------------------------
00614 !     PART 1 : INITIALISATION
00615 !------------------------------------------------------------------
00616 !
00617       PERCOU = PERICOU !UHM!!!!!!!
00618 !
00619       IF(PART==0.OR.PART==-1) THEN
00620         IF(DEBUG.GT.0) WRITE(LU,*) 'INITIALIZATION'
00621 !
00622         WRITE(LU,*) 'PART 0 : INITIALISING SISYPHE'
00623 !
00624 !       INITIALISES THE CONSTANTS
00625 !
00626         CALL INIT_CONSTANT(KARIM_HOLLY_YANG,KARMAN,PI)
00627 !
00628 !      READS THE WAVE DATA IN THE HYDRODYNAMIC FILE
00629 !
00630         IF(HOULE.AND.SIS_FILES(SISCOU)%NAME(1:1).EQ.' ') THEN
00631           ALIRE(12)=1
00632           ALIRE(13)=1
00633           ALIRE(14)=1
00634           ALIRE(22)=1
00635         ENDIF
00636 !
00637 !       READS THE SEDIMENTOLOGICAL DATA IN THE CONTINUATION FILE
00638 !
00639         IF(DEBU) THEN
00640           ALIRE(15)=1
00641           ALIRE(16)=1
00642           ALIRE(17)=1
00643           ALIRE(18)=1
00644           ALIRE(19)=1
00645           ALIRE(20)=1
00646           ALIRE(21)=1
00647 !         READS AVAI FROM THE PREVIOUS COMPUTATION FILE
00648           DO I=1,NSICLA*NOMBLAY
00649             ALIRE(22+I)=1
00650           ENDDO
00651 !         READS CS (CONCENTRATION) FROM THE PREVIOUS COMPUTATION FILE
00652           IF(SUSP) THEN
00653             DO I=1,NSICLA
00654              ALIRE(22+(NOMBLAY+1)*NSICLA+I)=1
00655             ENDDO
00656           ENDIF
00657 !         READS THE LAYER THICKNESSES
00658           DO I=1,NOMBLAY
00659             ALIRE(28+(NOMBLAY+4)*NSICLA+I)=1
00660           ENDDO
00661         ENDIF
00662 !       V6P2 CV lecture concentration des couches
00663         IF(TASS) THEN
00664           DO I=1,NOMBLAY
00665             ALIRE(28+(NOMBLAY+4)*NSICLA+I+NOMBLAY)=1
00666           ENDDO
00667         ENDIF
00668 !
00669 !       INITIALISING DEPOSITS ON THE BOTTOM
00670 !       THEY MUST BE USED IN BILAN_SISYPHE EVEN IF SUSP=.FALSE.
00671 !
00672         DO I=1,NSICLA
00673           MASDEP(I)=0.D0
00674           MASDEPT(I)=0.D0
00675         ENDDO
00676 !
00677 ! --------  INITIALISES (SETS TO 0) THE ARRAYS
00678 !
00679         CALL INIT_ZERO
00680 !
00681 ! --------  END OF INITIALISATION
00682 !
00683 !       DISCRETISATION : LINEAR FOR THE TIME BEING
00684 !
00685 !       IELMT HARD-CODED IN LECDON
00686 !
00687         IELMX  = MAX(IELMT,IELMH_SIS,IELMU_SIS)
00688         NELMAX = NELEM
00689 !
00690 !=======================================================================
00691 !
00692 ! : 1          READS, PREPARES AND CONTROLS THE DATA
00693 !
00694 !=======================================================================
00695 !
00696 !       READS THE BOUNDARY CONDITIONS AND INDICES FOR THE BOUNDARY NODES
00697 !       EBOR IS READ HERE FOR THE FIRST CLASS ONLY
00698 !       SEE CONLIT FOR OTHER CLASSES
00699 !
00700         IF(DEBUG.GT.0) WRITE(LU,*) 'LECLIS'
00701         CALL LECLIS(LIEBOR%I,LIHBOR%I,LIQBOR%I,EBOR%ADR(1)%P,Q2BOR,
00702      &              MESH%NPTFR,MESH%NBOR%I,3,
00703      &              SIS_FILES(SISCLI)%LU,KENT,
00704      &              MESH%ISEG%I,MESH%XSEG%R,MESH%YSEG%R,MESH%NACHB%I,
00705      &              NUMLIQ%I,NSICLA,AFBOR%R,BFBOR%R,BOUNDARY_COLOUR%I,
00706      &              MESH)
00707         IF(DEBUG.GT.0) WRITE(LU,*) 'END_LECLIS'
00708 !
00709 !-----------------------------------------------------------------------
00710 !
00711 !       COMPLEMENTS THE DATA STRUCTURE FOR BIEF
00712 !
00713         IF(DEBUG.GT.0) WRITE(LU,*) 'INBIEF'
00714         CALL INBIEF(IT1%I,KLOG,IT2,IT3,IT4,LVMAC,IELMX,
00715      &                 0.D0,SPHERI,MESH,T1,T2,OPTASS,PRODUC,EQUA)
00716 !
00717         IF(DEBUG.GT.0) WRITE(LU,*) 'END_INBIEF'
00718 !
00719 !-----------------------------------------------------------------------
00720 !
00721 !       LOCATES THE BOUNDARIES
00722 !
00723         IF(NCSIZE.GT.1) THEN
00724           NFRLIQ=0
00725           DO I=1,NPTFR
00726             NFRLIQ=MAX(NFRLIQ,NUMLIQ%I(I))
00727           ENDDO
00728           NFRLIQ=P_IMAX(NFRLIQ)
00729           WRITE(LU,*) ' '
00730           IF(LNG.EQ.1) WRITE(LU,*) 'FRONTIERES LIQUIDES :',NFRLIQ
00731           IF(LNG.EQ.2) WRITE(LU,*) 'LIQUID BOUNDARIES:',NFRLIQ
00732         ELSE
00733           IF(DEBUG.GT.0) WRITE(LU,*) 'APPEL DE FRONT2'
00734           CALL FRONT2(NFRLIQ,NFRSOL,DEBLIQ,FINLIQ,DEBSOL,FINSOL,
00735      &                LIEBOR%I,LIEBOR%I,
00736      &                MESH%X%R,MESH%Y%R,MESH%NBOR%I,MESH%KP1BOR%I,
00737      &                IT1%I,NPOIN,NPTFR,KLOG,.TRUE.,NUMLIQ%I,MAXFRO)
00738           IF(DEBUG.GT.0) WRITE(LU,*) 'RETOUR DE FRONT2'
00739         ENDIF
00740 !
00741 !-----------------------------------------------------------------------
00742 !       LOOKS FOR BOTTOM AND BOTTOM FRICTION IN THE GEOMETRY FILE :
00743 !-----------------------------------------------------------------------
00744 !
00745         IF(DEBUG.GT.0) WRITE(LU,*) 'FONSTR'
00746         CALL FONSTR(T1,ZF,T2,CHESTR,SIS_FILES(SISGEO)%LU,
00747      &              SIS_FILES(SISGEO)%FMT,
00748      &              SIS_FILES(SISFON)%LU,SIS_FILES(SISFON)%NAME,
00749      &              MESH,SFON,.TRUE.)
00750         IF(DEBUG.GT.0) WRITE(LU,*) 'END_FONSTR'
00751 !
00752 !-----------------------------------------------------------------------
00753 !
00754 !       PREPARES THE RESULTS FILE (OPTIONAL)
00755 !
00756 !       STANDARD SELAFIN FORMAT
00757 !
00758         IF(DEBUG.GT.0) WRITE(LU,*) 'ECRGEO'
00759 !       CREATES DATA FILE USING A GIVEN FILE FORMAT : FORMAT_RES
00760 !       THE DATA ARE CREATED IN THE FILE: SISRES, AND ARE
00761 !       CHARACTERISED BY A TITLE AND NAME OF OUTPUT VARIABLES
00762 !       CONTAINED IN THE FILE.
00763         CALL CREATE_DATASET(SIS_FILES(SISRES)%FMT, ! RESULTS FILE FORMAT
00764      &                      SIS_FILES(SISRES)%LU,  ! LU FOR RESULTS FILE
00765      &                      TITCA,      ! TITLE
00766      &                      MAXVAR,     ! MAX NUMBER OF OUTPUT VARIABLES
00767      &                      TEXTE,      ! NAMES OF OUTPUT VARIABLES
00768      &                      SORLEO)     ! PRINT TO FILE OR NOT
00769 !       WRITES THE MESH IN THE OUTPUT FILE :
00770 !       IN PARALLEL, REQUIRES NCSIZE AND NPTIR.
00771 !       THE REST OF THE INFORMATION IS IN MESH.
00772 !       ALSO WRITES : START DATE/TIME AND COORDINATES OF THE
00773 !       ORIGIN.
00774         CALL WRITE_MESH(SIS_FILES(SISRES)%FMT, ! RESULTS FILE FORMAT
00775      &                  SIS_FILES(SISRES)%LU,  ! LU FOR RESULTS FILE
00776      &                  MESH,          ! CHARACTERISES MESH
00777      &                  1,             ! NUMBER OF PLANES /NA/
00778      &                  MARDAT,        ! START DATE
00779      &                  MARTIM,        ! START TIME
00780      &                  I_ORIG,J_ORIG) ! COORDINATES OF THE ORIGIN.
00781         IF(DEBUG.GT.0) WRITE(LU,*) 'END_ECRGEO'
00782 !
00783 !   --- FILLS IN MASKEL BY DEFAULT
00784 !
00785         IF(MSK) CALL OS ('X=C     ', X=MASKEL, C=1.D0)
00786 !
00787 !       BUILDING THE MASK FOR LIQUID BOUNDARIES
00788 !       A SEGMENT IS LIQUID IF BOTH ENDS ARE NOT SOLID
00789 !
00790         DO IELEB = 1, MESH%NELEB
00791           K=MESH%IKLBOR%I(IELEB)
00792           KP1=MESH%IKLBOR%I(IELEB+MESH%NELEBX)
00793           IF(LIEBOR%I(K).NE.2.AND.LIEBOR%I(KP1).NE.2) THEN
00794             MASK%R(IELEB) = 1.D0
00795           ELSE
00796             MASK%R(IELEB) = 0.D0
00797           ENDIF
00798         ENDDO
00799 !
00800 !=======================================================================
00801 !
00802 ! : 2                  INITIALISES
00803 !
00804 !=======================================================================
00805 !
00806         PASS      = .TRUE.
00807 !
00808         PASS_SUSP = .TRUE.
00809         VCUMU     = 0.D0
00810         MASS_GF   = 0.D0
00811 !
00812 !
00813 !----   DETERMINES THE NUMBER OF EVENTS (1ST LOOP)       : NCALCU
00814 !                      NUMBER OF TIMESTEPS (2ND LOOP)    : NIDT
00815 !                      TOTAL NUMBER OF TIMESTEPS         : NIT
00816 !                      INITIAL TIME                      : AT0
00817 !                      TIMESTEP                          : DT
00818 !
00819 !
00820         IF(PART.EQ.0) THEN
00821           AT0=T_TEL
00822           DT = DT_TEL
00823           NCALCU = 1
00824           NIDT   = 1
00825           NIT=TELNIT
00826         ELSE
00827           AT0=MAX(TPREC,0.D0)
00828           DT=DELT
00829           IF(PERMA) THEN
00830             NCALCU=1
00831             NIDT=NPAS
00832             NIT=NIDT
00833             NSOUS=1
00834           ELSE
00835             NCALCU = NMAREE
00836 ! COMPUTES DT AFTER READING THE HYDRO FILE
00837 !           NIDT =  NINT ( PMAREE / DT + 0.1D0 )
00838 !           NIT=NIDT*NCALCU
00839 ! ELSE
00840             NIDT=NPAS
00841             NIT=NIDT*NCALCU
00842           ENDIF
00843         ENDIF
00844 !
00845 !       GETTING NUMEN: TOTAL NUMBER OF RECORDS, AND THE TIME STEP
00846 !
00847         IF(SIS_FILES(SISHYD)%NAME(1:1).NE.' ')  THEN
00848           IF(PERMA) THEN
00849 !           STEADY MODE
00850 !           TO GET NUMEN ONLY
00851             IF(DEBUG.GT.0) WRITE(LU,*) 'BIEF_SUITE'
00852             CALL BIEF_SUITE(VARSOR,VARCL,NUMEN,SIS_FILES(SISHYD)%LU,
00853      &                      SIS_FILES(SISHYD)%FMT,HIST,0,NPOIN,AT,
00854      &                      TEXTPR,VARCLA,
00855      &                      0,TROUVE,ALIR0,.TRUE.,.TRUE.,MAXVAR)
00856             DT=DELT
00857             IF(DEBUG.GT.0) WRITE(LU,*) 'SORTIE DE BIEF_SUITE'
00858           ELSE
00859 !           UNSTEADY MODE : DT IS TAKEN FROM THE HYDRO FILE
00860 !           TO GET NUMEN AND DT
00861             IF(DEBUG.GT.0) WRITE(LU,*) 'BIEF_SUITE'
00862             CALL BIEF_SUITE(VARSOR,VARCL,NUMEN,SIS_FILES(SISHYD)%LU,
00863      &                      SIS_FILES(SISHYD)%FMT,HIST,0,NPOIN,AT,
00864      &                      TEXTPR,VARCLA,
00865      &                      0,TROUVE,ALIR0,.TRUE.,.TRUE.,MAXVAR,DT=DT)
00866             IF(DEBUG.GT.0) WRITE(LU,*) 'SORTIE DE BIEF_SUITE'
00867             NIDT =  NINT ( PMAREE / DT + 0.1D0 )
00868             IF(ABS(NIDT*DT-PMAREE) > 1.D-3) THEN
00869               IF(LNG.EQ.1) WRITE(LU,101) NIDT*DT
00870               IF(LNG.EQ.2) WRITE(LU,102) NIDT*DT
00871             ENDIF
00872             NIT  = NCALCU * NIDT
00873           ENDIF
00874         ENDIF
00875 !
00876 !       VALIDATES AGAINST THE LAST TIMESTEP
00877 !
00878         VALNIT = NIT
00879 !
00880 101     FORMAT(/,'ATTENTION : LA PERIODE DE CALCUL NE CORRESPOND PAS A',
00881      &       /,'UN MULTIPLE DE LA PERIODE DE SORTIE HYDRODYNAMIQUE.',
00882      &       /,'LE CALCUL S''EFFECTUERA DONC SUR ',G16.7,'SECONDES')
00883 102     FORMAT(/,
00884      &         'CAUTION : THE PERIOD OF COMPUTATION IS NOT A MULTIPLE',
00885      &         /,'OF THE HYDRODYNAMIC FILE PRINTOUT PERIOD.',/,
00886      &       'THE LENGTH OF COMPUTATION WILL THEREFORE BE',G16.7,/,
00887      &       'SECONDS')
00888 !
00889 !  SISYPHE ONLY
00890 !  -----------------------------------------------------------------------
00891 !
00892 !  NUMEN : NUMBER OF RECORDS IN THE HYDRODYNAMIC FILE
00893 !  DT    : TIMESTEP OF THE HYDRODYNAMIC RECORDS
00894 !  NUMEN0: 1ST RECORD TO READ FROM HYDRODYNAMIC FILE
00895 !  TPREC : START TIME
00896 !
00897         IF(PART.EQ.-1) THEN
00898           IF(.NOT.PERMA) THEN
00899             IF(TPREC.GE.0.D0) THEN
00900               NUMEN0 = NINT((TPREC-AT0)/DT)+1
00901             ELSE
00902               NUMEN0 = NUMEN-INT(PMAREE/DT+1.1D0)
00903             ENDIF
00904           ELSE
00905             IF(TPREC.GE.0.D0) THEN
00906               NUMEN0 = NINT((TPREC-AT0)/DT)+1
00907             ELSE
00908               NUMEN0 = NUMEN
00909             ENDIF
00910           ENDIF
00911 !
00912           IF(SIS_FILES(SISHYD)%NAME(1:1).NE.' ')  THEN
00913             IF(DEBUG.GT.0) WRITE(LU,*) 'BIEF_SUITE'
00914             CALL BIEF_SUITE(VARSOR,VARCL,NUMEN0,SIS_FILES(SISHYD)%LU,
00915      &                      SIS_FILES(SISHYD)%FMT,HIST,0,NPOIN,AT,
00916      &                      TEXTPR,VARCLA,
00917      &                      0,TROUVE,ALIRE,.TRUE.,PERMA,MAXVAR)
00918             IF(DEBUG.GT.0) WRITE(LU,*) 'RETOUR DE BIEF_SUITE'
00919 !
00920 !           TRACES IF WAVE DATA HAVE BEEN FOUND
00921 !
00922             IF(HOULE) THEN
00923 !             NOTE: BIEF_ALLVEC SETS %TYPR TO '?'
00924               IF(TROUVE(12).EQ.1) HW%TYPR='Q'
00925               IF(TROUVE(13).EQ.1) TW%TYPR='Q'
00926               IF(TROUVE(14).EQ.1) THETAW%TYPR='Q'
00927               IF(TROUVE(22).EQ.1) UW%TYPR='Q'
00928               IF(UW%TYPR=='Q') THEN
00929                 WRITE(LU,*)
00930                 WRITE(LU,*) 'WAVE RESULTS IN :',SIS_FILES(SISHYD)%NAME
00931                 WRITE(LU,*)
00932                 WRITE(LU,*) 'BOTTOM ORBITAL VELOCITY FOUND'
00933                 WRITE(LU,*) 'THESE WILL BE USED DIRECTLY'
00934                 WRITE(LU,*)
00935               ELSEIF(HW%TYPR=='Q'.AND.TW%TYPR=='Q') THEN
00936                 WRITE(LU,*)
00937                 WRITE(LU,*) 'WAVE RESULTS IN :',SIS_FILES(SISHYD)%NAME
00938                 WRITE(LU,*)
00939                 WRITE(LU,*) 'WAVE HEIGHT AND PERIOD FOUND'
00940                 WRITE(LU,*) 'BOTTOM VELOCITY WILL BE COMPUTED IN CALCUW'
00941                 WRITE(LU,*)
00942               ENDIF
00943             ENDIF
00944             IF(DEBUG.GT.0) WRITE(LU,*) 'END_BIEF_SUITE'
00945             IF(DEBUG.GT.0) WRITE(LU,*) 'RESCUE_SISYPHE'
00946             CALL RESCUE_SISYPHE(QU%R,QV%R,Q%R,U2D%R,V2D%R,HN%R,Z%R,
00947      &                        ZF%R,HW%R,TW%R,THETAW%R,NPOIN,
00948      &                        TROUVE,ALIRE,PASS,ICF,.TRUE.,MAXVAR)
00949             IF(DEBUG.GT.0) WRITE(LU,*) 'END_RESCUE_SISYPHE'
00950           ENDIF
00951 !
00952         ENDIF
00953 !
00954 !---- RESUMES SISYPHE COMPUTATION
00955 !
00956         YAZR=.FALSE.
00957         IF(SIS_FILES(SISPRE)%NAME(1:1).NE.' ')  THEN
00958 !
00959 !         READS THE HYDRO AND SEDIMENTOLOGICAL VARIABLES
00960 !
00961           IF(DEBUG.GT.0) WRITE(LU,*) 'BIEF_SUITE'
00962           CALL BIEF_SUITE(VARSOR,VARCL,NUMENX,SIS_FILES(SISPRE)%LU,
00963      &                    SIS_FILES(SISPRE)%FMT,
00964      &                    HIST,0,NPOIN,AT0,TEXTPR,VARCLA,0,
00965      &                    TROUVE,ALIRE,.TRUE.,.TRUE.,MAXVAR)
00966           IF(TROUVE(9).EQ.1) YAZR=.TRUE.
00967           IF(DEBUG.GT.0) WRITE(LU,*) 'END_BIEF_SUITE'
00968 !
00969           IF(DEBUG.GT.0) WRITE(LU,*) 'RESCUE_SISYPHE'
00970           CALL RESCUE_SISYPHE(QU%R,QV%R,Q%R,U2D%R,V2D%R,HN%R,Z%R,ZF%R,
00971      &                        HW%R,TW%R,THETAW%R,NPOIN,TROUVE,ALIRE,
00972      &                        PASS,ICF,.TRUE.,MAXVAR)
00973           IF(DEBUG.GT.0) WRITE(LU,*) 'SORTIE DE BIEF_SUITE'
00974 !
00975 !         CHANGES THE UNITS OF CONCENTRATIONS
00976 !
00977           IF(SUSP.AND.UNIT) THEN
00978             DO I=1,NSICLA
00979               IF(TROUVE(21+(NOMBLAY+1)*NSICLA+I).EQ.1) THEN
00980                 CALL OS('X=CX    ',X=CS%ADR(I)%P,C=1.D0/XMVS)
00981               ENDIF
00982             ENDDO
00983           ENDIF
00984 !
00985         ENDIF
00986 !
00987 !----   READS THE LAST RECORD : WAVE FILE
00988 !
00989 !       NOTE : SIS_FILES(SISCOU)%NAME SET TO ' ' IF HOULE=NO
00990 !
00991         IF(SIS_FILES(SISCOU)%NAME(1:1).NE.' ')  THEN
00992 !
00993           IF(DEBUG.GT.0) WRITE(LU,*) 'SUITE_HOULE'
00994           WRITE(LU,*) ' LECTURE HOULE :',SIS_FILES(SISCOU)%NAME
00995           CALL BIEF_SUITE(VARSOR,VARCL,NUMENX,SIS_FILES(SISCOU)%LU,
00996      &                    SIS_FILES(SISCOU)%FMT,HIST,0,NPOIN,AT,
00997      &                    TEXTPR,VARCLA,0,
00998      &                    TROUVE,ALIRH,.TRUE.,.TRUE.,MAXVAR)
00999           IF(DEBUG.GT.0) WRITE(LU,*) 'END_SUITE_HOULE'
01000 !         TRACES IF WAVE DATA HAVE BEEN FOUND
01001           IF(TROUVE(12).EQ.1) HW%TYPR='Q'
01002           IF(TROUVE(13).EQ.1) TW%TYPR='Q'
01003           IF(TROUVE(14).EQ.1) THETAW%TYPR='Q'
01004 !
01005         ENDIF
01006 !
01007         IF(CODE(1:7) == 'TELEMAC'.AND.PART==0) THEN
01008 !
01009           AT0=T_TEL
01010 !
01011           WRITE(LU,*) 'INITIALISATION EN CAS DE COUPLAGE : PART=',PART
01012 !         INFORMATION ON SUSPENSION SENT BACK
01013           CHARR_TEL = CHARR
01014           SUSP_TEL = SUSP
01015 !
01016 !         OV INSTEAD OF OS IN ORDER TO AVOID PROBLEMS WITH QUASI-BUBBLE ELEMENTS
01017 !         OPERATES ONLY ON THE (1:NPOIN) RANGE OF THE TELEMAC FIELDS
01018 !         IT IS A HIDDEN DISCRETISATION CHANGE
01019 !
01020           CALL OV( 'X=Y     ', U2D%R, U_TEL%R, U_TEL%R, 0.D0, NPOIN)
01021           CALL OV( 'X=Y     ', V2D%R, V_TEL%R, V_TEL%R, 0.D0, NPOIN)
01022           CALL OV( 'X=Y     ', HN%R , H_TEL%R, H_TEL%R, 0.D0, NPOIN)
01023 !         BOTTOM GIVEN BY CALLING PROGRAMME
01024           CALL OS( 'X=Y     ', X=ZF,Y=ZF_TEL)
01025 !
01026 !         CLIPS NEGATIVE DEPTHS
01027 !
01028           IF(OPTBAN.GT.0) THEN
01029             DO I = 1,NPOIN
01030               IF(HN%R(I).LT.HMIN) THEN
01031                 U2D%R(I)=0.D0
01032                 V2D%R(I)=0.D0
01033                 HN%R(I) = HMIN
01034               ENDIF
01035             ENDDO
01036           ENDIF
01037 !
01038 !         CASE OF TRIPLE COUPLING
01039 !
01040           IF(INCLUS(COUPLING,'TOMAWAC')) THEN
01041 !           incident wave direction
01042             CALL OS( 'X=Y     ',THETAW,THETAW_TEL)
01043 !           Wave period
01044             CALL OS( 'X=Y     ', TW, TW_TEL)
01045 !           significant wave height
01046             CALL OS( 'X=Y     ', HW , HW_TEL)
01047 !           Bottom orbital velocity
01048             CALL OS( 'X=Y     ', UW , UW_TEL)
01049             HW%TYPR='Q'
01050             TW%TYPR='Q'
01051             THETAW%TYPR='Q'
01052             UW%TYPR='Q'
01053           ENDIF
01054 !
01055         ENDIF
01056 !
01057 !  ---- END COUPLING  -------------
01058 !
01059         IF(DEBUG.GT.0) WRITE(LU,*) 'CONDIM_SISYPHE'
01060         IF(.NOT.DEBU) THEN
01061         CALL CONDIM_SISYPHE
01062      &        (U2D%R,V2D%R,QU%R,QV%R,HN%R,ZF%R,Z%R,ESOMT%R,THETAW%R,
01063      &         Q%R,HW%R,TW%R,MESH%X%R,MESH%Y%R,NPOIN,AT0,PMAREE)
01064         ENDIF
01065         IF(DEBUG.GT.0) WRITE(LU,*) 'END_CONDIM_SISYPHE'
01066 !
01067 !       AT THIS LEVEL U2D,V2D,H AND ZF MUST HAVE BEEN DEFINED
01068 !       EITHER BY BIEF_SUITE, CONDIM_SISYPHE OR CALLING PROGRAM
01069 !
01070 !       NOW COMPUTES FUNCTIONS OF U2D,V2D,HN AND ZF
01071 !
01072 !       FREE SURFACE
01073         IF(DEBUG.GT.0) WRITE(LU,*) 'FREE SURFACE'
01074         CALL OS('X=Y+Z   ', X=Z, Y=ZF, Z=HN)
01075         IF(DEBUG.GT.0) WRITE(LU,*) 'END FREE SURFACE'
01076 !
01077         IF(CODE(1:7).NE.'TELEMAC') THEN
01078 !         PRODUCT H*
01079           CALL OS('X=YZ    ', X=QU, Y=U2D, Z=HN)
01080 !         PRODUCT H*V
01081           CALL OS('X=YZ    ', X=QV, Y=V2D, Z=HN)
01082 !         DISCHARGE
01083           CALL OS('X=N(Y,Z)', X=Q, Y=QU, Z=QV)
01084         ENDIF
01085 !
01086 !       CHECKS THE WAVE DATA
01087 !
01088         IF(HOULE) THEN
01089           IF(UW%TYPR.EQ.'Q') THEN
01090 !           IF(HW%TYPR    .NE.'Q'.OR.
01091           ELSEIF(HW%TYPR    .NE.'Q'.OR.
01092      &      TW%TYPR    .NE.'Q'.OR.
01093      &      THETAW%TYPR.NE.'Q') THEN
01094             WRITE(LU,*) ' '
01095             WRITE(LU,*) ' '
01096             IF(LNG.EQ.1) THEN
01097               WRITE(LU,*) 'DONNEES DE HOULE MANQUANTES'
01098               IF(HW%TYPR.NE.'Q') WRITE(LU,*) 'HAUTEUR HM0'
01099               IF(TW%TYPR.NE.'Q') WRITE(LU,*) 'PERIODE PIC TPR5'
01100               IF(THETAW%TYPR.NE.'Q') WRITE(LU,*) 'DIRECTION MOY'
01101             ENDIF
01102             IF(LNG.EQ.2) THEN
01103               WRITE(LU,*) 'MISSING WAVE DATA'
01104               IF(HW%TYPR.NE.'Q') WRITE(LU,*) 'WAVE HEIGHT HM0'
01105               IF(TW%TYPR.NE.'Q') WRITE(LU,*) 'PEAK PERIOD TPR5'
01106               IF(THETAW%TYPR.NE.'Q') WRITE(LU,*) 'MEAN DIRECTION'
01107             ENDIF
01108             CALL PLANTE(1)
01109             STOP
01110           ENDIF
01111         ENDIF
01112 !
01113 ! END OF HYDRODYNAMIC INITIALISATION
01114 !
01115 !
01116 !        COMPUTES AREAS (WITHOUT MASKING)
01117 !
01118         IF(DEBUG.GT.0) WRITE(LU,*) 'VECTOR FOR VOLU2D'
01119         CALL VECTOR(VOLU2D,'=','MASBAS          ',
01120      &              IELMH_SIS,1.D0,
01121      &              T1,T1,T1,T1,T1,T1,MESH,.FALSE.,MASKEL)
01122         IF(DEBUG.GT.0) WRITE(LU,*) 'END VECTOR FOR VOLU2D'
01123 !       V2DPAR : LIKE VOLU2D BUT IN PARALLEL VALUES COMPLETED AT
01124 !                INTERFACES BETWEEN SUBDOMAINS
01125         CALL OS('X=Y     ',X=V2DPAR,Y=VOLU2D)
01126         IF(NCSIZE.GT.1) CALL PARCOM(V2DPAR,2,MESH)
01127 !       INVERSE OF VOLUMES (DONE WITHOUT MASKING)
01128         CALL OS('X=1/Y   ',X=UNSV2D,Y=V2DPAR,
01129      &          IOPT=2,INFINI=0.D0,ZERO=1.D-12)
01130 !
01131 !       START OF MODIFICATIONS FOR MIXED SEDIMENTS
01132 !
01133 !       SETTING THE NON-ERODABLE BED (IT CAN BE SET BEFORE
01134 !                                     IF COMPUTATION CONTINUED, I.E. DEBU)
01135 !
01136         IF(.NOT.DEBU.OR..NOT.YAZR) THEN
01137           IF(DEBUG.GT.0) WRITE(LU,*) 'NOEROD'
01138           CALL NOEROD(HN%R,ZF%R,ZR%R,Z%R,MESH%X%R,
01139      &                MESH%Y%R,NPOIN,CHOIX,NLISS)
01140           IF(DEBUG.GT.0) WRITE(LU,*) 'END NOEROD'
01141         ENDIF
01142 !
01143 !       INITIALISATION FOR SEDIMENT
01144 !
01145         IF(DEBUG.GT.0) WRITE(LU,*) 'INIT_SEDIMENT'
01146         CALL INIT_SEDIMENT(NSICLA,ELAY,ZF,ZR,NPOIN,
01147      &                     AVAIL,FRACSED_GF,AVA0,LGRAFED,CALWC,
01148      &                     XMVS,XMVE,GRAV,VCE,XWC,FDM,CALAC,AC,
01149      &                     SEDCO,ES,ES_SABLE,ES_VASE,
01150      &                     NCOUCH_TASS,CONC_VASE,
01151      &                     MS_SABLE, MS_VASE,ACLADM,
01152      &                     UNLADM,TOCE_SABLE,
01153      &                     CONC,NLAYER,DEBU,MIXTE)
01154         IF(DEBUG.GT.0) WRITE(LU,*) 'END INIT_SEDIMENT'
01155 !
01156 !
01157 ! END OF MODIFICATIONS CV
01158 !
01159 !
01160 ! MEAN VELOCITY
01161 !======================================================================
01162         CALL OS('X=N(Y,Z)',X=UNORM,Y=U2D,Z=V2D)
01163 ! =====================================================================
01164 !  WAVE ORBITAL VELOCITY
01165 ! =====================================================================
01166         IF(HOULE) THEN
01167 !         JWI 31/05/2012 - added lines to use wave orbital velocities
01168 !         directly if found in hydro file; otherwise compute with CALCUW
01169           IF(UW%TYPR.NE.'Q') THEN
01170             CALL CALCUW(UW%R,HN%R,HW%R,TW%R,GRAV,NPOIN)
01171           ENDIF
01172         ENDIF
01173 !
01174 ! =====================================================================
01175 !
01176         IF(DEBUG.GT.0) WRITE(LU,*) 'TOB_SISYPHE'
01177         CALL TOB_SISYPHE(TOB,TOBW,MU,KS,KSP,KSR,CF,FW,
01178      &                   CHESTR,UETCAR,CF_TEL,KS_TEL,CODE ,
01179      &                   KFROT,ICR,KSPRATIO,HOULE,
01180      &                   GRAV,XMVE,XMVS,VCE,KARMAN,ZERO,
01181      &                   HMIN,HN,ACLADM,UNORM,UW,TW,NPOIN,KSPRED,IKS)
01182         IF(DEBUG.GT.0) WRITE(LU,*) 'END TOB_SISYPHE'
01183 !
01184 !       INITIALISATION FOR TRANSPORT
01185 !
01186         IF(DEBUG.GT.0) WRITE(LU,*) 'INIT_TRANSPORT'
01187         CALL INIT_TRANSPORT(TROUVE,DEBU,HIDING,NSICLA,NPOIN,
01188      &     T1,T2,T3,T4,T5,T6,T7,T8,T9,T10,T11,T12,T14,
01189      &     CHARR,QS_C,QSXC,QSYC,CALFA,SALFA,COEFPN,SLOPEFF,
01190      &     SUSP,QS_S,QS,QSCL,QSCL_C,QSCL_S,QSCLXS,QSCLYS,
01191      &     UNORM,U2D,V2D,HN,CF,MU,TOB,TOBW,UW,TW,THETAW,FW,HOULE,
01192      &     AVAIL,ACLADM,UNLADM,KSP,KSR,KS,
01193      &     ICF,HIDFAC,XMVS,XMVE,GRAV,VCE,HMIN,KARMAN,
01194      &     ZERO,PI,AC,IMP_INFLOW_C,ZREF,ICQ,CSTAEQ,CSRATIO,
01195      &     CMAX,CS,CS0,UCONV,VCONV,CORR_CONV,SECCURRENT,BIJK,
01196      &     IELMT,MESH,FDM,XWC,FD90,SEDCO,VITCE,PARTHENIADES,VITCD,
01197      &     U3D,V3D,CODE)
01198         IF(DEBUG.GT.0) WRITE(LU,*) 'END INIT_TRANSPORT'
01199 !
01200 ! ---------- DEBUT IMPRESSIION INITIALISATION =================
01201 !
01202         CALL ENTETE_SISYPHE(1,AT0,0)
01203 !       PREPARES RESULTS
01204 !
01205 ! CONCENTRATION OUTPUT IN G/L
01206 !
01207         IF(UNIT) CALL OS('X=CX    ',X=CS,C=XMVS)
01208         CALL PREDES(0,AT0)
01209 !
01210 !       PRINTS OUT THE RESULTS
01211 !
01212         IF(DEBUG.GT.0) WRITE(LU,*) 'BIEF_DESIMP'
01213         CALL BIEF_DESIMP(SIS_FILES(SISRES)%FMT,VARSOR,
01214      &                   HIST,0,NPOIN,SIS_FILES(SISRES)%LU,'STD',
01215      &                   AT0,0,LISPR,LEOPR,SORLEO,SORIMP,MAXVAR,
01216      &                   TEXTE,PTINIG,PTINIL)
01217         IF(DEBUG.GT.0) WRITE(LU,*) 'END BIEF_DESIMP'
01218 !
01219         IF(UNIT) CALL OS('X=CX    ',X=CS,C=1.D0/XMVS)
01220 !
01221 !===============FIN IMPRESSION CONDITIONS INITIALES =================
01222 !
01223 !      COUPLING
01224 !
01225         IF(DREDGESIM) THEN
01226           CALL DREDGESIM_INTERFACE(1)
01227           IF(LNG.EQ.1) WRITE(LU,*) 'SISYPHE COUPLE AVEC DREDGESIM'
01228           IF(LNG.EQ.2) WRITE(LU,*) 'SISYPHE COUPLED WITH DREDGESIM'
01229         ENDIF
01230 !
01231         IF(CODE(1:7).NE.'SISYPHE') THEN
01232           IF(LNG.EQ.1) WRITE(LU,*) 'SISYPHE COUPLE AVEC : ',CODE
01233           IF(LNG.EQ.2) WRITE(LU,*) 'SISYPHE COUPLED WITH: ',CODE
01234         ENDIF
01235 !
01236 !      COUPLING WITH TELEMAC-2D OR 3D
01237 !
01238         IF(CODE(1:7).EQ.'TELEMAC') NCALCU = 1
01239 !
01240 !=======================================================================
01241 !
01242 !     INITIAL CONDITION FOR CONSTANT FLOW DISCHARGE
01243 !
01244         IF(LCONDIS) THEN
01245           SISYPHE_CFD = LCONDIS
01246           NSIS_CFD    = NCONDIS
01247           CONSTFLOW   = .FALSE.
01248         ELSE
01249           SISYPHE_CFD = .FALSE.
01250           NSIS_CFD    = 1
01251           CONSTFLOW   = .FALSE.
01252         ENDIF
01253 !
01254 !=======================================================================
01255 !
01256 !     END OF INITIALISATIONS
01257         IF(DEBUG.GT.0) WRITE(LU,*) 'END_INITIALIZATION'
01258       ENDIF ! IF (PART==0 OR PART = -1)
01259 !
01260 !=======================================================================
01261 !
01262       IF(PART==1.OR.PART==-1) THEN
01263 !
01264         IF(DEBUG.GT.0) WRITE(LU,*) 'TIME_LOOP'
01265 !
01266 !=======================================================================
01267 !
01268 ! : 3                    /* LOOP ON TIME */
01269 !
01270 !=======================================================================
01271 !
01272 !----   STOPS THE COMPUTATION WHEN THE REQUIRED NUMBER OF ITERATIONS IS 0
01273 !
01274         IF(NIT == 0) THEN
01275           IF (LNG == 1) WRITE(LU,200)
01276           IF (LNG == 2) WRITE(LU,201)
01277 200       FORMAT(' ARRET DANS SISYPHE, NOMBRE D''ITERATIONS',/,
01278      &         ' DEMANDE NUL')
01279 201       FORMAT(' STOP IN SISYPHE, NUMBER OF ITERATIONS EQ.0')
01280           CALL PLANTE(1)
01281           STOP
01282         ENDIF
01283 !
01284 !---------------------------------------------------------------------
01285 !       STARTS THE COMPUTATIONS
01286 !---------------------------------------------------------------------
01287 !       LOOP ON THE NUMBER OF EVENTS
01288 !       (IN STEADY STATE: LOOP ON THE TIMESTEPS)
01289 !---------------------------------------------------------------------
01290 !
01291         IF(CODE(1:7) == 'TELEMAC') THEN
01292 !         VALNIT WILL BE USED FOR CALLING VALIDA
01293           VALNIT = (TELNIT/PERICOU)*PERICOU-PERICOU+1
01294 !         MODIFICATION JMH + CV: TO AVOID 2 SUCCESSIVE CALLS TO VALIDA
01295 !         WHEN BEDLOAD AND SUSPENSION
01296           IF(GRAFCOUNT.GT.TELNIT) VALNIT=NIT+1
01297 !         CHARR, SUSP AND TIME STEP MONITORED BY CALLING PROGRAM
01298           CHARR = CHARR_TEL
01299           SUSP= SUSP_TEL
01300           AT0=T_TEL
01301           DT=MOFAC*DT
01302 !
01303         ENDIF
01304 !
01305         DO MN = 1, NCALCU
01306 !
01307           IF(.NOT.PERMA.AND.CODE(1:7).NE.'TELEMAC') THEN
01308 !           DETERMINES THE FIRST RECORD TO BE READ :
01309 !           NUMDEB IS THE FIRST RECORD TO BE READ FROM THE HYDRO
01310 !           FILE
01311             NUMDEB=NUMEN0
01312             IF(NUMDEB+NIDT > NUMEN) THEN
01313               IF (LNG == 1) WRITE(LU,202)
01314               IF (LNG == 2) WRITE(LU,203)
01315 202           FORMAT(1X,'FICHIER HYDRODYNAMIQUE PAS ASSEZ LONG')
01316 203           FORMAT(1X,'THE HYDRODYNAMIC FILE IS NOT LONG ENOUGH')
01317               CALL PLANTE(1)
01318             ENDIF
01319           ENDIF
01320 !
01321 !         LOOP ON THE RECORDS (IF PERMA NIDT=1)
01322 !         ------------------------------
01323 
01324           DO MT = 1, NIDT
01325 !
01326 !  ----   DETERMINES THE TIMESTEP NUMBER :
01327 !
01328             LT = (MN-1)*NIDT +  MT
01329 !
01330             IF(CODE(1:7) == 'TELEMAC') THEN
01331               DT=DT_TEL*MOFAC
01332               LT    = LOOPCOUNT
01333               LEOPR = GRAFCOUNT
01334               LISPR = LISTCOUNT
01335               NSOUS=1
01336             ENDIF
01337 !
01338 !  ----     PRINTOUTS TO LISTING :
01339 !
01340             ENTETS = .FALSE.
01341             IF(LISPR*(LT/LISPR).EQ.LT) THEN
01342               ENTET = .TRUE.
01343             ELSE
01344               ENTET = .FALSE.
01345             ENDIF
01346 !
01347 !  ----     PRINTOUTS TO C-VSP !!! UHM
01348 !
01349             CVSM_OUT  = .FALSE.
01350             IF(CVSMPPERIOD*(LT/CVSMPPERIOD).EQ.LT) CVSM_OUT = .TRUE.
01351 !
01352 !----       READS AND UPDATES H AND ZF
01353 !----       IF 1ST PASS OR UNSTEADY AND NO COUPLING
01354 !
01355             DTS = DT / NSOUS
01356 !
01357             ISOUS = 0
01358 !
01359             IF(.NOT.PERMA.OR.PASS) THEN
01360 !
01361               IF(DEBUG.GT.0) WRITE(LU,*) 'CONDIM_SISYPHE'
01362               CALL CONDIM_SISYPHE
01363      &        (U2D%R,V2D%R,QU%R,QV%R,HN%R,ZF%R,Z%R,ESOMT%R,THETAW%R,
01364      &        Q%R,HW%R,TW%R,MESH%X%R,MESH%Y%R,NPOIN,AT0,PMAREE)
01365               IF(DEBUG.GT.0) WRITE(LU,*) 'END_CONDIM_SISYPHE'
01366 !
01367 !             BEWARE : THE VALUE FOR ESOMT IS NOT READ FROM THE FILE SISHYD
01368 !             NOTE : NAME FOR SISHYD SET TO ' ' IF COUPLING
01369 !
01370               IF(SIS_FILES(SISHYD)%NAME(1:1).NE.' ')  THEN
01371 !
01372 !               WORK ON ZF,QU,QV,Z WILL BE IN FACT DONE ON:
01373 !               T4,DEL_QU,DEL_QV AND DEL_Z
01374 !               BY PLAYING WITH POINTERS
01375 !               THEN THE INCREMENT IN ONE TIME STEP WILL BE COMPUTED
01376 !               AND THE INCREMENT FOR ONE SUB TIME-STEP DEDUCED
01377 !               THEN AT EVERY SUB TIME-STEP THIS INCREMENt WILL BE ADDED
01378 !               AND WILL PROGRESSIVELy UPDATE THE VARIABLES FROM TIME N
01379 !               TO TIME N+1
01380                 SAVEZF=>ZF%R
01381                 SAVEQU=>QU%R
01382                 SAVEQV=>QV%R
01383                 SAVEZ =>Z%R
01384                 IF(HOULE.AND.UW%TYPR=='Q') SAVEUW=>UW%R
01385                 ZF%R  =>T4%R
01386                 QU%R  =>DEL_QU%R
01387                 QV%R  =>DEL_QV%R
01388                 Z%R   =>DEL_Z%R
01389                 IF(HOULE.AND.UW%TYPR=='Q') UW%R  =>DEL_UW%R
01390 !
01391                 NUMDEB=NUMDEB+1
01392 !
01393                 IF(ENTET) WRITE(LU,*) 'DEFINITION INITIALE DES VITESSES'
01394 !
01395                 CALL BIEF_SUITE(VARSOR,VARCL,NUMDEB,
01396      &             SIS_FILES(SISHYD)%LU,SIS_FILES(SISHYD)%FMT,
01397      &             HIST,0,NPOIN,BID,TEXTPR,VARCLA,0,
01398      &             TROUVE,ALIRE,ENTET,PERMA,MAXVAR)
01399 !
01400                 IF(DEBUG.GT.0) WRITE(LU,*) 'RESCUE_SISYPHE_NOTPERMA'
01401                 CALL RESCUE_SISYPHE_NOTPERMA
01402      &               (QU%R,QV%R,Q%R,U2D%R,V2D%R,HN%R,Z%R,T4%R,
01403      &                HW%R,TW%R,THETAW%R,NPOIN,TROUVE,ALIRE,ICF,ENTET,
01404      &                MAXVAR)
01405                 IF(DEBUG.GT.0) WRITE(LU,*) 'END_RESCUE_SISYPHE_NOTPERMA'
01406 !
01407 !               BACK TO ORIGINAL ADDRESSES
01408                 ZF%R=>SAVEZF
01409                 QU%R=>SAVEQU
01410                 QV%R=>SAVEQV
01411                 Z%R=>SAVEZ
01412                 IF(HOULE.AND.UW%TYPR=='Q') UW%R=>SAVEUW
01413 !
01414 !               INCREMENT OF QU, QV AND Z PER SUB-TIME-STEP
01415                 DO I = 1,NPOIN
01416                   DEL_QU%R(I) = (DEL_QU%R(I)-QU%R(I))/NSOUS
01417                   DEL_QV%R(I) = (DEL_QV%R(I)-QV%R(I))/NSOUS
01418                   DEL_Z%R(I)  = (DEL_Z%R(I) -Z%R(I)) /NSOUS
01419                   IF(HOULE.AND.UW%TYPR=='Q') THEN
01420                     DEL_UW%R(I) = (DEL_UW%R(I)-UW%R(I))/NSOUS
01421                   ENDIF
01422                 ENDDO
01423 !
01424 !               UPDATES UNSTEADY HYDRO
01425 !               (TO BE MOVED TO RESCUE_SISYPHE_NOTPERMA)
01426 !               -----------------------------------
01427 !               CLIPS NEGATIVE DEPTHS
01428 !               COMPUTES U2D AND V2D
01429 !
01430                 CALL OS('X=Y-Z   ', X=HN, Y=Z, Z=ZF)
01431 !
01432                 IF(OPTBAN.GT.0) THEN
01433                   DO I = 1,NPOIN
01434                     IF(HN%R(I).LT.HMIN) THEN
01435                       U2D%R(I)=0.D0
01436                       V2D%R(I)=0.D0
01437                       HN%R(I) = MAX(HN%R(I),HMIN)
01438                     ELSE
01439                       U2D%R(I)=QU%R(I)/HN%R(I)
01440                       V2D%R(I)=QV%R(I)/HN%R(I)
01441                     ENDIF
01442                   ENDDO
01443                 ELSE
01444                   CALL OS('X=Y/Z   ', X=U2D, Y=QU,   Z=HN)
01445                   CALL OS('X=Y/Z   ', X=V2D, Y=QV,   Z=HN)
01446                 ENDIF
01447 !
01448               ENDIF ! (SIS_FILES(SISHYD)%NAME(1:1) /=' ')
01449             ENDIF ! (NOT.PERMA.OR.PASS)
01450 !
01451         IF(PASS) THEN
01452 !         IN STEADY STATE LOGICAL FOR READING SET TO FALSE
01453           IF (PERMA) PASS = .FALSE.
01454         ELSE
01455 !         COMPUTES THE WATER DEPTH
01456           CALL OS('X=Y-Z   ', X=HN, Y=Z, Z=ZF)
01457         ENDIF
01458 !
01459 !       COUPLING
01460 !
01461         IF(CODE(1:7).EQ.'TELEMAC') THEN
01462 !
01463 !         OV INSTEAD OF OS IN ORDER TO AVOID PROBLEMS WITH QUASI-BUBBLE ELEMENTS
01464 !         OPERATES ONLY ON THE (1:NPOIN) RANGE OF THE TELEMAC FIELDS
01465 !         IT IS A HIDDEN DISCRETISATION CHANGE
01466 !
01467           CALL OV( 'X=Y     ',U2D%R, U_TEL%R, U_TEL%R, 0.D0, NPOIN)
01468           CALL OV( 'X=Y     ',V2D%R, V_TEL%R, V_TEL%R, 0.D0, NPOIN)
01469           CALL OV( 'X=Y     ', HN%R, H_TEL%R, H_TEL%R, 0.D0, NPOIN)
01470 !         ADDED BY JMH 01/07/2004 (ZF MAY BE MODIFIED BY CALLING PROGRAM)
01471           CALL OS('X=Y     ', X=ZF, Y=ZF_TEL)
01472 !         ELAY MUST BE UPDATED CONSEQUENTLY (ADDED BY JMH 02/05/2014)
01473           CALL OS('X=Y-Z   ',X=ELAY,Y=ZF,Z=ZR)
01474 !         CLIPS NEGATIVE DEPTHS
01475           IF(OPTBAN.GT.0) THEN
01476             DO I = 1,HN%DIM1
01477               IF(HN%R(I).LT.HMIN) THEN
01478                 U2D%R(I)=0.D0
01479                 V2D%R(I)=0.D0
01480                 HN%R(I)=HMIN
01481               ENDIF
01482             ENDDO
01483           ENDIF
01484 !         FREE SURFACE
01485           CALL OS('X=Y+Z   ', X=Z, Y=ZF, Z=HN)
01486 !
01487 !         CV: COPY OF TOMAWAC VARIABLES
01488 !
01489           IF(INCLUS(COUPLING,'TOMAWAC')) THEN
01490 !           INCIDENT WAVE DIRECTION
01491             CALL OS( 'X=Y     ',THETAW,THETAW_TEL)
01492 !           Wave period
01493             CALL OS( 'X=Y     ', TW, TW_TEL)
01494 !           SIGNIFICANT HEIGHT
01495             CALL OS( 'X=Y     ', HW , HW_TEL)
01496 !           SIGNIFICANT HEIGHT
01497             CALL OS( 'X=Y     ', UW , UW_TEL)
01498             HW%TYPR='Q'
01499             TW%TYPR='Q'
01500             THETAW%TYPR='Q'
01501             UW%TYPR='Q'
01502           ENDIF
01503 !
01504         ENDIF
01505 !
01506 !       END OF COUPLING
01507 ! =========================================================================
01508 ! TREATMENT OF TIDAL FLATS, DEFINITION OF THE MASKS
01509 ! =====================================================================!
01510 !
01511         IF(OPTBAN.EQ.2) THEN
01512 !
01513 ! ----    BUILDS MASKING BY ELEMENTS
01514 !
01515           CALL OS ('X=Y     ', X=MSKTMP, Y=MASKEL)
01516           CALL OS ('X=C     ', X=MASKEL, C=1.D0)
01517           IF(CODE(1:7) == 'TELEMAC') THEN
01518 !           MASKS ARE DERIVED FROM THE NON-CLIPPED VALUES OF H
01519 !           PROVIDED BY TELEMAC
01520             CALL MASKTF(MASKEL%R,H_TEL%R,HMIN,MESH%IKLE%I,
01521      &                  NELEM,NPOIN)
01522           ELSE
01523             CALL MASKTF(MASKEL%R,HN%R,HMIN,MESH%IKLE%I,
01524      &                  NELEM,NPOIN)
01525           ENDIF
01526 !
01527 !        JMH 17/12/2009
01528 !
01529 !        ELSEIF(OPTBAN.EQ.1) THEN
01530 !
01531 !          CANCELS Q QU AND QV IF HN.LE.0.D0
01532 !          CALL MASKAB(HN%R,Q%R,QU%R,QV%R,NPOIN)
01533 !
01534         ENDIF
01535 !
01536 ! ----   BUILDS THE MASK OF THE POINTS FROM THE MASK OF THE ELEMENTS
01537 ! ----   AND CHANGES IFAMAS (IFABOR WITH MASKING)
01538 !
01539         IF(MSK) CALL MASKTO(MASKEL%R,MASKPT,IFAMAS%I,
01540      &                      MESH%IKLE%I,
01541      &                      MESH%IFABOR%I,MESH%ELTSEG%I,MESH%NSEG,
01542      &                      NELEM,NPOIN,IELMT,MESH)
01543 !
01544 ! ------------------------------------------------------------------
01545 !  START OF SUB-ITERATIONS IN UNSTEADY STATE
01546 !
01547 ! ------------------------------------------------------------------
01548 !
01549 702     CONTINUE
01550 !
01551         ISOUS = ISOUS + 1
01552         IF(CODE(1:7).NE.'TELEMAC') AT0=AT0+DTS
01553         IF(ENTET.AND.ISOUS.EQ.1) CALL ENTETE_SISYPHE(2,AT0,LT)
01554         IF(ENTET.AND.ISOUS.EQ.NSOUS) ENTETS=.TRUE.
01555 !
01556 !---------------------------------------------------------------------
01557 !        FRICTION COEFFICIENT VARIABLE IN TIME
01558 !---------------------------------------------------------------------
01559 !
01560         CALL CORSTR_SISYPHE
01561 !
01562 ! ----  TREATS THE BOUNDARY CONDITIONS
01563 !
01564         IF(DEBUG.GT.0) WRITE(LU,*) 'CONLIT'
01565         CALL CONLIT(MESH%NBOR%I,AT0)
01566         IF(DEBUG.GT.0) WRITE(LU,*) 'END CONLIT'
01567 !
01568 ! =======================================================================
01569 !
01570 !       IF 'VARIABLE TIME-STEP = YES' NSOUS WILL BE COMPUTED FURTHER DOWN
01571 !       THE CONPUTATION OF THE TIMESTEP SIS HAS BEEN MOVED BEFORE READING
01572 !       THE HYDRO CONDITIONS
01573 !
01574 !  ---  MEAN DIAMETER FOR THE ACTIVE-LAYER AND UNDER-LAYER
01575 !
01576         IF(.NOT.MIXTE.AND.NSICLA.GT.1) CALL MEAN_GRAIN_SIZE
01577 !
01578 !  ---  MEAN VELOCITY UNORM
01579 !
01580         CALL OS('X=N(Y,Z)',X=UNORM,Y=U2D,Z=V2D)
01581 !
01582 !  ---  WAVE ORBITAL VELOCITY --> UW
01583 !
01584         IF(HOULE) THEN
01585           IF(UW%TYPR.NE.'Q') THEN
01586             CALL CALCUW(UW%R,HN%R,HW%R,TW%R,GRAV,NPOIN)
01587           ENDIF
01588         ENDIF
01589 !
01590         CALL TOB_SISYPHE
01591      &    (TOB,TOBW, MU, KS, KSP,KSR,CF, FW,
01592      &     CHESTR, UETCAR, CF_TEL,KS_TEL, CODE ,
01593      &     KFROT, ICR, KSPRATIO,HOULE,
01594      &     GRAV,XMVE,  XMVS, VCE, KARMAN,ZERO,
01595      &     HMIN,HN, ACLADM, UNORM,UW, TW, NPOIN,KSPRED,IKS)
01596 !
01597 !  END OF INITIALISATION
01598 !
01599       ! ******************** !
01600       ! BEDLOAD COMPUTATION  !
01601       ! ******************** !
01602         IF(CHARR) THEN
01603           IF(DEBUG.GT.0) WRITE(LU,*) 'BEDLOAD_MAIN'
01604           CALL BEDLOAD_MAIN
01605      &        (ACLADM,KSP,KSR,VOLU2D,UNSV2D,
01606      &         CF,EBOR,FW,HN,LIQBOR,MASK,MASKEL,
01607      &         MASKPT,Q,QBOR,U2D,V2D,S,UNLADM,UW,THETAW,
01608      &         MU,TOB,TOBW,TW,ZF,DEBUG,HIDFAC,ICF,
01609      &         IELMT,ISOUS,KDDL,KDIR,KENT,KINC,KLOG,KNEU,KSORT,
01610      &         LOADMETH,LT,NPOIN,NPTFR,NSICLA,
01611      &         OPTBAN,LS0,BETA,FD90,FDM,GRAV,HIDI,HMIN,
01612      &         VCE,CSF_SABLE,XMVE,XMVS,XWC,PI,KARMAN,ZERO,
01613      &         KARIM_HOLLY_YANG,MSK,SUSP,VF,ENTET,
01614      &         CONST_ALAYER,LCONDIS,LGRAFED,MESH,
01615      &         ELAY,LIEBOR,LIMTEC,MASKTR,
01616      &         IT1,T1,T2,T3,T4,T5,T6,T7,T8,T9,
01617      &         T10,T11,T12,T13,UNORM,AC,AT0,DTS,ELAY0,FRACSED_GF,
01618      &         AVAIL,BREACH,CALFA,COEFPN,DZF_GF,HIDING,
01619      &         QSCL_C,QSCL_S,QS_C,QSCLXC,QSXC,QSCLYC,
01620      &         QSYC,SALFA,ZF_C,ZFCL_C,NSOUS,ENTETS,
01621      &         SECCURRENT,SLOPEFF,PHISED,DEVIA,BETA2,BIJK,
01622      &         SEDCO,HOULE,U3D,V3D,CODE,FLBCLA,MAXADV)
01623 !
01624           IF(DEBUG.GT.0) WRITE(LU,*) 'END_BEDLOAD_MAIN'
01625 !
01626 !         UPDATES THE BOTTOM
01627 !
01628           IF(.NOT.STAT_MODE) THEN
01629             CALL OS('X=X+Y   ',X=ZF,Y=ZF_C)
01630           ENDIF
01631 !
01632 !         UPDATES THE LAYERS  --> ELAY
01633 !
01634           IF(.NOT.MIXTE.AND.NSICLA.GT.1) THEN
01635             IF(DEBUG.GT.0) WRITE(LU,*) 'LAYER AFTER BEDLOAD'
01636             IF(VSMTYPE.EQ.0) THEN
01637               CALL LAYER(ZFCL_C,NLAYER,ZR,ZF,ESTRAT,ELAY,VOLU2D,
01638      &                   ACLADM,NSICLA,NPOIN,ELAY0,VOLTOT,ES,
01639      &                   AVAIL,CONST_ALAYER,DTS,T2%R,IT1%I)
01640             ELSE
01641               CALL CVSP_MAIN(ZFCL_C,NLAYER,ZR,ZF,ESTRAT,ELAY,VOLU2D,
01642      &                       ACLADM,NSICLA,NPOIN,ELAY0,VOLTOT,ES,
01643      &                       AVAIL,CONST_ALAYER,DTS,T2%R,IT1%I)
01644             ENDIF
01645             IF(DEBUG.GT.0) WRITE(LU,*) 'END_LAYER'
01646           ELSE
01647             CALL OS('X=Y-Z   ',X=ELAY,Y=ZF,Z=ZR)
01648           ENDIF
01649 ! END OF BEDLOAD
01650         ENDIF
01651       ! ********************** !
01652       ! SUSPENSION COMPUTATION !
01653       ! ********************** !
01654         IF(SUSP) THEN
01655 !
01656           IF(DEBUG.GT.0) WRITE(LU,*) 'SUSPENSION_MAIN'
01657           CALL SUSPENSION_MAIN
01658      &(SLVTRA,HN,HN_TEL,MU,TOB,FDM,FD90,KSP,KSR,KS,
01659      & VOLU2D,V2DPAR,UNSV2D,AFBOR,BFBOR,ZF,LICBOR,
01660      & IFAMAS,MASKEL,MASKPT,U2D,V2D,NSICLA,
01661      & NPOIN,NPTFR,IELMT,OPTDIF,RESOL,LT,NIT,OPTBAN,OPTADV,
01662      & OPDTRA,KENT,KSORT,KLOG,KINC,KNEU,KDIR,KDDL,
01663      & DEBUG,DTS,CSF_SABLE,ZERO,GRAV,XKX,XKY,
01664      & KARMAN,XMVE,XMVS,VCE,HMIN,XWC,VITCD,VITCE,PARTHENIADES,
01665      & BILMA,MSK,CHARR,IMP_INFLOW_C,MESH,ZF_S,CS,
01666      & CST,CTILD,CBOR,DISP,IT1,IT2,IT3,IT4,TB,T1,T2,T3,T4,T5,T6,
01667      & T7,T8,T9,T10,T11,T12,T14,W1,TE1,CLT,TE2,TE3,S,AM1_S,AM2_S,MBOR,
01668      & ELAY,LIMDIF,MASKTR,TETA_SUSP,AC,MASED0,MASINI,MASTEN,MASTOU,
01669      & ES,ES_SABLE, ES_VASE,AVAIL,ENTETS,PASS_SUSP,
01670      & ZFCL_S,HPROP,FLUDPT,FLUDP,FLUER,DISP_C,KX,KY,KZ,UCONV,
01671      & VCONV,QSXS,QSYS,QSCLXS,QSCLYS,QSCL_S,QS_S,QS_C,
01672      & CSTAEQ,CSRATIO,ICQ,MASTCP,MASFIN,MASDEPT,MASDEP,MASSOU,CORR_CONV,
01673      & ZREF,SEDCO,VISC_TEL,CODE,DIFT,DM1,UCONV_TEL,VCONV_TEL,
01674      & ZCONV,SOLSYS,FLBOR_TEL,FLBOR_SIS,FLBORTRA,NUMLIQ%I,NFRLIQ,
01675      & MIXTE,NCOUCH_TASS,CONC,TOCE_VASE,TOCE_SABLE,
01676      & FLUER_VASE,TOCE_MIXTE,MS_SABLE,MS_VASE,TASS,DIRFLU,MAXADV)
01677           IF(DEBUG.GT.0) WRITE(LU,*) 'END_SUSPENSION_MAIN'
01678 !
01679 !      UPDATES THE BOTTOM
01680 !
01681 ! mak:
01682           IF(.NOT.STAT_MODE) THEN
01683             CALL OS('X=X+Y   ',X=ZF,Y=ZF_S)
01684           ENDIF
01685 ! end mak
01686 !        CALL OS('X=X+Y   ',X=ZF,Y=ZF_S)
01687 !
01688 !      UPDATES THE LAYERS
01689 !      REDEFINES THE LAYER OF ERODABLE SEDIMENT
01690 !      EXTENDED GRANULOMETRY (TO BE REPLACED WITH NOMBLAY>1
01691 !
01692           IF(.NOT.MIXTE.AND.NSICLA.GT.1) THEN
01693             IF(DEBUG.GT.0) WRITE(LU,*) 'LAYER AFTER SUSPENSION'
01694             IF(VSMTYPE.EQ.0) THEN
01695               CALL LAYER(ZFCL_S,NLAYER,ZR,ZF,ESTRAT,ELAY,VOLU2D,
01696      &                   ACLADM,NSICLA,NPOIN,ELAY0,VOLTOT,ES,
01697      &                   AVAIL,CONST_ALAYER,DTS,T2%R,IT1%I)
01698             ELSE
01699               CALL CVSP_MAIN(ZFCL_S,NLAYER,ZR,ZF,ESTRAT,ELAY,VOLU2D,
01700      &                       ACLADM,NSICLA,NPOIN,ELAY0,VOLTOT,ES,
01701      &                       AVAIL,CONST_ALAYER,DTS,T2%R,IT1%I)
01702             ENDIF
01703             IF(DEBUG.GT.0) WRITE(LU,*) 'END_LAYER'
01704           ELSE
01705             CALL OS('X=Y-Z   ',X=ELAY,Y=ZF,Z=ZR)
01706           ENDIF
01707 ! END OF SUSPENSION
01708         ENDIF
01709 !
01710 ! RECONSTITUTES THE BEDLOAD AND/OR SUSPENSION DATA
01711 ! -----------------------------------------------------
01712 !
01713         IF(DEBUG.GT.0) WRITE(LU,*) 'QS_RESULT'
01714 !
01715         DO I = 1, NSICLA
01716           CALL OS('X=Y+Z   ',X=T1,Y=QSCLXC%ADR(I)%P,Z=QSCLXS%ADR(I)%P)
01717           CALL OS('X=Y+Z   ',X=T2,Y=QSCLYC%ADR(I)%P,Z=QSCLYS%ADR(I)%P)
01718           CALL OS('X=N(Y,Z)', X=QSCL%ADR(I)%P,Y=T1,Z=T2)
01719           IF(I.EQ.1) THEN
01720             CALL OS('X=Y     ', X=QSX, Y=T1)
01721             CALL OS('X=Y     ', X=QSY, Y=T2)
01722           ELSE
01723             CALL OS('X=X+Y   ', X=QSX, Y=T1)
01724             CALL OS('X=X+Y   ', X=QSY, Y=T2)
01725           ENDIF
01726         ENDDO
01727         CALL OS('X=N(Y,Z)', X=QS, Y=QSX, Z=QSY)
01728         IF(DEBUG.GT.0) WRITE(LU,*) 'END_QS_RESULT'
01729 !
01730 !=======================================================================
01731 !
01732 !     MAXIMUM BOTTOM SLOPE : EVOL IN T1
01733 !
01734       IF(SLIDE) THEN
01735 !
01736         IF(ENTET) CALL ENTETE_SISYPHE(14,AT0,LT)
01737 !
01738         IF(DEBUG.GT.0) WRITE(LU,*) 'APPEL DE MAXSLOPE'
01739         CALL MAXSLOPE(PHISED,ZF%R,ZR%R,MESH%XEL%R,MESH%YEL%R,
01740      &                MESH%NELEM,
01741      &                MESH%NELMAX,NPOIN,MESH%IKLE%I,T1,UNSV2D,MESH,
01742      &                ZFCL_MS,AVAIL,NOMBLAY,NSICLA)
01743         IF(DEBUG.GT.0) WRITE(LU,*) 'RETOUR DE MAXSLOPE'
01744         CALL OS('X=X+Y   ',X=ZF,Y=T1)
01745 !
01746         IF(.NOT.MIXTE.AND.NSICLA.GT.1) THEN
01747           IF(DEBUG.GT.0) WRITE(LU,*) 'LAYER AFTER SLIDE'
01748           IF(VSMTYPE.EQ.0) THEN
01749             CALL LAYER(ZFCL_MS,NLAYER,ZR,ZF,ESTRAT,ELAY,VOLU2D,
01750      &                 ACLADM,NSICLA,NPOIN,ELAY0,VOLTOT,ES,
01751      &                 AVAIL,CONST_ALAYER,DTS,T2%R,IT1%I)
01752           ELSE
01753             CALL CVSP_MAIN(ZFCL_MS,NLAYER,ZR,ZF,ESTRAT,ELAY,VOLU2D,
01754      &                     ACLADM,NSICLA,NPOIN,ELAY0,VOLTOT,ES,
01755      &                     AVAIL,CONST_ALAYER,DTS,T2%R,IT1%I)
01756           ENDIF
01757           IF(DEBUG.GT.0) WRITE(LU,*) 'END_LAYER'
01758         ELSE
01759           CALL OS('X=Y-Z   ',X=ELAY,Y=ZF,Z=ZR)
01760         ENDIF
01761 !
01762       ENDIF
01763 !
01764 !========================================================================
01765 !
01766 !     SETTLING: EVOLUTION COMPUTED IN T3
01767 !
01768       IF(TASS) THEN
01769 !
01770         IF(ENTET) THEN
01771           IF(.NOT.CHARR.AND..NOT.SUSP.AND..NOT.SLIDE) THEN
01772             CALL ENTETE_SISYPHE(2,AT0,LT)
01773           ENDIF
01774           CALL ENTETE_SISYPHE(15,AT0,LT)
01775         ENDIF
01776 !
01777         IF(DEBUG.GT.0) WRITE(LU,*) 'APPEL DE TASSEMENT'
01778         IF(ITASS.EQ.1) THEN
01779           CALL TASSEMENT(ZF,NPOIN,DTS,ELAY,T3,T2,LT,AVAIL,NSICLA,
01780      &                 ES,XMVS,XKV,TRANS_MASS,CONC_VASE,NCOUCH_TASS,
01781      &                 MS_SABLE%R,MS_VASE%R)
01782         ELSEIF(ITASS.EQ.2.OR.ITASS.EQ.3) THEN
01783 !!!! ONLY FOR ONE CLASS
01784           CALL TASSEMENT_2(ZF,NPOIN,DTS,ELAY,T3,T2,LT,XMVS,XMVE,GRAV,
01785      &        NOMBLAY,ES,CONC_VASE,MS_VASE%R,XWC(1),
01786      &        COEF_N,CONC_GEL,CONC_MAX)
01787 !        ELSEIF(ITASS.EQ.3) THEN
01788 !!!! ONLY FOR ONE CLASS
01789 !          CALL TASSEMENT_3(ZF,NPOIN,DTS,ELAY,
01790 !     &               T3,T2,LT,XMVS,XMVE,GRAV,NOMBLAY,
01791 !     &               ES,CONC_VASE,CONC,IVIDE,MS_VASE%R,XWC(1),
01792 !     &                 TRA01,TRA02,TRA03,CONC_GEL,COEF_N,CONC_MAX)
01793         ENDIF
01794         IF(DEBUG.GT.0) WRITE(LU,*) 'RETOUR DE TASSEMENT'
01795 !
01796 !       UPDATES ZF (ELAY HAS BEEN UPDATED IN TASSEMENT)
01797 !
01798 ! mak:
01799         IF(.NOT.STAT_MODE) THEN
01800            CALL OS('X=X+Y   ',X=ZF,Y=T3)
01801         ENDIF
01802 ! end mak
01803 !        CALL OS('X=X+Y   ',X=ZF,Y=T3)
01804 !
01805       ENDIF
01806 !
01807 !=======================================================================
01808 ! : 5        COMPUTES THE EVOLUTIONS FOR THIS CYCLE OF TIMESTEP
01809 !            AND UPDATES AFTER THIS COMPUTATION
01810 !=======================================================================
01811 !
01812 ! ----  COMPUTES  THE EVOLUTIONS FOR THIS (SUB) TIMESTEP
01813 !
01814       IF(CHARR) THEN
01815         CALL OS('X=Y     ',X=E,Y=ZF_C)
01816       ELSE
01817         CALL OS('X=0     ',X=E)
01818       ENDIF
01819       IF(SUSP)  CALL OS('X=X+Y   ',X=E,Y=ZF_S)
01820       IF(SLIDE) CALL OS('X=X+Y   ',X=E,Y=T1)
01821       IF(TASS)  CALL OS('X=X+Y   ',X=E,Y=T3)
01822 !
01823       CALL OS('X=X+Y   ', X=ESOMT, Y=E)
01824 !
01825 !  UPDATES
01826 !
01827       IF(PART.EQ.-1) THEN
01828 !
01829 ! mak:
01830         IF(.NOT.STAT_MODE) THEN
01831           CALL OS('X=X-Y   ',X=HN,Y=E)
01832         ENDIF
01833 ! end mak
01834 !      CALL OS('X=X-Y   ',X=HN,Y=E)
01835         IF(OPTBAN.GT.0) THEN
01836           DO I = 1,HN%DIM1
01837             IF(HN%R(I).LT.HMIN) THEN
01838               U2D%R(I)=0.D0
01839               V2D%R(I)=0.D0
01840               HN%R(I) =HMIN
01841             ELSE
01842               U2D%R(I)= QU%R(I)/HN%R(I)
01843               V2D%R(I)= QV%R(I)/HN%R(I)
01844             ENDIF
01845           ENDDO
01846         ELSE
01847           CALL OS('X=Y/Z   ', X=U2D, Y=QU,   Z=HN)
01848           CALL OS('X=Y/Z   ', X=V2D, Y=QV,   Z=HN)
01849         ENDIF
01850 !
01851 !=======================================================================
01852 ! : 6     STOPS IF EVOLUTIONS GREATER THAN EMAX = RC*(INITIAL DEPTH)
01853 !=======================================================================
01854 !
01855 !       DETERMINES THE MAXIMUM EVOLUTION THRESHOLD
01856         DO I = 1, NPOIN
01857           EMAX%R(I) = RC*MAX(HN%R(I),HMIN)
01858         ENDDO
01859 !
01860 ! ----  STOPS WHEN THE EVOLUTIONS ARE GREATER THAN A CERTAIN THRESHOLD
01861 !      THIS TEST IS ONLY CALLED IN 'SISYPHE ONLY' MODE
01862 !
01863         IF(DEBUG.GT.0) WRITE(LU,*) 'ARRET'
01864           CALL SIS_ARRET(ESOMT,EMAX,HN,VARSOR,NPOIN,MN,
01865      &                   SIS_FILES(SISRES)%LU,SIS_FILES(SISRES)%FMT,
01866      &                   MAXVAR,AT0,RC,HIST,BINRESSIS,TEXTE,
01867      &                   SORLEO,SORIMP,T1,T2)
01868           IF(DEBUG.GT.0) WRITE(LU,*) 'END_ARRET'
01869         ENDIF
01870 !
01871 ! ----     CONSTANT FLOW DISCHARGE
01872 !
01873         IF(LCONDIS) THEN
01874           CALL CONDIS_SISYPHE(CONSTFLOW)
01875         ELSE
01876           CONSTFLOW =.FALSE.
01877         ENDIF
01878 !
01879 !=======================================================================
01880 ! : 8     MASS BALANCE
01881 !=======================================================================
01882 !       COMPUTES THE COMPONENTS OF SAND TRANSPORT FOR THE MASS BALANCE,
01883 !       GRAPHIC OUTPUTS AND VALIDATION STAGE
01884 !
01885         IF(BILMA) THEN
01886 ! CVL Bilan in volume
01887 ! not valid for cohesif or mixte
01888           IF(.NOT.MIXTE.AND.(.NOT.SEDCO(1))) THEN
01889           IF(DEBUG.GT.0) WRITE(LU,*) 'BILAN_SISYPHE'
01890           CALL BILAN_SISYPHE(E,ESOMT,MESH,MSK,MASKEL,T1,T2,S,
01891      &                       IELMU_SIS,VCUMU,DTS,NPTFR,ENTETS,
01892      &                       ZFCL_C,ZFCL_S,ZFCL_MS,
01893      &                       QSCLXC,QSCLYC,NSICLA,VOLTOT,DZF_GF,MASS_GF,
01894      &                       LGRAFED,NUMLIQ%I,NFRLIQ,FLBCLA,VF,LT,NIT,
01895      &                       NPOIN,VOLU2D,CSF_SABLE,MASDEP,MASDEPT,
01896      &                       CHARR,SUSP,SLIDE)
01897           ENDIF
01898           IF(DEBUG.GT.0) WRITE(LU,*) 'END_BILAN_SISYPHE'
01899         ENDIF
01900 !
01901 !       CONTROL SECTIONS
01902 !
01903         IF(NCP.GT.0) THEN
01904           IF(DEBUG.GT.0) WRITE(LU,*) 'FLUSEC_SISYPHE'
01905           CALL FLUSEC_SISYPHE(U2D,V2D,HN,
01906      &                        QSXC,QSYC,CHARR,QSXS,QSYS,SUSP,
01907      &                        MESH%IKLE%I,
01908      &                        MESH%NELMAX,MESH%NELEM,
01909      &                        MESH%X%R,MESH%Y%R,
01910      &                        DTS,NCP,CTRLSC,ENTETS,AT0)
01911           IF(DEBUG.GT.0) WRITE(LU,*) 'END_FLUSEC_SISYPHE'
01912         ENDIF
01913 !
01914 !-----------------------------------------------------------------------
01915 !
01916         IF(.NOT.PERMA.AND.SIS_FILES(SISHYD)%NAME(1:1).NE.' ') THEN
01917 !
01918 !         UPDATES THE HYDRO
01919 !
01920 !         IF READING ON HYDRODYNAMIC FILE, INCREMENTS QU, QV AND Z
01921           IF(SIS_FILES(SISHYD)%NAME(1:1).NE.' ')  THEN
01922             CALL OS('X=X+Y   ', X=QU, Y=DEL_QU)
01923             CALL OS('X=X+Y   ', X=QV, Y=DEL_QV)
01924             CALL OS('X=X+Y   ', X=Z , Y=DEL_Z)
01925 ! JWI 31/05/2012 - added line to use wave orbital velocities directly if found in hydro file
01926             IF(HOULE.AND.UW%TYPR=='Q') CALL OS('X=X+Y   ',X=UW,Y=DEL_UW)
01927 ! JWI END
01928           ENDIF
01929           CALL OS('X=Y-Z   ', X=HN, Y=Z, Z=ZF)
01930 !         CLIPS NEGATIVE DEPTHS
01931           IF(OPTBAN.GT.0) THEN
01932             DO I = 1, NPOIN
01933               IF(HN%R(I).LT.HMIN) THEN
01934                 U2D%R(I)=0.D0
01935                 V2D%R(I)=0.D0
01936                 HN%R(I) = MAX(HN%R(I),HMIN)
01937               ELSE
01938                 U2D%R(I)= QU%R(I)/HN%R(I)
01939                 V2D%R(I)= QV%R(I)/HN%R(I)
01940               ENDIF
01941             ENDDO
01942           ELSE
01943             CALL OS('X=Y/Z   ', X=U2D, Y=QU,   Z=HN)
01944             CALL OS('X=Y/Z   ', X=V2D, Y=QV,   Z=HN)
01945           ENDIF
01946 !
01947         ENDIF
01948 !
01949 !       END OF THE LOOP ON SUB-TIMESTEPS NSOUS
01950 !
01951         IF(ISOUS.LT.NSOUS) GO TO 702
01952 !
01953 !=======================================================================
01954 ! : 9        PRINTS OUT EXTREME VALUES
01955 !=======================================================================
01956 !
01957         IF(ENTET.AND.CHARR) THEN
01958           WRITE(LU,*)
01959           CALL MAXI(XMAX,IMAX,E%R,NPOIN)
01960           IF(NCSIZE.GT.1) THEN
01961             XMA=P_DMAX(XMAX)
01962             IF(XMAX.EQ.XMA) THEN
01963               IMA=MESH%KNOLG%I(IMAX)
01964             ELSE
01965               IMA=0
01966             ENDIF
01967             IMA=P_IMAX(IMA)
01968           ELSE
01969             IMA=IMAX
01970             XMA=XMAX
01971           ENDIF
01972           IF(LNG.EQ.1) WRITE(LU,371) XMA,IMA
01973           IF(LNG.EQ.2) WRITE(LU,372) XMA,IMA
01974 371       FORMAT(' EVOLUTION MAXIMUM        : ',G16.7,' NOEUD : ',I6)
01975 372       FORMAT(' MAXIMAL EVOLUTION        : ',G16.7,' NODE  : ',I6)
01976           CALL MINI(XMIN,IMIN,E%R,NPOIN)
01977           IF(NCSIZE.GT.1) THEN
01978             XMI=P_DMIN(XMIN)
01979             IF(XMIN.EQ.XMI) THEN
01980               IMI=MESH%KNOLG%I(IMIN)
01981             ELSE
01982               IMI=0
01983             ENDIF
01984             IMI=P_IMAX(IMI)
01985           ELSE
01986             IMI=IMIN
01987             XMI=XMIN
01988           ENDIF
01989           IF(LNG.EQ.1) WRITE(LU,373) XMI,IMI
01990           IF(LNG.EQ.2) WRITE(LU,374) XMI,IMI
01991 373       FORMAT(' EVOLUTION MINIMUM        : ',G16.7,' NOEUD : ',I6)
01992 374       FORMAT(' MINIMAL EVOLUTION        : ',G16.7,' NODE  : ',I6)
01993 !
01994           IF(CONST_ALAYER) THEN
01995             IF(NSICLA.GT.1.AND.XMI.LT.-0.5D0*ELAY0) THEN
01996               IF(LNG.EQ.1) WRITE(LU,885)
01997               IF(LNG.EQ.2) WRITE(LU,886)
01998 885           FORMAT(' EROSION SUPERIEURE A EPAISSEUR DE COUCHE !')
01999 886           FORMAT(' EROSION GREATER THAN ONE LAYER THICKNESS !')
02000             ENDIF
02001             IF(NSICLA.GT.1.AND.XMA.GT.ELAY0) THEN
02002               IF(LNG.EQ.1) WRITE(LU,887)
02003               IF(LNG.EQ.2) WRITE(LU,888)
02004 887           FORMAT(' DEPOT SUPERIEUR A EPAISSEUR DE COUCHE !')
02005 888           FORMAT(' DEPOSITION MORE THAN ONE LAYER THICKNESS !')
02006             ENDIF
02007           ELSE
02008             DO J=1,NPOIN
02009               ELAY0 = 3.D0*ACLADM%R(J)
02010               IF(NSICLA.GT.1.AND.E%R(J).LT.-0.5D0*ELAY0) THEN
02011                 IF(LNG.EQ.1) WRITE(LU,885)
02012                 IF(LNG.EQ.2) WRITE(LU,886)
02013               ENDIF
02014               IF(NSICLA.GT.1.AND.E%R(J).GT.ELAY0) THEN
02015                 IF(LNG.EQ.1) WRITE(LU,887)
02016                 IF(LNG.EQ.2) WRITE(LU,888)
02017               ENDIF
02018             ENDDO
02019           ENDIF
02020         ENDIF
02021         IF(ENTET) THEN
02022           CALL MAXI(XMAX,IMAX,ESOMT%R,NPOIN)
02023           IF(NCSIZE.GT.1) THEN
02024             XMA=P_DMAX(XMAX)
02025             IF(XMAX.EQ.XMA) THEN
02026               IMA=MESH%KNOLG%I(IMAX)
02027             ELSE
02028               IMA=0
02029             ENDIF
02030             IMA=P_IMAX(IMA)
02031           ELSE
02032             IMA=IMAX
02033             XMA=XMAX
02034           ENDIF
02035           IF (LNG.EQ.1) WRITE(LU,881) XMA,IMA
02036           IF (LNG.EQ.2) WRITE(LU,882) XMA,IMA
02037 881       FORMAT(' EVOLUTION MAXIMUM TOTALE : ',G16.7,' NOEUD : ',I6)
02038 882       FORMAT(' TOTAL MAXIMAL EVOLUTION  : ',G16.7,' NODE  : ',I6)
02039           CALL MINI(XMIN,IMIN,ESOMT%R,NPOIN)
02040           IF(NCSIZE.GT.1) THEN
02041             XMI=P_DMIN(XMIN)
02042             IF(XMIN.EQ.XMI) THEN
02043               IMI=MESH%KNOLG%I(IMIN)
02044             ELSE
02045               IMI=0
02046             ENDIF
02047             IMI=P_IMAX(IMI)
02048           ELSE
02049             IMI=IMIN
02050             XMI=XMIN
02051           ENDIF
02052           IF (LNG.EQ.1) WRITE(LU,883) XMI,IMI
02053           IF (LNG.EQ.2) WRITE(LU,884) XMI,IMI
02054 883       FORMAT(' EVOLUTION MINIMUM TOTALE : ',G16.7,' NOEUD : ',I6)
02055 884       FORMAT(' TOTAL MINIMAL EVOLUTION  : ',G16.7,' NODE  : ',I6)
02056         ENDIF
02057 !
02058 !=======================================================================
02059 ! : 10         PRINTS OUT RESULTS AT THIS TIMESTEP
02060 !              AND COMPARES AGAINST A REFERENCE FILE
02061 !=======================================================================
02062 !
02063         IF(UNIT) CALL OS('X=CX    ', X=CS, C= XMVS)
02064 !
02065         IF(DEBUG.GT.0) WRITE(LU,*) 'APPEL DE PREDES'
02066         CALL PREDES(LT,AT0)
02067         IF(DEBUG.GT.0) WRITE(LU,*) 'RETOUR DE PREDES'
02068         IF(DEBUG.GT.0) WRITE(LU,*) 'APPEL DE BIEF_DESIMP'
02069         CALL BIEF_DESIMP(SIS_FILES(SISRES)%FMT,VARSOR,
02070      &                   HIST,0,NPOIN,SIS_FILES(SISRES)%LU,
02071      &                   'STD',AT0,LT,LISPR,LEOPR,
02072      &                   SORLEO,SORIMP,MAXVAR,TEXTE,PTINIG,PTINIL)
02073         IF(DEBUG.GT.0) WRITE(LU,*) 'RETOUR DE BIEF_DESIMP'
02074 !
02075         IF(UNIT) CALL OS('X=CX    ', X=CS,  C= 1.D0/XMVS)
02076 !
02077 !       SENDS THE NEW ZF TO TELEMAC-2D OR 3D
02078 !
02079         IF(CODE(1:7).EQ.'TELEMAC') THEN
02080           CALL OV ('X=Y     ', ZF_TEL%R, ZF%R, ZF%R, 0.D0, NPOIN)
02081         ENDIF
02082 !
02083 !       THE SUBROUTINE VALIDA FROM THE LIBRARY IS STANDARD
02084 !       IT CAN BE MODIFIED FOR EACH PARTICULAR CASE
02085 !       BUT ITS CALL MUST BE LEFT IN THE LOOP ON TIME
02086 !
02087         IF(VALID) THEN
02088           IF(DEBUG.GT.0) WRITE(LU,*) 'APPEL DE BIEF_VALIDA'
02089           CALL BIEF_VALIDA(TB,TEXTPR,SIS_FILES(SISREF)%LU,
02090      &                     SIS_FILES(SISREF)%FMT,
02091      &                     VARSOR,TEXTE,SIS_FILES(SISRES)%LU,
02092      &                     SIS_FILES(SISRES)%FMT,
02093      &                     MAXVAR,NPOIN,LT,VALNIT,ALIRV)
02094           IF(DEBUG.GT.0) WRITE(LU,*) 'RETOUR DE BIEF_VALIDA'
02095         ENDIF
02096 !
02097 !       END OF THE LOOP ON THE RECORDS
02098         ENDDO !MT
02099 !
02100 !=======================================================================
02101 !
02102 !       END OF THE LOOP ON THE NUMBER OF EVENTS
02103 !
02104         ENDDO !MN
02105 !
02106 !-----------------------------------------------------------------------
02107 !
02108         IF(DEBUG.GT.0) WRITE(LU,*) 'END_TIME_LOOP'
02109       ENDIF
02110 !
02111 !-----------------------------------------------------------------------
02112 !
02113       IF(DREDGESIM.AND.(LOOPCOUNT.EQ.TELNIT.AND.PART.EQ.1.
02114      &                                      .OR.PART.EQ.-1)) THEN
02115         CALL DREDGESIM_INTERFACE(3)
02116       ENDIF
02117 !
02118 !-----------------------------------------------------------------------
02119 !
02120       RETURN
02121       END

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