ROMS
Loading...
Searching...
No Matches
ana_stflux.h
Go to the documentation of this file.
1!!
2 SUBROUTINE ana_stflux (ng, tile, model, itrc)
3!
4!! git $Id$
5!!======================================================================
6!! Copyright (c) 2002-2025 The ROMS Group !
7!! Licensed under a MIT/X style license !
8!! See License_ROMS.md !
9!=======================================================================
10! !
11! Sets surface flux of tracer type variables stflux(:,:,itrc) using !
12! analytical expressions (TracerUnits m/s). The surface fluxes are !
13! processed and loaded to state variable "stflx" in "set_vbc". !
14! !
15!=======================================================================
16!
17 USE mod_param
18 USE mod_forces
19 USE mod_ncparam
20!
21! Imported variable declarations.
22!
23 integer, intent(in) :: ng, tile, model, itrc
24!
25! Local variable declarations.
26!
27 character (len=*), parameter :: MyFile = &
28 & __FILE__
29!
30#include "tile.h"
31!
32 CALL ana_stflux_tile (ng, tile, model, itrc, &
33 & lbi, ubi, lbj, ubj, &
34 & imins, imaxs, jmins, jmaxs, &
35#ifdef SHORTWAVE
36 & forces(ng) % srflx, &
37#endif
38 & forces(ng) % stflux)
39!
40! Set analytical header file name used.
41!
42#ifdef DISTRIBUTE
43 IF (lanafile) THEN
44#else
45 IF (lanafile.and.(tile.eq.0)) THEN
46#endif
47 ananame(31)=myfile
48 END IF
49!
50 RETURN
51 END SUBROUTINE ana_stflux
52!
53!***********************************************************************
54 SUBROUTINE ana_stflux_tile (ng, tile, model, itrc, &
55 & LBi, UBi, LBj, UBj, &
56 & IminS, ImaxS, JminS, JmaxS, &
57#ifdef SHORTWAVE
58 & srflx, &
59#endif
60 & stflux)
61!***********************************************************************
62!
63 USE mod_param
64 USE mod_scalars
65!
67#ifdef DISTRIBUTE
69#endif
70!
71! Imported variable declarations.
72!
73 integer, intent(in) :: ng, tile, model, itrc
74 integer, intent(in) :: LBi, UBi, LBj, UBj
75 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
76!
77#ifdef ASSUMED_SHAPE
78# ifdef SHORTWAVE
79 real(r8), intent(in) :: srflx(LBi:,LBj:)
80# endif
81 real(r8), intent(inout) :: stflux(LBi:,LBj:,:)
82#else
83# ifdef SHORTWAVE
84 real(r8), intent(in) :: srflx(LBi:UBi,LBj:UBj)
85# endif
86 real(r8), intent(inout) :: stflux(LBi:UBi,LBj:UBj,NT(ng))
87#endif
88!
89! Local variable declarations.
90!
91 integer :: i, j
92
93#include "set_bounds.h"
94!
95!-----------------------------------------------------------------------
96! Set surface net heat flux (degC m/s) at horizontal RHO-points.
97!-----------------------------------------------------------------------
98!
99 IF (itrc.eq.itemp) THEN
100 DO j=jstrt,jendt
101 DO i=istrt,iendt
102#ifdef BL_TEST
103 stflux(i,j,itrc)=srflx(i,j)
104#else
105 stflux(i,j,itrc)=0.0_r8
106#endif
107 END DO
108 END DO
109!
110!-----------------------------------------------------------------------
111! Set surface freshwater flux (m/s) at horizontal RHO-points. The
112! scaling by surface salinity is done in "set_vbc".
113!-----------------------------------------------------------------------
114!
115 ELSE IF (itrc.eq.isalt) THEN
116 DO j=jstrt,jendt
117 DO i=istrt,iendt
118 stflux(i,j,itrc)=0.0_r8
119 END DO
120 END DO
121!
122!-----------------------------------------------------------------------
123! Set surface flux (Tunits m/s) of passive tracers at RHO-points,
124! if any.
125!-----------------------------------------------------------------------
126!
127 ELSE
128 DO j=jstrt,jendt
129 DO i=istrt,iendt
130 stflux(i,j,itrc)=0.0_r8
131 END DO
132 END DO
133 END IF
134!
135! Exchange boundary data.
136!
137 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
138 CALL exchange_r2d_tile (ng, tile, &
139 & lbi, ubi, lbj, ubj, &
140 & stflux(:,:,itrc))
141 END IF
142
143#ifdef DISTRIBUTE
144!
145 CALL mp_exchange2d (ng, tile, model, 1, &
146 & lbi, ubi, lbj, ubj, &
147 & nghostpoints, &
148 & ewperiodic(ng), nsperiodic(ng), &
149 & stflux(:,:,itrc))
150#endif
151!
152 RETURN
153 END SUBROUTINE ana_stflux_tile
subroutine ana_stflux_tile(ng, tile, model, itrc, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, srflx, stflux)
Definition ana_stflux.h:61
subroutine ana_stflux(ng, tile, model, itrc)
Definition ana_stflux.h:3
subroutine exchange_r2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
type(t_forces), dimension(:), allocatable forces
Definition mod_forces.F:554
logical lanafile
character(len=256), dimension(39) ananame
integer nghostpoints
Definition mod_param.F:710
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic
integer isalt
integer itemp
subroutine mp_exchange2d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, nghost, ew_periodic, ns_periodic, a, b, c, d)