ROMS
Loading...
Searching...
No Matches
tl_frc_adjust.F
Go to the documentation of this file.
1#include "cppdefs.h"
3
4#if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
5!
6!git $Id$
7!================================================== Hernan G. Arango ===
8! Copyright (c) 2002-2025 The ROMS Group !
9! Licensed under a MIT/X style license !
10! See License_ROMS.md !
11!=======================================================================
12! !
13! This routine time-interpolates 4DVar tangent linear model surface !
14! forcing increments. The increments can be constant (Nfrec=1) or !
15! time interpolated between snapshots (Nfrec>1). !
16! !
17! On Input: !
18! !
19! ng Nested grid number. !
20! tile Domain partition. !
21! Linp Tangent linear state time index to process. !
22! !
23!=======================================================================
24!
25 implicit none
26
27 PRIVATE
28 PUBLIC :: tl_frc_adjust
29
30 CONTAINS
31!
32!***********************************************************************
33 SUBROUTINE tl_frc_adjust (ng, tile, Linp)
34!***********************************************************************
35!
36 USE mod_param
37 USE mod_forces
38!
39! Imported variable declarations.
40!
41 integer, intent(in) :: ng, tile, linp
42!
43! Local variable declarations.
44!
45 character (len=*), parameter :: myfile = &
46 & __FILE__
47!
48# include "tile.h"
49!
50# ifdef PROFILE
51 CALL wclock_on (ng, itlm, 7, __line__, myfile)
52# endif
53 CALL tl_frc_adjust_tile (ng, tile, &
54 & lbi, ubi, lbj, ubj, &
55 & imins, imaxs, jmins, jmaxs, &
56# ifdef ADJUST_WSTRESS
57 & forces(ng) % tl_ustr, &
58 & forces(ng) % tl_vstr, &
59 & forces(ng) % tl_sustr, &
60 & forces(ng) % tl_svstr, &
61# endif
62# if defined ADJUST_STFLUX && defined SOLVE3D
63 & forces(ng) % tl_tflux, &
64 & forces(ng) % tl_stflx, &
65# endif
66 & linp)
67# ifdef PROFILE
68 CALL wclock_off (ng, itlm, 7, __line__, myfile)
69# endif
70!
71 RETURN
72 END SUBROUTINE tl_frc_adjust
73!
74!***********************************************************************
75 SUBROUTINE tl_frc_adjust_tile (ng, tile, &
76 & LBi, UBi, LBj, UBj, &
77 & IminS, ImaxS, JminS, JmaxS, &
78# ifdef ADJUST_WSTRESS
79 & tl_ustr, tl_vstr, &
80 & tl_sustr, tl_svstr, &
81# endif
82# if defined ADJUST_STFLUX && defined SOLVE3D
83 & tl_tflux, tl_stflx, &
84# endif
85 & Linp)
86!***********************************************************************
87!
88 USE mod_param
89 USE mod_scalars
90!
91! Imported variable declarations.
92!
93 integer, intent(in) :: ng, tile
94 integer, intent(in) :: LBi, UBi, LBj, UBj
95 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
96 integer, intent(in) :: Linp
97!
98# ifdef ASSUMED_SHAPE
99# ifdef ADJUST_WSTRESS
100 real(r8), intent(in) :: tl_ustr(LBi:,LBj:,:,:)
101 real(r8), intent(in) :: tl_vstr(LBi:,LBj:,:,:)
102# endif
103# if defined ADJUST_STFLUX && defined SOLVE3D
104 real(r8), intent(inout) :: tl_tflux(LBi:,LBj:,:,:,:)
105# endif
106# ifdef ADJUST_WSTRESS
107 real(r8), intent(inout) :: tl_sustr(LBi:,LBj:)
108 real(r8), intent(inout) :: tl_svstr(LBi:,LBj:)
109# endif
110# if defined ADJUST_STFLUX && defined SOLVE3D
111 real(r8), intent(inout) :: tl_stflx(LBi:,LBj:,:)
112# endif
113# else
114# ifdef ADJUST_WSTRESS
115 real(r8), intent(in) :: tl_ustr(LBi:UBi,LBj:UBj,Nfrec(ng),2)
116 real(r8), intent(in) :: tl_vstr(LBi:UBi,LBj:UBj,Nfrec(ng),2)
117# endif
118# if defined ADJUST_STFLUX && defined SOLVE3D
119 real(r8), intent(in) :: tl_tflux(LBi:UBi,LBj:UBj, &
120 & Nfrec(ng),2,NT(ng))
121# endif
122# ifdef ADJUST_WSTRESS
123 real(r8), intent(inout) :: tl_sustr(LBi:UBi,LBj:UBj)
124 real(r8), intent(inout) :: tl_svstr(LBi:UBi,LBj:UBj)
125# endif
126# if defined ADJUST_STFLUX && defined SOLVE3D
127 real(r8), intent(inout) :: tl_stflx(LBi:UBi,LBj:UBj,NT(ng))
128# endif
129# endif
130!
131! Local variable declarations.
132!
133 integer :: i, it1, it2, j
134# ifdef SOLVE3D
135 integer :: itrc
136# endif
137 real(r8) :: fac, fac1, fac2
138
139# include "set_bounds.h"
140!
141!-----------------------------------------------------------------------
142! Adjust tangent linear surface forcing fields with 4DVar increments.
143!-----------------------------------------------------------------------
144!
145! Set time records and interpolation factor, if any.
146!
147 IF (nfrec(ng).eq.1) THEN
148 it1=1
149 it2=1
150 fac1=1.0_r8
151 fac2=0.0_r8
152 ELSE
153# ifdef GENERIC_DSTART
154 it1=max(0,(iic(ng)-ntstart(ng))/nsff(ng))+1
155# else
156 it1=max(0,(iic(ng)-1)/nsff(ng))+1
157# endif
158 it2=min(it1+1,nfrec(ng))
159 fac1=sf_time(it2,ng)-(time(ng)+dt(ng))
160 fac2=(time(ng)+dt(ng))-sf_time(it1,ng)
161 fac=1.0_r8/(fac1+fac2)
162 fac1=fac*fac1
163 fac2=fac*fac2
164 END IF
165
166# ifdef ADJUST_WSTRESS
167!
168! Adjust surface wind stress. Interpolate between surface forcing
169! increments, if appropriate.
170!
171 DO j=jstrr,jendr
172 DO i=istr,iendr
173 tl_sustr(i,j)=fac1*tl_ustr(i,j,it1,linp)+ &
174 & fac2*tl_ustr(i,j,it2,linp)
175 END DO
176 END DO
177 DO j=jstr,jendr
178 DO i=istrr,iendr
179 tl_svstr(i,j)=fac1*tl_vstr(i,j,it1,linp)+ &
180 & fac2*tl_vstr(i,j,it2,linp)
181 END DO
182 END DO
183# endif
184# if defined ADJUST_STFLUX && defined SOLVE3D
185!
186! Adjust surface tracer fluxes. Interpolate between surface forcing
187! increments, if appropriate.
188!
189 IF (lstflux(itemp,ng)) THEN
190# ifdef QCORRECTION
191 DO j=jstrr,jendr
192 DO i=istrr,iendr
193 tl_stflx(i,j,itemp)=tl_stflx(i,j,itemp)+ &
194 & fac1*tl_tflux(i,j,it1,linp,itemp)+ &
195 & fac2*tl_tflux(i,j,it2,linp,itemp)
196 END DO
197 END DO
198# else
199 DO j=jstrr,jendr
200 DO i=istrr,iendr
201 tl_stflx(i,j,itemp)=fac1*tl_tflux(i,j,it1,linp,itemp)+ &
202 & fac2*tl_tflux(i,j,it2,linp,itemp)
203 END DO
204 END DO
205# endif
206 END IF
207
208# ifdef SALINITY
209 IF (lstflux(isalt,ng)) THEN
210 DO j=jstrr,jendr
211 DO i=istrr,iendr
212 tl_stflx(i,j,isalt)=tl_stflx(i,j,isalt)+ &
213 & fac1*tl_tflux(i,j,it1,linp,isalt)+ &
214 & fac2*tl_tflux(i,j,it2,linp,isalt)
215 END DO
216 END DO
217 END IF
218# endif
219
220 DO itrc=nat+1,nt(ng)
221 IF (lstflux(itrc,ng)) THEN
222 DO j=jstrr,jendr
223 DO i=istrr,iendr
224 tl_stflx(i,j,itrc)=fac1*tl_tflux(i,j,it1,linp,itrc)+ &
225 & fac2*tl_tflux(i,j,it2,linp,itrc)
226 END DO
227 END DO
228 END IF
229 END DO
230# endif
231!
232 RETURN
233 END SUBROUTINE tl_frc_adjust_tile
234#endif
235 END MODULE tl_frc_adjust_mod
type(t_forces), dimension(:), allocatable forces
Definition mod_forces.F:554
integer nat
Definition mod_param.F:499
integer, parameter itlm
Definition mod_param.F:663
integer, dimension(:), allocatable iic
real(dp), dimension(:), allocatable dt
logical, dimension(:,:), allocatable lstflux
integer isalt
integer itemp
real(dp), dimension(:), allocatable time
integer, dimension(:), allocatable nsff
integer, dimension(:), allocatable ntstart
real(dp), dimension(:,:), allocatable sf_time
subroutine tl_frc_adjust_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, tl_ustr, tl_vstr, tl_sustr, tl_svstr, tl_tflux, tl_stflx, linp)
subroutine, public tl_frc_adjust(ng, tile, linp)
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