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

Functions/Subroutines

subroutine, public frc_adjust (ng, tile, linp)
 
subroutine frc_adjust_tile (ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, tl_ustr, tl_vstr, ustr, vstr, sustr, svstr, tl_tflux, tflux, stflx, linp)
 
subroutine, public load_frc (ng, tile, lout)
 
subroutine load_frc_tile (ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, lout)
 

Function/Subroutine Documentation

◆ frc_adjust()

subroutine, public frc_adjust_mod::frc_adjust ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) linp )

Definition at line 35 of file frc_adjust.F.

36!***********************************************************************
37!
38 USE mod_param
39 USE mod_forces
40!
41! Imported variable declarations.
42!
43 integer, intent(in) :: ng, tile, Linp
44!
45! Local variable declarations.
46!
47 character (len=*), parameter :: MyFile = &
48 & __FILE__
49!
50# include "tile.h"
51!
52# ifdef PROFILE
53 CALL wclock_on (ng, inlm, 7, __line__, myfile)
54# endif
55 CALL frc_adjust_tile (ng, tile, &
56 & lbi, ubi, lbj, ubj, &
57 & imins, imaxs, jmins, jmaxs, &
58# ifdef ADJUST_WSTRESS
59 & forces(ng) % tl_ustr, &
60 & forces(ng) % tl_vstr, &
61 & forces(ng) % ustr, &
62 & forces(ng) % vstr, &
63 & forces(ng) % sustr, &
64 & forces(ng) % svstr, &
65# endif
66# if defined ADJUST_STFLUX && defined SOLVE3D
67 & forces(ng) % tl_tflux, &
68 & forces(ng) % tflux, &
69 & forces(ng) % stflx, &
70# endif
71 & linp)
72# ifdef PROFILE
73 CALL wclock_off (ng, inlm, 7, __line__, myfile)
74# endif
75!
76 RETURN
type(t_forces), dimension(:), allocatable forces
Definition mod_forces.F:554
integer, parameter inlm
Definition mod_param.F:662
recursive subroutine wclock_off(ng, model, region, line, routine)
Definition timers.F:148
recursive subroutine wclock_on(ng, model, region, line, routine)
Definition timers.F:3

References mod_forces::forces, frc_adjust_tile(), mod_param::inlm, wclock_off(), and wclock_on().

Referenced by main3d().

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

◆ frc_adjust_tile()

subroutine frc_adjust_mod::frc_adjust_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,
real(r8), dimension(lbi:,lbj:,:,:), intent(in) tl_ustr,
real(r8), dimension(lbi:,lbj:,:,:), intent(in) tl_vstr,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) ustr,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) vstr,
real(r8), dimension(lbi:,lbj:), intent(inout) sustr,
real(r8), dimension(lbi:,lbj:), intent(inout) svstr,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(inout) tl_tflux,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(inout) tflux,
real(r8), dimension(lbi:,lbj:,:), intent(inout) stflx,
integer, intent(in) linp )
private

Definition at line 80 of file frc_adjust.F.

