var/da/da_transfer_model/da_transfer_xatowrftl_adj_lbc.inc

   1 subroutine da_transfer_xatowrftl_adj_lbc(grid, config_flags, filnam)
   2
   3    !---------------------------------------------------------------------------
   4    ! Purpose: Convert WRFTL variables to analysis increments
   5    !           (inverse of the incremental part of xatowrf)
   6    !---------------------------------------------------------------------------
   7
   8    implicit none
   9
  10
  11    type(domain), intent(inout)               :: grid
  12    type(grid_config_rec_type), intent(inout) :: config_flags
  13
  14    character*4, intent(in) :: filnam
  15
  16 #ifdef VAR4D
  17
  18    ! Local variables
  19
  20    integer :: i, j, k, ii, jj, spec_bdy_width
  21    real    :: sdmd, s1md
  22    real    :: rho_cgrid
  23
  24    real, dimension(ims:ime,jms:jme,kms:kme) :: a_press
  25 #ifndef A2C
  26    real, dimension(ims:ime,jms:jme,kms:kme) :: utmp
  27    real, dimension(ims:ime,jms:jme,kms:kme) :: vtmp
  28 #endif
  29
  30    integer ndynopt
  31
  32    if (trace_use) call da_trace_entry("da_transfer_xatowrftl_adj_lbc")
  33
  34    spec_bdy_width = grid%spec_bdy_width
  35
  36    ALLOCATE ( ubdy3dtemp2(ims:ime,jms:jme,kms:kme) )
  37    ALLOCATE ( vbdy3dtemp2(ims:ime,jms:jme,kms:kme) )
  38    ALLOCATE ( tbdy3dtemp2(ims:ime,jms:jme,kms:kme) )
  39    ALLOCATE ( pbdy3dtemp2(ims:ime,jms:jme,kms:kme) )
  40    ALLOCATE ( qbdy3dtemp2(ims:ime,jms:jme,kms:kme) )
  41    ALLOCATE ( mbdy2dtemp2(ims:ime,1:1,    jms:jme) )
  42 !  ALLOCATE ( wbdy3dtemp2(ims:ime,jms:jme,kms:kme) )
  43
  44    ubdy3dtemp2 = 0.0
  45    vbdy3dtemp2 = 0.0
  46    tbdy3dtemp2 = 0.0
  47    pbdy3dtemp2 = 0.0
  48    qbdy3dtemp2 = 0.0
  49    mbdy2dtemp2 = 0.0
  50 !  wbdy3dtemp2 = 0.0
  51
  52    !---------------------------------------------------------------------------
  53    ! [6.5] Add adjoint from LBC
  54    !---------------------------------------------------------------------------
  55
  56 !  CALL a_stuff_bdytend_new ( wbdy3dtemp2 , REAL(interval_seconds) ,                 &
  57 !                             model_grid%a_w_btxs, model_grid%a_w_btxe,     &
  58 !                             model_grid%a_w_btys, model_grid%a_w_btye,     &
  59 !                             'W' , &
  60 !                             spec_bdy_width      , &
  61 !                             grid%sd31, grid%ed31, grid%sd32, grid%ed32, grid%sd33, grid%ed33, &
  62 !                             grid%sm31, grid%em31, grid%sm32, grid%em32, grid%sm33, grid%em33, &
  63 !                             grid%sp31, grid%ep31, grid%sp32, grid%ep32, grid%sp33, grid%ep33 )
  64
  65    CALL a_stuff_bdytend_new ( ubdy3dtemp2 , REAL(interval_seconds) ,                 &
  66                               model_grid%a_u_btxs, model_grid%a_u_btxe,     &
  67                               model_grid%a_u_btys, model_grid%a_u_btye,     &
  68                               'U' , &
  69                               spec_bdy_width      , &
  70                               grid%sd31, grid%ed31, grid%sd32, grid%ed32, grid%sd33, grid%ed33, &
  71                               grid%sm31, grid%em31, grid%sm32, grid%em32, grid%sm33, grid%em33, &
  72                               grid%sp31, grid%ep31, grid%sp32, grid%ep32, grid%sp33, grid%ep33 )
  73
  74    CALL a_stuff_bdytend_new ( vbdy3dtemp2 , REAL(interval_seconds) ,                 &
  75                               model_grid%a_v_btxs, model_grid%a_v_btxe,     &
  76                               model_grid%a_v_btys, model_grid%a_v_btye,     &
  77                               'V' , &
  78                               spec_bdy_width      , &
  79                               grid%sd31, grid%ed31, grid%sd32, grid%ed32, grid%sd33, grid%ed33, &
  80                               grid%sm31, grid%em31, grid%sm32, grid%em32, grid%sm33, grid%em33, &
  81                               grid%sp31, grid%ep31, grid%sp32, grid%ep32, grid%sp33, grid%ep33 )
  82
  83    CALL a_stuff_bdytend_new ( tbdy3dtemp2 , REAL(interval_seconds) ,                 &
  84                               model_grid%a_t_btxs, model_grid%a_t_btxe,     &
  85                               model_grid%a_t_btys, model_grid%a_t_btye,     &
  86                               'T' , &
  87                               spec_bdy_width      , &
  88                               grid%sd31, grid%ed31, grid%sd32, grid%ed32, grid%sd33, grid%ed33, &
  89                               grid%sm31, grid%em31, grid%sm32, grid%em32, grid%sm33, grid%em33, &
  90                               grid%sp31, grid%ep31, grid%sp32, grid%ep32, grid%sp33, grid%ep33 )
  91
  92    CALL a_stuff_bdytend_new ( pbdy3dtemp2 , REAL(interval_seconds) ,                 &
  93                               model_grid%a_ph_btxs, model_grid%a_ph_btxe,     &
  94                               model_grid%a_ph_btys, model_grid%a_ph_btye,     &
  95                               'W' , &
  96                               spec_bdy_width      , &
  97                               grid%sd31, grid%ed31, grid%sd32, grid%ed32, grid%sd33, grid%ed33, &
  98                               grid%sm31, grid%em31, grid%sm32, grid%em32, grid%sm33, grid%em33, &
  99                               grid%sp31, grid%ep31, grid%sp32, grid%ep32, grid%sp33, grid%ep33 )
 100
 101    CALL a_stuff_bdytend_new ( qbdy3dtemp2 , REAL(interval_seconds) ,                 &
 102                               model_grid%a_moist_btxs(:,:,:,P_A_QV), model_grid%a_moist_btxe(:,:,:,P_A_QV),     &
 103                               model_grid%a_moist_btys(:,:,:,P_A_QV), model_grid%a_moist_btye(:,:,:,P_A_QV),     &
 104                               'T' , &
 105                               spec_bdy_width      , &
 106                               grid%sd31, grid%ed31, grid%sd32, grid%ed32, grid%sd33, grid%ed33, &
 107                               grid%sm31, grid%em31, grid%sm32, grid%em32, grid%sm33, grid%em33, &
 108                               grid%sp31, grid%ep31, grid%sp32, grid%ep32, grid%sp33, grid%ep33 )
 109
 110    CALL a_stuff_bdytend_new ( mbdy2dtemp2 , REAL(interval_seconds) ,                 &
 111                               model_grid%a_mu_btxs, model_grid%a_mu_btxe,     &
 112                               model_grid%a_mu_btys, model_grid%a_mu_btye,     &
 113                               'M' , &
 114                               spec_bdy_width      , &
 115                               grid%sd31, grid%ed31, grid%sd32, grid%ed32, 1, 1, &
 116                               grid%sm31, grid%em31, grid%sm32, grid%em32, 1, 1, &
 117                               grid%sp31, grid%ep31, grid%sp32, grid%ep32, 1, 1 )
 118
 119 #ifdef DM_PARALLEL
 120 #include "HALO_EM_E.inc"
 121 #endif
 122
 123 !  CALL a_couple ( model_config_flags, grid%mu_2 , grid%a_mu_2, grid%mub , wbdy3dtemp2 , grid%w_2, grid%a_w_2, 'h' , grid%msfty , &
 124 !                  grid%sd31, grid%ed31, grid%sd32, grid%ed32, grid%sd33, grid%ed33, &
 125 !                  grid%sm31, grid%em31, grid%sm32, grid%em32, grid%sm33, grid%em33, &
 126 !                  grid%sp31, grid%ep31, grid%sp32, grid%ep32, grid%sp33, grid%ep33 )
 127
 128    CALL a_couple ( model_config_flags, grid%mu_2 , grid%a_mu_2, grid%mub , ubdy3dtemp2 , grid%u_2 , grid%a_u_2 , 'u' , grid%msfuy , &
 129                    grid%sd31, grid%ed31, grid%sd32, grid%ed32, grid%sd33, grid%ed33, &
 130                    grid%sm31, grid%em31, grid%sm32, grid%em32, grid%sm33, grid%em33, &
 131                    grid%sp31, grid%ep31, grid%sp32, grid%ep32, grid%sp33, grid%ep33 )
 132
 133    CALL a_couple ( model_config_flags, grid%mu_2 , grid%a_mu_2, grid%mub , vbdy3dtemp2 , grid%v_2 , grid%a_v_2 , 'v' , grid%msfvx , &
 134                    grid%sd31, grid%ed31, grid%sd32, grid%ed32, grid%sd33, grid%ed33, &
 135                    grid%sm31, grid%em31, grid%sm32, grid%em32, grid%sm33, grid%em33, &
 136                    grid%sp31, grid%ep31, grid%sp32, grid%ep32, grid%sp33, grid%ep33 )
 137
 138    CALL a_couple ( model_config_flags, grid%mu_2 , grid%a_mu_2, grid%mub , tbdy3dtemp2 , grid%t_2 , grid%a_t_2 , 't' , grid%msfty , &
 139                    grid%sd31, grid%ed31, grid%sd32, grid%ed32, grid%sd33, grid%ed33, &
 140                    grid%sm31, grid%em31, grid%sm32, grid%em32, grid%sm33, grid%em33, &
 141                    grid%sp31, grid%ep31, grid%sp32, grid%ep32, grid%sp33, grid%ep33 )
 142
 143    CALL a_couple ( model_config_flags, grid%mu_2 , grid%a_mu_2, grid%mub , pbdy3dtemp2 , grid%ph_2, grid%a_ph_2, 'h' , grid%msfty , &
 144                    grid%sd31, grid%ed31, grid%sd32, grid%ed32, grid%sd33, grid%ed33, &
 145                    grid%sm31, grid%em31, grid%sm32, grid%em32, grid%sm33, grid%em33, &
 146                    grid%sp31, grid%ep31, grid%sp32, grid%ep32, grid%sp33, grid%ep33 )
 147
 148    CALL a_couple ( model_config_flags, grid%mu_2 , grid%a_mu_2, grid%mub , qbdy3dtemp2 , grid%moist(:,:,:,P_A_QV), grid%a_moist(:,:,:,P_A_QV),  't' , grid%msfty , &
 149                    grid%sd31, grid%ed31, grid%sd32, grid%ed32, grid%sd33, grid%ed33, &
 150                    grid%sm31, grid%em31, grid%sm32, grid%em32, grid%sm33, grid%em33, &
 151                    grid%sp31, grid%ep31, grid%sp32, grid%ep32, grid%sp33, grid%ep33 )
 152
 153    DO j = grid%sp32 , MIN(grid%ed32-1,grid%ep32)
 154       DO i = grid%sp31 , MIN(grid%ed31-1,grid%ep31)
 155          grid%a_mu_2(i,j) = grid%a_mu_2(i,j) + mbdy2dtemp2(i,1,j)
 156          mbdy2dtemp2(i,1,j) = 0.0
 157       END DO
 158    END DO
 159
 160 #ifdef DM_PARALLEL
 161    call da_halo_em_e_ad ( grid%a_mu_2 )
 162 #endif
 163
 164    model_grid%a_w_btxs = 0.0
 165    model_grid%a_w_btxe = 0.0
 166    model_grid%a_w_btys = 0.0
 167    model_grid%a_w_btye = 0.0
 168    model_grid%a_w_bxs = 0.0
 169    model_grid%a_w_bxe = 0.0
 170    model_grid%a_w_bys = 0.0
 171    model_grid%a_w_bye = 0.0
 172
 173    DEALLOCATE ( ubdy3dtemp2 )
 174    DEALLOCATE ( vbdy3dtemp2 )
 175    DEALLOCATE ( tbdy3dtemp2 )
 176    DEALLOCATE ( pbdy3dtemp2 )
 177    DEALLOCATE ( qbdy3dtemp2 )
 178    DEALLOCATE ( mbdy2dtemp2 )
 179 !  DEALLOCATE ( wbdy3dtemp2 )
 180
 181 #ifndef A2C
 182    !---------------------------------------------------------------------------
 183    ! [5.0] Adjoint of CONVERT FROM A-GRID TO C-GRID
 184    !---------------------------------------------------------------------------
 185
 186    ! Fill the halo region for a_u and a_v.
 187    utmp=grid%x6a%u
 188    vtmp=grid%x6a%v
 189 #endif
 190    grid%x6a%u=grid%a_u_2
 191    grid%x6a%v=grid%a_v_2
 192
 193 #ifdef A2C
 194   if ((fg_format==fg_format_wrf_arw_regional  .or. &
 195        fg_format==fg_format_wrf_arw_global  ) .and. ide == ipe ) then
 196      ipe = ipe + 1
 197      ide = ide + 1
 198   end if
 199
 200   if ((fg_format==fg_format_wrf_arw_regional  .or. &
 201        fg_format==fg_format_wrf_arw_global  ) .and. jde == jpe ) then
 202      jpe = jpe + 1
 203      jde = jde + 1
 204   end if
 205 #endif
 206
 207 #ifdef DM_PARALLEL
 208 #include "HALO_X6A_A.inc"
 209 #endif
 210
 211 #ifdef A2C
 212   if ((fg_format==fg_format_wrf_arw_regional  .or. &
 213        fg_format==fg_format_wrf_arw_global  ) .and. ide == ipe ) then
 214      ipe = ipe - 1
 215      ide = ide - 1
 216   end if
 217
 218   if ((fg_format==fg_format_wrf_arw_regional  .or. &
 219        fg_format==fg_format_wrf_arw_global  ) .and. jde == jpe ) then
 220      jpe = jpe - 1
 221      jde = jde - 1
 222   end if
 223 #else
 224    grid%a_u_2=grid%x6a%u
 225    grid%a_v_2=grid%x6a%v
 226    grid%x6a%u=utmp
 227    grid%x6a%v=vtmp
 228    utmp=0.0
 229    vtmp=0.0
 230
 231    do k=kts,kte
 232       do j=jts,jte
 233          do i=its,ite
 234             utmp(i,j,k)=utmp(i,j,k)+0.5*(grid%a_u_2(i+1,j  ,k)+grid%a_u_2(i,j,k))
 235             vtmp(i,j,k)=vtmp(i,j,k)+0.5*(grid%a_v_2(i  ,j+1,k)+grid%a_v_2(i,j,k))
 236          end do
 237       end do
 238    end do
 239
 240    utmp(its-1,jts:jte,kts:kte)=utmp(its-1,jts:jte,kts:kte)+0.5*grid%a_u_2(its,jts:jte,kts:kte)
 241    utmp(ite+1,jts:jte,kts:kte)=utmp(ite+1,jts:jte,kts:kte)+0.5*grid%a_u_2(ite+1,jts:jte,kts:kte)
 242    vtmp(its:ite,jts-1,kts:kte)=vtmp(its:ite,jts-1,kts:kte)+0.5*grid%a_v_2(its:ite,jts,kts:kte)
 243    vtmp(its:ite,jte+1,kts:kte)=vtmp(its:ite,jte+1,kts:kte)+0.5*grid%a_v_2(its:ite,jte+1,kts:kte)
 244
 245    ! The western boundary
 246    if (its == grid%xp%ids) then
 247       grid%x6a%u(its  ,jts:jte,kts:kte)=grid%x6a%u(its  ,jts:jte,kts:kte)+2.0*utmp(its-1,jts:jte,kts:kte)
 248       grid%x6a%u(its+1,jts:jte,kts:kte)=grid%x6a%u(its+1,jts:jte,kts:kte)-utmp(its-1,jts:jte,kts:kte)
 249    end if
 250
 251    ! The eastern boundary
 252    if (ite == grid%xp%ide) then
 253       grid%x6a%u(ite  ,jts:jte,kts:kte)=grid%x6a%u(ite  ,jts:jte,kts:kte)+2.0*utmp(ite+1,jts:jte,kts:kte)
 254       grid%x6a%u(ite-1,jts:jte,kts:kte)=grid%x6a%u(ite-1,jts:jte,kts:kte)-utmp(ite+1,jts:jte,kts:kte)
 255    end if
 256
 257    grid%x6a%u=grid%x6a%u+utmp
 258
 259    ! The southern boundary
 260    if (jts == grid%xp%jds) then
 261       grid%x6a%v(its:ite,jts  ,kts:kte)=grid%x6a%v(its:ite,jts  ,kts:kte)+2.0*vtmp(its:ite,jts-1,kts:kte)
 262       grid%x6a%v(its:ite,jts+1,kts:kte)=grid%x6a%v(its:ite,jts+1,kts:kte)-vtmp(its:ite,jts-1,kts:kte)
 263    end if
 264
 265    ! The northern boundary
 266    if (jte == grid%xp%jde) then
 267       grid%x6a%v(its:ite,jte  ,kts:kte)=grid%x6a%v(its:ite,jte  ,kts:kte)+2.0*vtmp(its:ite,jte+1,kts:kte)
 268       grid%x6a%v(its:ite,jte-1,kts:kte)=grid%x6a%v(its:ite,jte-1,kts:kte)-vtmp(its:ite,jte+1,kts:kte)
 269    end if
 270
 271    grid%x6a%v=grid%x6a%v+vtmp
 272 #endif
 273
 274    grid%a_u_2 = 0.0
 275    grid%a_v_2 = 0.0
 276
 277    !---------------------------------------------------------------------------
 278    ! [4.0] Adjoint of CONVERT TEMPERATURE inCREMENTS inTO THETA inCREMENTS
 279    !       EVALUATE ALSO THE inCREMENTS OF (1/rho) AND GEOPOTENTIAL
 280    !---------------------------------------------------------------------------
 281
 282    a_press(its:ite,jts:jte,kts:kte+1)=0.0
 283 !  do k=kte,kts,-1
 284 !     do j=jts,jte
 285 !        do i=its,ite
 286 !           rho_cgrid=-(grid%ph_2(i,j,k+1)-grid%ph_2(i,j,k))*grid%a_ph_2(i,j,k+1)/grid%xb%rho(i,j,k)
 287 !           a_press(i,j,k )=a_press(i,j,k )+grid%a_ph_2(i,j,k+1)/grid%xb%rho(i,j,k)
 288 !           a_press(i,j,k+1)=a_press(i,j,k+1)-grid%a_ph_2(i,j,k+1)/grid%xb%rho(i,j,k)
 289 !           grid%a_ph_2(i,j,k ) =grid%a_ph_2(i,j,k)   +grid%a_ph_2(i,j,k+1)
 290 !           grid%x6a%q(i,j,k)=grid%x6a%q(i,j,k)-grid%xb%rho(i,j,k)*0.61*rho_cgrid/(1.0+0.61*grid%xb%q(i,j,k))
 291 !           grid%x6a%t(i,j,k)=grid%x6a%t(i,j,k)-grid%xb%rho(i,j,k)*rho_cgrid/grid%xb%t(i,j,k)
 292 !           grid%x6a%p(i,j,k)= grid%x6a%p(i,j,k)+grid%xb%rho(i,j,k)*rho_cgrid/grid%xb%p(i,j,k)
 293 !        end do
 294 !     end do
 295 !  end do
 296
 297    do k=kts,kte
 298       do j=jts,jte
 299          do i=its,ite
 300             grid%x6a%p(i,j,k)=grid%x6a%p(i,j,k)-(t0+grid%t_2(i,j,k))*kappa*grid%a_t_2(i,j,k)/grid%xb%p(i,j,k)
 301             grid%x6a%t(i,j,k)=grid%x6a%t(i,j,k)+(t0+grid%t_2(i,j,k))*grid%a_t_2(i,j,k)/grid%xb%t(i,j,k)
 302          end do
 303       end do
 304    end do
 305
 306    grid%a_t_2 = 0.0
 307    grid%a_ph_2 = 0.0
 308
 309    !---------------------------------------------------------------------------
 310    ! [3.0] Adjoint of COMPUTE pressure increments (for computing theta increments)
 311    !---------------------------------------------------------------------------
 312
 313    do k=kts,kte
 314       do j=jts,jte
 315          do i=its,ite
 316             a_press(i,j,k+1)=a_press(i,j,k+1)+0.5*grid%x6a%p(i,j,k)
 317             a_press(i,j,k )=a_press(i,j,k )+0.5*grid%x6a%p(i,j,k)
 318             grid%x6a%p(i,j,k)=0.0
 319             grid%a_moist(i,j,k,P_A_QV)=grid%a_moist(i,j,k,P_A_QV)-(grid%mu_2(i,j)+grid%mub(i,j))*a_press(i,j,k)*grid%dn(k)
 320             grid%a_mu_2(i,j)=grid%a_mu_2(i,j)-a_press(i,j,k)*(1.0+grid%moist(i,j,k,P_QV))*grid%dn(k)
 321             a_press(i,j,k+1)=a_press(i,j,k+1)+a_press(i,j,k)
 322          end do
 323       end do
 324    end do
 325
 326    !---------------------------------------------------------------------------
 327    ! [2.0] Adjoint of COMPUTE increments of dry-column air mass per unit area
 328    !---------------------------------------------------------------------------
 329
 330    do j=jts,jte
 331       do i=its,ite
 332          sdmd=0.0
 333          s1md=0.0
 334          do k=kts,kte
 335             s1md=s1md+(1.0+grid%moist(i,j,k,P_QV))*grid%dnw(k)
 336          end do
 337          sdmd=sdmd-grid%xb%psac(i,j)*grid%a_mu_2(i,j)/s1md
 338          grid%x6a%psfc(i,j)=grid%x6a%psfc(i,j)-grid%a_mu_2(i,j)/s1md
 339          do k=kts,kte
 340             grid%a_moist(i,j,k,P_A_QV)=grid%a_moist(i,j,k,P_A_QV)+sdmd*grid%dnw(k)
 341          end do
 342       end do
 343    end do
 344
 345    grid%a_mu_2 = 0.0
 346    !---------------------------------------------------------------------------
 347    ! [1.0] Adjoint of Get the mixing ratio of moisture
 348    !---------------------------------------------------------------------------
 349    do k=kts,kte
 350       do j=jts,jte
 351          do i=its,ite
 352             grid%x6a%q(i,j,k)=grid%x6a%q(i,j,k)+grid%a_moist(i,j,k,P_A_QV)/(1.0-grid%xb%q(i,j,k))**2
 353          end do
 354       end do
 355    end do
 356
 357    grid%a_moist = 0.0
 358
 359    do k=kts,kte+1
 360       do j=jts,jte
 361          do i=its,ite
 362             grid%x6a%w(i,j,k)=grid%x6a%w(i,j,k)+grid%a_w_2(i,j,k)
 363          end do
 364       end do
 365    end do
 366
 367    grid%a_w_2 = 0.0
 368
 369    if (trace_use) call da_trace_exit("da_transfer_xatowrftl_adj_lbc")
 370
 371 #endif
 372 end subroutine da_transfer_xatowrftl_adj_lbc
 373
 374