ROMS
Loading...
Searching...
No Matches
ad_force_dual_mod Module Reference

Functions/Subroutines

subroutine, public ad_force_dual (ng, tile, kfrc, nfrc)
 
subroutine ad_force_dual_tile (ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, kfrc, nfrc, f_t, f_u, f_v, f_zeta, ad_t, ad_u, ad_v, ad_zeta)
 

Function/Subroutine Documentation

◆ ad_force_dual()

subroutine, public ad_force_dual_mod::ad_force_dual ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) kfrc,
integer, intent(in) nfrc )

Definition at line 25 of file ad_force_dual.F.

26!***********************************************************************
27!
28 USE mod_param
29 USE mod_ocean
30!
31! Imported variable declarations.
32!
33 integer, intent(in) :: ng, tile, Kfrc, Nfrc
34!
35! Local variable declarations.
36!
37# include "tile.h"
38!
39 CALL ad_force_dual_tile (ng, tile, &
40 & lbi, ubi, lbj, ubj, &
41 & imins, imaxs, jmins, jmaxs, &
42 & kfrc, nfrc, &
43# ifdef SOLVE3D
44 & ocean(ng) % f_t, &
45 & ocean(ng) % f_u, &
46 & ocean(ng) % f_v, &
47# else
48 & ocean(ng) % f_ubar, &
49 & ocean(ng) % f_vbar, &
50# endif
51 & ocean(ng) % f_zeta, &
52# ifdef SOLVE3D
53 & ocean(ng) % ad_t, &
54 & ocean(ng) % ad_u, &
55 & ocean(ng) % ad_v, &
56# else
57 & ocean(ng) % ad_ubar, &
58 & ocean(ng) % ad_vbar, &
59# endif
60 & ocean(ng) % ad_zeta)
61
62 RETURN
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351

References ad_force_dual_tile(), and mod_ocean::ocean.

Referenced by ad_main3d().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ ad_force_dual_tile()

subroutine ad_force_dual_mod::ad_force_dual_tile ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) lbi,
integer, intent(in) ubi,
integer, intent(in) lbj,
integer, intent(in) ubj,
integer, intent(in) imins,
integer, intent(in) imaxs,
integer, intent(in) jmins,
integer, intent(in) jmaxs,
integer, intent(in) kfrc,
integer, intent(in) nfrc,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) f_t,
real(r8), dimension(lbi:,lbj:,:), intent(inout) f_u,
real(r8), dimension(lbi:,lbj:,:), intent(inout) f_v,
real(r8), dimension(lbi:,lbj:), intent(inout) f_zeta,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(inout) ad_t,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) ad_u,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) ad_v,
real(r8), dimension(lbi:,lbj:,:), intent(inout) ad_zeta )
private

Definition at line 66 of file ad_force_dual.F.

