bedload_solvs_fe.f

Go to the documentation of this file.

!                   ***************************
                    SUBROUTINE bedload_solvs_fe
!                   ***************************
!
     &(mesh,s,ebor,maskel,mask,qsx,qsy,ielmt,npoin,nptfr,kent,kdir,kddl,
     & limtec,dt,msk,entet,t1,t2,t3,t4,t8,zfcl,hz,hzn,gloseg,dimglo,
     & flodel,flulim,nseg,unsv2d,csf_sable,icla,flbcla,ava,liqbor,qbor,
     & maxadv)
!
!***********************************************************************
! SISYPHE   V6P2                                   21/07/2011
!***********************************************************************
!
!brief    SOLVES:
!code
!+     D(HZ)
!+     ---- + DIV(QS) = 0
!+      DT
!
!warning
!+     LIMTEC is used here instead of LIEBOR. The difference is that
!+     LIMTEC is LIEBOR corrected in view of sign of u.n at boundaries
!+     then KENT, KSORT apply to LIEBOR, while KDIR and KDDL apply on
!+     LIMTEC, see bedload_diffin.f
!
!history  E. PELTIER; C. LENORMANT; J.-M. HERVOUET
!+        11/09/1995
!+        V5P1
!+
!
!history  B. MINH DUC
!+        **/**/2002
!+        V5P3
!+
!
!history  F. HUVELIN
!+        14/09/2004
!+        V5P5
!+
!
!history  J.-M. HERVOUET
!+        29/10/2007
!+        V5P8
!+
!
!history  J.-M. HERVOUET
!+        05/09/2009
!+        V6P0
!+   NEW METHOD
!
!history  N.DURAND (HRW), S.E.BOURBAN (HRW)
!+        13/07/2010
!+        V6P0
!+   Translation of French comments within the FORTRAN sources into
!+   English comments
!
!history  N.DURAND (HRW), S.E.BOURBAN (HRW)
!+        21/08/2010
!+        V6P0
!+   Creation of DOXYGEN tags for automated documentation and
!+   cross-referencing of the FORTRAN sources
!
!history  C.VILLARET (EDF-LNHE), P.TASSI (EDF-LNHE)
!+        19/07/2011
!+        V6P1
!+  Name of variables
!+
!history  J-M HERVOUET (EDF-LNHE)
!+        27/01/2012
!+        V6P2
!+  Argument ICLA added
!
!history  J-M HERVOUET (EDF-LNHE)
!+        14/02/2012
!+        V6P2
!+  Optimisation, and FLBCLA built and kept for use in bilan_sisyphe
!+  Treatment of QBOR added
!
!history  J-M HERVOUET (EDF-LNHE)
!+        25/12/2012
!+        V6P3
!+  3 arguments added to VECTOS.
!
!history  L. STADLER (BAW)
!+        17/03/2016
!+        V7P2
!+  Call to flusec_sis added for new computation of discharges through
!+  cross-sections.
!
!history  J-M HERVOUET (EDF-LNHE)
!+        10/06/2016
!+        V7P2
!+  Cancelling sediment fluxes to and from dry nodes.
!
!history  J-M HERVOUET (EDF-LNHE)
!+        09/09/2016
!+        V7P2
!+  Adding fake arguments PLUIE and RAIN in the call to positive_depths.
!
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!| DIMGLO         |-->| FIRST DIMENSION OF GLOSEG
!| DT             |-->| TIME STEP
!| EBOR           |<->| BOUNDARY CONDITION FOR BED EVOLUTION (DIRICHLET)
!| ENTET          |-->| LOGICAL, IF YES INFORMATION IS GIVEN ON MASS CONSERVATION
!| FLBCLA         |<->| FLUXES AT BOUNDARY FOR THE CLASS
!| FLODEL         |<--| FLUXES BETWEEN POINTS (PER SEGMENT)
!| FLULIM         |<--| LIMITATION OF FLUXES
!| GLOSEG         |-->| CONNECTIVITY TABLE FOR SEGMENTS
!| HZ             |<--| NEW AVAILABLE LAYER OF SEDIMENT
!| HZN            |-->| OLD AVAILABLE LAYER OF SEDIMENT
!| ICLA           |-->| CLASS NUMBER
!| IELMT          |-->| NUMBER OF ELEMENTS
!| KDDL           |-->| CONVENTION FOR DEGREE OF FREEDOM
!| KDIR           |-->| CONVENTION FOR DIRICHLET POINT
!| KENT           |-->| CONVENTION FOR LIQUID INPUT WITH PRESCRIBED VALUE
!| LIMTEC         |-->| TYPE OF BOUNDARY CONDITION
!| LIQBOR         |-->| TYPE OF BOUNDARY CONDITION ON BEDLOAD DISCHARGE
!| MASK           |-->| BLOCK OF MASKS, EVERY ONE FOR A TYPE OF BOUNDARY
!|                |   | SEE DIFFIN.F IN LIBRARY BIEF.
!| MASKEL         |-->| MASKING OF ELEMENTS
!| MAXADV         |-->| MAXIMUM NUMBER OF ITERATIONS (IN POSITIVE_DEPTH)
!| MESH           |<->| MESH STRUCTURE
!| MSK            |-->| IF YES, THERE IS MASKED ELEMENTS
!| NPOIN          |-->| NUMBER OF POINTS
!| NPTFR          |-->| NUMBER OF BOUNDARY POINTS
!| NSEG           |-->| NUMBER OF SEGMENTS PER CONTROL SECTION
!| QBOR           |-->| PRESCRIBED BEDLOAD DISCHARGES
!| QSX            |-->| SOLID DISCHARGE X
!| QSY            |-->| SOLID DISCHARGE Y
!| S              |-->| VOID STRUCTURE
!| T1             |<->| WORK BIEF_OBJ STRUCTURE
!| T2             |<->| WORK BIEF_OBJ STRUCTURE
!| T3             |<->| WORK BIEF_OBJ STRUCTURE
!| T4             |<->| WORK BIEF_OBJ STRUCTURE
!| T8             |<->| WORK BIEF_OBJ STRUCTURE
!| UNSV2D         |-->| INVERSE OF INTEGRALS OF TEST FUNCTIONS
!| ZFCL           |<--| ZFCL=HZ-HZN
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!
      USE bief
      USE interface_sisyphe, ex_bedload_solvs_fe => bedload_solvs_fe
      USE declarations_sisyphe, ONLY : doflux,hn,hmin_bedload