92!***********************************************************************
93!
94 USE mod_param
95 USE mod_scalars
96!
97! Imported variable declarations.
98!
99 integer, intent(in) :: ng, tile
100 integer, intent(in) :: LBi, UBi, LBj, UBj
101 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
102 integer, intent(in) :: Linp
103!
104# ifdef ASSUMED_SHAPE
105# ifdef ADJUST_WSTRESS
106 real(r8), intent(in) :: tl_ustr(LBi:,LBj:,:,:)
107 real(r8), intent(in) :: tl_vstr(LBi:,LBj:,:,:)
108# endif
109# if defined ADJUST_STFLUX && defined SOLVE3D
110 real(r8), intent(inout) :: tflux(LBi:,LBj:,:,:,:)
111 real(r8), intent(inout) :: tl_tflux(LBi:,LBj:,:,:,:)
112# endif
113# ifdef ADJUST_WSTRESS
114 real(r8), intent(inout) :: ustr(LBi:,LBj:,:,:)
115 real(r8), intent(inout) :: vstr(LBi:,LBj:,:,:)
116 real(r8), intent(inout) :: sustr(LBi:,LBj:)
117 real(r8), intent(inout) :: svstr(LBi:,LBj:)
118# endif
119# if defined ADJUST_STFLUX && defined SOLVE3D
120 real(r8), intent(inout) :: stflx(LBi:,LBj:,:)
121# endif
122# else
123# ifdef ADJUST_WSTRESS
124 real(r8), intent(in) :: tl_ustr(LBi:UBi,LBj:UBj,Nfrec(ng),2)
125 real(r8), intent(in) :: tl_vstr(LBi:UBi,LBj:UBj,Nfrec(ng),2)
126# endif
127# if defined ADJUST_STFLUX && defined SOLVE3D
128 real(r8), intent(in) :: tflux(LBi:UBi,LBj:UBj, &
129 & Nfrec(ng),2,NT(ng))
130 real(r8), intent(in) :: tl_tflux(LBi:UBi,LBj:UBj, &
131 & Nfrec(ng),2,NT(ng))
132# endif
133# ifdef ADJUST_WSTRESS
134 real(r8), intent(inout) :: ustr(LBi:UBi,LBj:UBj,Nfrec(ng),2)
135 real(r8), intent(inout) :: vstr(LBi:UBi,LBj:UBj,Nfrec(ng),2)
136 real(r8), intent(inout) :: sustr(LBi:UBi,LBj:UBj)
137 real(r8), intent(inout) :: svstr(LBi:UBi,LBj:UBj)
138# endif
139# if defined ADJUST_STFLUX && defined SOLVE3D
140 real(r8), intent(inout) :: stflx(LBi:UBi,LBj:UBj,NT(ng))
141# endif
142# endif
143!
144! Local variable declarations.
145!
146 integer :: i, it1, it2, j
147# ifdef SOLVE3D
148 integer :: itrc
149# endif
150 real(r8) :: fac, fac1, fac2
151
152# include "set_bounds.h"
153!
154!-----------------------------------------------------------------------
155! Adjust nonlinear surface forcing fields with 4DVar increments.
156!-----------------------------------------------------------------------
157!
158! Set time records and interpolation factor, if any.
159!
160 IF (nfrec(ng).eq.1) THEN
161 it1=1
162 it2=1
163 fac1=1.0_r8
164 fac2=0.0_r8
165 ELSE
166# ifdef GENERIC_DSTART
167 it1=max(0,(iic(ng)-ntstart(ng))/nsff(ng))+1
168# else
169 it1=max(0,(iic(ng)-1)/nsff(ng))+1
170# endif
171 it2=min(it1+1,nfrec(ng))
172 fac1=sf_time(it2,ng)-(time(ng)+dt(ng))
173 fac2=(time(ng)+dt(ng))-sf_time(it1,ng)
174 fac=1.0_r8/(fac1+fac2)
175 fac1=fac*fac1
176 fac2=fac*fac2
177 END IF
178
179# ifdef ADJUST_WSTRESS
180!
181! Adjust surface momentum stress. Interpolate between surface forcing
182! increments, if appropriate.
183!
184 DO j=jstrr,jendr
185 DO i=istr,iendr
186 sustr(i,j)=sustr(i,j)+ &
187 & fac1*tl_ustr(i,j,it1,linp)+ &
188 & fac2*tl_ustr(i,j,it2,linp)
189 END DO
190 END DO
191 DO j=jstr,jendr
192 DO i=istrr,iendr
193 svstr(i,j)=svstr(i,j)+ &
194 & fac1*tl_vstr(i,j,it1,linp)+ &
195 & fac2*tl_vstr(i,j,it2,linp)
196 END DO
197 END DO
198# endif
199# if defined ADJUST_STFLUX && defined SOLVE3D
200!
201! Adjust surface tracer fluxes. Interpolate between surface forcing
202! increments, if appropriate.
203!
204 DO itrc=1,nt(ng)
205 IF (lstflux(itrc,ng)) THEN
206 DO j=jstrr,jendr
207 DO i=istrr,iendr
208 stflx(i,j,itrc)=stflx(i,j,itrc)+ &
209 & fac1*tl_tflux(i,j,it1,linp,itrc)+ &
210 & fac2*tl_tflux(i,j,it2,linp,itrc)
211 END DO
212 END DO
213 END IF
214 END DO
215# endif
216!
217 RETURN
integer, dimension(:), allocatable nt
Definition mod_param.F:489
integer, dimension(:), allocatable iic
real(dp), dimension(:), allocatable dt
logical, dimension(:,:), allocatable lstflux
integer, dimension(:), allocatable nfrec
real(dp), dimension(:), allocatable time
integer, dimension(:), allocatable nsff
integer, dimension(:), allocatable ntstart
real(dp), dimension(:,:), allocatable sf_time

