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

Functions/Subroutines

subroutine, public wrt_aug_imp (ng, tile, model, iinp, iout)
 
subroutine, private wrt_aug_imp_nf90 (ng, tile, model, iinp, iout, lbi, ubi, lbj, ubj)
 
subroutine, private wrt_aug_imp_pio (ng, tile, model, iinp, iout, lbi, ubi, lbj, ubj)
 

Function/Subroutine Documentation

◆ wrt_aug_imp()

subroutine, public wrt_aug_imp_mod::wrt_aug_imp ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) model,
integer, intent(in) iinp,
integer, intent(in) iout )

Definition at line 42 of file wrt_aug_imp.F.

43!***********************************************************************
44!
45! Imported variable declarations.
46!
47 integer, intent(in) :: ng, tile, model, Iinp, Iout
48!
49! Local variable declarations.
50!
51 integer :: LBi, UBi, LBj, UBj
52!
53 character (len=*), parameter :: MyFile = &
54 & __FILE__
55!
56!-----------------------------------------------------------------------
57! Write out time-averaged fields according to IO type.
58!-----------------------------------------------------------------------
59!
60 lbi=bounds(ng)%LBi(tile)
61 ubi=bounds(ng)%UBi(tile)
62 lbj=bounds(ng)%LBj(tile)
63 ubj=bounds(ng)%UBj(tile)
64!
65 SELECT CASE (tlf(ng)%IOtype)
66 CASE (io_nf90)
67 CALL wrt_aug_imp_nf90 (ng, tile, model, iinp, iout, &
68 & lbi, ubi, lbj, ubj)
69
70# if defined PIO_LIB && defined DISTRIBUTE
71 CASE (io_pio)
72 CALL wrt_aug_imp_pio (ng, tile, model, iinp, iout, &
73 & lbi, ubi, lbj, ubj)
74# endif
75 CASE DEFAULT
76 IF (master) WRITE (stdout,10) tlf(ng)%IOtype
77 exit_flag=3
78 END SELECT
79 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
80!
81 10 FORMAT (' WRT_AUG_IMP - Illegal output tile type, io_type = ',i0, &
82 & /,15x,'Check KeyWord ''OUT_LIB'' in ''roms.in''.' )
83!
84 RETURN

References mod_param::bounds, mod_scalars::exit_flag, strings_mod::founderror(), mod_ncparam::io_nf90, mod_ncparam::io_pio, mod_parallel::master, mod_scalars::noerror, mod_iounits::stdout, mod_iounits::tlf, wrt_aug_imp_nf90(), and wrt_aug_imp_pio().

Referenced by rbl4dvar_mod::increment().

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

◆ wrt_aug_imp_nf90()

subroutine, private wrt_aug_imp_mod::wrt_aug_imp_nf90 ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) model,
integer, intent(in) iinp,
integer, intent(in) iout,
integer, intent(in) lbi,
integer, intent(in) ubi,
integer, intent(in) lbj,
integer, intent(in) ubj )
private

Definition at line 88 of file wrt_aug_imp.F.

