da_transform_vtoy_adj.inc
References to this file elsewhere.
1 subroutine da_transform_vtoy_adj(iter,cv_size, be, ep, cv, iv, vp, vv, xbx, y, &
2 grid, config_flags)
3
4 !-------------------------------------------------------------------------
5 ! Purpose: Does Adjoint of control variable (V) transform to Obs-space(Y)
6 !-------------------------------------------------------------------------
7
8 implicit none
9
10 integer, intent(in) :: iter
11 integer, intent(in) :: cv_size ! Size of cv array.
12 type(be_type), intent(in) :: be ! background error structure.
13 type(ep_type), intent(in) :: ep ! ensemble perturbation structure.
14 real, intent(out) :: cv(1:cv_size) ! control variables.
15 type(iv_type), intent(inout) :: iv ! innovation vector (o-b).
16 type(vp_type), intent(inout) :: vp ! Grdipt/level CV.
17 type(vp_type), intent(inout) :: vv ! Grdipt/EOF CV.
18 type(xbx_type), intent(in) :: xbx ! For header & non-grid arrays.
19 type(y_type), intent(inout) :: y ! y = H(x_inc).
20 type(domain), intent(inout) :: grid
21 type(grid_config_rec_type), intent(inout) :: config_flags
22 type (xbx_type) :: xbx_tmp
23
24 integer :: nobwin,ndynopt
25 #ifdef DM_PARALLEL
26 integer :: wrf_done_unit
27 #endif
28
29 character(len=4) :: filnam
30
31 call da_trace_entry("da_transform_vtoy_adj")
32
33 if (var4d) then
34 ndynopt = grid%dyn_opt
35 grid%dyn_opt = DYN_EM_TL
36 call nl_set_dyn_opt (1 , DYN_EM_TL)
37
38 if (jcdfi_use .AND. iter > 0) then
39 call da_med_initialdata_input(grid , config_flags, 'tldf')
40 grid%g_u_2 = - grid%g_u_2
41 grid%g_v_2 = - grid%g_v_2
42 grid%g_w_2 = - grid%g_w_2
43 grid%g_t_2 = - grid%g_t_2
44 grid%g_ph_2 = - grid%g_ph_2
45 grid%g_mu_2 = - grid%g_mu_2
46 grid%g_moist = - grid%g_moist
47 else
48 grid%g_u_2 = 0.0
49 grid%g_v_2 = 0.0
50 grid%g_w_2 = 0.0
51 grid%g_t_2 = 0.0
52 grid%g_ph_2 = 0.0
53 grid%g_mu_2 = 0.0
54 grid%g_moist = 0.0
55 end if
56
57 call med_hist_out(grid , 3 , config_flags)
58
59 grid%dyn_opt = ndynopt
60 call nl_set_dyn_opt (1 , DYN_EM)
61
62 do nobwin=num_fgat_time,1,-1
63 write(unit=filnam, fmt='(a, i2.2)') 'fg', nobwin
64 call da_med_initialdata_input( grid , config_flags, &
65 filnam)
66 call da_setup_firstguess( xbx_tmp, grid)
67
68 iv%time = nobwin
69 iv%info(:)%n1 = iv%info(:)%plocal(iv%time-1) + 1
70 iv%info(:)%n2 = iv%info(:)%plocal(iv%time)
71 call da_zero_x(grid%xa)
72 call da_transform_xtoy_adj(grid, iv, y, grid%xa)
73 write(unit=filnam,fmt='(a2,i2.2)') 'af',nobwin
74
75 call da_transfer_wrftltoxa_adj(grid, config_flags, filnam)
76 end do
77
78 call da_trace("da_transform_vtoy_adj","Starting da_run_wrfplus_ad.ksh")
79
80 #ifdef DM_PARALLEL
81 if (var4d_coupling == var4d_coupling_disk_simul) then
82
83 if (rootproc) then
84 call da_system("da_run_wrfplus_ad.ksh pre")
85 call da_system("rm -rf wrf_done")
86 call da_system("touch wrf_go_ahead")
87 call da_get_unit(wrf_done_unit)
88 do while (.true.)
89 open(wrf_done_unit,file="wrf_done",status="old",err=303)
90 close(wrf_done_unit)
91 exit
92 303 continue
93 call da_system("sleep 1")
94 end do
95 call da_free_unit(wrf_done_unit)
96 call da_system("da_run_wrfplus_ad.ksh post")
97 end if
98 ! Wait until PE 0 agrees that AD finished
99 call wrf_get_dm_communicator ( comm )
100 call mpi_barrier(comm, ierr)
101 end if
102 #else
103 call da_system("da_run_wrfplus_ad.ksh")
104 #endif
105 call da_trace("da_transform_vtoy_adj","Finished da_run_wrfplus_ad.ksh")
106
107 call da_transfer_xatowrftl_adj(grid, config_flags, 'gr01')
108 else ! var4d
109 iv%time = 1
110 iv%info(:)%n1 = iv%info(:)%plocal(iv%time-1) + 1
111 iv%info(:)%n2 = iv%info(:)%plocal(iv%time)
112 call da_zero_x(grid%xa)
113 call da_transform_xtoy_adj(grid, iv, y,grid%xa)
114 end if ! var4d
115
116 cv = 0.0
117 call da_transform_vtox_adj(grid, cv_size, xbx, be, ep, vp, vv, cv)
118
119 call da_trace_exit("da_transform_vtoy_adj")
120
121 end subroutine da_transform_vtoy_adj
122
123