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

Data Types

type  t_sedbed
 

Functions/Subroutines

subroutine, public allocate_sedbed (ng, lbi, ubi, lbj, ubj)
 
subroutine, public deallocate_sedbed (ng)
 
subroutine, public initialize_sedbed (ng, tile, model)
 

Variables

type(t_sedbed), dimension(:), allocatable sedbed
 

Function/Subroutine Documentation

◆ allocate_sedbed()

subroutine, public mod_sedbed::allocate_sedbed ( integer, intent(in) ng,
integer, intent(in) lbi,
integer, intent(in) ubi,
integer, intent(in) lbj,
integer, intent(in) ubj )

Definition at line 161 of file sedbed_mod.h.

162!
163!=======================================================================
164! !
165! This routine allocates all variables in the module for all nested !
166! grids. !
167! !
168!=======================================================================
169!
170 USE mod_param
171 USE mod_ncparam
172 USE mod_sediment
173!
174! Imported variable declarations.
175!
176 integer, intent(in) :: ng, LBi, UBi, LBj, UBj
177!
178!-----------------------------------------------------------------------
179! Allocate structure variables.
180!-----------------------------------------------------------------------
181!
182 IF (ng.eq.1) allocate ( sedbed(ngrids) )
183!
184! Nonlinear model state.
185!
186#if defined BEDLOAD && \
187 defined averages || \
188 (defined ad_averages && defined adjoint) || \
189 (defined rp_averages && defined tl_ioms) || \
190 (defined tl_averages && defined tangent)
191 IF (any(aout(idubld(:),ng))) THEN
192 allocate ( sedbed(ng) % avgbedldu(lbi:ubi,lbj:ubj,nst) )
193 END IF
194 IF (any(aout(idvbld(:),ng))) THEN
195 allocate ( sedbed(ng) % avgbedldv(lbi:ubi,lbj:ubj,nst) )
196 END IF
197#endif
198#if defined SEDIMENT
199 allocate ( sedbed(ng) % bed(lbi:ubi,lbj:ubj,nbed,mbedp) )
200 allocate ( sedbed(ng) % bed_frac(lbi:ubi,lbj:ubj,nbed,nst) )
201 allocate ( sedbed(ng) % bed_mass(lbi:ubi,lbj:ubj,nbed,2,nst) )
202#endif
203#if defined SEDIMENT && defined SED_MORPH
204 allocate ( sedbed(ng) % bed_thick0(lbi:ubi,lbj:ubj) )
205 allocate ( sedbed(ng) % bed_thick(lbi:ubi,lbj:ubj,3) )
206#endif
207#ifdef BEDLOAD
208 allocate ( sedbed(ng) % bedldu(lbi:ubi,lbj:ubj,nst) )
209 allocate ( sedbed(ng) % bedldv(lbi:ubi,lbj:ubj,nst) )
210#endif
211 allocate ( sedbed(ng) % bottom(lbi:ubi,lbj:ubj,mbotp) )
212#if defined SEDIMENT && defined SUSPLOAD
213 allocate ( sedbed(ng) % ero_flux(lbi:ubi,lbj:ubj,nst) )
214 allocate ( sedbed(ng) % settling_flux(lbi:ubi,lbj:ubj,nst) )
215#endif
216
217#if defined TANGENT || defined TL_IOMS
218!
219! Tangent linear model state.
220!
221# if defined SEDIMENT
222 allocate ( sedbed(ng) % tl_bed(lbi:ubi,lbj:ubj,nbed,mbedp) )
223 allocate ( sedbed(ng) % tl_bed_frac(lbi:ubi,lbj:ubj,nbed,nst) )
224 allocate ( sedbed(ng) % tl_bed_mass(lbi:ubi,lbj:ubj,nbed,2,nst) )
225# endif
226# if defined SEDIMENT && defined SED_MORPH
227 allocate ( sedbed(ng) % tl_bed_thick0(lbi:ubi,lbj:ubj) )
228 allocate ( sedbed(ng) % tl_bed_thick(lbi:ubi,lbj:ubj,3) )
229# endif
230# ifdef BEDLOAD
231 allocate ( sedbed(ng) % tl_bedldu(lbi:ubi,lbj:ubj,nst) )
232 allocate ( sedbed(ng) % tl_bedldv(lbi:ubi,lbj:ubj,nst) )
233# endif
234 allocate ( sedbed(ng) % tl_bottom(lbi:ubi,lbj:ubj,mbotp) )
235# if defined SEDIMENT && defined SUSPLOAD
236 allocate ( sedbed(ng) % tl_ero_flux(lbi:ubi,lbj:ubj,nst) )
237 allocate ( sedbed(ng) % tl_settling_flux(lbi:ubi,lbj:ubj,nst) )
238# endif
239#endif
240
241#ifdef ADJOINT
242!
243! Adjoint model state.
244!
245# if defined SEDIMENT
246 allocate ( sedbed(ng) % ad_bed(lbi:ubi,lbj:ubj,nbed,mbedp) )
247 allocate ( sedbed(ng) % ad_bed_frac(lbi:ubi,lbj:ubj,nbed,nst) )
248 allocate ( sedbed(ng) % ad_bed_mass(lbi:ubi,lbj:ubj,nbed,2,nst) )
249# endif
250# if defined SEDIMENT && defined SED_MORPH
251 allocate ( sedbed(ng) % ad_bed_thick0(lbi:ubi,lbj:ubj) )
252 allocate ( sedbed(ng) % ad_bed_thick(lbi:ubi,lbj:ubj,3) )
253# endif
254# ifdef BEDLOAD
255 allocate ( sedbed(ng) % ad_bedldu(lbi:ubi,lbj:ubj,nst) )
256 allocate ( sedbed(ng) % ad_bedldv(lbi:ubi,lbj:ubj,nst) )
257# endif
258 allocate ( sedbed(ng) % ad_bottom(lbi:ubi,lbj:ubj,mbotp) )
259# if defined SEDIMENT && defined SUSPLOAD
260 allocate ( sedbed(ng) % ad_ero_flux(lbi:ubi,lbj:ubj,nst) )
261 allocate ( sedbed(ng) % ad_settling_flux(lbi:ubi,lbj:ubj,nst) )
262# endif
263#endif
264!
265 RETURN
logical, dimension(:,:), allocatable aout
integer nbed
Definition mod_param.F:517
integer ngrids
Definition mod_param.F:113
integer nst
Definition mod_param.F:521
integer, dimension(:), allocatable idubld
integer, parameter mbotp
integer, dimension(:), allocatable idvbld
integer, parameter mbedp

