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

Data Types

type  t_sources
 

Functions/Subroutines

subroutine, public allocate_sources (ng)
 
subroutine, public deallocate_sources (ng)
 
subroutine check_sources (ng, ncname, npsrc)
 

Variables

type(t_sources), dimension(:), allocatable sources
 
integer, dimension(:), allocatable msrc
 
integer, dimension(:), allocatable nsrc
 

Function/Subroutine Documentation

◆ allocate_sources()

subroutine, public mod_sources::allocate_sources ( integer ng)

Definition at line 101 of file mod_sources.F.

102!
103!=======================================================================
104! !
105! This routine allocates and initializes all variables in the module !
106! for all nested grids. !
107! !
108!=======================================================================
109!
110 USE mod_param
111#ifndef ANA_PSOURCE
112 USE mod_parallel
113 USE mod_iounits
114 USE mod_ncparam
115 USE mod_netcdf
116# if defined PIO_LIB && defined DISTRIBUTE
118# endif
119 USE mod_scalars
120#endif
121!
122 USE strings_mod, ONLY : founderror
123!
124! Imported variable declarations.
125!
126 integer :: ng
127!
128! Local variable declarations.
129!
130#ifndef ANA_PSOURCE
131 logical :: foundit
132!
133 integer :: Vid, ifile, nvatt, nvdim
134#endif
135 integer :: is, itrc, k, mg
136
137 real(r8), parameter :: IniVal = 0.0_r8
138!
139 character (len=*), parameter :: MyFile = &
140 & __FILE__//", allocate_sources"
141
142#ifndef ANA_PSOURCE
143# if defined PIO_LIB && defined DISTRIBUTE
144!
145 TYPE (Var_desc_t) :: my_pioVar
146# endif
147#endif
148!
149!-----------------------------------------------------------------------
150! Allocate module variables.
151!-----------------------------------------------------------------------
152!
153 IF (.not.allocated(msrc)) THEN
154 allocate ( msrc(ngrids) )
155 END IF
156
157 IF (.not.allocated(nsrc)) THEN
158 allocate ( nsrc(ngrids) )
159 END IF
160
161#ifndef ANA_PSOURCE
162!
163! Inquire about the number of point Sources/Sinks.
164!
165 IF (ng.eq.1) THEN
166 DO mg=1,ngrids
167 foundit=.false.
168 IF (luvsrc(mg).or.lwsrc(mg).or.any(ltracersrc(:,mg))) THEN
169 SELECT CASE (ssf(ng)%IOtype)
170 CASE (io_nf90)
171 CALL netcdf_inq_var (ng, inlm, ssf(mg)%name, &
172 & myvarname = vname(1,idrxpo), &
173 & searchvar = foundit, &
174 & varid = vid, &
175 & nvardim = nvdim, &
176 & nvaratt = nvatt)
177
178# if defined PIO_LIB && defined DISTRIBUTE
179 CASE (io_pio)
180 CALL pio_netcdf_inq_var (ng, inlm, ssf(mg)%name, &
181 & myvarname = vname(1,idrxpo), &
182 & searchvar = foundit, &
183 & piovar = my_piovar, &
184 & nvardim = nvdim, &
185 & nvaratt = nvatt)
186# endif
187 END SELECT
188 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
189 IF (foundit) THEN
190 nsrc(mg)=var_dsize(1) ! first dimension
191 msrc(mg)=nsrc(mg)
192 CALL check_sources (mg, ssf(mg)%name, nsrc(mg))
193 IF (founderror(exit_flag, noerror, __line__, &
194 & myfile)) RETURN
195
196 END IF
197 END IF
198 END DO
199 END IF
200#else
201!
202! Set number of point sources to maximum number of analytical sources.
203! Notice that a maximum of 200 analytical sources are set-up here.
204!
205 msrc(ng)=200
206 nsrc(ng)=msrc(ng)
207#endif
208!
209! Allocate structure.
210!
211 IF (ng.eq.1) allocate ( sources(ngrids) )
212!
213! Allocate point Sources/Sinks variables.
214!
215 allocate ( sources(ng) % Isrc(nsrc(ng)) )
216 dmem(ng)=dmem(ng)+real(nsrc(ng),r8)
217
218 allocate ( sources(ng) % Jsrc(nsrc(ng)) )
219 dmem(ng)=dmem(ng)+real(nsrc(ng),r8)
220
221 allocate ( sources(ng) % Dsrc(nsrc(ng)) )
222 dmem(ng)=dmem(ng)+real(nsrc(ng),r8)
223
224 allocate ( sources(ng) % Fsrc(nsrc(ng)) )
225 dmem(ng)=dmem(ng)+real(nsrc(ng),r8)
226
227 allocate ( sources(ng) % Qbar(nsrc(ng)) )
228 dmem(ng)=dmem(ng)+real(nsrc(ng),r8)
229
230 allocate ( sources(ng) % Qshape(nsrc(ng),n(ng)) )
231 dmem(ng)=dmem(ng)+real(nsrc(ng)*n(ng),r8)
232
233 allocate ( sources(ng) % Qsrc(nsrc(ng),n(ng)) )
234 dmem(ng)=dmem(ng)+real(nsrc(ng)*n(ng),r8)
235
236 allocate ( sources(ng) % Tsrc(nsrc(ng),n(ng),nt(ng)) )
237 dmem(ng)=dmem(ng)+real(nsrc(ng)*n(ng)*nt(ng),r8)
238
239 allocate ( sources(ng) % Xsrc(nsrc(ng)) )
240 dmem(ng)=dmem(ng)+real(nsrc(ng),r8)
241
242 allocate ( sources(ng) % Ysrc(nsrc(ng)) )
243 dmem(ng)=dmem(ng)+real(nsrc(ng),r8)
244
245#ifndef ANA_PSOURCE
246 allocate ( sources(ng) % QbarG(nsrc(ng),2) )
247 dmem(ng)=dmem(ng)+2.0_r8*real(nsrc(ng),r8)
248
249 allocate ( sources(ng) % TsrcG(nsrc(ng),n(ng),2,nt(ng)) )
250 dmem(ng)=dmem(ng)+2.0_r8*real(nsrc(ng)*n(ng)*nt(ng),r8)
251#endif
252
253#ifdef ADJOINT
254 allocate ( sources(ng) % ad_Qbar(nsrc(ng)) )
255 dmem(ng)=dmem(ng)+real(nsrc(ng),r8)
256
257 allocate ( sources(ng) % ad_Qsrc(nsrc(ng),n(ng)) )
258 dmem(ng)=dmem(ng)+real(nsrc(ng)*n(ng),r8)
259
260 allocate ( sources(ng) % ad_Tsrc(nsrc(ng),n(ng),nt(ng)) )
261 dmem(ng)=dmem(ng)+real(nsrc(ng)*n(ng)*nt(ng),r8)
262#endif
263
264#ifdef TANGENT
265 allocate ( sources(ng) % tl_Qbar(nsrc(ng)) )
266 dmem(ng)=dmem(ng)+real(nsrc(ng),r8)
267
268 allocate ( sources(ng) % tl_Qsrc(nsrc(ng),n(ng)) )
269 dmem(ng)=dmem(ng)+real(nsrc(ng)*n(ng),r8)
270
271 allocate ( sources(ng) % tl_Tsrc(nsrc(ng),n(ng),nt(ng)) )
272 dmem(ng)=dmem(ng)+real(nsrc(ng)*n(ng)*nt(ng),r8)
273#endif
274
275!
276!-----------------------------------------------------------------------
277! Initialize module variables.
278!-----------------------------------------------------------------------
279!
280 DO is=1,nsrc(ng)
281 sources(ng) % Isrc(is) = 0
282 sources(ng) % Jsrc(is) = 0
283 sources(ng) % Dsrc(is) = inival
284 sources(ng) % Fsrc(is) = inival
285 sources(ng) % Xsrc(is) = inival
286 sources(ng) % Ysrc(is) = inival
287 sources(ng) % Qbar(is) = inival
288#ifndef ANA_PSOURCE
289 sources(ng) % QbarG(is,1) = inival
290 sources(ng) % QbarG(is,2) = inival
291#endif
292#ifdef ADJOINT
293 sources(ng) % ad_Qbar(is) = inival
294#endif
295#ifdef TANGENT
296 sources(ng) % tl_Qbar(is) = inival
297#endif
298 END DO
299 DO k=1,n(ng)
300 DO is=1,nsrc(ng)
301 sources(ng) % Qshape(is,k) = inival
302 sources(ng) % Qsrc(is,k) = inival
303#ifdef ADJOINT
304 sources(ng) % ad_Qsrc(is,k) = inival
305#endif
306#ifdef TANGENT
307 sources(ng) % tl_Qsrc(is,k) = inival
308#endif
309 END DO
310 END DO
311 DO itrc=1,nt(ng)
312 DO k=1,n(ng)
313 DO is=1,nsrc(ng)
314 sources(ng) % Tsrc(is,k,itrc) = inival
315#ifdef ADJOINT
316 sources(ng) % ad_Tsrc(is,k,itrc) = inival
317#endif
318#ifdef TANGENT
319 sources(ng) % tl_Tsrc(is,k,itrc) = inival
320#endif
321#ifndef ANA_PSOURCE
322 sources(ng) % TsrcG(is,k,1,itrc) = inival
323 sources(ng) % TsrcG(is,k,2,itrc) = inival
324#endif
325 END DO
326 END DO
327 END DO
328!
329 RETURN
type(t_io), dimension(:), allocatable ssf
integer, parameter io_nf90
Definition mod_ncparam.F:95
integer, parameter io_pio
Definition mod_ncparam.F:96
character(len=maxlen), dimension(6, 0:nv) vname
integer idrxpo
integer, dimension(nvard) var_dsize
Definition mod_netcdf.F:177
subroutine, public netcdf_inq_var(ng, model, ncname, ncid, myvarname, searchvar, varid, nvardim, nvaratt)
integer, parameter inlm
Definition mod_param.F:662
integer, dimension(:), allocatable n
Definition mod_param.F:479
real(r8), dimension(:), allocatable dmem
Definition mod_param.F:137
integer ngrids
Definition mod_param.F:113
integer, dimension(:), allocatable nt
Definition mod_param.F:489
subroutine, public pio_netcdf_inq_var(ng, model, ncname, piofile, myvarname, searchvar, piovar, nvardim, nvaratt)
logical, dimension(:), allocatable luvsrc
logical, dimension(:,:), allocatable ltracersrc
logical, dimension(:), allocatable lwsrc
integer exit_flag
integer noerror
logical function, public founderror(flag, noerr, line, routine)
Definition strings.F:52

