CALDRG.f

Go to the documentation of this file.

!     
      SUBROUTINE caldrg(DRAGCO)

!     
!
      use vrbls3d, only: uh, vh
      use vrbls2d, only: uz0, vz0, ustar, u10, v10
      use masks, only: lmh
      use params_mod, only: d00, d50, d25
      use ctlblk_mod, only: jsta, jend, jsta_m, jend_m, modelname, spval, im, jm,  &
                            jsta_2l, jend_2u
      use gridspec_mod, only: gridtype
!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
      implicit none
!
!     INCLUDE/SET PARAMETERS.
!     
!     DECLARE VARIABLES.
      REAL,intent(inout) ::  dragco(im,jsta_2l:jend_2u)
      INTEGER ihe(jm),ihw(jm)
      integer i,j,lhmk,ie,iw,lmhk
      real ubar,vbar,wspdsq,ustrsq,sumu,sumv,ulmh,vlmh,uz0h,vz0h
!     
!********************************************************************
!     START CALDRG HERE.
!     
!     INITIALIZE DRAG COEFFICIENT ARRAY TO ZERO.
!     
!$omp parallel do private(i,j)
      DO j=jsta,jend
        DO i=1,im
!          DRAGCO(I,J) = D00
          dragco(i,j) = 0.0 

        ENDDO
      ENDDO
!

      IF(gridtype=='A')THEN 
       DO j=jsta,jend
       DO i=1,im
!     

       IF (ustar(i,j) /= spval) THEN

!        LMHK=NINT(LMH(I,J))
!     
!        COMPUTE A MEAN MASS POINT WIND SPEED BETWEEN THE
!        FIRST ATMOSPHERIC ETA LAYER AND Z0.  ACCORDING TO
!        NETCDF OUTPUT, UZ0 AND VZ0 ARE AT MASS POINTS. (MEB 6/11/02)
!
!        UBAR=D50*(UH(I,J,LMHK)+UZ0(I,J))
!        VBAR=D50*(VH(I,J,LMHK)+VZ0(I,J))
!        WSPDSQ=UBAR*UBAR+VBAR*VBAR

! dong use 10m wind
         wspdsq=u10(i,j)*u10(i,j)+v10(i,j)*v10(i,j)
!     
!        COMPUTE A DRAG COEFFICIENT.
!
        ustrsq=ustar(i,j)*ustar(i,j)
        IF(wspdsq > 1.0) dragco(i,j)=ustrsq/wspdsq

       END IF
       ENDDO
       ENDDO
      ELSE IF(gridtype=='E')THEN
      
       DO j=jsta_m,jend_m
        ihe(j)=mod(j+1,2)
        ihw(j)=ihe(j)-1
       ENDDO
       
       DO j=jsta_m,jend_m
       DO i=2,im-1
!
!        COMPUTE A MEAN MASS POINT WIND IN THE
!        FIRST ATMOSPHERIC ETA LAYER.
!
        lmhk=nint(lmh(i,j))
        ie=i+ihe(j)
        iw=i+ihw(j)
        sumu=uh(ie,j,lmhk)+uh(iw,j,lmhk)+uh(i,j-1,lmhk)          &
          +uh(i,j+1,lmhk)
        sumv=vh(ie,j,lmhk)+vh(iw,j,lmhk)+vh(i,j-1,lmhk)          &
          +vh(i,j+1,lmhk)
        ulmh=d25*sumu
        vlmh=d25*sumv
!
!        COMPUTE A MEAN MASS POINT WIND AT HEIGHT Z0.
!
        uz0h=d25*(uz0(ie,j)+uz0(iw,j)+uz0(i,j-1)+uz0(i,j+1))
        vz0h=d25*(vz0(ie,j)+vz0(iw,j)+vz0(i,j-1)+vz0(i,j+1))
!
!        COMPUTE A MEAN MASS POINT WIND SPEED BETWEEN THE
!        FIRST ATMOSPHERIC ETA LAYER AND Z0.
!
        ubar=d50*(ulmh+uz0h)
        vbar=d50*(vlmh+vz0h)
        wspdsq=ubar*ubar+vbar*vbar
!jjt  WSPDSQ=MIN(WSPDSQ,0.1)
!
!        COMPUTE A DRAG COEFFICIENT.
!
        ustrsq=ustar(i,j)*ustar(i,j)
        IF(wspdsq > 1.0e-6)dragco(i,j)=ustrsq/wspdsq
!
       END DO
       END DO
      ELSE IF(gridtype=='B')THEN 
       DO j=jsta_m,jend_m
       DO i=2,im-1
!
!        COMPUTE A MEAN MASS POINT WIND IN THE
!        FIRST ATMOSPHERIC ETA LAYER.
!
        lmhk=nint(lmh(i,j))
        ie=i
        iw=i-1
        sumu=uh(ie,j,lmhk)+uh(iw,j,lmhk)+uh(i,j-1,lmhk)          &
          +uh(iw,j-1,lmhk)
        sumv=vh(ie,j,lmhk)+vh(iw,j,lmhk)+vh(i,j-1,lmhk)          &
          +vh(iw,j-1,lmhk)
        ulmh=d25*sumu
        vlmh=d25*sumv
!
!        COMPUTE A MEAN MASS POINT WIND AT HEIGHT Z0.
!
! NEMS-NMMB is now putting uz0 and vz0 on mass points to save double interpolation time in
! the model
        IF(modelname == 'NMM')THEN
          uz0h=uz0(i,j)
          vz0h=vz0(i,j)
        ELSE   
          uz0h=d25*(uz0(ie,j)+uz0(iw,j)+uz0(i,j-1)+uz0(iw,j-1))
          vz0h=d25*(vz0(ie,j)+vz0(iw,j)+vz0(i,j-1)+vz0(iw,j-1))
        END IF  
!
!        COMPUTE A MEAN MASS POINT WIND SPEED BETWEEN THE
!        FIRST ATMOSPHERIC ETA LAYER AND Z0.
!
        ubar=d50*(ulmh+uz0h)
        vbar=d50*(vlmh+vz0h)
        wspdsq=ubar*ubar+vbar*vbar
!jjt  WSPDSQ=MIN(WSPDSQ,0.1)
!
!        COMPUTE A DRAG COEFFICIENT.
!
        ustrsq=ustar(i,j)*ustar(i,j)
        IF(wspdsq > 1.0e-6)dragco(i,j)=ustrsq/wspdsq
!
       END DO
       END DO
      ELSE 
      
!$omp parallel do private(i,j)
        DO j=jsta,jend
          DO i=1,im
            dragco(i,j) = spval
          ENDDO
        ENDDO

      END IF 
!     
!     END OF ROUTINE.
!     
      RETURN
      END