subroutine da_balance_cycloterm_adj (rho, ub, vb, u, v, coefx, coefy, term_x, term_y) !--------------------------------------------------------------------------- ! Purpose: Adjoint of da_balance_cycloterm !--------------------------------------------------------------------------- implicit none real, intent(in) :: rho(ims:ime,jms:jme) ! Density. real, intent(in) :: ub(ims:ime,jms:jme) ! Background u wind real, intent(in) :: vb(ims:ime,jms:jme) ! Background u wind real, intent(in) :: term_x(ims:ime,jms:jme) ! x component of term real, intent(in) :: term_y(ims:ime,jms:jme) ! y component of term real, intent(in) :: coefx(ims:ime,jms:jme) real, intent(in) :: coefy(ims:ime,jms:jme) ! Mulplicative coeff. real, intent(inout) :: u(ims:ime,jms:jme) ! u wind increment real, intent(inout) :: v(ims:ime,jms:jme) ! v wind increment integer :: i, j ! Loop counters. integer :: is, ie ! 1st dim. end points. integer :: js, je ! 2nd dim. end points. real :: data(ims:ime,jms:jme) ! Work array. real :: var, varb, uar if (trace_use) call da_trace_entry("da_balance_cycloterm_adj") !--------------------------------------------------------------------------- ! [1.0] Initialise: !--------------------------------------------------------------------------- ! Computation to check for edge of domain: is = its; ie = ite; js = jts; je = jte if (.not. global .and. its == ids) is = ids+1 if (.not. global .and. ite == ide) ie = ide-1 if (jts == jds) js = jds+1 if (jte == jde) je = jde-1 !--------------------------------------------------------------------------- ! [3.0] Calculate term_y = rho M ( u'dv/dx + v'dv/dy + udv'/dx + vdv'/dy ): !--------------------------------------------------------------------------- ! [3.7] Multiply by rho and add to term_y data(its:ite,jts:jte) = rho(its:ite,jts:jte) * term_y(its:ite,jts:jte) if (.NOT. global) then ! [3.6] Corner points: if (its == ids .AND. jts == jds ) then data(its,jts+1) = data(its,jts+1) + 0.5 * data(its,jts) data(its+1,jts) = data(its+1,jts) + 0.5 * data(its,jts) end if if (ite == ide .AND. jts == jds ) then data(ite-1,jts) = data(ite-1,jts) + 0.5 * data(ite,jts) data(ite,jts+1) = data(ite,jts+1) + 0.5 * data(ite,jts) end if if (its == ids .AND. jte == jde ) then data(its,jde-1) = data(its,jde-1) + 0.5 * data(its,jde) data(its+1,jde) = data(its+1,jde) + 0.5 * data(its,jde) end if if (ite == ide .AND. jte == jde ) then data(ite-1,jte) = data(ite-1,jte) + 0.5 * data(ite,jte) data(ite,jte-1) = data(ite,jte-1) + 0.5 * data(ite,jte) end if ! [3.5] Right boundaries: if (jte == jde ) then j = jte do i = is, ie varb = 3.0*vb(i,j)-4.0*vb(i,j-1)+vb(i,j-2) var = coefy(i,j)* vb(i,j) * data(i,j) uar = coefx(i,j)* data(i,j) * ub(i,j) u(i,j) = u(i,j) + coefx(i,j)*data(i,j) * ( vb(i+1,j) - vb(i-1,j) ) v(i,j) = v(i,j) + coefy(i,j)*data(i,j) * varb v(i+1,j) = v(i+1,j) + uar v(i-1,j) = v(i-1,j) - uar v(i,j) = v(i,j) + 3.0*var v(i,j-1) = v(i,j-1) -4.0*var v(i,j-2) = v(i,j-2) + var end do end if ! [3.4] Left boundaries: if (jts == jds ) then j = jts do i = is, ie varb = -3.0*vb(i,j)+4.0*vb(i,j+1)-vb(i,j+2) var = coefy(i,j)*vb(i,j) * data(i,j) uar = coefx(i,j)*ub(i,j) * data(i,j) v(i,j) = v(i,j) + coefy(i,j)*data(i,j) * varb v(i+1,j) = v(i+1,j) + uar v(i-1,j) = v(i-1,j) - uar v(i,j) = v(i,j) - 3.0*var v(i,j+1) = v(i,j+1) +4.0*var v(i,j+2) = v(i,j+2) - var end do end if ! [3.3] Top boundaries: if (ite == ide ) then i = ite do j = js, je varb = 3.0*vb(i,j)-4.0*vb(i-1,j)+vb(i-2,j) var = coefx(i,j)* ub(i,j) * data(i,j) uar = coefy(i,j)* vb(i,j) * data(i,j) u(i,j) = u(i,j) + coefx(i,j)*data(i,j) * varb v(i,j) = v(i,j) + coefy(i,j)*data(i,j) * ( vb(i,j+1) - vb(i,j-1) ) v(i,j+1) = v(i,j+1) + uar v(i,j-1) = v(i,j-1) - uar v(i,j) = v(i,j) + 3.0*var v(i-1,j) = v(i-1,j) -4.0**var v(i-2,j) = v(i-2,j) + var end do end if ! [3.2] Bottom boundaries: if (its == ids ) then i = its do j = js, je varb = -3.0*vb(i,j)+4.0*vb(i+1,j)-vb(i+2,j) var = coefx(i,j)* ub(i,j) * data(i,j) uar = coefy(i,j)* vb(i,j) * data(i,j) u(i,j) = u(i,j) + coefx(i,j)*data(i,j) * varb v(i,j) = v(i,j) + coefy(i,j)*data(i,j) * ( vb(i,j+1) - vb(i,j-1) ) v(i,j+1) = v(i,j+1) + uar v(i,j-1) = v(i,j-1) - uar v(i,j) = v(i,j) - 3.0*var v(i+1,j) = v(i+1,j) +4.0**var v(i+2,j) = v(i+2,j) - var end do end if end if ! not global ! [3.1] Interior points: do j = je, js, -1 do i = ie, is, -1 uar = coefx(i,j) * ub(i,j) * data(i,j) var = coefy(i,j) * vb(i,j) * data(i,j) u(i,j) = u(i,j) + coefx(i,j)*data(i,j)*( vb(i+1,j) - vb(i-1,j) ) v(i,j) = v(i,j) + coefy(i,j)*data(i,j)*( vb(i,j+1) - vb(i,j-1) ) v(i+1,j) = v(i+1,j) + uar v(i-1,j) = v(i-1,j) - uar v(i,j+1) = v(i,j+1) + var v(i,j-1) = v(i,j-1) - var end do end do !--------------------------------------------------------------------------- ! [2.0] Calculate term_x = rho M ( u'du/dx + v'du/dy + udu'/dx + vdu'/dy ): !--------------------------------------------------------------------------- ! [2.7] Multiply by rho and add to term_x: data(its:ite,jts:jte) = rho(its:ite,jts:jte) * term_x(its:ite,jts:jte) if( .NOT. global) then ! [2.6] Corner points: if (its == ids .AND. jts == jds ) then data(its,jts+1) = data(its,jts+1) + 0.5 * data(its,jts) data(its+1,jts) = data(its+1,jts) + 0.5 * data(its,jts) end if if (ite == ide .AND. jts == jds ) then data(ite-1,jts) = data(ite-1,jts) + 0.5 * data(ite,jts) data(ite,jts+1) = data(ite,jts+1) + 0.5 * data(ite,jts) end if if (its == ids .AND. jte == jde ) then data(its,jde-1) = data(its,jde-1) + 0.5 * data(its,jde) data(its+1,jde) = data(its+1,jde) + 0.5 * data(its,jde) end if if (ite == ide .AND. jte == jde ) then data(ite-1,jte) = data(ite-1,jte) + 0.5 * data(ite,jte) data(ite,jte-1) = data(ite,jte-1) + 0.5 * data(ite,jte) end if ! [2.5] Right boundaries: if (jte == jde ) then j = jte do i = is, ie varb = 3.0*ub(i,j)-4.0*ub(i,j-1)+ub(i,j-2) var = coefy(i,j) * vb(i,j) * data(i,j) uar = coefx(i,j) * ub(i,j) * data(i,j) u(i+1,j) = u(i+1,j) + uar u(i-1,j) = u(i-1,j) - uar v(i,j) = v(i,j) + coefy(i,j)*data(i,j) * varb u(i,j) = u(i,j) + coefx(i,j)*data(i,j) * ( ub(i+1,j) - ub(i-1,j) ) u(i,j) = u(i,j) + 3.0*var u(i,j-1) = u(i,j-1) -4.0*var u(i,j-2) = u(i,j-2) + var end do end if ! [2.4] Left boundaries: if (jts == jds ) then j = jts do i = is, ie varb = -3.0*ub(i,j)+4.0*ub(i,j+1)-ub(i,j+2) var = coefy(i,j)*vb(i,j) * data(i,j) uar = coefx(i,j)*ub(i,j) * data(i,j) u(i+1,j) = u(i+1,j) + uar u(i-1,j) = u(i-1,j) - uar v(i,j) = v(i,j) + coefy(i,j)*data(i,j) * varb u(i,j) = u(i,j) + coefx(i,j)*data(i,j) * ( ub(i+1,j) - ub(i-1,j) ) u(i,j) = u(i,j) - 3.0*var u(i,j+1) = u(i,j+1) +4.0*var u(i,j+2) = u(i,j+2) - var end do end if ! [2.3] Top boundaries: if (ite == ide ) then i = ite do j = js, je varb = 3.0*ub(i,j)-4.0*ub(i-1,j)+ub(i-2,j) var = coefx(i,j)*ub(i,j) * data(i,j) uar = coefy(i,j)*vb(i,j) * data(i,j) u(i,j+1) = u(i,j+1) + uar u(i,j-1) = u(i,j-1) - uar v(i,j) = v(i,j) + coefy(i,j)*data(i,j) * ( ub(i,j+1) - ub(i,j-1) ) u(i,j) = u(i,j) + coefx(i,j)*data(i,j) * varb u(i,j) = u(i,j) + 3.0*var u(i-1,j) = u(i-1,j) - 4.0*var u(i-2,j) = u(i-2,j) + var end do end if ! [2.2] Bottom boundaries: if (its == ids ) then i = its do j = js, je varb = -3.0*ub(i,j)+4.0*ub(i+1,j)-ub(i+2,j) var = coefy(i,j)*ub(i,j) * data(i,j) uar = coefy(i,j)*vb(i,j) * data(i,j) u(i,j+1) = u(i,j+1) + uar u(i,j-1) = u(i,j-1) - uar v(i,j) = v(i,j) + coefy(i,j)*data(i,j) * ( ub(i,j+1) - ub(i,j-1) ) u(i,j) = u(i,j) + coefx(i,j)*data(i,j) * varb u(i,j) = u(i,j) - 3.0*var u(i+1,j) = u(i+1,j) + 4.0*var u(i+2,j) = u(i+2,j) - var end do end if end if ! not global ! [2.1] Interior points: do j = je, js, -1 do i = ie, is, -1 uar = coefx(i,j) * ub(i,j) * data(i,j) var = coefy(i,j) * vb(i,j) * data(i,j) u(i,j) = u(i,j) + coefx(i,j)*( ub(i+1,j) - ub(i-1,j) ) * data(i,j) v(i,j) = v(i,j) + coefy(i,j)*( ub(i,j+1) - ub(i,j-1) ) * data(i,j) u(i+1,j) = u(i+1,j) + uar u(i-1,j) = u(i-1,j) - uar u(i,j+1) = u(i,j+1) + var u(i,j-1) = u(i,j-1) - var end do end do if (trace_use) call da_trace_exit("da_balance_cycloterm_adj") end subroutine da_balance_cycloterm_adj