82!***********************************************************************
83!
84 USE mod_param
85!
86! Imported variable declarations.
87!
88 integer, intent(in) :: ng, tile
89 integer, intent(in) :: LBi, UBi, LBj, UBj
90 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
91 integer, intent(in) :: Kfrc
92 integer, intent(in) :: Nfrc
93!
94# ifdef ASSUMED_SHAPE
95# ifdef SOLVE3D
96 real(r8), intent(inout) :: f_t(LBi:,LBj:,:,:)
97 real(r8), intent(inout) :: f_u(LBi:,LBj:,:)
98 real(r8), intent(inout) :: f_v(LBi:,LBj:,:)
99# else
100 real(r8), intent(inout) :: f_ubar(LBi:,LBj:)
101 real(r8), intent(inout) :: f_vbar(LBi:,LBj:)
102# endif
103 real(r8), intent(inout) :: f_zeta(LBi:,LBj:)
104# ifdef SOLVE3D
105 real(r8), intent(inout) :: ad_t(LBi:,LBj:,:,:,:)
106 real(r8), intent(inout) :: ad_u(LBi:,LBj:,:,:)
107 real(r8), intent(inout) :: ad_v(LBi:,LBj:,:,:)
108# else
109 real(r8), intent(inout) :: ad_ubar(LBi:,LBj:,:)
110 real(r8), intent(inout) :: ad_vbar(LBi:,LBj:,:)
111# endif
112 real(r8), intent(inout) :: ad_zeta(LBi:,LBj:,:)
113# else
114# ifdef SOLVE3D
115 real(r8), intent(inout) :: f_t(LBi:UBi,LBj:UBj,N(ng),NT(ng))
116 real(r8), intent(inout) :: f_u(LBi:UBi,LBj:UBj,N(ng))
117 real(r8), intent(inout) :: f_v(LBi:UBi,LBj:UBj,N(ng))
118# else
119 real(r8), intent(inout) :: f_ubar(LBi:UBi,LBj:UBj)
120 real(r8), intent(inout) :: f_vbar(LBi:UBi,LBj:UBj)
121# endif
122 real(r8), intent(inout) :: f_zeta(LBi:UBi,LBj:UBj)
123# ifdef SOLVE3D
124 real(r8), intent(inout) :: ad_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
125 real(r8), intent(inout) :: ad_u(LBi:UBi,LBj:UBj,N(ng),2)
126 real(r8), intent(inout) :: ad_v(LBi:UBi,LBj:UBj,N(ng),2)
127# else
128 real(r8), intent(inout) :: ad_ubar(LBi:UBi,LBj:UBj,:)
129 real(r8), intent(inout) :: ad_vbar(LBi:UBi,LBj:UBj,:)
130# endif
131 real(r8), intent(inout) :: ad_zeta(LBi:UBi,LBj:UBj,:)
132# endif
133!
134! Local variable declarations.
135!
136 integer :: i, j
137# ifdef SOLVE3D
138 integer :: itrc, k
139# endif
140
141# include "set_bounds.h"
142!
143!-----------------------------------------------------------------------
144! Add given forcing to adjoint state.
145!-----------------------------------------------------------------------
146!
147! Adjoint free-surface.
148!
149 DO j=jstrr,jendr
150 DO i=istrr,iendr
151 ad_zeta(i,j,kfrc)=ad_zeta(i,j,kfrc)+f_zeta(i,j)
152 f_zeta(i,j)=0.0_r8
153 END DO
154 END DO
155
156# ifndef SOLVE3D
157!
158! Adjoint 2D momentum.
159!
160 DO j=jstrr,jendr
161 DO i=istr,iendr
162 ad_ubar(i,j,kfrc)=ad_ubar(i,j,kfrc)+f_ubar(i,j)
163 f_ubar(i,j)=0.0_r8
164 END DO
165 END DO
166!
167 DO j=jstr,jendr
168 DO i=istrr,iendr
169 ad_vbar(i,j,kfrc)=ad_vbar(i,j,kfrc)+f_vbar(i,j)
170 f_vbar(i,j)=0.0_r8
171 END DO
172 END DO
173# else
174!
175! Adjoint 3D momentum.
176!
177 DO k=1,n(ng)
178 DO j=jstrr,jendr
179 DO i=istr,iendr
180 ad_u(i,j,k,nfrc)=ad_u(i,j,k,nfrc)+f_u(i,j,k)
181 f_u(i,j,k)=0.0_r8
182 END DO
183 END DO
184 DO j=jstr,jendr
185 DO i=istrr,iendr
186 ad_v(i,j,k,nfrc)=ad_v(i,j,k,nfrc)+f_v(i,j,k)
187 f_v(i,j,k)=0.0_r8
188 END DO
189 END DO
190 END DO
191!
192! Adjoint tracers.
193!
194 DO itrc=1,nt(ng)
195 DO k=1,n(ng)
196 DO j=jstrr,jendr
197 DO i=istrr,iendr
198 ad_t(i,j,k,nfrc,itrc)=ad_t(i,j,k,nfrc,itrc)+ &
199 & f_t(i,j,k,itrc)
200 f_t(i,j,k,itrc)=0.0_r8
201 END DO
202 END DO
203 END DO
204 END DO
205# endif
206
207 RETURN
integer, dimension(:), allocatable n
Definition mod_param.F:479
integer, dimension(:), allocatable nt
Definition mod_param.F:489

Referenced by ad_force_dual().

Here is the caller graph for this function: