da_check_vtox_adjoint.inc
References to this file elsewhere.
1 subroutine da_check_vtox_adjoint(cv_size, xb, xbx, be, ep, cv1, vv, vp, xp, &
2 xa)
3
4 !---------------------------------------------------------------------------
5 ! Purpose: Test V to X routine and adjoint for compatibility.
6 !
7 ! Method: Standard adjoint test: < x, x > = < v_adj, v >.
8 !---------------------------------------------------------------------------
9
10 implicit none
11
12 integer, intent(in) :: cv_size
13 type (xb_type), intent(in) :: xb ! first guess (local).
14 type (xbx_type),intent(in) :: xbx ! For header & non-grid arrays.
15 type (be_type), intent(in) :: be ! background error structure.
16 type (ep_type), intent(in) :: ep ! Ensemble perturbation structure.
17 real, intent(in) :: cv1(1:cv_size) ! control variable (local).
18 type (vp_type), intent(inout) :: vv ! Grdipt/EOF CV.
19 type (vp_type), intent(inout) :: vp ! Grdipt/level CV.
20 type (xpose_type), intent(inout):: xp ! Domain decomposition vars.
21 type (x_type) , intent(inout) :: xa ! gridded analy. incs. (local)
22
23 real :: cv2(1:cv_size) ! control variable (local).
24 real :: adj_par_lhs ! < x, x >
25 real :: adj_par_rhs ! < x, x >
26 real :: adj_sum_lhs ! < x, x >
27 real :: adj_sum_rhs ! < v_adj, v >
28
29 if (trace_use) call da_trace_entry("da_check_vtox_adjoint")
30
31 write(unit=stdout, fmt='(/a/)') &
32 'da_check_vtox_adjoint: Adjoint Test Results:'
33
34 !----------------------------------------------------------------------
35 ! [1.0] Initialise:
36 !----------------------------------------------------------------------
37
38 cv2(:) = 0.0
39
40 !----------------------------------------------------------------------
41 ! [2.0] Perform x = U v transform:
42 !----------------------------------------------------------------------
43
44 call da_zero_x (xa)
45
46 call da_transform_vtox(cv_size, xb, xbx, be, ep, cv1, vv, vp, xp, xa)
47
48 !----------------------------------------------------------------------
49 ! [3.0] Calculate LHS of adjoint test equation:
50 !----------------------------------------------------------------------
51
52 adj_par_lhs = sum(xa % u(its:ite, jts:jte, kts:kte)**2) / typical_u_rms**2 &
53 + sum(xa % v(its:ite, jts:jte, kts:kte)**2) / typical_v_rms**2 &
54 + sum(xa % p(its:ite, jts:jte, kts:kte)**2) / typical_p_rms**2 &
55 + sum(xa % t(its:ite, jts:jte, kts:kte)**2) / typical_t_rms**2 &
56 + sum(xa % q(its:ite, jts:jte, kts:kte)**2) / typical_q_rms**2 &
57 + sum(xa % rho(its:ite,jts:jte,kts:kte)**2)/ typical_rho_rms**2 &
58 + sum(xa % psfc(its:ite, jts:jte)**2) / typical_p_rms**2
59
60 if (cv_options_hum == 2) then
61 adj_par_lhs = adj_par_lhs &
62 + sum(xa % rh(its:ite, jts:jte, kts:kte)**2) / typical_rh_rms**2
63 end if
64
65 if (cv_options_hum == 3) then
66 adj_par_lhs = adj_par_lhs &
67 + sum(xa % qcw(its:ite, jts:jte, kts:kte)**2)/typical_qcw_rms**2 &
68 + sum(xa % qrn(its:ite, jts:jte, kts:kte)**2)/typical_qrn_rms**2
69 adj_par_lhs = adj_par_lhs &
70 + sum(xa % qt (its:ite, jts:jte, kts:kte)**2)/typical_q_rms**2
71 end if
72
73 if (Use_RadarObs) then
74 adj_par_lhs = adj_par_lhs &
75 + sum(xa % wh (its:ite, jts:jte, kts:kte)**2)/typical_w_rms**2
76 else
77 adj_par_lhs = adj_par_lhs &
78 + sum(xa % w (its:ite, jts:jte, kts:kte)**2)/typical_w_rms**2
79 end if
80
81 !-------------------------------------------------------------------------
82 ! [4.0] Rescale input to adjoint routine:
83 !-------------------------------------------------------------------------
84
85 xa % u(:,:,:) = xa % u(:,:,:) / typical_u_rms**2
86 xa % v(:,:,:) = xa % v(:,:,:) / typical_v_rms**2
87 xa % p(:,:,:) = xa % p(:,:,:) / typical_p_rms**2
88 xa % t(:,:,:) = xa % t(:,:,:) / typical_t_rms**2
89 xa % q(:,:,:) = xa % q(:,:,:) / typical_q_rms**2
90 xa % rho(:,:,:) = xa % rho(:,:,:) / typical_rho_rms**2
91
92 xa % psfc(:,:) = xa % psfc(:,:) / typical_p_rms**2
93
94 if (cv_options_hum == 2) then
95 xa % rh(:,:,:) = xa % rh(:,:,:) / typical_rh_rms**2
96 end if
97
98 if (cv_options_hum == 3) then
99 xa % qcw(:,:,:) = xa % qcw(:,:,:) / typical_qcw_rms**2
100 xa % qrn(:,:,:) = xa % qrn(:,:,:) / typical_qrn_rms**2
101 xa % qt (:,:,:) = xa % qt (:,:,:) / typical_q_rms**2
102 end if
103
104 if (Use_RadarObs) then
105 xa %wh(:,:,:) = xa %wh(:,:,:) / typical_w_rms**2
106 xa % w(:,:,:) = 0.0
107 else
108 xa %w (:,:,:) = xa %w (:,:,:) / typical_w_rms**2
109 end if
110
111 !-------------------------------------------------------------------------
112 ! [5.0] Perform adjoint operation:
113 !-------------------------------------------------------------------------
114
115 call da_transform_vtox_adj(cv_size, xb, xbx, be, ep, xa, xp, vp, vv, cv2)
116
117 !-------------------------------------------------------------------------
118 ! [6.0] Calculate RHS of adjoint test equation:
119 !-------------------------------------------------------------------------
120
121 adj_par_rhs = sum(cv1(1:cv_size) * cv2(1:cv_size))
122
123 !-------------------------------------------------------------------------
124 ! [7.0] Print output:
125 !-------------------------------------------------------------------------
126
127 adj_sum_lhs = wrf_dm_sum_real(adj_par_lhs)
128
129 if (global) then
130 adj_sum_rhs = adj_par_rhs
131 else
132 adj_sum_rhs = wrf_dm_sum_real(adj_par_rhs)
133 end if
134 write (unit=stdout,fmt='(A,2F10.2)') &
135 'TEST_COVERAGE_da_check_vtox_adjoint: adj_sum_lhs,adj_sum_rhs = ', &
136 adj_sum_lhs,adj_sum_rhs
137
138 if (rootproc) then
139 write(unit=stdout, fmt='(/)')
140 write(unit=stdout, fmt='(a,1pe22.14)') &
141 'Whole Domain: < x, x > = ', adj_sum_lhs, &
142 'Whole Domain: < v_adj, v > = ', adj_sum_rhs
143 end if
144
145 write(unit=stdout, fmt='(/a/)') 'da_check_vtox_adjoint: Finished'
146
147 if (trace_use) call da_trace_exit("da_check_vtox_adjoint")
148
149 end subroutine da_check_vtox_adjoint
150
151