References check_sources(), mod_param::dmem, mod_scalars::exit_flag, strings_mod::founderror(), mod_ncparam::idrxpo, mod_param::inlm, mod_ncparam::io_nf90, mod_ncparam::io_pio, mod_scalars::ltracersrc, mod_scalars::luvsrc, mod_scalars::lwsrc, msrc, mod_param::n, mod_netcdf::netcdf_inq_var(), mod_param::ngrids, mod_scalars::noerror, nsrc, mod_param::nt, mod_pio_netcdf::pio_netcdf_inq_var(), mod_kinds::r8, sources, mod_iounits::ssf, mod_netcdf::var_dsize, and mod_ncparam::vname.

Referenced by mod_arrays::roms_allocate_arrays().

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

◆ check_sources()

subroutine mod_sources::check_sources ( integer, intent(in) ng,
character (len=*), intent(in) ncname,
integer, intent(in) npsrc )

Definition at line 446 of file mod_sources.F.

447!
448!=======================================================================
449! !
450! It checks input NetCDF file for correctness in the spefication of !
451! the point Source/Sink grid-cell face flag values (0, 1, 2) and the !
452! implementation methodology using "LuvSrc" and or "LwSrc". !
453! !
454!=======================================================================
455!
456 USE mod_param
457 USE mod_parallel
458 USE mod_iounits
459 USE mod_netcdf
460#if defined PIO_LIB && defined DISTRIBUTE
462#endif
463 USE mod_scalars
464!
465 USE strings_mod, ONLY : founderror
466
467! Imported variable declarations.
468!
469 integer, intent(in) :: ng, Npsrc
470!
471 character (len=*), intent(in) :: ncname
472!
473! Local variable declarations.
474!
475 integer :: i, ic_u, ic_v, ic_w
476!
477 real(r8) :: Dsrc_min, Dsrc_max
478!
479 real(r8) :: Dsrc(Npsrc)
480!
481 character (len=*), parameter :: MyFile = &
482 & __FILE__//", check_sources"
483!
484!-----------------------------------------------------------------------
485! Read in Point Source/Sink cell-face flag and check for correct values
486! according to the implementation methodology:
487!
488! If LuvSrc = T, needs Dsrc = 0 (flow across grid cell u-face) or
489! Dsrc = 1 (flow across grid cell v-face)
490!
491! If LwSrc = T, needs Dsrc = 2 (flow across grid cell w-face)
492!-----------------------------------------------------------------------
493!
494! Read in Source/Sink cell-face flag.
495!
496 SELECT CASE (ssf(ng)%IOtype)
497 CASE (io_nf90)
498 CALL netcdf_get_fvar (ng, inlm, ncname, &
499 & vname(1,idrdir), dsrc, &
500 & min_val = dsrc_min, &
501 & max_val = dsrc_max)
502
503# if defined PIO_LIB && defined DISTRIBUTE
504 CASE (io_pio)
505 CALL pio_netcdf_get_fvar (ng, inlm, ncname, &
506 & vname(1,idrdir), dsrc, &
507 & min_val = dsrc_min, &
508 & max_val = dsrc_max)
509# endif
510 END SELECT
511 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
512!
513! Count and report Source/Sink for each cell-face flag value.
514!
515 ic_u=0
516 ic_v=0
517 ic_w=0
518!
519 DO i=1,npsrc
520 IF (int(dsrc(i)).eq.0) ic_u=ic_u+1
521 IF (int(dsrc(i)).eq.1) ic_v=ic_v+1
522 IF (int(dsrc(i)).eq.2) ic_w=ic_w+1
523 END DO
524!
525 IF (master) THEN
526 IF (ng.eq.1) WRITE (stdout,10)
527 WRITE (stdout,20) ng, trim(ncname), &
528 luvsrc(ng), ic_u, &
529 luvsrc(ng), ic_v, &
530 lwsrc(ng), ic_w
531 END IF
532!
533! Stop if illegal configuration.
534!
535 IF (.not.lwsrc(ng).and.luvsrc(ng).and. &
536 (ic_u.eq.0).and.(ic_v.eq.0)) THEN
537 IF (master) WRITE (stdout,30) 'LuvSrc'
538 exit_flag=5
539 END IF
540 IF (.not.luvsrc(ng).and.lwsrc(ng).and.(ic_w.eq.0)) THEN
541 IF (master) WRITE (stdout,30) 'LwSrc'
542 exit_flag=5
543 END IF
544 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
545!
546 10 FORMAT (/,1x,'Point Sources/Sinks grid-cell flag locations ', &
547 & 'counter:',/)
548 20 FORMAT (4x,'Grid: ',i0,', file: ',a,/, &
549 & 19x,'LuvSrc = ',l1,2x,'u-face = ',i0,/, &
550 & 19x,'LuvSrc = ',l1,2x,'v-face = ',i0,/, &
551 & 19x,'LwSrc = ',l1,2x,'w-face = ',i0)
552 30 FORMAT (/,' CHECK_SOURCES - Cannot find point Souces/Sinks ', &
553 & "the '",a,"' method.")
554!
555 RETURN
integer stdout
logical master