References mod_ncparam::aout, mod_sediment::idubld, mod_sediment::idvbld, mod_sediment::mbedp, mod_sediment::mbotp, mod_param::nbed, mod_param::ngrids, mod_param::nst, and sedbed.

Referenced by mod_arrays::roms_allocate_arrays().

Here is the caller graph for this function:

◆ deallocate_sedbed()

subroutine, public mod_sedbed::deallocate_sedbed ( integer, intent(in) ng)

Definition at line 268 of file sedbed_mod.h.

269!
270!=======================================================================
271! !
272! This routine deallocates all variables in the module for all nested !
273! grids. !
274! !
275!=======================================================================
276!
277 USE mod_param, ONLY : ngrids
278#ifdef SUBOBJECT_DEALLOCATION
279 USE destroy_mod, ONLY : destroy
280#endif
281
282!
283! Imported variable declarations.
284!
285 integer, intent(in) :: ng
286
287#ifdef SUBOBJECT_DEALLOCATION
288!
289!-----------------------------------------------------------------------
290! Deallocate each variable in the derived-type T_BBL structure
291! separately.
292!-----------------------------------------------------------------------
293!
294! Nonlinear model state.
295!
296# if defined BEDLOAD && \
297 defined averages || \
298 (defined ad_averages && defined adjoint) || \
299 (defined rp_averages && defined tl_ioms) || \
300 (defined tl_averages && defined tangent)
301 IF (.not.destroy(ng, sedbed(ng)%avgbedldu, myfile, &
302 & __line__, 'SEDBED(ng)%avgbedldu')) RETURN
303
304 IF (.not.destroy(ng, sedbed(ng)%avgbedldv, myfile, &
305 & __line__, 'SEDBED(ng)%avgbedldv')) RETURN
306# endif
307
308# if defined SEDIMENT
309 IF (.not.destroy(ng, sedbed(ng)%bed, myfile, &
310 & __line__, 'SEDBED(ng)%bed')) RETURN
311
312 IF (.not.destroy(ng, sedbed(ng)%bed_frac, myfile, &
313 & __line__, 'SEDBED(ng)%bed_frac')) RETURN
314
315 IF (.not.destroy(ng, sedbed(ng)%bed_mass, myfile, &
316 & __line__, 'SEDBED(ng)%bed_mass')) RETURN
317# endif
318
319# if defined SEDIMENT && defined SED_MORPH
320 IF (.not.destroy(ng, sedbed(ng)%bed_thick0, myfile, &
321 & __line__, 'SEDBED(ng)%bed_thick0')) RETURN
322
323 IF (.not.destroy(ng, sedbed(ng)%bed_thick, myfile, &
324 & __line__, 'SEDBED(ng)%bed_thick')) RETURN
325# endif
326
327# ifdef BEDLOAD
328 IF (.not.destroy(ng, sedbed(ng)%bedldu, myfile, &
329 & __line__, 'SEDBED(ng)%bedldu')) RETURN
330
331 IF (.not.destroy(ng, sedbed(ng)%bedldv, myfile, &
332 & __line__, 'SEDBED(ng)%bedldv')) RETURN
333# endif
334
335 IF (.not.destroy(ng, sedbed(ng)%bottom, myfile, &
336 & __line__, 'SEDBED(ng)%bottom')) RETURN
337
338# if defined SEDIMENT && defined SUSPLOAD
339 IF (.not.destroy(ng, sedbed(ng)%ero_flux, myfile, &
340 & __line__, 'SEDBED(ng)%ero_flux')) RETURN
341
342 IF (.not.destroy(ng, sedbed(ng)%settling_flux, myfile, &
343 & __line__, 'SEDBED(ng)%settling_flux')) RETURN
344# endif
345
346# if defined TANGENT || defined TL_IOMS
347!
348! Tangent linear model state.
349!
350# if defined SEDIMENT
351 IF (.not.destroy(ng, sedbed(ng)%tl_bed, myfile, &
352 & __line__, 'SEDBED(ng)%tl_bed')) RETURN
353
354 IF (.not.destroy(ng, sedbed(ng)%tl_bed_frac, myfile, &
355 & __line__, 'SEDBED(ng)%tl_bed_frac')) RETURN
356
357 IF (.not.destroy(ng, sedbed(ng)%tl_bed_mass, myfile, &
358 & __line__, 'SEDBED(ng)%tl_bed_mass')) RETURN
359# endif
360
361# if defined SEDIMENT && defined SED_MORPH
362 IF (.not.destroy(ng, sedbed(ng)%tl_bed_thick0, myfile, &
363 & __line__, 'SEDBED(ng)%tl_bed_thick0')) RETURN
364
365 IF (.not.destroy(ng, sedbed(ng)%tl_bed_thick, myfile, &
366 & __line__, 'SEDBED(ng)%tl_bed_thick')) RETURN
367# endif
368
369# ifdef BEDLOAD
370 IF (.not.destroy(ng, sedbed(ng)%tl_bedldu, myfile, &
371 & __line__, 'SEDBED(ng)%tl_bedldu')) RETURN
372
373 IF (.not.destroy(ng, sedbed(ng)%tl_bedldv, myfile, &
374 & __line__, 'SEDBED(ng)%tl_bedldv')) RETURN
375# endif
376
377 IF (.not.destroy(ng, sedbed(ng)%bottom, myfile, &
378 & __line__, 'SEDBED(ng)%bottom')) RETURN
379
380# if defined SEDIMENT && defined SUSPLOAD
381 IF (.not.destroy(ng, sedbed(ng)%tl_ero_flux, myfile, &
382 & __line__, 'SEDBED(ng)%tl_ero_flux')) RETURN
383
384 IF (.not.destroy(ng, sedbed(ng)%tl_settling_flux, myfile, &
385 & __line__, 'SEDBED(ng)%tl_settling_flux')) RETURN
386# endif
387# endif
388
389# ifdef ADJOINT
390!
391! Adjoint model state.
392!
393# if defined SEDIMENT
394 IF (.not.destroy(ng, sedbed(ng)%ad_bed, myfile, &
395 & __line__, 'SEDBED(ng)%ad_bed')) RETURN
396
397 IF (.not.destroy(ng, sedbed(ng)%ad_bed_frac, myfile, &
398 & __line__, 'SEDBED(ng)%ad_bed_frac')) RETURN
399
400 IF (.not.destroy(ng, sedbed(ng)%ad_bed_mass, myfile, &
401 & __line__, 'SEDBED(ng)%ad_bed_mass')) RETURN
402# endif
403
404# if defined SEDIMENT && defined SED_MORPH
405 IF (.not.destroy(ng, sedbed(ng)%ad_bed_thick0, myfile, &
406 & __line__, 'SEDBED(ng)%ad_bed_thick0')) RETURN
407
408 IF (.not.destroy(ng, sedbed(ng)%ad_bed_thick, myfile, &
409 & __line__, 'SEDBED(ng)%ad_bed_thick')) RETURN
410# endif
411
412# ifdef BEDLOAD
413 IF (.not.destroy(ng, sedbed(ng)%ad_bedldu, myfile, &
414 & __line__, 'SEDBED(ng)%ad_bedldu')) RETURN
415
416 IF (.not.destroy(ng, sedbed(ng)%ad_bedldv, myfile, &
417 & __line__, 'SEDBED(ng)%ad_bedldv')) RETURN
418# endif
419
420 IF (.not.destroy(ng, sedbed(ng)%bottom, myfile, &
421 & __line__, 'SEDBED(ng)%bottom')) RETURN
422
423# if defined SEDIMENT && defined SUSPLOAD
424 IF (.not.destroy(ng, sedbed(ng)%ad_ero_flux, myfile, &
425 & __line__, 'SEDBED(ng)%ad_ero_flux')) RETURN
426
427 IF (.not.destroy(ng, sedbed(ng)%ad_settling_flux, myfile, &
428 & __line__, 'SEDBED(ng)%ad_settling_flux')) RETURN
429# endif
430# endif
431#endif
432!
433!-----------------------------------------------------------------------
434! Deallocate derived-type SEDBED structure.
435!-----------------------------------------------------------------------
436!
437 IF (ng.eq.ngrids) THEN
438 IF (allocated(sedbed)) deallocate ( sedbed )
439 END IF
440!
441 RETURN

