Ignore:
Timestamp:
Dec 17, 2008, 2:30:13 PM (16 years ago)
Author:
Laurent Fairhead
Message:
  • Modifications lie au premier niveau du modele pour la diffusion turbulent

du vent.

  • Preparation pour un couplage des courrant oceaniques.

JG

File:
1 edited

Legend:

Unmodified
Added
Removed

  • LMDZ4/trunk/libf/phylmd/climb_wind_mod.F90

    r793 r1067  
    1 !
    2 ! $Header$
    31!
    42MODULE climb_wind_mod
    … …  
    2119  REAL, DIMENSION(:,:), ALLOCATABLE  :: Ccoef_V, Dcoef_V
    2220  !$OMP THREADPRIVATE(Ccoef_V, Dcoef_V)
     21  REAL, DIMENSION(:), ALLOCATABLE   :: Acoef_U, Bcoef_U
     22  !$OMP THREADPRIVATE(Acoef_U, Bcoef_U)
     23  REAL, DIMENSION(:), ALLOCATABLE   :: Acoef_V, Bcoef_V
     24  !$OMP THREADPRIVATE(Acoef_V, Bcoef_V)
    2325  LOGICAL                            :: firstcall=.TRUE.
    2426  !$OMP THREADPRIVATE(firstcall)
    2527
    2628 
    27   PUBLIC :: climb_wind_down, calcul_wind_flux, climb_wind_up
     29  PUBLIC :: climb_wind_down, climb_wind_up
    2830
    2931CONTAINS
    … …  
    6264    IF (ierr /= 0) CALL abort_gcm(modname,'Pb in allocation Dcoef_V',1)
    6365
     66    ALLOCATE(Acoef_U(klon), Bcoef_U(klon), Acoef_V(klon), Bcoef_V(klon), STAT=ierr)
     67    IF ( ierr /= 0 )  PRINT*,' pb in allloc Acoef_U and Bcoef_U, ierr=', ierr
     68
    6469    firstcall=.FALSE.
    6570
    … …  
    6873!****************************************************************************************
    6974!
    70   SUBROUTINE climb_wind_down(knon, dtime, coef_in, pplay, paprs, temp, delp, u_old, v_old)
     75  SUBROUTINE climb_wind_down(knon, dtime, coef_in, pplay, paprs, temp, delp, u_old, v_old, &
     76       Acoef_U_out, Acoef_V_out, Bcoef_U_out, Bcoef_V_out)
    7177!
    7278! This routine calculates for the wind components u and v,
    … …  
    8894    REAL, DIMENSION(klon,klev), INTENT(IN)   :: v_old
    8995
     96! Output arguments
     97!****************************************************************************************
     98    REAL, DIMENSION(klon), INTENT(OUT)       :: Acoef_U_out
     99    REAL, DIMENSION(klon), INTENT(OUT)       :: Acoef_V_out
     100    REAL, DIMENSION(klon), INTENT(OUT)       :: Bcoef_U_out
     101    REAL, DIMENSION(klon), INTENT(OUT)       :: Bcoef_V_out
     102
    90103! Local variables
    91104!****************************************************************************************
    … …  
    102115!
    103116!****************************************************************************************
    104 
    105117! - Define alpha (alf1 and alf2)
    106118    alf1(:) = 1.0
    107119    alf2(:) = 1.0 - alf1(:)
    108120
    109 ! - Calculte the wind components for the first layer
    110     u1lay(1:knon) = u_old(1:knon,1)*alf1(1:knon) + u_old(1:knon,2)*alf2(1:knon)
    111     v1lay(1:knon) = v_old(1:knon,1)*alf1(1:knon) + v_old(1:knon,2)*alf2(1:knon)   
    112 
    113 ! - Calculate K
     121! - Calculate the coefficients K
    114122    Kcoefm(:,:) = 0.0
    115     DO i = 1, knon
    116        Kcoefm(i,1) = coef_in(i,1) * (1.0+SQRT(u1lay(i)**2+v1lay(i)**2))* &
    117             pplay(i,1)/(RD*temp(i,1))
    118        Kcoefm(i,1) = Kcoefm(i,1) * dtime*RG
    119     END DO
    120 
    121123    DO k = 2, klev
    122124       DO i=1,knon
    123           Kcoefm(i,k) = coef_in(i,k)*RG/(pplay(i,k-1)-pplay(i,k)) &
     125          Kcoefm(i,k) = coef_in(i,k)*RG*RG*dtime/(pplay(i,k-1)-pplay(i,k)) &
    124126               *(paprs(i,k)*2/(temp(i,k)+temp(i,k-1))/RD)**2
    125           Kcoefm(i,k) = Kcoefm(i,k) * dtime*RG
    126127       END DO
    127128    END DO
    … …  
    130131    CALL calc_coef(knon, Kcoefm(:,:), delp(:,:), &
    131132         u_old(:,:), alf1(:), alf2(:),  &
    132          Ccoef_U(:,:), Dcoef_U(:,:))
     133         Ccoef_U(:,:), Dcoef_U(:,:), Acoef_U(:), Bcoef_U(:))
    133134
    134135! - Calculate the coefficients C and D, component "v"
    135136    CALL calc_coef(knon, Kcoefm(:,:), delp(:,:), &
    136137         v_old(:,:), alf1(:), alf2(:),  &
    137          Ccoef_V(:,:), Dcoef_V(:,:))
     138         Ccoef_V(:,:), Dcoef_V(:,:), Acoef_V(:), Bcoef_V(:))
     139
     140!****************************************************************************************
     141! 6)
     142! Return the first layer in output variables
     143!
     144!****************************************************************************************
     145    Acoef_U_out = Acoef_U
     146    Bcoef_U_out = Bcoef_U
     147    Acoef_V_out = Acoef_V
     148    Bcoef_V_out = Bcoef_V
    138149
    139150  END SUBROUTINE climb_wind_down
    … …  
    141152!****************************************************************************************
    142153!
    143   SUBROUTINE calc_coef(knon, Kcoef, dels, X, alfa1, alfa2, Ccoef, Dcoef)
    144 !
    145 ! Find the coefficients C and D in fonction of alfa, K and dels
     154  SUBROUTINE calc_coef(knon, Kcoef, delp, X, alfa1, alfa2, Ccoef, Dcoef, Acoef, Bcoef)
     155!
     156! Find the coefficients C and D in fonction of alfa, K and delp
    146157!
    147158! Input arguments
    148159!****************************************************************************************
    149160    INTEGER, INTENT(IN)                      :: knon
    150     REAL, DIMENSION(klon,klev), INTENT(IN)   :: Kcoef, dels
     161    REAL, DIMENSION(klon,klev), INTENT(IN)   :: Kcoef, delp
    151162    REAL, DIMENSION(klon,klev), INTENT(IN)   :: X
    152163    REAL, DIMENSION(klon), INTENT(IN)        :: alfa1, alfa2
    … …  
    154165! Output arguments
    155166!****************************************************************************************
     167    REAL, DIMENSION(klon), INTENT(OUT)       :: Acoef, Bcoef
    156168    REAL, DIMENSION(klon,klev), INTENT(OUT)  :: Ccoef, Dcoef
    157169 
    … …  
    161173    REAL                                     :: buf
    162174
     175    INCLUDE "YOMCST.h"
    163176!****************************************************************************************
    164177!
    165 ! Niveau au sommet, k=klev
     178
     179! Calculate coefficients C and D at top level, k=klev
    166180!
    167181    Ccoef(:,:) = 0.0
    … …  
    169183
    170184    DO i = 1, knon
    171        buf = dels(i,klev) + Kcoef(i,klev)
    172 
    173        Ccoef(i,klev) = X(i,klev)*dels(i,klev)/buf
     185       buf = delp(i,klev) + Kcoef(i,klev)
     186
     187       Ccoef(i,klev) = X(i,klev)*delp(i,klev)/buf
    174188       Dcoef(i,klev) = Kcoef(i,klev)/buf
    175189    END DO
    176190   
    177191!
    178 ! Niveau (klev-1) <= k <= 2
     192! Calculate coefficients C and D at top level (klev-1) <= k <= 2
    179193!
    180194    DO k=(klev-1),2,-1
    181195       DO i = 1, knon
    182           buf = dels(i,k) + Kcoef(i,k) + &
    183                Kcoef(i,k+1)*(1.-Dcoef(i,k+1))
     196          buf = delp(i,k) + Kcoef(i,k) + Kcoef(i,k+1)*(1.-Dcoef(i,k+1))
    184197         
    185           Ccoef(i,k) = X(i,k)*dels(i,k) + &
    186                Kcoef(i,k+1)*Ccoef(i,k+1)
    187           Ccoef(i,k) = Ccoef(i,k)/buf         
     198          Ccoef(i,k) = (X(i,k)*delp(i,k) + Kcoef(i,k+1)*Ccoef(i,k+1))/buf
    188199          Dcoef(i,k) = Kcoef(i,k)/buf
    189200       END DO
    190201    END DO
    191202
    192 ! 
    193 ! Niveau k=1
     203!
     204! Calculate coeffiecent A and B at surface
    194205!
    195206    DO i = 1, knon
    196        buf = dels(i,1) + &
    197             (alfa1(i) + alfa2(i)*Dcoef(i,2)) * Kcoef(i,1) + &
    198             (1.-Dcoef(i,2))*Kcoef(i,2)
    199        
    200        Ccoef(i,1) = X(i,1)*dels(i,1) + &
    201             (Kcoef(i,2)-Kcoef(i,1)*alfa2(i)) * Ccoef(i,2)
    202        Ccoef(i,1) = Ccoef(i,1)/buf
    203        Dcoef(i,1) = Kcoef(i,1)/buf
     207       buf = delp(i,1) + Kcoef(i,2)*(1-Dcoef(i,2))
     208       Acoef(i) = (X(i,1)*delp(i,1) + Kcoef(i,2)*Ccoef(i,2))/buf
     209       Bcoef(i) = -RG/buf
    204210    END DO
    205211
    … …  
    208214!****************************************************************************************
    209215!
    210   SUBROUTINE calcul_wind_flux(knon, dtime, &
    211        flux_u, flux_v)
    212 
    213     INCLUDE "YOMCST.h"
    214 
    215 ! Input arguments
    216 !****************************************************************************************
    217     INTEGER, INTENT(IN)                  :: knon
    218     REAL, INTENT(IN)                     :: dtime
    219 
    220 ! Output arguments
    221 !****************************************************************************************
    222     REAL, DIMENSION(klon), INTENT(OUT)   :: flux_u
    223     REAL, DIMENSION(klon), INTENT(OUT)   :: flux_v
    224 
    225 ! Local variables
    226 !****************************************************************************************
    227     INTEGER                              :: i
    228     REAL, DIMENSION(klon)                :: u0, v0  ! u and v at niveau 0
    229     REAL, DIMENSION(klon)                :: u1, v1  ! u and v at niveau 1
    230     REAL, DIMENSION(klon)                :: u2, v2  ! u and v at niveau 2
    231    
    232 
    233 !****************************************************************************************
    234 ! Les vents de surface sont supposes nuls
    235 !
    236 !****************************************************************************************
    237     u0(:) = 0.0
    238     v0(:) = 0.0
    239 
    240 !****************************************************************************************
    241 ! On calcule les vents du couhes 1 et 2 recurviement
    242 !
    243 !****************************************************************************************
    244     DO i = 1, knon
    245        u1(i) = Ccoef_U(i,1) + Dcoef_U(i,1)*u0(i)
    246        v1(i) = Ccoef_V(i,1) + Dcoef_V(i,1)*v0(i)
    247        u2(i) = Ccoef_U(i,2) + Dcoef_U(i,2)*u1(i)
    248        v2(i) = Ccoef_V(i,2) + Dcoef_V(i,2)*v1(i)
    249     END DO
    250 
    251 !****************************************************************************************
    252 ! On calcule le flux
    253 !
    254 !****************************************************************************************
    255     flux_u(:) = 0.0
    256     flux_v(:) = 0.0
    257 
    258     DO i=1,knon
    259        flux_u(i) = Kcoefm(i,1)/RG/dtime * (u1(i)*alf1(i) + u2(i)*alf2(i) - u0(i))
    260        flux_v(i) = Kcoefm(i,1)/RG/dtime * (v1(i)*alf1(i) + v2(i)*alf2(i) - v0(i))
    261     END DO
    262 
    263   END SUBROUTINE calcul_wind_flux
    264 !
    265 !****************************************************************************************
    266 !
    267   SUBROUTINE climb_wind_up(knon, dtime, u_old, v_old, &
     216
     217  SUBROUTINE climb_wind_up(knon, dtime, u_old, v_old, flx_u1, flx_v1,  &
    268218       flx_u_new, flx_v_new, d_u_new, d_v_new)
    269219!
    270 ! Diffuse the wind components from the surface and up to the top layer. Coefficents
    271 ! C and D are known from before. The values for U and V at surface are supposed to be
    272 ! zero (this could be modified).
    273 !
    274 ! u(k) = Cu(k) + Du(k)*u(k-1)
    275 ! v(k) = Cv(k) + Dv(k)*v(k-1)
    276 ! [1 <= k <= klev]
     220! Diffuse the wind components from the surface layer and up to the top layer.
     221! Coefficents A, B, C and D are known from before. Start values for the diffusion are the
     222! momentum fluxes at surface.
     223!
     224! u(k=1) = A + B*flx*dtime
     225! u(k)   = C(k) + D(k)*u(k-1)  [2 <= k <= klev]
    277226!
    278227!****************************************************************************************
    … …  
    285234    REAL, DIMENSION(klon,klev), INTENT(IN)  :: u_old
    286235    REAL, DIMENSION(klon,klev), INTENT(IN)  :: v_old
     236    REAL, DIMENSION(klon), INTENT(IN)       :: flx_u1, flx_v1 ! momentum flux
    287237
    288238! Output arguments
    … …  
    294244!****************************************************************************************
    295245    REAL, DIMENSION(klon,klev)              :: u_new, v_new
    296     REAL, DIMENSION(klon)                   :: u0, v0
    297246    INTEGER                                 :: k, i
    298247   
    … …  
    300249!****************************************************************************************
    301250
    302 ! Niveau 0
    303     u0(1:knon) = 0.0
    304     v0(1:knon) = 0.0
    305 
    306251! Niveau 1
    307252    DO i = 1, knon
    308        u_new(i,1) = Ccoef_U(i,1) + Dcoef_U(i,1) * u0(i)
    309        v_new(i,1) = Ccoef_V(i,1) + Dcoef_V(i,1) * v0(i)
     253       u_new(i,1) = Acoef_U(i) + Bcoef_U(i)*flx_u1(i)*dtime
     254       v_new(i,1) = Acoef_V(i) + Bcoef_V(i)*flx_v1(i)*dtime
    310255    END DO
    311256
    … …  
    329274    flx_v_new(:,:) = 0.0
    330275
    331 ! Niveau 1
    332     DO i = 1, knon
    333        flx_u_new(i,1) = Kcoefm(i,1)/RG/dtime * &
    334             (u_new(i,1)*alf1(i)+u_new(i,2)*alf2(i) - u0(i))
    335        flx_v_new(i,1) = Kcoefm(i,1)/RG/dtime * &
    336             (v_new(i,1)*alf1(i)+v_new(i,2)*alf2(i) - v0(i))
    337     END DO
     276    flx_u_new(1:knon,1)=flx_u1(1:knon)
     277    flx_v_new(1:knon,1)=flx_v1(1:knon)
    338278
    339279! Niveau 2->klev
Note: See TracChangeset for help on using the changeset viewer.