source_waq.F

Go to the documentation of this file.

!                 **********************
                  SUBROUTINE source_waq
!                 **********************
!
     &(npoin3,npoin2,texp,timp,tn,hprop,u,v,cf,
     & t1,t2,t3,t4,t5,t6,t7,t8,t9,t10,t11,t12, t2_1,t2_2,t2_3,
     & patmos,dimm,nplan,
     & latit,longit,at,mardat,martim,zprop,rho)
!
!***********************************************************************
! WAQTEL      V8P2
!***********************************************************************
!
!brief    GIVES CONTRIBUTION OF WAQ PROCESSES TO SOURCE TERMS
!+                FOR THE TRACER.
!
!history  R. ATA
!+        21/09/2014
!+        V7P0
!+   CREATION
!
!history  M. JODEAU (LNHE)
!+        18/05/2016
!+        V7P3
!+   + AED2 source terms
!
!history  S.E. BOURBAN (HRW)
!+        21/09/2017
!+        V7P3
!+        WAQPROCESS is now a prime number, so that multiple processes
!+        can be called by multiplication of the prime numbers.
!
!history  C.-T. PHAM (LNHE)
!+        31/07/2019
!+        V8P1
!+        Density RHO computed in drsurr in TELEMAC-3D and given to
!+        WAQTEL as optional (none in 2D) rather than computed again
!
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!| AT             |-->| TIME IN SECONDS
!| DIMM           |-->| 2D OR 3D
!| HPROP          |-->| PROPAGATION DEPTH (2D)
!| LONGIT         |-->| LONGITUTE OF ORIGIN POINT
!| LATIT          |-->| LATITUDE OF ORIGIN POINT
!| MAXSCE         |-->| MAXIMUM NUMBER OF SOURCES
!| MAXTRA         |-->| MAXIMUM NUMBER OF TRACERS
!| NPLAN          |-->| NUMBER OF VERTICAL PLANES
!| RHO            |-->| WATER DENSITY
!| T1,..,T12      |<->| WORKING STRUCTURES
!| T2_1,T2_2      |<->| 2D WORKING STRUCTURES
!| TETAT          |-->| COEFFICIENT OF IMPLICITATION FOR TRACERS.
!| TEXP           |-->| EXPLICIT SOURCE TERM.
!| TIMP           |-->| IMPLICIT SOURCE TERM.
!| TN             |-->| TRACERS AT TIME N
!| TSCE           |-->| PRESCRIBED VALUES OF TRACERS AT POINT SOURCES
!| TSCEXP         |<--| EXPLICIT SOURCE TERM OF POINT SOURCES
!|                |   | IN TRACER EQUATION, EQUAL TO:
!|                |   | TSCE - ( 1 - TETAT ) TN
!| YASMI          |<--| IF YES, THERE ARE IMPLICIT SOURCE TERMS
!| ZPROP          |-->| Z COORDINATES FOR 3D NODES
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!
      USE bief
      USE declarations_special
      USE declarations_waqtel
      USE interface_waqtel, ex_source_waq => source_waq
!
      IMPLICIT NONE
!
!+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
!
      INTEGER       , INTENT(IN)    :: NPOIN3,DIMM
      INTEGER       , INTENT(IN)    :: NPLAN,NPOIN2
      INTEGER       , INTENT(IN)    :: MARDAT(3),MARTIM(3)
      DOUBLE PRECISION,INTENT(IN  ) :: LATIT,LONGIT,AT
      TYPE(bief_obj), INTENT(IN)    :: TN,HPROP,CF,U,V
      TYPE(bief_obj), INTENT(IN)    :: PATMOS,ZPROP
      TYPE(bief_obj), INTENT(INOUT) :: TEXP,TIMP
      TYPE(bief_obj), INTENT(INOUT) :: T1,T2,T3,T4,T5,T6
      TYPE(bief_obj), INTENT(INOUT) :: T7,T8,T9,T10,T11,T12
      TYPE(bief_obj), INTENT(INOUT) :: T2_1,T2_2,T2_3
