SURFEX v7.3
General documentation of Surfex
 All Classes Files Functions Variables Typedefs
/home/dasprezs/EXPORT_v7_3/src/SURFEX/roof_impl_coef.F90
Go to the documentation of this file.
00001 !   ###############################################
00002 SUBROUTINE ROOF_IMPL_COEF(PTSTEP,PTDEEP_A,PTDEEP_B)
00003 !   ###############################################
00004 !
00005 !!
00006 !!    PURPOSE
00007 !!    -------
00008 !
00009 !     Computes the corfficients for implicitation of upper
00010 !     roof layer with what is above
00011 !         
00012 !     
00013 !!**  METHOD
00014 !     ------
00015 !
00016 !  One computes a guess assuming a zero flux condition at the base
00017 !  of the roof. One solves the half part of the tridiagonal matrix
00018 !  fromm bottom to top.
00019 !
00020 !!     The classical tridiagonal algorithm is used to invert the 
00021 !!     implicit operator (from bottom to top only). Its matrix is given by:
00022 !!
00023 !!     (  b(1)      c(1)      0        0        0         0        0        0  )
00024 !!     (  a(2)      b(2)     c(2)      0  ...    0        0        0        0  ) 
00025 !!     (   0        a(3)     b(3)     c(3)       0        0        0        0  ) 
00026 !!      .......................................................................
00027 !!     (   0   ...   0      a(k)      b(k)     c(k)       0   ...  0        0  ) 
00028 !!      .......................................................................
00029 !!     (   0         0        0        0        0 ...  a(n-1)   b(n-1)   c(n-1))
00030 !!     (   0         0        0        0        0 ...     0      a(n)     b(n) )
00031 !!
00032 !!
00033 !!    AUTHOR
00034 !!    ------
00035 !!
00036 !!      V. Masson           * Meteo-France *
00037 !!
00038 !!    MODIFICATIONS
00039 !!    -------------
00040 !!      Original    01/2013
00041 !!
00042 !-------------------------------------------------------------------------------
00043 !
00044 !*       0.     DECLARATIONS
00045 !               ------------
00046 !
00047 USE YOMHOOK   ,ONLY : LHOOK,   DR_HOOK
00048 USE PARKIND1  ,ONLY : JPRB
00049 !
00050 USE MODD_TEB_n, ONLY : XD_ROOF, XTC_ROOF, XHC_ROOF, XT_ROOF, NROOF_LAYER
00051 !
00052 USE MODI_LAYER_E_BUDGET_GET_COEF
00053 !
00054 IMPLICIT NONE
00055 !
00056 !
00057 !*      0.1    Declarations of arguments
00058 !
00059 REAL              , INTENT(IN)  :: PTSTEP             ! time step
00060 REAL, DIMENSION(:), INTENT(OUT) :: PTDEEP_A, PTDEEP_B     
00061                                               ! Deep soil temperature (prescribed)
00062 !                                      PTDEEP_A = Deep soil temperature
00063 !                                                 coefficient depending on flux
00064 !                                      PTDEEP_B = Deep soil temperature (prescribed)
00065 !                                               which models heating/cooling from
00066 !                                               below the diurnal wave penetration
00067 !                                               (surface temperature) depth. If it
00068 !                                               is FLAGGED as undefined, then the zero
00069 !                                               flux lower BC is applied.
00070 !                                      Tdeep = PTDEEP_B + PTDEEP_A * PDEEP_FLUX
00071 !                                              (with PDEEP_FLUX in W/m2)
00072 !
00073 !*      0.2    Local variables
00074 !
00075 REAL :: ZIMPL = 1.0        ! implicit coefficient
00076 INTEGER :: JK              ! loop counter
00077 !
00078 REAL, DIMENSION(SIZE(PTDEEP_A),NROOF_LAYER) :: ZA, ! lower diag.
00079                                                ZB, ! main  diag.
00080                                                ZC, ! upper diag.
00081                                                ZY   ! r.h.s.
00082 
00083 REAL, DIMENSION(SIZE(PTDEEP_A))             :: ZDET ! work array
00084 REAL, DIMENSION(SIZE(PTDEEP_A),NROOF_LAYER) :: ZW   ! work array
00085 REAL, DIMENSION(SIZE(PTDEEP_A),NROOF_LAYER) :: ZT   ! guess of T
00086 !
00087 REAL(KIND=JPRB) :: ZHOOK_HANDLE
00088 !-------------------------------------------------------------------------------
00089 IF (LHOOK) CALL DR_HOOK('ROOF_IMPL_COEF',0,ZHOOK_HANDLE)
00090 !
00091 !*      1.0    Coefficients of the tridioagonal matrix for heat conduction eq.
00092 !              ---------------------------------------------------------------
00093 !
00094  CALL LAYER_E_BUDGET_GET_COEF( XT_ROOF, PTSTEP, ZIMPL, XHC_ROOF, XTC_ROOF, XD_ROOF, &
00095                               ZA, ZB, ZC, ZY )
00096 !
00097 !-------------------------------------------------------------------------------
00098 !
00099 !*      2.0    Solving of the equation from bottom to top
00100 !              ------------------------------------------
00101 !
00102 ! layer at bottom of roof
00103 !
00104 ZDET(:)   = ZB(:,NROOF_LAYER)
00105 !
00106 ZT  (:,NROOF_LAYER) = ZY(:,NROOF_LAYER) / ZDET(:)
00107 !
00108 ! internal layers & top layer (but without the external heat flux term)
00109 !
00110 DO JK=NROOF_LAYER-1,1,-1
00111   ZW  (:,JK)  = ZA(:,JK+1)/ZDET(:)
00112   ZDET(:)     = ZB(:,JK  ) - ZC(:,JK)*ZW(:,JK)
00113   ZT  (:,JK)  = ( ZY(:,JK) - ZC(:,JK)*ZT(:,JK+1) ) / ZDET(:) ! + FLUX / ZDET
00114   !                                                          ! for layer 1
00115   !                                                          ! because the external
00116   !                                                          ! flux would be
00117   !                                                          ! included in the Y
00118   !                                                          ! term
00119 END DO
00120 !
00121 ! Implicit coefficients for the heat flux
00122 !
00123 PTDEEP_B = ZT  (:,1)
00124 PTDEEP_A = 1. / ZDET(:) 
00125 !
00126 !* The following lines are here if you want to test the explicit coupling
00127 !PTDEEP_B = XT_ROOF(:,1)
00128 !PTDEEP_A = 0.
00129 !-------------------------------------------------------------------------------
00130 IF (LHOOK) CALL DR_HOOK('ROOF_IMPL_COEF',1,ZHOOK_HANDLE)
00131 !-------------------------------------------------------------------------------
00132 !
00133 END SUBROUTINE ROOF_IMPL_COEF