!
      USE declarations_special
      IMPLICIT NONE
!
!+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
!
      TYPE(bief_mesh), INTENT(INOUT)  :: MESH
      TYPE(bief_obj),   INTENT(IN)    :: S,MASKEL,MASK,QSX,QSY
      INTEGER,          INTENT(IN)    :: IELMT,NPOIN,NPTFR,KENT,KDIR
      INTEGER,          INTENT(IN)    :: DIMGLO,NSEG,ICLA,KDDL,MAXADV
      INTEGER,          INTENT(IN)    :: GLOSEG(dimglo,2)
      DOUBLE PRECISION, INTENT(IN)    :: DT,CSF_SABLE,AVA(npoin)
      DOUBLE PRECISION, INTENT(INOUT) :: FLULIM(nseg)
      LOGICAL,          INTENT(IN)    :: MSK,ENTET
      TYPE(bief_obj),   INTENT(INOUT) :: FLODEL,T1,T2,T3,T4,T8
      TYPE(bief_obj),   INTENT(INOUT) :: HZ,EBOR,LIMTEC
      TYPE(bief_obj),   INTENT(INOUT) :: ZFCL,FLBCLA
      TYPE(bief_obj),   INTENT(IN)    :: HZN,UNSV2D,LIQBOR,QBOR
!
!+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
!
      INTEGER K,N,I1,I2,ISEG
!
!-----------------------------------------------------------------------
!
!     BOUNDARY FLUXES
!
      CALL vector(flbcla,'=','FLUBOR          ',ielbor(ielmt,1),1.d0,
     &            s,s,s,qsx,qsy,s,mesh,.true.,mask)