References mod_scalars::dt, mod_scalars::iic, mod_scalars::lstflux, mod_scalars::nsff, mod_scalars::ntstart, mod_scalars::sf_time, and mod_scalars::time.

Referenced by frc_adjust().

Here is the caller graph for this function:

◆ load_frc()

subroutine, public frc_adjust_mod::load_frc ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) lout )

Definition at line 220 of file frc_adjust.F.

221!
222!=======================================================================
223! !
224! This routine surface forcing into nonlinear storage arrays. In !
225! 4DVAR, the surface forcing adjustment are stored in arrays with !
226! extra dimensions to aid minimization at other times in addition !
227! to initialization time. !
228! !
229! On Input: !
230! !
231! ng Nested grid number. !
232! tile Domain partition. !
233! Lout Time index to process in storage arrays. !
234! !
235!=======================================================================
236!
237 USE mod_param
238!
239! Imported variable declarations.
240!
241 integer, intent(in) :: ng, tile, Lout
242!
243! Local variable declarations.
244!
245 character (len=*), parameter :: MyFile = &
246 & __FILE__//", load_frc"
247!
248# include "tile.h"
249!
250# ifdef PROFILE
251 CALL wclock_on (ng, inlm, 8, __line__, myfile)
252# endif
253 CALL load_frc_tile (ng, tile, &
254 & lbi, ubi, lbj, ubj, &
255 & imins, imaxs, jmins, jmaxs, &
256 & lout)
257# ifdef PROFILE
258 CALL wclock_off (ng, inlm, 8, __line__, myfile)
259# endif
260!
261 RETURN

References mod_param::inlm, load_frc_tile(), wclock_off(), and wclock_on().

Referenced by main3d().

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

◆ load_frc_tile()

subroutine frc_adjust_mod::load_frc_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) lout )
private

Definition at line 265 of file frc_adjust.F.

269!***********************************************************************
270!
271 USE mod_param
272 USE mod_forces
273 USE mod_ncparam
274 USE mod_scalars
275!
276! Imported variable declarations.
277!
278 integer, intent(in) :: ng, tile
279 integer, intent(in) :: LBi, UBi, LBj, UBj
280 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
281 integer, intent(in) :: Lout
282!
283! Local variable declarations.
284!
285 integer :: i, ic, j
286# ifdef SOLVE3D
287 integer :: itrc, k
288# endif
289
290# include "set_bounds.h"
291!
292!-----------------------------------------------------------------------
293! Load nonlinear open boundary fields into storage arrays.
294!-----------------------------------------------------------------------
295!
296 load_sf : IF (mod(iic(ng)-1,nsff(ng)).eq.0) THEN
297 sfcount(ng)=sfcount(ng)+1
298 ic=sfcount(ng)
299
300# ifdef ADJUST_WSTRESS
301!
302! Adjust surface momentum stress.
303!
304 DO j=jstrr,jendr
305 DO i=istr,iendr
306 forces(ng)%ustr(i,j,ic,lout)=forces(ng)%sustr(i,j)
307 END DO
308 END DO
309 DO j=jstr,jendr
310 DO i=istrr,iendr
311 forces(ng)%vstr(i,j,ic,lout)=forces(ng)%svstr(i,j)
312 END DO
313 END DO
314# endif
315
316# if defined ADJUST_STFLUX && defined SOLVE3D
317!
318! Adjust surface tracer fluxes.
319!
320 DO itrc=1,nt(ng)
321 IF (lstflux(itrc,ng)) THEN
322 DO j=jstrr,jendr
323 DO i=istrr,iendr
324 forces(ng)%tflux(i,j,ic,lout,itrc)= &
325 & forces(ng)%stflx(i,j,itrc)
326 END DO
327 END DO
328 END IF
329 END DO
330# endif
331 END IF load_sf
332!
333 RETURN
integer, dimension(:), allocatable sfcount

References mod_forces::forces, mod_scalars::iic, mod_scalars::lstflux, mod_scalars::nsff, mod_param::nt, and mod_scalars::sfcount.

Referenced by load_frc().

Here is the caller graph for this function: