da_balance_equation_adj.inc
References to this file elsewhere.
1 subroutine da_balance_equation_adj( xb, xbx, xp, phi_b, u, v)
2
3 !---------------------------------------------------------------------------
4 ! Purpose: Adjoint of da_balance_equation
5 !---------------------------------------------------------------------------
6
7 implicit none
8
9 type (xb_type), intent(in) :: xb ! First guess structure.
10 type (xbx_type),intent(in) :: xbx ! Header & non-gridded vars.
11 type (xpose_type),intent(inout) :: xp ! Dimensions and xpose buffers.
12 real, intent(in) :: phi_b(ims:ime,jms:jme,kms:kme) ! Balanced mass increment.
13 real, intent(inout) :: u(ims:ime,jms:jme,kms:kme) ! u wind comp. (dot pts)
14 real, intent(inout) :: v(ims:ime,jms:jme,kms:kme) ! v wind comp. (dot pts)
15
16 integer :: i, j, k ! Loop counters.
17 integer :: is, ie ! 1st dim. end points.
18 integer :: js, je ! 2nd dim. end points.
19
20 real, dimension(ims:ime,jms:jme) :: coefx, & ! Multiplicative coefficient.
21 coefy, & ! Multiplicative coefficient.
22 term_x, & ! Balance eqn x term
23 term_y ! Balance eqn y term
24
25 real, dimension(ims:ime,jms:jme,kms:kme) :: del2phi_b ! Del**2 phi_b/M**2
26 real :: coeff1, coeff2 ! Multiplicative coefficient.
27
28 !---------------------------------------------------------------------------
29 ! [1.0] Initialise:
30 !---------------------------------------------------------------------------
31
32 ! Computation to check for edge of domain:
33 is = its-1; ie = ite+1; js = jts-1; je = jte+1
34 if (.not. global .and. its == ids ) is = ids+1
35 if (.not. global .and. ite == ide ) ie = ide-1
36 if (jts == jds ) js = jds+1
37 if (jte == jde ) je = jde-1
38
39 if( fg_format == fg_format_kma_global) then
40 coefx(is:ie,js:je) = xb%coefx(is:ie,js:je)
41 coefy(is:ie,js:je) = xb%coefy(is:ie,js:je)
42 else if( fg_format == fg_format_wrf) then
43 coefx(is:ie,js:je) = xb%coefz(is:ie,js:je)
44 coefy(is:ie,js:je) = coefx(is:ie,js:je)
45 else
46 write (unit=message(1),fmt='(A,I3)') ' Wrong FG_FORMAT = ',fg_format
47 call da_error(__FILE__,__LINE__,message(1:1))
48 end if
49
50 ! [1.1] Multiplicative coefficent for conversion RHS->Del**2 phi_b/M**2:
51
52 del2phi_b(:,:,:) = 0.0
53
54 !---------------------------------------------------------------------------
55 ! [3.0] Solve Del**2 phi_b = RHS for phi_b:
56 !---------------------------------------------------------------------------
57
58 call da_solve_poissoneqn_fst_adj(xbx, xp, phi_b, &
59 del2phi_b)
60
61 !---------------------------------------------------------------------------
62 ! [2.0] Calculate RHS of balance equation in gridpt space:
63 !---------------------------------------------------------------------------
64
65 do k = kts, kte
66 ! [2.4] Del**2 Phi_b boundary conditions (null as zero boundary conditions):
67
68 ! [2.3] Take divergence to get Del**2 phi_b/M**2:
69
70 term_x(ims:ime,jms:jme) = 0.0
71 term_y(ims:ime,jms:jme) = 0.0
72
73 do j = je, js, -1
74 do i = ie, is, -1
75 coeff1 = coefx(i,j) * del2phi_b(i,j,k)
76 coeff2 = coefy(i,j) * del2phi_b(i,j,k)
77
78 term_x(i+1,j) = term_x(i+1,j) - coeff1
79 term_x(i-1,j) = term_x(i-1,j) + coeff1
80 term_y(i,j+1) = term_y(i,j+1) - coeff2
81 term_y(i,j-1) = term_y(i,j-1) + coeff2
82 end do
83 end do
84
85 ! [2.2] Include cyclostrophic terms in balance eqn if requested:
86
87 if (balance_type == balance_cyc .OR. &
88 balance_type == balance_geocyc ) then
89 call da_balance_cycloterm_adj( xb % rho(:,:,k),xb % u(:,:,k),&
90 xb % v(:,:,k), u(:,:,k), v(:,:,k),&
91 xb % coefx(:,:), xb % coefy(:,:),&
92 term_x(:,:), term_y(:,:))
93 end if
94
95
96 ! [2.1] Calculate geostrophic terms in balance eqn:
97
98 if (balance_type == balance_geo .OR. &
99 balance_type == balance_geocyc ) then
100 call da_balance_geoterm_adj( xb % cori, xb % rho(:,:,k), &
101 term_x, term_y, &
102 u(:,:,k), v(:,:,k))
103 end if
104 end do
105
106 end subroutine da_balance_equation_adj
107
108