References mod_scalars::exit_flag, strings_mod::founderror(), mod_param::inlm, mod_scalars::luvsrc, mod_scalars::lwsrc, mod_parallel::master, mod_scalars::noerror, mod_iounits::ssf, and mod_iounits::stdout.

Referenced by allocate_sources().

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

◆ deallocate_sources()

subroutine, public mod_sources::deallocate_sources ( integer ng)

Definition at line 332 of file mod_sources.F.

333!
334!=======================================================================
335! !
336! This routine deallocates all variables in the module for all nested !
337! grids. !
338! !
339!=======================================================================
340
341#ifdef SUBOBJECT_DEALLOCATION
342!
343 USE destroy_mod, ONLY : destroy
344#endif
345!
346! Imported variable declarations.
347!
348 integer :: ng
349!
350! Local variable declarations.
351!
352 character (len=*), parameter :: MyFile = &
353 & __FILE__//", deallocate_mixing"
354
355#ifdef SUBOBJECT_DEALLOCATION
356!
357!-----------------------------------------------------------------------
358! Deallocate each variable in the derived-type T_SOURCES structure
359! separately.
360!-----------------------------------------------------------------------
361!
362 IF (.not.destroy(ng, sources(ng)%Isrc, myfile, &
363 & __line__, 'SOURCES(ng)%Isrc')) RETURN
364
365 IF (.not.destroy(ng, sources(ng)%Jsrc, myfile, &
366 & __line__, 'SOURCES(ng)%Jsrc')) RETURN
367
368 IF (.not.destroy(ng, sources(ng)%Dsrc, myfile, &
369 & __line__, 'SOURCES(ng)%Dsrc')) RETURN
370
371 IF (.not.destroy(ng, sources(ng)%Fsrc, myfile, &
372 & __line__, 'SOURCES(ng)%Fsrc')) RETURN
373
374 IF (.not.destroy(ng, sources(ng)%Qbar, myfile, &
375 & __line__, 'SOURCES(ng)%Qbar')) RETURN
376
377 IF (.not.destroy(ng, sources(ng)%Qshape, myfile, &
378 & __line__, 'SOURCES(ng)%Qshape')) RETURN
379
380 IF (.not.destroy(ng, sources(ng)%Qsrc, myfile, &
381 & __line__, 'SOURCES(ng)%Qsrc')) RETURN
382
383 IF (.not.destroy(ng, sources(ng)%Tsrc, myfile, &
384 & __line__, 'SOURCES(ng)%Tsrc')) RETURN
385
386 IF (.not.destroy(ng, sources(ng)%Xsrc, myfile, &
387 & __line__, 'SOURCES(ng)%Xsrc')) RETURN
388
389 IF (.not.destroy(ng, sources(ng)%Ysrc, myfile, &
390 & __line__, 'SOURCES(ng)%Ysrc')) RETURN
391
392# ifndef ANA_PSOURCE
393 IF (.not.destroy(ng, sources(ng)%QbarG, myfile, &
394 & __line__, 'SOURCES(ng)%QbarG')) RETURN
395
396 IF (.not.destroy(ng, sources(ng)%TsrcG, myfile, &
397 & __line__, 'SOURCES(ng)%TsrcG')) RETURN
398# endif
399
400# ifdef ADJOINT
401 IF (.not.destroy(ng, sources(ng)%ad_Qbar, myfile, &
402 & __line__, 'SOURCES(ng)%ad_Qbar')) RETURN
403
404 IF (.not.destroy(ng, sources(ng)%ad_Qsrc, myfile, &
405 & __line__, 'SOURCES(ng)%ad_Qsrc')) RETURN
406
407 IF (.not.destroy(ng, sources(ng)%ad_Tsrc, myfile, &
408 & __line__, 'SOURCES(ng)%ad_Tsrc')) RETURN
409# endif
410
411# ifdef TANGENT
412 IF (.not.destroy(ng, sources(ng)%tl_Qbar, myfile, &
413 & __line__, 'SOURCES(ng)%tl_Qbar')) RETURN
414
415 IF (.not.destroy(ng, sources(ng)%tl_Qsrc, myfile, &
416 & __line__, 'SOURCES(ng)%tl_Qsrc')) RETURN
417
418 IF (.not.destroy(ng, sources(ng)%tl_Tsrc, myfile, &
419 & __line__, 'SOURCES(ng)%tl_Tsrc')) RETURN
420# endif
421#endif
422!
423!-----------------------------------------------------------------------
424! Deallocate derived-type SOURCES structure.
425!-----------------------------------------------------------------------
426!
427 IF (ng.eq.ngrids) THEN
428 IF (allocated(sources)) deallocate ( sources )
429 END IF
430!
431!-----------------------------------------------------------------------
432! Deallocate other variables in module.
433!-----------------------------------------------------------------------
434!
435 IF (allocated(msrc)) THEN
436 deallocate ( msrc )
437 END IF
438
439 IF (allocated(nsrc)) THEN
440 deallocate ( nsrc )
441 END IF
442!
443 RETURN

References msrc, nsrc, and sources.

Referenced by mod_arrays::roms_deallocate_arrays().

Here is the caller graph for this function:

Variable Documentation

◆ msrc

integer, dimension(:), allocatable mod_sources::msrc

Definition at line 96 of file mod_sources.F.

96 integer, allocatable :: Msrc(:)

Referenced by allocate_sources(), analytical_mod::ana_psource_tile(), and deallocate_sources().

◆ nsrc

◆ sources