90!***********************************************************************
91!
92 USE mod_netcdf
93!
94! Imported variable declarations.
95!
96 integer, intent(in) :: ng, tile, model, Iinp, Iout
97 integer, intent(in) :: LBi, UBi, LBj, UBj
98!
99! Local variable declarations.
100!
101 integer :: itrc, gfactor, gtype, status
102!
103 real(dp) :: scale
104!
105 character (len=*), parameter :: MyFile = &
106 & __FILE__//", wrt_aug_imp_nf90"
107
108# include "set_bounds.h"
109!
110 sourcefile=myfile
111!
112!-----------------------------------------------------------------------
113! Determine variables to read and their availability.
114!-----------------------------------------------------------------------
115!
116 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
117!
118! Report.
119!
120 IF (master) WRITE (stdout,10) iout, trim(tlf(ng)%name)
121!
122! Set grid type factor to write full (gfactor=1) fields or water
123! points (gfactor=-1) fields only.
124!
125# if defined WRITE_WATER && defined MASKING
126 gfactor=-1
127# else
128 gfactor=1
129# endif
130!
131! Process free-surface weak-constraint impulse forcing.
132!
133 scale=1.0_dp
134 gtype=gfactor*r2dvar
135 status=nf_fwrite2d(ng, model, tlf(ng)%ncid, idztlf, &
136 & tlf(ng)%Vid(idztlf), &
137 & iout, gtype, &
138 & lbi, ubi, lbj, ubj, scale, &
139# ifdef MASKING
140 & grid(ng) % rmask, &
141# endif
142 & ocean(ng) % tl_zeta(:,:,iinp))
143 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
144 IF (master) THEN
145 WRITE (stdout,20) trim(vname(1,idztlf)), iout, &
146 & trim(tlf(ng)%name)
147 END IF
148 exit_flag=3
149 ioerror=status
150 RETURN
151 END IF
152
153# ifndef SOLVE3D
154!
155! Process 2D U-momentum weak-constraint impulse forcing.
156!
157 scale=1.0_dp
158 gtype=gfactor*u2dvar
159 status=nf_fwrite2d(ng, model, tlf(ng)%ncid, idubtf, &
160 & tlf(ng)%Vid(idubtf), &
161 & iout, gype, &
162 & lbi, ubi, lbj, ubj, scale, &
163# ifdef MASKING
164 & grid(ng) % umask_full, &
165# endif
166 & ocean(ng) % tl_ubar(:,:,iinp))
167 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
168 IF (master) THEN
169 WRITE (stdout,20) trim(vname(1,idubtf)), iout, &
170 & trim(tlf(ng)%name)
171 END IF
172 exit_flag=3
173 ioerror=status
174 RETURN
175 END IF
176!
177! Process 2D V-momentum weak-constraint impulse forcing.
178!
179 scale=1.0_dp
180 gtype=gfactor*v2dvar
181 status=nf_fwrite2d(ng, model, tlf(ng)%ncid, idvbtf, &
182 & tlf(ng)%Vid(idvbtf), &
183 & iout, gtype, &
184 & lbi, ubi, lbj, ubj, scale, &
185# ifdef MASKING
186 & grid(ng) % vmask_full, &
187# endif
188 & ocean(ng) % tl_vbar(:,:,iinp))
189 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
190 IF (master) THEN
191 WRITE (stdout,20) trim(vname(1,idvbtf)), iout, &
192 & trim(tlf(ng)%name)
193 END IF
194 exit_flag=3
195 ioerror=status
196 RETURN
197 END IF
198# endif
199# ifdef SOLVE3D
200!
201! Process 3D U-momentum weak-constraint impulse forcing.
202!
203 scale=1.0_dp
204 gtype=gfactor*u3dvar
205 status=nf_fwrite3d(ng, model, tlf(ng)%ncid, idutlf, &
206 & tlf(ng)%Vid(idutlf), &
207 & iout, gtype, &
208 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
209# ifdef MASKING
210 & grid(ng) % umask_full, &
211# endif
212 & ocean(ng) % tl_u(:,:,:,iinp))
213 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
214 IF (master) THEN
215 WRITE (stdout,20) trim(vname(1,idutlf)), iout, &
216 & trim(tlf(ng)%name)
217 END IF
218 exit_flag=3
219 ioerror=status
220 RETURN
221 END IF
222!
223! Process 3D V-momentum weak-constraint impulse forcing.
224!
225 scale=1.0_dp
226 gtype=gfactor*v3dvar
227 status=nf_fwrite3d(ng, model, tlf(ng)%ncid, idvtlf, &
228 & tlf(ng)%Vid(idvtlf), &
229 & iout, gtype, &
230 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
231# ifdef MASKING
232 & grid(ng) % vmask_full, &
233# endif
234 & ocean(ng) % tl_v(:,:,:,iinp))
235 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
236 IF (master) THEN
237 WRITE (stdout,20) trim(vname(1,idvtlf)), iout, &
238 & trim(tlf(ng)%name)
239 END IF
240 exit_flag=3
241 ioerror=status
242 RETURN
243 END IF
244!
245! Process tracer type variables impulses.
246!
247 scale=1.0_dp
248 gtype=gfactor*r3dvar
249 DO itrc=1,nt(ng)
250 status=nf_fwrite3d(ng, model, tlf(ng)%ncid, idttlf(itrc), &
251 & tlf(ng)%Tid(itrc), &
252 & iout, gtype, &
253 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
254# ifdef MASKING
255 & grid(ng) % rmask_full, &
256# endif
257 & ocean(ng) % tl_t(:,:,:,iinp,itrc))
258 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
259 IF (master) THEN
260 WRITE (stdout,20) trim(vname(1,idttlf(itrc))), iout, &
261 & trim(tlf(ng)%name)
262 END IF
263 exit_flag=3
264 ioerror=status
265 RETURN
266 END IF
267 END DO
268# endif
269!
270!-----------------------------------------------------------------------
271! Synchronize impulse NetCDF file to disk to allow other processes
272! to access data immediately after it is written.
273!-----------------------------------------------------------------------
274!
275 CALL netcdf_sync (ng, model, tlf(ng)%name, tlf(ng)%ncid)
276 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
277!
278 10 FORMAT (2x,'WRT_AUG_IMP_NF90 - writing augmented adjoint', &
279 & ' impulses, record: ',i0,/,22x,'file: ',a)
280 20 FORMAT (/,' WRT_AUG_IMP_NF90 - error while writing variable: ',a, &
281 & 2x,'at time record = ',i0, &
282 & /,20x,'into NetCDF file: ',a)
283!
284 RETURN
subroutine, public netcdf_sync(ng, model, ncname, ncid)

