ROMS
Loading...
Searching...
No Matches
propagator_fte.h
Go to the documentation of this file.
1 MODULE propagator_mod
2!
3!git $Id$
4!================================================== Hernan G. Arango ===
5! Copyright (c) 2002-2025 The ROMS Group Andrew M. Moore !
6! Licensed under a MIT/X style license !
7! See License_ROMS.md !
8!=======================================================================
9! !
10! Finite Time Eigenvalues Propagator: !
11! !
12! This routine it is used to compute the Finite Time Eigenmodes !
13! (FTEs) of the propagator R(0,t) linearized about a time evolving !
14! circulation. It involves a single integration of an arbitrary !
15! perturbation state vector forward in time over the interval [0,t] !
16! by the tangent linear model. !
17! !
18! Reference: !
19! !
20! Moore, A.M. et al., 2004: A comprehensive ocean prediction and !
21! analysis system based on the tangent linear and adjoint of a !
22! regional ocean model, Ocean Modelling, 7, 227-258. !
23! !
24!=======================================================================
25!
26 USE mod_kinds
27!
28 implicit none
29!
30 PRIVATE
31 PUBLIC :: propagator_fte
32!
33 CONTAINS
34!
35!***********************************************************************
36 SUBROUTINE propagator_fte (RunInterval, state, tl_state)
37!***********************************************************************
38!
39 USE mod_param
40 USE mod_parallel
41#ifdef SOLVE3D
42 USE mod_coupling
43#endif
44 USE mod_iounits
45 USE mod_ocean
46 USE mod_scalars
47 USE mod_stepping
48!
49 USE close_io_mod, ONLY : close_inp
51 USE packing_mod, ONLY : tl_unpack, tl_pack
52#ifdef SOLVE3D
53 USE set_depth_mod, ONLY : set_depth
54#endif
55 USE strings_mod, ONLY : founderror
56!
57! Imported variable declarations.
58!
59 real(dp), intent(in) :: runinterval
60!
61 TYPE (t_gst), intent(in) :: state(ngrids)
62 TYPE (t_gst), intent(inout) :: tl_state(ngrids)
63!
64! Local variable declarations.
65!
66#ifdef SOLVE3D
67 logical :: firstpass = .true.
68!
69#endif
70 integer :: ng, tile
71!
72 real(r8) :: statenorm(ngrids)
73!
74 character (len=*), parameter :: myfile = &
75 & __FILE__
76!
77!=======================================================================
78! Forward integration of the tangent linear model.
79!=======================================================================
80!
81!$OMP MASTER
82 nrun=nrun+1
83 IF (master) THEN
84 DO ng=1,ngrids
85 WRITE (stdout,10) ' PROPAGATOR - Grid: ', ng, &
86 & ', Iteration: ', nrun, &
87 & ', number converged RITZ values: ', &
88 & nconv(ng)
89 END DO
90 END IF
91!$OMP END MASTER
92!
93! Initialize time stepping indices and counters.
94!
95 DO ng=1,ngrids
96 iif(ng)=1
97 indx1(ng)=1
98 kstp(ng)=1
99 krhs(ng)=1
100 knew(ng)=1
101 predictor_2d_step(ng)=.false.
102!
103 iic(ng)=0
104 nstp(ng)=1
105 nrhs(ng)=1
106 nnew(ng)=1
107!
108 synchro_flag(ng)=.true.
109 tdays(ng)=dstart
110 time(ng)=tdays(ng)*day2sec
111!$OMP MASTER
112 ntstart(ng)=int((time(ng)-dstart*day2sec)/dt(ng))+1
113 ntend(ng)=ntimes(ng)
114 ntfirst(ng)=ntstart(ng)
115!$OMP END MASTER
116 END DO
117!$OMP BARRIER
118!
119!-----------------------------------------------------------------------
120! Clear tangent linear state variables. There is not need to clean
121! the basic state arrays since they were zeroth out at initialization
122! and bottom of previous iteration.
123!-----------------------------------------------------------------------
124!
125 DO ng=1,ngrids
126 DO tile=first_tile(ng),last_tile(ng),+1
127 CALL initialize_ocean (ng, tile, itlm)
128 END DO
129!$OMP BARRIER
130 END DO
131
132#ifdef SOLVE3D
133!
134!-----------------------------------------------------------------------
135! Compute basic state initial level thicknesses used for state norm
136! scaling. It uses zero time averaged free-surface (rest state).
137! Therefore, the norm scaling is time invariant.
138!-----------------------------------------------------------------------
139!
140 DO ng=1,ngrids
141 DO tile=last_tile(ng),first_tile(ng),-1
142 CALL set_depth (ng, tile, itlm)
143 END DO
144!$OMP BARRIER
145 END DO
146#endif
147!
148!-----------------------------------------------------------------------
149! Unpack tangent linear initial conditions from state vector.
150!-----------------------------------------------------------------------
151!
152 DO ng=1,ngrids
153 DO tile=first_tile(ng),last_tile(ng),+1
154 CALL tl_unpack (ng, tile, nstr(ng), nend(ng), &
155 & state(ng)%vector)
156 END DO
157!$OMP BARRIER
158 END DO
159!
160!-----------------------------------------------------------------------
161! Compute initial tangent linear state dot product norm.
162!-----------------------------------------------------------------------
163!
164 DO ng=1,ngrids
165 DO tile=last_tile(ng),first_tile(ng),-1
166 CALL tl_statenorm (ng, tile, kstp(ng), nstp(ng), &
167 & statenorm(ng))
168 END DO
169!$OMP BARRIER
170
171!$OMP MASTER
172 IF (master) THEN
173 WRITE (stdout,20) ' PROPAGATOR - Grid: ', ng, &
174 & ', Tangent Initial Norm: ', statenorm(ng)
175 END IF
176!$OMP END MASTER
177 END DO
178!
179!-----------------------------------------------------------------------
180! Read in initial forcing, climatology and assimilation data from
181! input NetCDF files. It loads the first relevant data record for
182! the time-interpolation between snapshots.
183!-----------------------------------------------------------------------
184!
185 DO ng=1,ngrids
186!$OMP MASTER
187 CALL close_inp (ng, itlm)
188 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
189
190 CALL tl_get_idata (ng)
191 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
192
193 CALL tl_get_data (ng)
194!$OMP END MASTER
195!$OMP BARRIER
196 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
197 END DO
198!
199!-----------------------------------------------------------------------
200! Time-step the tangent linear model.
201!-----------------------------------------------------------------------
202!
203 DO ng=1,ngrids
204!$OMP MASTER
205 IF (master) THEN
206 WRITE (stdout,30) 'TL', ng, ntstart(ng), ntend(ng)
207 END IF
208 time(ng)=time(ng)-dt(ng)
209!$OMP END MASTER
210 iic(ng)=ntstart(ng)-1
211 END DO
212!$OMP BARRIER
213
214#ifdef SOLVE3D
215 CALL tl_main3d (runinterval)
216#else
217 CALL tl_main2d (runinterval)
218#endif
219!$OMP BARRIER
220 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
221!
222!-----------------------------------------------------------------------
223! Clear nonlinear state (basic state) variables and insure that the
224! time averaged free-surface is zero for scaling below and next
225! iteration.
226!-----------------------------------------------------------------------
227!
228 DO ng=1,ngrids
229 DO tile=first_tile(ng),last_tile(ng),+1
230 CALL initialize_ocean (ng, tile, inlm)
231#ifdef SOLVE3D
232 CALL initialize_coupling (ng, tile, 0)
233#endif
234 END DO
235!$OMP BARRIER
236 END DO
237
238#ifdef SOLVE3D
239!
240!-----------------------------------------------------------------------
241! Compute basic state final level thicknesses used for state norm
242! scaling. It uses zero time averaged free-surface (rest state).
243! Therefore, the norm scaling is time invariant.
244!-----------------------------------------------------------------------
245!
246 DO ng=1,ngrids
247 DO tile=last_tile(ng),first_tile(ng),-1
248 CALL set_depth (ng, tile, itlm)
249 END DO
250!$OMP BARRIER
251 END DO
252#endif
253!
254!-----------------------------------------------------------------------
255! Compute final tangent linear state dot product norm.
256!-----------------------------------------------------------------------
257!
258 DO ng=1,ngrids
259 DO tile=first_tile(ng),last_tile(ng),+1
260 CALL tl_statenorm (ng, tile, knew(ng), nstp(ng), &
261 & statenorm(ng))
262 END DO
263!$OMP BARRIER
264
265!$OMP MASTER
266 IF (master) THEN
267 WRITE (stdout,20) ' PROPAGATOR - Grid: ', ng, &
268 & ', Tangent Final Norm: ', statenorm(ng)
269 END IF
270!$OMP END MASTER
271 END DO
272!
273!-----------------------------------------------------------------------
274! Pack final tangent linear solution into tangent state vector.
275!-----------------------------------------------------------------------
276!
277 DO ng=1,ngrids
278 DO tile=last_tile(ng),first_tile(ng),-1
279 CALL tl_pack (ng, tile, nstr(ng), nend(ng), &
280 & tl_state(ng)%vector)
281 END DO
282!$OMP BARRIER
283 END DO
284!
285 10 FORMAT (/,a,i2.2,a,i3.3,a,i3.3/)
286 20 FORMAT (/,a,i2.2,a,1p,e15.6,/)
287 30 FORMAT (/,1x,a,1x,'ROMS: started time-stepping:', &
288 & ' (Grid: ',i2.2,' TimeSteps: ',i8.8,' - ',i8.8,')')
289
290 RETURN
291 END SUBROUTINE propagator_fte
292
293 END MODULE propagator_mod
subroutine, public close_inp(ng, model)
Definition close_io.F:92
subroutine, public tl_statenorm(ng, tile, kinp, ninp, statenorm)
Definition dotproduct.F:699
subroutine, public initialize_coupling(ng, tile, model)
integer stdout
integer, parameter r8
Definition mod_kinds.F:28
integer, parameter dp
Definition mod_kinds.F:25
subroutine, public initialize_ocean(ng, tile, model)
Definition mod_ocean.F:1526
integer, dimension(:), allocatable first_tile
logical master
integer, dimension(:), allocatable last_tile
integer, parameter inlm
Definition mod_param.F:662
integer, dimension(:), allocatable nstr
Definition mod_param.F:646
integer ngrids
Definition mod_param.F:113
integer, parameter itlm
Definition mod_param.F:663
integer, dimension(:), allocatable nend
Definition mod_param.F:647
real(dp), parameter day2sec
integer, dimension(:), allocatable ntimes
integer, dimension(:), allocatable iic
real(dp), dimension(:), allocatable dt
integer, dimension(:), allocatable nconv
logical, dimension(:), allocatable synchro_flag
logical, dimension(:), allocatable predictor_2d_step
real(dp), dimension(:), allocatable tdays
real(dp) dstart
integer, dimension(:), allocatable ntend
integer exit_flag
integer, dimension(:), allocatable indx1
integer, dimension(:), allocatable ntfirst
real(dp), dimension(:), allocatable time
integer, dimension(:), allocatable ntstart
integer nrun
integer, dimension(:), allocatable iif
integer noerror
integer, dimension(:), allocatable kstp
integer, dimension(:), allocatable knew
integer, dimension(:), allocatable nrhs
integer, dimension(:), allocatable nnew
integer, dimension(:), allocatable krhs
integer, dimension(:), allocatable nstp
subroutine, public propagator_fte(runinterval, state, tl_state)
subroutine, public set_depth(ng, tile, model)
Definition set_depth.F:34
logical function, public founderror(flag, noerr, line, routine)
Definition strings.F:52
subroutine tl_pack(ng, tile, mstr, mend, tl_state)
Definition packing.F:6015
subroutine tl_unpack(ng, tile, mstr, mend, state)
Definition packing.F:6561
subroutine tl_get_data(ng)
Definition tl_get_data.F:4
subroutine tl_get_idata(ng)
Definition tl_get_idata.F:4
subroutine tl_main2d
Definition tl_main2d.F:429
subroutine tl_main3d(runinterval)
Definition tl_main3d.F:4