!     OPTIONAL
      TYPE(bief_obj), INTENT(IN),OPTIONAL :: RHO
!
!+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
!
      INTEGER :: J,ITRAC,IERR
      LOGICAL :: FOUND
!
      found = .false.
!
!-----------------------------------------------------------------------
!
!     O2 MODULE
!
      IF( 2*int(waqprocess/2).EQ.waqprocess ) THEN
        IF(debug.GT.0) WRITE(lu,*) 'CALL OF CALCS~D_O2'
!
        IF(dimm.EQ.2)THEN
          CALL calcs2d_o2
     &         (npoin2,wattemp,o2satu,demben,formk2,k1,k44,k22,
     &          photo,resp,tn,texp,timp,t2,t3,t4,hprop,u,v,debug)
!
!          YASMI(IND_O2) = YASMI(IND_O2) .OR. .FALSE.
!          YASMI(IND_OL) = .TRUE.
!          YASMI(IND_NH4) = .TRUE.
!
        ELSEIF(dimm.EQ.3)THEN
          CALL calcs3d_o2
     &         (npoin3,npoin2,nplan,wattemp,o2satu,demben,
     &          formk2,k1,k44,k22,photo,resp,tn,texp,timp,
     &          t1,t2,t2_1,hprop,zprop,u,v)
!
          texp%ADR(ind_o2)%P%TYPR='Q'
          timp%ADR(ind_ol)%P%TYPR='Q'
          timp%ADR(ind_nh4)%P%TYPR='Q'
!
        ENDIF
        found = .true.
        IF(debug.GT.0) WRITE(lu,*) 'BACK FROMCALCS~D_O2'
      ENDIF
!
!-----------------------------------------------------------------------
!
!     BIOMASS MODULE
!
      IF( 3*int(waqprocess/3).EQ.waqprocess ) THEN
        IF(debug.GT.0) WRITE(lu,*) 'CALL OF CALCS~D_BIOMASS'
!
        IF(dimm.EQ.2)THEN
          CALL calcs2d_biomass(npoin2,wattemp,tn,texp,rayeff,
     &                         hprop,t1,t2,t3,t4,t5,t6,debug)
!
!          YASMI(IND_PHY) = YASMI(IND_PHY) .OR. .FALSE.
!          YASMI(IND_PO4) = YASMI(IND_PO4) .OR. .FALSE.
!          YASMI(IND_POR) = YASMI(IND_POR) .OR. .FALSE.
!          YASMI(IND_NO3) = YASMI(IND_NO3) .OR. .FALSE.
!          YASMI(IND_NOR) = YASMI(IND_NOR) .OR. .FALSE.
!
        ELSEIF(dimm.EQ.3)THEN
          CALL calcs3d_biomass(npoin3,npoin2,nplan,wattemp,tn,texp,
     &                         rayeff,zprop,t1,t2,t3,t4,t5,t6)
!
          texp%ADR(ind_phy)%P%TYPR='Q'
          texp%ADR(ind_po4)%P%TYPR='Q'
          texp%ADR(ind_por)%P%TYPR='Q'
          texp%ADR(ind_no3)%P%TYPR='Q'
          texp%ADR(ind_nor)%P%TYPR='Q'
!
        ENDIF
        found = .true.
        IF(debug.GT.0) WRITE(lu,*) 'BACK FROM CALCS~D_BIOMASS'
      ENDIF
!
!-----------------------------------------------------------------------
!
!     EUTRO MODULE
!
      IF( 5*int(waqprocess/5).EQ.waqprocess ) THEN
        IF(debug.GT.0) WRITE(lu,*) 'CALL OF CALCS~D_EUTRO'
!
        IF(dimm.EQ.2)THEN
          CALL calcs2d_eutro(npoin2,wattemp,tn,texp,timp,rayeff,
     &                       hprop,t1,t2,t3,t4,t5,t6,t7,t8,t9,t10,
     &                       t11,t12,debug,u,v)
