da_test_vtoy_transform.inc
References to this file elsewhere.
1 subroutine da_test_vtoy_transform(grid, config_flags, vp, vv, xbx, be, iv, y)
2
3 !---------------------------------------------------------------------------
4 ! Purpose: Perform V To Y Adjoint transform test
5 !
6 ! Method: Perform adjoint test on complete transform: <y,y> = <v_adj,v>.
7 !---------------------------------------------------------------------------
8
9 implicit none
10
11 type(grid_config_rec_type), intent(inout) :: config_flags
12 type(domain), intent(inout) :: grid
13 type(vp_type), intent(inout) :: vv ! Grdipt/EOF CV.
14 type(vp_type), intent(inout) :: vp ! Grdipt/level CV.
15 type(xbx_type), intent(in) :: xbx ! Header & non-gridded vars.
16 type(be_type), intent(in) :: be ! background error structure.
17 type(iv_type), intent(inout) :: iv ! ob. increment vector.
18 type(y_type), intent(inout) :: y ! y = h (xa)
19
20
21 real :: cv(1:cv_size) ! Test control variable.
22 real :: cv_2(1:cv_size)
23 real :: field(ims:ime, jms:jme)
24 real :: adj_sum_lhs ! < y, y >
25 real :: adj_rhs,adj_sum_rhs ! < v, v_adj >
26 real :: partial_lhs ! < y, y >
27 real :: pertile_lhs ! < y, y >
28
29 if (trace_use) call da_trace_entry("da_test_vtoy_transform")
30
31 write(unit=stdout, fmt='(a)') ' da_test_vtoy_transform:'
32 write(unit=stdout, fmt='(a)') '---------------------------------------'
33
34 call random_number(cv(:))
35 cv(:) = cv(:) - 0.5
36
37 cv_2(1:cv_size) = cv(1:cv_size)
38
39 call da_zero_x(grid%xa)
40 call da_zero_vp_type(vp)
41 call da_zero_vp_type(vv)
42
43 call da_transform_vtoy(be, cv, iv, vp, vv, xbx, y, grid, config_flags )
44
45 !-------------------------------------------------------------------------
46 ! [3.0] Calculate LHS of adjoint test equation and
47 ! Rescale input to adjoint routine :
48 !-------------------------------------------------------------------------
49
50 call da_get_y_lhs_value( iv, y, partial_lhs, pertile_lhs, adj_sum_lhs )
51
52 cv = 0.0
53
54 ! call da_zero_vp_type(vp)
55 ! call da_zero_vp_type(vv)
56 ! call da_zero_x(grid%xa)
57
58 call da_transform_vtoy_adj(be, cv, iv, vp, vv, xbx, y, grid, config_flags, .true. )
59
60 adj_rhs = sum( cv(1:cv_size) * cv_2(1:cv_size) )
61 !-------------------------------------------------------------------------
62 ! Print output:
63 !-------------------------------------------------------------------------
64 ! FIX? with wrf_dm_sum_real
65
66 #ifdef DM_PARALLEL
67 if( global ) then
68 adj_sum_rhs = adj_rhs
69 else
70 call mpi_allreduce( adj_rhs, adj_sum_rhs, 1, true_mpi_real, mpi_sum, comm, ierr)
71 end if
72 #else
73 adj_sum_rhs = adj_rhs
74 adj_sum_lhs = partial_lhs
75 #endif
76
77 #ifdef DM_PARALLEL
78 if ( rootproc ) then
79 write(unit=stdout, fmt='(A,1pe22.14)') &
80 'Whole Domain < y, y > = ', adj_sum_lhs
81 write(unit=stdout, fmt='(A,1pe22.14)') &
82 'Whole Domain < v, v_adj > = ', adj_sum_rhs
83 end if
84 #else
85 write(unit=stdout, fmt='(A,1pe22.14)') &
86 'Whole Domain < y, y > = ', adj_sum_lhs
87 write(unit=stdout, fmt='(A,1pe22.14)') &
88 'Whole Domain < v, v_adj > = ', adj_sum_rhs
89 #endif
90
91 if (trace_use) call da_trace_exit("da_test_vtoy_transform")
92
93 end subroutine da_test_vtoy_transform
94
95