ROMS
Loading...
Searching...
No Matches
r4dvar_roms.h
Go to the documentation of this file.
1 MODULE roms_kernel_mod
2!
3!git $Id$
4!=================================================== Andrew M. Moore ===
5! Copyright (c) 2002-2025 The ROMS Group Hernan G. Arango !
6! Licensed under a MIT/X style license !
7! See License_ROMS.md !
8!=======================================================================
9! !
10! ROMS Strong/Weak Constraint 4-Dimensional Variational Data !
11! Assimilation Driver: Indirect Representer Approach !
12! (R4D-Var). !
13! !
14! This driver is used for the dual formulation (observation space), !
15! strong or weak constraint 4D-Var where errors may be considered !
16! in both model and observations. !
17! !
18! The routines in this driver control the initialization, time- !
19! stepping, and finalization of ROMS model following ESMF/NUOPC !
20! conventions: !
21! !
22! ROMS_initialize !
23! ROMS_run !
24! ROMS_finalize !
25! !
26! References: !
27! !
28! Moore, A.M., H.G. Arango, G. Broquet, B.S. Powell, A.T. Weaver, !
29! and J. Zavala-Garay, 2011: The Regional Ocean Modeling System !
30! (ROMS) 4-dimensional variational data assimilations systems, !
31! Part I - System overview and formulation, Prog. Oceanogr., 91, !
32! 34-49, doi:10.1016/j.pocean.2011.05.004. !
33! !
34! Moore, A.M., H.G. Arango, G. Broquet, C. Edward, M. Veneziani, !
35! B. Powell, D. Foley, J.D. Doyle, D. Costa, and P. Robinson, !
36! 2011: The Regional Ocean Modeling System (ROMS) 4-dimensional !
37! variational data assimilations systems, Part II - Performance !
38! and application to the California Current System, Prog. !
39! Oceanogr., 91, 50-73, doi:10.1016/j.pocean.2011.05.003. !
40! !
41!=======================================================================
42!
43 USE mod_param
44 USE mod_parallel
45 USE mod_arrays
46 USE mod_fourdvar
47 USE mod_iounits
48 USE mod_ncparam
49 USE mod_netcdf
50#if defined PIO_LIB && defined DISTRIBUTE
52#endif
53 USE mod_scalars
54 USE mod_stepping
55!
56 USE r4dvar_mod
57!
59 USE def_dai_mod, ONLY : def_dai
60 USE get_state_mod, ONLY : get_state
61 USE inp_par_mod, ONLY : inp_par
62#ifdef MCT_LIB
63# ifdef ATM_COUPLING
64 USE mct_coupler_mod, ONLY : initialize_ocn2atm_coupling
65# endif
66# ifdef WAV_COUPLING
67 USE mct_coupler_mod, ONLY : initialize_ocn2wav_coupling
68# endif
69#endif
72 USE strings_mod, ONLY : founderror
73 USE wrt_dai_mod, ONLY : wrt_dai
74 USE wrt_rst_mod, ONLY : wrt_rst
75!
76 implicit none
77!
78 PUBLIC :: roms_initialize
79 PUBLIC :: roms_run
80 PUBLIC :: roms_finalize
81!
82 CONTAINS
83!
84 SUBROUTINE roms_initialize (first, mpiCOMM)
85!
86!=======================================================================
87! !
88! This routine allocates and initializes ROMS state variables and !
89! internal parameters. It reads standard input parameters. !
90! !
91!=======================================================================
92!
93! Imported variable declarations.
94!
95 logical, intent(inout) :: first
96!
97 integer, intent(in), optional :: mpiCOMM
98!
99! Local variable declarations.
100!
101 logical :: allocate_vars = .true.
102!
103#ifdef DISTRIBUTE
104 integer :: MyError, MySize
105#endif
106 integer :: chunk_size, ng, thread
107#ifdef _OPENMP
108 integer :: my_threadnum
109#endif
110!
111 character (len=*), parameter :: MyFile = &
112 & __FILE__//", ROMS_initialize"
113
114#ifdef DISTRIBUTE
115!
116!-----------------------------------------------------------------------
117! Set distribute-memory (mpi) world communictor.
118!-----------------------------------------------------------------------
119!
120 IF (PRESENT(mpicomm)) THEN
121 ocn_comm_world=mpicomm
122 ELSE
123 ocn_comm_world=mpi_comm_world
124 END IF
125 CALL mpi_comm_rank (ocn_comm_world, myrank, myerror)
126 CALL mpi_comm_size (ocn_comm_world, mysize, myerror)
127#endif
128!
129!-----------------------------------------------------------------------
130! On first pass, initialize model parameters a variables for all
131! nested/composed grids. Notice that the logical switch "first"
132! is used to allow multiple calls to this routine during ensemble
133! configurations.
134!-----------------------------------------------------------------------
135!
136 IF (first) THEN
137 first=.false.
138!
139! Initialize parallel control switches. These scalars switches are
140! independent from standard input parameters.
141!
143!
144! Set the ROMS standard output unit to write verbose execution info.
145! Notice that the default standard out unit in Fortran is 6.
146!
147! In some applications like coupling or disjointed mpi-communications,
148! it is advantageous to write standard output to a specific filename
149! instead of the default Fortran standard output unit 6. If that is
150! the case, it opens such formatted file for writing.
151!
152 IF (set_stdoutunit) THEN
154 set_stdoutunit=.false.
155 END IF
156!
157! Read in model tunable parameters from standard input. Allocate and
158! initialize variables in several modules after the number of nested
159! grids and dimension parameters are known.
160!
161 CALL inp_par (inlm)
162 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
163!
164! Initialize counters. The 'Nrun' counter will be recomputed in the
165! RBL4D-Var phases to process the obervation operator correctly.
166!
167 nrun=1 ! run counter
168 erstr=1 ! ensemble start counter
169 erend=nouter ! ensemble end counter
170!
171! Set domain decomposition tile partition range. This range is
172! computed only once since the "first_tile" and "last_tile" values
173! are private for each parallel thread/node.
174!
175#if defined _OPENMP
176 mythread=my_threadnum()
177#elif defined DISTRIBUTE
179#else
180 mythread=0
181#endif
182 DO ng=1,ngrids
183 chunk_size=(ntilex(ng)*ntilee(ng)+numthreads-1)/numthreads
184 first_tile(ng)=mythread*chunk_size
185 last_tile(ng)=first_tile(ng)+chunk_size-1
186 END DO
187!
188! Initialize internal wall clocks. Notice that the timings does not
189! includes processing standard input because several parameters are
190! needed to allocate clock variables.
191!
192 IF (master) THEN
193 WRITE (stdout,10)
194 10 FORMAT (/,' Process Information:',/)
195 END IF
196!
197 DO ng=1,ngrids
198 DO thread=thread_range
199 CALL wclock_on (ng, inlm, 0, __line__, myfile)
200 END DO
201 END DO
202!
203! Allocate and initialize modules variables.
204!
205 CALL roms_allocate_arrays (allocate_vars)
207 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
208
209 END IF
210
211#if defined MCT_LIB && (defined ATM_COUPLING || defined WAV_COUPLING)
212!
213!-----------------------------------------------------------------------
214! Initialize coupling streams between model(s).
215!-----------------------------------------------------------------------
216!
217 DO ng=1,ngrids
218# ifdef ATM_COUPLING
219 CALL initialize_ocn2atm_coupling (ng, myrank)
220# endif
221# ifdef WAV_COUPLING
222 CALL initialize_ocn2wav_coupling (ng, myrank)
223# endif
224 END DO
225#endif
226!
227!-----------------------------------------------------------------------
228! Set application grid, metrics, and associated variables. Then,
229! Proccess background and model prior error covariance standard
230! deviations and normalization coefficients.
231!-----------------------------------------------------------------------
232!
233 DO ng=1,ngrids
234#ifdef STD_MODEL
235 ldefstd(ng)=.true.
236 lwrtstd(ng)=.true.
237#else
238 lreadstd(ng)=.true.
239#endif
240 CALL prior_error (ng)
241 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
242 setgridconfig(ng)=.false.
243 END DO
244!
245 RETURN
246 END SUBROUTINE roms_initialize
247!
248 SUBROUTINE roms_run (RunInterval)
249!
250!=======================================================================
251! !
252! This subroutine runs the Strong or Weak constraint, Indirect !
253! Representers 4D-Var data assimilation (R4D-Var) algorithm. It !
254! time-steps ROMS nonlinear, representer, tangent linear, and !
255! adjoint kernels. !
256! !
257! On Input: !
258! !
259! RunInterval Execution time stepping window (seconds) !
260! !
261!=======================================================================
262!
263! Imported variable declarations
264!
265 real(dp), intent(in) :: RunInterval
266!
267! Local variable declarations.
268!
269 integer :: my_outer, ng
270!
271 character (len=*), parameter :: MyFile = &
272 & __FILE__//", ROMS_run"
273!
274!=======================================================================
275! Run R4D-Var Data Assimilation algorithm.
276!=======================================================================
277!
278! Initialize several global parameters.
279!
280 DO ng=1,ngrids
281#if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
282 lfinp(ng)=1 ! forcing index for input
283 lfout(ng)=1 ! forcing index for output history files
284#endif
285#ifdef ADJUST_BOUNDARY
286 lbinp(ng)=1 ! boundary index for input
287 lbout(ng)=1 ! boundary index for output history files
288#endif
289 lold(ng)=1 ! old minimization time index
290 lnew(ng)=2 ! new minimization time index
291 END DO
292!
293 outer=0 ! outer-loop counter
294 inner=0 ! inner-loop counter
295!
296! Compute nonlinear background state trajectory, Xb(t)|n-1. Interpolate
297! the background at the observation locations, and compute the quality
298! control accept/reject flag, ObsScale. The background state is used
299! to linearize the tangent linear and adjoint models during the
300! minimization.
301!
302 CALL background (outer, runinterval)
303 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
304!
305! Start outer loop iterations.
306!
307 outer_loop : DO my_outer=1,nouter
308 outer=my_outer
309 inner=0
310!
311! Compute 4D-Var data assimilation increment, dXa, by iterating over
312! the inner loops, and minimizing the cost function.
313!
314 CALL increment (my_outer, runinterval)
315 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
316!
317! Compute 4D-Var data assimilation analysis, Xa = Xb + dXa. Set
318! nonlinear model initial conditions for next outer loop.
319!
320 CALL analysis (my_outer, runinterval)
321 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
322
323 END DO outer_loop
324
325#if defined POSTERIOR_ERROR_I || defined POSTERIOR_ERROR_F || \
326 defined posterior_eofs
327!
328! Compute full (diagonal) posterior analysis error covariance matrix.
329! (NOTE: Currently, this code only works for a single outer-loop).
330!
331 CALL posterior_error (runinterval)
332 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
333#endif
334!
335 RETURN
336 END SUBROUTINE roms_run
337!
338 SUBROUTINE roms_finalize
339!
340!=======================================================================
341! !
342! This routine terminates ROMS R4D-Var execution. !
343! !
344!=======================================================================
345!
346! Local variable declarations.
347!
348 integer :: Fcount, InpRec, Nfiles, Tindex
349 integer :: ifile, lstr, ng, thread
350!
351 character (len=10) :: suffix
352
353 character (len=*), parameter :: MyFile = &
354 & __FILE__//", ROMS_finalize"
355!
356!-----------------------------------------------------------------------
357! Write out 4D-Var analysis fields that can be used as the initial
358! conditions for the next data assimilation cycle. Here, use the
359! last record of the RPM for the final outer loop.
360!-----------------------------------------------------------------------
361!
362 IF (exit_flag.eq.noerror) THEN
363 DO ng=1,ngrids
364 ldefdai(ng)=.true.
365 CALL def_dai (ng)
366 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
367!
368 WRITE (tlm(ng)%name,10) trim(fwd(ng)%head), nouter
369 10 FORMAT (a,'_outer',i0,'.nc')
370 lstr=len_trim(tlm(ng)%name)
371 tlm(ng)%base=tlm(ng)%name(1:lstr-3)
372 IF (tlm(ng)%Nfiles.gt.1) THEN
373 nfiles=tlm(ng)%Nfiles
374 DO ifile=1,nfiles
375 WRITE (suffix,"('_',i4.4,'.nc')") ifile
376 tlm(ng)%files(ifile)=trim(tlm(ng)%base)//trim(suffix)
377 END DO
378 tlm(ng)%name=trim(tlm(ng)%files(nfiles))
379 ELSE
380 tlm(ng)%files(1)=trim(tlm(ng)%name)
381 END IF
382!
383 SELECT CASE (tlm(ng)%IOtype)
384 CASE (io_nf90)
385 CALL netcdf_get_dim (ng, irpm, tlm(ng)%name, &
386 & dimname = 'ocean_time', &
387 & dimsize = inprec)
388
389#if defined PIO_LIB && defined DISTRIBUTE
390 CASE (io_pio)
391 CALL pio_netcdf_get_dim (ng, irpm, tlm(ng)%name, &
392 & dimname = 'ocean_time', &
393 & dimsize = inprec)
394#endif
395 END SELECT
396 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
397!
398 tindex=1
399 CALL get_state (ng, irpm, 1, tlm(ng), inprec, tindex)
400 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
401!
402 kout=tindex
403 nout=tindex
404#ifdef DISTRIBUTE
405 CALL wrt_dai (ng, myrank)
406#else
407 CALL wrt_dai (ng, -1)
408#endif
409 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
410 END DO
411 END IF
412!
413!-----------------------------------------------------------------------
414! Compute and report model-observation comparison statistics.
415!-----------------------------------------------------------------------
416!
417 DO ng=1,ngrids
418#ifdef DISTRIBUTE
419 CALL stats_modobs (ng, myrank)
420#else
421 CALL stats_modobs (ng, -1)
422#endif
423 END DO
424!
425!-----------------------------------------------------------------------
426! If blowing-up, save latest model state into RESTART NetCDF file.
427!-----------------------------------------------------------------------
428!
429! If cycling restart records, write solution into record 3.
430!
431 IF (exit_flag.eq.1) THEN
432 DO ng=1,ngrids
433 IF (lwrtrst(ng)) THEN
434 IF (master) WRITE (stdout,20)
435 20 FORMAT (/,' Blowing-up: Saving latest model state into ', &
436 & ' RESTART file',/)
437 fcount=rst(ng)%load
438 IF (lcyclerst(ng).and.(rst(ng)%Nrec(fcount).ge.2)) THEN
439 rst(ng)%Rindex=2
440 lcyclerst(ng)=.false.
441 END IF
444#ifdef DISTRIBUTE
445 CALL wrt_rst (ng, myrank)
446#else
447 CALL wrt_rst (ng, -1)
448#endif
449 END IF
450 END DO
451 END IF
452!
453!-----------------------------------------------------------------------
454! Stop model and time profiling clocks, report memory requirements,
455! and close output NetCDF files.
456!-----------------------------------------------------------------------
457!
458! Stop time clocks.
459!
460 IF (master) THEN
461 WRITE (stdout,30)
462 30 FORMAT (/,'Elapsed wall CPU time for each process (seconds):',/)
463 END IF
464!
465 DO ng=1,ngrids
466 DO thread=thread_range
467 CALL wclock_off (ng, inlm, 0, __line__, myfile)
468 END DO
469 END DO
470!
471! Report dynamic memory and automatic memory requirements.
472!
473 CALL memory
474!
475! Close IO files.
476!
477 DO ng=1,ngrids
478 CALL close_inp (ng, inlm)
479 END DO
480 CALL close_out
481!
482 RETURN
483 END SUBROUTINE roms_finalize
484
485 END MODULE roms_kernel_mod
subroutine memory
Definition memory.F:3
subroutine, public close_out
Definition close_io.F:175
subroutine, public close_inp(ng, model)
Definition close_io.F:92
subroutine, public def_dai(ng)
Definition def_dai.F:52
subroutine, public get_state(ng, model, msg, s, inirec, tindex)
Definition get_state.F:90
subroutine, public inp_par(model)
Definition inp_par.F:56
subroutine, public roms_initialize_arrays
Definition mod_arrays.F:351
subroutine, public roms_allocate_arrays(allocate_vars)
Definition mod_arrays.F:114
type(t_io), dimension(:), allocatable tlm
type(t_io), dimension(:), allocatable fwd
type(t_io), dimension(:), allocatable rst
integer stdout
integer, parameter io_nf90
Definition mod_ncparam.F:95
integer, parameter io_pio
Definition mod_ncparam.F:96
subroutine, public netcdf_get_dim(ng, model, ncname, ncid, dimname, dimsize, dimid)
Definition mod_netcdf.F:330
subroutine, public initialize_parallel
integer numthreads
integer, dimension(:), allocatable first_tile
integer mythread
logical master
integer, dimension(:), allocatable last_tile
integer ocn_comm_world
integer, parameter inlm
Definition mod_param.F:662
integer, parameter irpm
Definition mod_param.F:664
integer, dimension(:), allocatable ntilex
Definition mod_param.F:685
integer ngrids
Definition mod_param.F:113
integer, dimension(:), allocatable ntilee
Definition mod_param.F:686
subroutine, public pio_netcdf_get_dim(ng, model, ncname, piofile, dimname, dimsize, dimid)
integer nouter
logical, dimension(:), allocatable lreadstd
logical, dimension(:), allocatable setgridconfig
integer blowup
integer erend
integer exit_flag
integer erstr
logical, dimension(:), allocatable lwrtrst
integer nrun
logical, dimension(:), allocatable ldefdai
integer inner
integer noerror
logical, dimension(:), allocatable lcyclerst
integer outer
integer, dimension(:), allocatable lold
integer, dimension(:), allocatable lbout
integer, dimension(:), allocatable lfinp
integer, dimension(:), allocatable lbinp
integer, dimension(:), allocatable lnew
integer, dimension(:), allocatable lfout
subroutine, public posterior_error(runinterval)
Definition r4dvar.F:1946
subroutine, public increment(my_outer, runinterval)
Definition r4dvar.F:380
subroutine, public analysis(my_outer, runinterval)
Definition r4dvar.F:1384
subroutine, public prior_error(ng)
Definition r4dvar.F:1782
subroutine, public background(my_outer, runinterval)
Definition r4dvar.F:158
subroutine, public roms_finalize
Definition ad_roms.h:283
subroutine, public roms_run(runinterval)
Definition ad_roms.h:239
subroutine, public roms_initialize(first, mpicomm)
Definition ad_roms.h:52
subroutine, public stats_modobs(ng, tile)
integer function, public stdout_unit(mymaster)
Definition stdout_mod.F:48
logical, save set_stdoutunit
Definition stdout_mod.F:41
logical function, public founderror(flag, noerr, line, routine)
Definition strings.F:52
subroutine, public wrt_dai(ng, tile)
Definition wrt_dai.F:46
subroutine, public wrt_rst(ng, tile)
Definition wrt_rst.F:63
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