!
!          YASMI(IND_PHY) = YASMI(IND_PHY) .OR. .FALSE.
!          YASMI(IND_PO4) = YASMI(IND_PO4) .OR. .FALSE.
!          YASMI(IND_POR) = YASMI(IND_POR) .OR. .FALSE.
!          YASMI(IND_NO3) = YASMI(IND_NO3) .OR. .FALSE.
!          YASMI(IND_NOR) = YASMI(IND_NOR) .OR. .FALSE.
!          YASMI(IND_NH4) = .TRUE.
!          YASMI(IND_OL) = .TRUE.
!          YASMI(IND_O2) = YASMI(IND_O2) .OR. .FALSE.
!
        ELSEIF(dimm.EQ.3)THEN
          CALL calcs3d_eutro(npoin3,npoin2,nplan,wattemp,tn,texp,
     &                       timp,rayeff,hprop,zprop,t1,t2_1,t3,t4,
     &                       t5,t6,t7,t8,t9,t10,t11,t12,debug,u,v)
!
          texp%ADR(ind_phy)%P%TYPR='Q'
          texp%ADR(ind_po4)%P%TYPR='Q'
          texp%ADR(ind_por)%P%TYPR='Q'
          texp%ADR(ind_no3)%P%TYPR='Q'
          texp%ADR(ind_nor)%P%TYPR='Q'
          texp%ADR(ind_nh4)%P%TYPR='Q'
          timp%ADR(ind_nh4)%P%TYPR='Q'
          texp%ADR(ind_ol)%P%TYPR='Q'
          timp%ADR(ind_ol)%P%TYPR='Q'
          texp%ADR(ind_o2)%P%TYPR='Q'
!
        ENDIF
        found = .true.
        IF(debug.GT.0) WRITE(lu,*) 'BACK FROM CALCS~D_EUTRO'
      ENDIF
!
!-----------------------------------------------------------------------
!
!     MICROPOL MODULE
!
      IF( 7*int(waqprocess/7).EQ.waqprocess ) THEN
        IF(debug.GT.0) WRITE(lu,*) 'CALL OF CALCS~D_MICROPOL'
!
        IF(dimm.EQ.2)THEN
          CALL calcs2d_micropol(npoin2,tn,texp,timp,hprop,
     &                          cf,u,v,t1,t2,t3,t4)
!
!          YASMI(IND_SS) = YASMI(IND_SS) .OR. .FALSE.
!          YASMI(IND_SF) = YASMI(IND_SF) .OR. .FALSE.
!          YASMI(IND_C) = .TRUE.
!          YASMI(IND_CSS) = .TRUE.
!          YASMI(IND_CSF) = YASMI(IND_CSF) .OR. .FALSE.
!
        ELSEIF(dimm.EQ.3)THEN
          CALL calcs3d_micropol(npoin2,tn,texp,
     &                          timp,zprop,cf,u,v,t2_1,t2_2,t2_3,t1,t2,
     &                          debug)
!
          texp%ADR(ind_ss)%P%TYPR='Q'
          texp%ADR(ind_sf)%P%TYPR='Q'
          timp%ADR(ind_c)%P%TYPR='Q'
          texp%ADR(ind_c)%P%TYPR='Q'
          timp%ADR(ind_css)%P%TYPR='Q'
          texp%ADR(ind_css)%P%TYPR='Q'
          texp%ADR(ind_csf)%P%TYPR='Q'
!
        ENDIF
        found = .true.
        IF(debug.GT.0) WRITE(lu,*) 'BACK FROM CALCS~D_MICROPOL'
      ENDIF
!
!-----------------------------------------------------------------------
!
!     THERMIC MODULE
!
      IF( 11*int(waqprocess/11).EQ.waqprocess ) THEN
!
        IF(dimm.EQ.2)THEN
          IF(debug.GT.0) WRITE(lu,*) 'CALL OF CALCS2D_THERMIC'
          CALL calcs2d_thermic(npoin2,tn,texp,hprop,patmos)
          IF(debug.GT.0) WRITE(lu,*) 'BACK FROM CALCS2D_THERMIC'
