da_check.inc
References to this file elsewhere.
1 subroutine da_check(grid, config_flags, cv_size, xbx, be, ep, iv, vv, vp, y)
2
3 !-----------------------------------------------------------------------
4 ! Purpose: TBD
5 !-----------------------------------------------------------------------
6
7 implicit none
8
9 type(domain), intent(inout) :: grid
10 type(grid_config_rec_type), intent(inout) :: config_flags
11
12 integer, intent(in) :: cv_size ! Size of cv array.
13 type (xbx_type), intent(in) :: xbx ! Header & non-gridded vars.
14 type (be_type), intent(in) :: be ! background error structure.
15 type (ep_type), intent(in) :: ep ! Ensemble perturbation structure.
16 type (iv_type), intent(in) :: iv ! ob. increment vector.
17 type (vp_type), intent(inout) :: vv ! Grdipt/EOF CV.
18 type (vp_type), intent(inout) :: vp ! Grdipt/level CV.
19 type (y_type), intent(inout) :: y ! y = h (xa)
20
21 integer :: sizec
22 real :: cvtest(1:cv_size) ! background error structure.
23 real :: field(its:ite,jts:jte) ! Field for spectral transform test.
24
25 call da_trace_entry("da_check")
26
27 !----------------------------------------------------------------------------
28 ! [1] Set up test data:
29 !----------------------------------------------------------------------------
30
31 ! Initialize cv values with random data:
32 call random_number(cvtest(:))
33 cvtest(:) = cvtest(:) - 0.5
34
35 ! vv arrays initialized already.
36 ! vp arrays initialized already.
37
38 !----------------------------------------------------------------------------
39 ! [2] Perform vtox adjoint tests:
40 !----------------------------------------------------------------------------
41
42 call da_message((/"Performing vtox adjoint tests"/))
43
44 ! v_to_vv adjoint test:
45
46 call da_check_cvtovv_adjoint(grid, cv_size, xbx, be, cvtest, vv)
47
48 !-------------------------------------------------------------------------
49 ! vv_to_vp adjoint test:
50 !-------------------------------------------------------------------------
51
52 call da_check_vvtovp_adjoint(be % ne, grid%xb, be, vv, vp)
53
54 !-------------------------------------------------------------------------
55 ! vptox adjoint test:
56 !-------------------------------------------------------------------------
57
58 call da_check_vptox_adjoint(grid,be % ne, be, ep, vp)
59
60 !-------------------------------------------------------------------------
61 ! vtox adjoint test: <x,x> = <v_adj,v>
62 !-------------------------------------------------------------------------
63
64 call da_check_vtox_adjoint(grid, cv_size, xbx, be, ep, cvtest, vv, vp)
65
66 !----------------------------------------------------------------------------
67 ! [2] Perform xtoy adjoint tests:
68 !----------------------------------------------------------------------------
69
70 call da_message((/"Performing xtoy adjoint tests"/))
71
72 call da_allocate_y(iv, y)
73 call da_zero_x(grid%xa)
74
75 call da_setup_testfield(grid)
76
77 ! WHY?
78 ! Make cv_array random.
79
80 ! call random_number(cvtest(1:cv_size))
81 ! cvtest(1:cv_size) = cvtest(1:cv_size) - 0.5
82
83 ! call da_transform_vtox(grid, cv_size, xbx, be, ep, cvtest, vv, vp)
84
85 call da_check_xtoy_adjoint(grid, config_flags, iv, y)
86
87 call da_deallocate_y(y)
88
89 !----------------------------------------------------------------------------
90 ! [4] Perform spectral test:
91 !----------------------------------------------------------------------------
92
93 if (global) then
94
95 call da_message((/"Performing spectral tests"/))
96
97 call random_number(field(:,:))
98 field(:,:) = field(:,:) - 0.5
99
100 sizec = (be % max_wave+1) * (be % max_wave+2)/2
101 call da_test_spectral(be % max_wave, sizec, xbx, field)
102
103 end if
104
105 call da_trace_exit("da_check")
106
107 end subroutine da_check
108
109