!
!     BOUNDARY CONDITIONS: EITHER EBOR OR QBOR PRESCRIBED (NOT THE 2)
!
      DO k=1,nptfr
        IF(liqbor%I(k).EQ.kent) THEN
!         QBOR IS GIVEN BY USER, AND POSITIVE IF ENTERING
!         HERE WE PUT THE INTERNAL USAGE <0 = ENTERING
          flbcla%R(k)=-qbor%R(k)
!         EVEN IF USER HAS SPECIFIED LIEBOR=KSORT, LIMTEC MAY HAVE BEEN
!         SET TO KDIR BY CHECKING IF VELOCITY IS ENTERING, THIS IS
!         UNWANTED HERE AS QBOR ONLY IS TAKEN INTO ACCOUNT,
!         SO DDL IS PUT TO AVOID A DIRICHLET TREATMENT IN POSITIVE_DEPTHS.
          limtec%I(k)=kddl
        ELSEIF(limtec%I(k).EQ.kdir) THEN
!         HERE THE VARIABLE WILL BE THE LAYER DEPTH OF THE SEDIMENT CLASS,
!         PUT IN T8, NOT THE EVOLUTION
          n=mesh%NBOR%I(k)
          t8%R(k)=ava(n)*ebor%R(k)*csf_sable+hzn%R(n)
        ENDIF
      ENDDO
!
!     HERE T1 MAY NOT BE ASSEMBLED, WE WORK DIRECTLY ON MESH%W%R AFTER,
!     FOR CALLING FLUX_EF_VF (IT IS T1 IN NON ASSEMBLED FORM)
!
      CALL vector(t1,'=','VGRADP          ',qsx%ELM,-1.d0,
     &            s,s,s,qsx,qsy,s,mesh,msk,maskel,lego=.false.)
!
!     !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!
!     FLODEL COMPUTED HERE, NOT IN POSITIVE_DEPTHS
!
      CALL flux_ef_vf(flodel%R,mesh%W%R,mesh%NSEG,mesh%NELEM,
     &                mesh%NELMAX,mesh%ELTSEG%I,mesh%ORISEG%I,
     &                mesh%IKLE%I,.true.,2)
!
!     CANCELLING THE FLUXES TO AND FROM DRY POINTS
!
      DO iseg=1,nseg
        i1=gloseg(iseg,1)
        i2=gloseg(iseg,2)
        IF(hn%R(i1).LT.hmin_bedload.OR.
     &     hn%R(i2).LT.hmin_bedload) flodel%R(iseg)=0.d0
      ENDDO
!
!     END OF NEW SECTION (ONLY TRUE CHANGED INTO FALSE AFTER FLODEL IN THE NEXT CALL)
!
!     !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!
      CALL positive_depths(t1,t2,t3,t4,hz,hzn,mesh,
!                                  !!!!!!
!    &                     FLODEL, .TRUE.,FLBCLA,DT,UNSV2D,NPOIN,
     &                     flodel,.false.,flbcla,dt,unsv2d,npoin,
     &                     gloseg(1:dimglo,1),gloseg(1:dimglo,2),
     &                     mesh%NBOR%I,nptfr,t8,.false.,t8,.false.,
!                                            VOID (SMH) VOID (PLUIE)
     &                     1,flulim,
     &                     limtec%I,t8%R  ,kdir,entet,mesh%W%R,
!                                   EBOR%R
     &                     'SISYPHE                 ',2,maxadv)
!                                                     2 : HARDCODED
!                             OPTION FOR POSITIVE DEPTHS ALGORITHMS
!                             HERE CHOICE OF OPTION INDEPENDENT OF
!                             SEGMENT NUMBERING
!
      CALL os('X=Y-Z   ' ,x=zfcl,y=hz,z=hzn)
!
!-----------------------------------------------------------------------
!
!     NEW FLUXES ACROSS CROSS-SECTIONS
!
      IF(doflux) THEN
        CALL flusec_sis(gloseg,dimglo,dt,mesh,
     &                  flodel,icla,entet)
      ENDIF
!
!-----------------------------------------------------------------------
!
      RETURN
      END