References mod_scalars::exit_flag, strings_mod::founderror(), mod_grid::grid, mod_ncparam::idttlf, mod_ncparam::idubtf, mod_ncparam::idutlf, mod_ncparam::idvbtf, mod_ncparam::idvtlf, mod_ncparam::idztlf, mod_iounits::ioerror, mod_parallel::master, mod_param::n, mod_netcdf::netcdf_sync(), mod_scalars::noerror, mod_param::nt, mod_ocean::ocean, mod_param::r2dvar, mod_param::r3dvar, mod_iounits::sourcefile, mod_iounits::stdout, mod_iounits::tlf, mod_param::u2dvar, mod_param::u3dvar, mod_param::v2dvar, mod_param::v3dvar, and mod_ncparam::vname.

Referenced by wrt_aug_imp().

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

◆ wrt_aug_imp_pio()

subroutine, private wrt_aug_imp_mod::wrt_aug_imp_pio ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) model,
integer, intent(in) iinp,
integer, intent(in) iout,
integer, intent(in) lbi,
integer, intent(in) ubi,
integer, intent(in) lbj,
integer, intent(in) ubj )
private

Definition at line 290 of file wrt_aug_imp.F.

292!***********************************************************************
293!
295!
296! Imported variable declarations.
297!
298 integer, intent(in) :: ng, tile, model, Iinp, Iout
299 integer, intent(in) :: LBi, UBi, LBj, UBj
300!
301! Local variable declarations.
302!
303 integer :: itrc, status
304!
305 real(dp) :: scale
306!
307 character (len=*), parameter :: MyFile = &
308 & __FILE__", wrt_aug_imp_pio"
309!
310 TYPE (IO_desc_t), pointer :: ioDesc
311
312# include "set_bounds.h"
313!
314 sourcefile=myfile
315!
316!-----------------------------------------------------------------------
317! Determine variables to read and their availability.
318!-----------------------------------------------------------------------
319!
320 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
321!
322! Report.
323!
324 IF (master) WRITE (stdout,20) iout, trim(tlf(ng)%name)
325!
326! Process free-surface weak-constraint impulse forcing.
327!
328 scale=1.0_dp
329 IF (tlf(ng)%pioVar(idztlf)%dkind.eq.pio_double) THEN
330 iodesc => iodesc_dp_r2dvar(ng)
331 ELSE
332 iodesc => iodesc_sp_r2dvar(ng)
333 END IF
334!
335 status=nf_fwrite2d(ng, model, tlf(ng)%pioFile, idztlf, &
336 & tlf(ng)%pioVar(idztlf), iout, &
337 & iodesc, &
338 & lbi, ubi, lbj, ubj, scale, &
339# ifdef MASKING
340 & grid(ng) % rmask, &
341# endif
342 & ocean(ng) % tl_zeta(:,:,iinp))
343 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
344 IF (master) THEN
345 WRITE (stdout,20) trim(vname(1,idztlf)), iout, &
346 & trim(tlf(ng)%name)
347 END IF
348 exit_flag=3
349 ioerror=status
350 RETURN
351 END IF
352
353# ifndef SOLVE3D
354!
355! Process 2D U-momentum weak-constraint impulse forcing.
356!
357 scale=1.0_dp
358 IF (tlf(ng)%pioVar(idubtf)%dkind.eq.pio_double) THEN
359 iodesc => iodesc_dp_u2dvar(ng)
360 ELSE
361 iodesc => iodesc_sp_u2dvar(ng)
362 END IF
363!
364 status=nf_fwrite2d(ng, model, tlf(ng)%pioFile, idubtf, &
365 & tlf(ng)%pioVar(idubtf), iout, &
366 & iodesc, &
367 & lbi, ubi, lbj, ubj, scale, &
368# ifdef MASKING
369 & grid(ng) % umask_full, &
370# endif
371 & ocean(ng) % tl_ubar(:,:,iinp))
372 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
373 IF (master) THEN
374 WRITE (stdout,20) trim(vname(1,idubtf)), iout, &
375 & trim(tlf(ng)%name)
376 END IF
377 exit_flag=3
378 ioerror=status
379 RETURN
380 END IF
381!
382! Process 2D V-momentum weak-constraint impulse forcing.
383!
384 scale=1.0_dp
385 IF (tlf(ng)%pioVar(idvbtf)%dkind.eq.pio_double) THEN
386 iodesc => iodesc_dp_v2dvar(ng)
387 ELSE
388 iodesc => iodesc_sp_v2dvar(ng)
389 END IF
390!
391 status=nf_fwrite2d(ng, model, tlf(ng)%pioFile, idvbtf, &
392 & tlf(ng)%pioVar(idvbtf), iout, &
393 & iodesc, &
394 & lbi, ubi, lbj, ubj, scale, &
395# ifdef MASKING
396 & grid(ng) % vmask_full, &
397# endif
398 & ocean(ng) % tl_vbar(:,:,iinp))
399 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
400 IF (master) THEN
401 WRITE (stdout,20) trim(vname(1,idvbtf)), iout, &
402 & trim(tlf(ng)%name)
403 END IF
404 exit_flag=3
405 ioerror=status
406 RETURN
407 END IF
408# endif
409# ifdef SOLVE3D
410!
411! Process 3D U-momentum weak-constraint impulse forcing.
412!
413 scale=1.0_dp
414 IF (tlf(ng)%pioVar(idutlf)%dkind.eq.pio_double) THEN
415 iodesc => iodesc_dp_u3dvar(ng)
416 ELSE
417 iodesc => iodesc_sp_u3dvar(ng)
418 END IF
419!
420 status=nf_fwrite3d(ng, model, tlf(ng)%pioFile, idutlf, &
421 & tlf(ng)%pioVar(idutlf), iout, &
422 & iodesc, &
423 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
424# ifdef MASKING
425 & grid(ng) % umask_full, &
426# endif
427 & ocean(ng) % tl_u(:,:,:,iinp))
428 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
429 IF (master) THEN
430 WRITE (stdout,20) trim(vname(1,idutlf)), iout, &
431 & trim(tlf(ng)%name)
432 END IF
433 exit_flag=3
434 ioerror=status
435 RETURN
436 END IF
437!
438! Process 3D V-momentum weak-constraint impulse forcing.
439!
440 scale=1.0_dp
441 IF (tlf(ng)%pioVar(idvtlf)%dkind.eq.pio_double) THEN
442 iodesc => iodesc_dp_v3dvar(ng)
443 ELSE
444 iodesc => iodesc_sp_v3dvar(ng)
445 END IF
446!
447 status=nf_fwrite3d(ng, model, tlf(ng)%pioFile, idvtlf, &
448 & tlf(ng)%pioVar(idvtlf), iout, &
449 & iodesc, &
450 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
451# ifdef MASKING
452 & grid(ng) % vmask_full, &
453# endif
454 & ocean(ng) % tl_v(:,:,:,iinp))
455 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
456 IF (master) THEN
457 WRITE (stdout,20) trim(vname(1,idvtlf)), iout, &
458 & trim(tlf(ng)%name)
459 END IF
460 exit_flag=3
461 ioerror=status
462 RETURN
463 END IF
464!
465! Process tracer type variables impulses.
466!
467 DO itrc=1,nt(ng)
468 scale=1.0_dp
469 IF (tlf(ng)%pioTrc(itrc)%dkind.eq.pio_double) THEN
470 iodesc => iodesc_dp_r3dvar(ng)
471 ELSE
472 iodesc => iodesc_sp_r3dvar(ng)
473 END IF
474!
475 status=nf_fwrite3d(ng, model, tlf(ng)%pioFile, idttlf(itrc), &
476 & tlf(ng)%pioTrc(itrc), iout, &
477 & iodesc, &
478 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
479# ifdef MASKING
480 & grid(ng) % rmask_full, &
481# endif
482 & ocean(ng) % tl_t(:,:,:,iinp,itrc))
483 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
484 IF (master) THEN
485 WRITE (stdout,20) trim(vname(1,idttlf(itrc))), iout, &
486 & trim(tlf(ng)%name)
487 END IF
488 exit_flag=3
489 ioerror=status
490 RETURN
491 END IF
492 END DO
493# endif
494!
495!-----------------------------------------------------------------------
496! Synchronize impulse NetCDF file to disk to allow other processes
497! to access data immediately after it is written.
498!-----------------------------------------------------------------------
499!
500 CALL pio_netcdf_sync (ng, model, tlf(ng)%name, tlf(ng)%pioFile)
501 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
502!
503 10 FORMAT (2x,'WRT_AUG_IMP_PIO - writing augmented adjoint', &
504 & ' impulses, record: ',i0,/,22x,'file: ',a)
505 20 FORMAT (/,' WRT_AUG_IMP_PIO - error while writing variable: ',a, &
506 & 2x,'at time record = ',i0, &
507 & /,20x,'into NetCDF file: ',a)
508!
509 RETURN
type(io_desc_t), dimension(:), pointer iodesc_dp_r3dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_u3dvar
subroutine, public pio_netcdf_sync(ng, model, ncname, piofile)
type(io_desc_t), dimension(:), pointer iodesc_sp_u2dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_u3dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_r2dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_u2dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_v3dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_v2dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_v3dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_r3dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_r2dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_v2dvar

