var/da/da_test/da_check_vchemtox_adjoint.inc

   1 subroutine da_check_vchemtox_adjoint(grid, vchem, be, vchem2)
   2
   3    !---------------------------------------------------------------------------
   4    ! Purpose: Test V to X routine and adjoint for compatibility.
   5    !
   6    ! Method:  Standard adjoint test: < x, x > = < v_adj, v >.
   7    !---------------------------------------------------------------------------
   8
   9    use module_state_description, only : num_chem, PARAM_FIRST_SCALAR
  10    use da_define_structures, only : da_zero_xchem_type
  11    use module_domain, only : xchem_type
  12
  13    implicit none
  14
  15    type(domain),  intent(inout) :: grid
  16    type(xchem_type), intent(in) :: vchem   ! Grid point/EOF equivalent.
  17    type(be_type), intent(in)    :: be   ! Background error structure.
  18
  19    real                   :: adj_par_lhs ! < x, x >
  20    real                   :: adj_par_rhs ! < x, x >
  21    real                   :: adj_sum_lhs ! < x, x >
  22    real                   :: adj_sum_rhs ! < v_adj, v >
  23    integer                :: ic
  24    type(xchem_type), intent(out) :: vchem2   ! Grid point/EOF equivalent.
  25
  26    if (trace_use) call da_trace_entry("da_check_vchemtox_adjoint")
  27
  28    write(unit=stdout, fmt='(/a/)') &
  29       'da_check_vchemtox_adjoint: Adjoint Test Results:'
  30
  31    !----------------------------------------------------------------------
  32    ! [1.0] Initialise:
  33    !----------------------------------------------------------------------
  34    vchem2%chem_ic(:,:,:,:) = 0.D0
  35    adj_par_rhs = 0.0
  36    adj_par_lhs = 0.0
  37    !----------------------------------------------------------------------
  38    ! [2.0] Perform x = U v transform:
  39    !----------------------------------------------------------------------
  40
  41    call da_zero_xchem_type (grid%xachem)
  42    call da_transform_vchemtox(grid, vchem, be)
  43
  44    !----------------------------------------------------------------------
  45    ! [3.0] Calculate LHS of adjoint test equation:
  46    !----------------------------------------------------------------------
  47
  48    do ic = PARAM_FIRST_SCALAR ,num_chem
  49      adj_par_lhs = adj_par_lhs + sum(grid%xachem%chem_ic(its:ite,jts:jte,kts:kte,ic)**2) / typical_t_rms**2
  50    end do
  51
  52    !-------------------------------------------------------------------------
  53    ! [4.0] Rescale input to adjoint routine:
  54    !-------------------------------------------------------------------------
  55
  56     do ic = PARAM_FIRST_SCALAR ,num_chem
  57 !      grid%xachem%chem_ic(:,:,:,ic) = grid%xachem%chem_ic(:,:,:,ic) / typical_t_rms**2
  58       grid%xachem%chem_ic(its:ite,jts:jte,kts:kte,ic) = grid%xachem%chem_ic(its:ite,jts:jte,kts:kte,ic) / typical_t_rms**2
  59     end do
  60
  61    !-------------------------------------------------------------------------
  62    ! [5.0] Perform adjoint operation:
  63    !-------------------------------------------------------------------------
  64
  65    call da_transform_vchemtox_adj(grid, vchem2, be)
  66
  67    !-------------------------------------------------------------------------
  68    ! [6.0] Calculate RHS of adjoint test equation:
  69    !-------------------------------------------------------------------------
  70
  71     do ic = PARAM_FIRST_SCALAR ,num_chem
  72       adj_par_rhs = adj_par_rhs+sum(vchem%chem_ic(its:ite,jts:jte,1:be%v12(ic-1)%mz,ic) * vchem2%chem_ic(its:ite,jts:jte,1:be%v12(ic-1)%mz,ic))
  73     end do
  74    !-------------------------------------------------------------------------
  75    ! [7.0] Print output:
  76    !-------------------------------------------------------------------------
  77
  78    if (.not. global ) then
  79     if( num_procs == 1) then
  80       write(unit=stdout, fmt='(/)')
  81       write(unit=stdout, fmt='(a,1pe22.14)') &
  82          'Single Domain: < x, x >     = ', adj_par_lhs, &
  83          'Single Domain: < v_adj, v > = ', adj_par_rhs
  84     else
  85       write(unit=stdout, fmt='(/a/,a/)')&
  86         'It is Multi Processor Run: ',&
  87         'For Single Domain: da_check_vchemtox_adjoint Test: Not Performed'
  88     endif
  89    end if
  90
  91    adj_sum_lhs = wrf_dm_sum_real(adj_par_lhs)
  92
  93    if (global) then
  94       adj_sum_rhs = adj_par_rhs
  95    else
  96       adj_sum_rhs = wrf_dm_sum_real(adj_par_rhs)
  97    end if
  98
  99    if (rootproc) then
 100       write(unit=stdout, fmt='(/)')
 101       write(unit=stdout, fmt='(a,1pe22.14)') &
 102          'Whole  Domain: < x, x >     = ', adj_sum_lhs, &
 103          'Whole  Domain: < vchem_adj, vchem > = ', adj_sum_rhs
 104    end if
 105
 106    write(unit=stdout, fmt='(/a/)') 'da_check_vchemtox_adjoint: Finished'
 107
 108    if (trace_use) call da_trace_exit("da_check_vchemtox_adjoint")
 109
 110 end subroutine da_check_vchemtox_adjoint
 111
 112