References mod_param::ngrids, and sedbed.

Referenced by mod_arrays::roms_deallocate_arrays().

Here is the caller graph for this function:

◆ initialize_sedbed()

subroutine, public mod_sedbed::initialize_sedbed ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) model )

Definition at line 444 of file sedbed_mod.h.

445!
446!=======================================================================
447! !
448! This routine initialize structure variables in the module using !
449! first touch distribution policy. In shared-memory configuration, !
450! this operation actually performs the propagation of the shared !
451! arrays across the cluster, unless another policy is specified !
452! to override the default. !
453! !
454!=======================================================================
455!
456 USE mod_param
457 USE mod_ncparam
458 USE mod_sediment
459!
460! Imported variable declarations.
461!
462 integer, intent(in) :: ng, tile, model
463!
464! Local variable declarations.
465!
466 integer :: Imin, Imax, Jmin, Jmax
467 integer :: i, itrc, j, k
468
469 real(r8), parameter :: IniVal = 0.0_r8
470
471#include "set_bounds.h"
472!
473! Set array initialization range.
474!
475#ifdef _OPENMP
476 IF (domain(ng)%Western_Edge(tile)) THEN
477 imin=bounds(ng)%LBi(tile)
478 ELSE
479 imin=istr
480 END IF
481 IF (domain(ng)%Eastern_Edge(tile)) THEN
482 imax=bounds(ng)%UBi(tile)
483 ELSE
484 imax=iend
485 END IF
486 IF (domain(ng)%Southern_Edge(tile)) THEN
487 jmin=bounds(ng)%LBj(tile)
488 ELSE
489 jmin=jstr
490 END IF
491 IF (domain(ng)%Northern_Edge(tile)) THEN
492 jmax=bounds(ng)%UBj(tile)
493 ELSE
494 jmax=jend
495 END IF
496#else
497 imin=bounds(ng)%LBi(tile)
498 imax=bounds(ng)%UBi(tile)
499 jmin=bounds(ng)%LBj(tile)
500 jmax=bounds(ng)%UBj(tile)
501#endif
502!
503!-----------------------------------------------------------------------
504! Initialize sediment structure variables.
505!-----------------------------------------------------------------------
506!
507! Nonlinear model state.
508!
509 IF ((model.eq.0).or.(model.eq.inlm)) THEN
510
511#if defined BEDLOAD && \
512 defined averages || \
513 (defined ad_averages && defined adjoint) || \
514 (defined rp_averages && defined tl_ioms) || \
515 (defined tl_averages && defined tangent)
516 IF (any(aout(idubld(:),ng))) THEN
517 DO itrc=1,nst
518 DO j=jmin, jmax
519 DO i=imin,imax
520 sedbed(ng) % avgbedldu(i,j,itrc) = inival
521 END DO
522 END DO
523 END DO
524 END IF
525 IF (any(aout(idvbld(:),ng))) THEN
526 DO itrc=1,nst
527 DO j=jmin, jmax
528 DO i=imin,imax
529 sedbed(ng) % avgbedldv(i,j,itrc) = inival
530 END DO
531 END DO
532 END DO
533 END IF
534#endif
535 DO j=jmin,jmax
536#ifdef SEDIMENT
537 DO itrc=1,mbedp
538 DO k=1,nbed
539 DO i=imin,imax
540 sedbed(ng) % bed(i,j,k,itrc) = inival
541 END DO
542 END DO
543 END DO
544 DO itrc=1,nst
545 DO k=1,nbed
546 DO i=imin,imax
547 sedbed(ng) % bed_frac(i,j,k,itrc) = inival
548 sedbed(ng) % bed_mass(i,j,k,1,itrc) = inival
549 sedbed(ng) % bed_mass(i,j,k,2,itrc) = inival
550 END DO
551 END DO
552 END DO
553#endif
554#if defined SEDIMENT && defined SED_MORPH
555 DO i=imin,imax
556 sedbed(ng) % bed_thick0(i,j) = inival
557 sedbed(ng) % bed_thick(i,j,1) = inival
558 sedbed(ng) % bed_thick(i,j,2) = inival
559 sedbed(ng) % bed_thick(i,j,3) = inival
560 END DO
561#endif
562#ifdef BEDLOAD
563 DO itrc=1,nst
564 DO i=imin,imax
565 sedbed(ng) % bedldu(i,j,itrc) = inival
566 sedbed(ng) % bedldv(i,j,itrc) = inival
567 END DO
568 END DO
569#endif
570 DO itrc=1,mbotp
571 DO i=imin,imax
572 sedbed(ng) % bottom(i,j,itrc) = inival
573 END DO
574 END DO
575#if defined SEDIMENT && defined SUSPLOAD
576 DO itrc=1,nst
577 DO i=imin,imax
578 sedbed(ng) % ero_flux(i,j,itrc) = inival
579 sedbed(ng) % settling_flux(i,j,itrc) = inival
580 END DO
581 END DO
582#endif
583 END DO
584 END IF
585
586#if defined TANGENT || defined TL_IOMS
587!
588! Tangent linear model state.
589!
590 IF ((model.eq.0).or.(model.eq.itlm).or.(model.eq.irpm)) THEN
591 DO j=jmin,jmax
592# ifdef SEDIMENT
593 DO itrc=1,mbedp
594 DO k=1,nbed
595 DO i=imin,imax
596 sedbed(ng) % tl_bed(i,j,k,itrc) = inival
597 END DO
598 END DO
599 END DO
600 DO itrc=1,nst
601 DO k=1,nbed
602 DO i=imin,imax
603 sedbed(ng) % tl_bed_frac(i,j,k,itrc) = inival
604 sedbed(ng) % tl_bed_mass(i,j,k,1,itrc) = inival
605 sedbed(ng) % tl_bed_mass(i,j,k,2,itrc) = inival
606 END DO
607 END DO
608 END DO
609# endif
610# if defined SEDIMENT && defined SED_MORPH
611 DO i=imin,imax
612 sedbed(ng) % tl_bed_thick0(i,j) = inival
613 sedbed(ng) % tl_bed_thick(i,j,1) = inival
614 sedbed(ng) % tl_bed_thick(i,j,2) = inival
615 sedbed(ng) % tl_bed_thick(i,j,3) = inival
616 END DO
617# endif
618# ifdef BEDLOAD
619 DO itrc=1,nst
620 DO i=imin,imax
621 sedbed(ng) % tl_bedldu(i,j,itrc) = inival
622 sedbed(ng) % tl_bedldv(i,j,itrc) = inival
623 END DO
624 END DO
625# endif
626 DO itrc=1,mbotp
627 DO i=imin,imax
628 sedbed(ng) % tl_bottom(i,j,itrc) = inival
629 END DO
630 END DO
631# if defined SEDIMENT && defined SUSPLOAD
632 DO itrc=1,nst
633 DO i=imin,imax
634 sedbed(ng) % tl_ero_flux(i,j,itrc) = inival
635 sedbed(ng) % tl_settling_flux(i,j,itrc) = inival
636 END DO
637 END DO
638# endif
639 END DO
640 END IF
641#endif
642
643#ifdef ADJOINT
644!
645! Adjoint model state.
646!
647 IF ((model.eq.0).or.(model.eq.iadm)) THEN
648 DO j=jmin,jmax
649# ifdef SEDIMENT
650 DO itrc=1,mbedp
651 DO k=1,nbed
652 DO i=imin,imax
653 sedbed(ng) % ad_bed(i,j,k,itrc) = inival
654 END DO
655 END DO
656 END DO
657 DO itrc=1,nst
658 DO k=1,nbed
659 DO i=imin,imax
660 sedbed(ng) % ad_bed_frac(i,j,k,itrc) = inival
661 sedbed(ng) % ad_bed_mass(i,j,k,1,itrc) = inival
662 sedbed(ng) % ad_bed_mass(i,j,k,2,itrc) = inival
663 END DO
664 END DO
665 END DO
666# endif
667# if defined SEDIMENT && defined SED_MORPH
668 DO i=imin,imax
669 sedbed(ng) % ad_bed_thick0(i,j) = inival
670 sedbed(ng) % ad_bed_thick(i,j,1) = inival
671 sedbed(ng) % ad_bed_thick(i,j,2) = inival
672 sedbed(ng) % ad_bed_thick(i,j,3) = inival
673 END DO
674# endif
675# ifdef BEDLOAD
676 DO itrc=1,nst
677 DO i=imin,imax
678 sedbed(ng) % ad_bedldu(i,j,itrc) = inival
679 sedbed(ng) % ad_bedldv(i,j,itrc) = inival
680 END DO
681 END DO
682# endif
683 DO itrc=1,mbotp
684 DO i=imin,imax
685 sedbed(ng) % ad_bottom(i,j,itrc) = inival
686 END DO
687 END DO
688# if defined SEDIMENT && defined SUSPLOAD
689 DO itrc=1,nst
690 DO i=imin,imax
691 sedbed(ng) % ad_ero_flux(i,j,itrc) = inival
692 sedbed(ng) % ad_settling_flux(i,j,itrc) = inival
693 END DO
694 END DO
695# endif
696 END DO
697 END IF
698#endif
699!
700 RETURN
integer, parameter inlm
Definition mod_param.F:662
integer, parameter irpm
Definition mod_param.F:664
type(t_bounds), dimension(:), allocatable bounds
Definition mod_param.F:232
integer, parameter iadm
Definition mod_param.F:665
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
integer, parameter itlm
Definition mod_param.F:663

References mod_ncparam::aout, mod_param::bounds, mod_param::domain, mod_param::iadm, mod_sediment::idubld, mod_sediment::idvbld, mod_param::inlm, mod_param::irpm, mod_param::itlm, mod_sediment::mbedp, mod_sediment::mbotp, mod_param::nbed, mod_param::nst, and sedbed.

Referenced by mod_arrays::roms_initialize_arrays().

Here is the caller graph for this function:

Variable Documentation

◆ sedbed