References mod_scalars::exit_flag, strings_mod::founderror(), mod_grid::grid, mod_ncparam::idttlf, mod_ncparam::idubtf, mod_ncparam::idutlf, mod_ncparam::idvbtf, mod_ncparam::idvtlf, mod_ncparam::idztlf, mod_pio_netcdf::iodesc_dp_r2dvar, mod_pio_netcdf::iodesc_dp_r3dvar, mod_pio_netcdf::iodesc_dp_u2dvar, mod_pio_netcdf::iodesc_dp_u3dvar, mod_pio_netcdf::iodesc_dp_v2dvar, mod_pio_netcdf::iodesc_dp_v3dvar, mod_pio_netcdf::iodesc_sp_r2dvar, mod_pio_netcdf::iodesc_sp_r3dvar, mod_pio_netcdf::iodesc_sp_u2dvar, mod_pio_netcdf::iodesc_sp_u3dvar, mod_pio_netcdf::iodesc_sp_v2dvar, mod_pio_netcdf::iodesc_sp_v3dvar, mod_iounits::ioerror, mod_parallel::master, mod_param::n, mod_scalars::noerror, mod_param::nt, mod_ocean::ocean, mod_pio_netcdf::pio_netcdf_sync(), mod_iounits::sourcefile, mod_iounits::stdout, mod_iounits::tlf, and mod_ncparam::vname.

Referenced by wrt_aug_imp().

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