!
!          YASMI(IND_T) = YASMI(IND_T) .OR. .FALSE.
!
        ELSEIF(dimm.EQ.3)THEN
!         SOURCE TERMS (VOLUME)
          IF(debug.GT.0) WRITE(lu,*) 'CALL OF CALCS3D_THERMICV'
          IF(atmosexch.EQ.2) THEN 
            CALL calcs3d_thermicv(npoin3,npoin2,nplan,zprop%R,rho%R,
     &                            tn,texp,longit,latit,
     &                            at,mardat,martim)
!
            texp%ADR(ind_t)%P%TYPR='Q'
          ENDIF
!
          IF(debug.GT.0) WRITE(lu,*) 'BACK FROM CALCS3D_THERMICV'
!
        ENDIF
        found = .true.
      ENDIF
!
!-----------------------------------------------------------------------
!
!     AED2 COUPLING
!
      IF( 13*int(waqprocess/13).EQ.waqprocess ) THEN
#if defined HAVE_AED2
        IF(.NOT.deja_sw) THEN
          ALLOCATE(extcaed2(npoin2,nplan),stat=ierr)
          deja_sw = .true.
        ENDIF
        IF(dimm.EQ.3)THEN
          IF(debug.GT.0) WRITE(lu,*) 'CALL OF CALCS3D_AED2'
          CALL calcs3d_aed2(npoin3,npoin2,nplan,zprop%R,rho%R,
     &                      tn,texp,timp,longit,latit,at,
     &                      mardat,martim,extcaed2,rayaed2)
          IF(debug.GT.0) WRITE(lu,*) 'BACK FROM CALCS3D_AED2'
        ENDIF
        found = .true.
#else
        WRITE(lu,*) 'ERROR: AED2 LIBRARY NOT COMPILED'
        CALL plante(1)
#endif
      ENDIF
!
!-----------------------------------------------------------------------
!
!     DEGRADATION LAW
!
      IF( 17*int(waqprocess/17).EQ.waqprocess ) THEN
!
        IF(dimm.EQ.2) THEN
          IF(debug.GT.0) WRITE(lu,*) 'CALL OF CALCS2D_DEGRADATION'
          CALL calcs2d_degradation(npoin2,tn,texp,timp,hprop,nwaq_degra,
     &                             rank_degra,loitrac,coef1trac)
          IF(debug.GT.0) WRITE(lu,*) 'BACK FROM CALCS2D_DEGRADATION'
!          DO J = 1,NWAQ_DEGRA
!            ITRAC = RANK_DEGRA(J)
!            YASMI(ITRAC) = .TRUE.
!          ENDDO
        ELSEIF(dimm.EQ.3) THEN
          IF(debug.GT.0) WRITE(lu,*) 'CALL OF CALCS3D_DEGRADATION'
          CALL calcs3d_degradation(npoin3,npoin2,nplan,tn,texp,timp,
     &                             zprop,nwaq_degra,rank_degra,
     &                             loitrac,coef1trac)
          IF(debug.GT.0) WRITE(lu,*) 'BACK FROM CALCS3D_DEGRADATION'
!
          DO j = 1,nwaq_degra
            itrac = rank_degra(j)
            timp%ADR(itrac)%P%TYPR='Q'
          ENDDO
        ENDIF
        found = .true.
!
      ENDIF
!
!-----------------------------------------------------------------------
!
!     GHOST PROCESS IN WAITING FOR THE MERGE WITH ICE PROCESS
!
      IF( 19*int(waqprocess/19).EQ.waqprocess ) THEN
        found = .true.
      ENDIF
!
!-----------------------------------------------------------------------
!
!      UNKNOWN PROCESS
!
      IF( .NOT.found ) THEN
        WRITE(lu,20) waqprocess
20    FORMAT(1x,'SOURCE_WAQ: UNKNOWN WAQ MODULE: ',i4)
        CALL plante(1)
        stop
      ENDIF
!
!
!-----------------------------------------------------------------------
!
      RETURN
      END