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

Functions/Subroutines

subroutine, public tl_t3drelax (ng, tile)
 
subroutine tl_t3drelax_tile (ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, nnew, umask, vmask, hz, pmon_u, pnom_v, pm, pn, tl_t)
 

Function/Subroutine Documentation

◆ tl_t3drelax()

subroutine, public tl_t3drelax_mod::tl_t3drelax ( integer, intent(in) ng,
integer, intent(in) tile )

Definition at line 26 of file tl_t3drelax.F.

27!***********************************************************************
28!
29 USE mod_param
30 USE mod_grid
31 USE mod_ocean
32 USE mod_stepping
33!
34! Imported variable declarations.
35!
36 integer, intent(in) :: ng, tile
37!
38! Local variable declarations.
39!
40 character (len=*), parameter :: MyFile = &
41 & __FILE__
42!
43# include "tile.h"
44!
45# ifdef PROFILE
46 CALL wclock_on (ng, itlm, 24, __line__, myfile)
47# endif
48 CALL tl_t3drelax_tile (ng, tile, &
49 & lbi, ubi, lbj, ubj, &
50 & imins, imaxs, jmins, jmaxs, &
51 & nrhs(ng), nnew(ng), &
52# ifdef MASKING
53 & grid(ng) % umask, &
54 & grid(ng) % vmask, &
55# endif
56 & grid(ng) % Hz, &
57 & grid(ng) % pmon_u, &
58 & grid(ng) % pnom_v, &
59 & grid(ng) % pm, &
60 & grid(ng) % pn, &
61 & ocean(ng) % tl_t)
62# ifdef PROFILE
63 CALL wclock_off (ng, itlm, 24, __line__, myfile)
64# endif
65!
66 RETURN
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
integer, parameter itlm
Definition mod_param.F:663
integer, dimension(:), allocatable nrhs
integer, dimension(:), allocatable nnew
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

References mod_grid::grid, mod_param::itlm, mod_stepping::nnew, mod_stepping::nrhs, mod_ocean::ocean, tl_t3drelax_tile(), wclock_off(), and wclock_on().

Referenced by tl_rhs3d_mod::tl_rhs3d().

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

◆ tl_t3drelax_tile()

subroutine tl_t3drelax_mod::tl_t3drelax_tile ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) lbi,
integer, intent(in) ubi,
integer, intent(in) lbj,
integer, intent(in) ubj,
integer, intent(in) imins,
integer, intent(in) imaxs,
integer, intent(in) jmins,
integer, intent(in) jmaxs,
integer, intent(in) nrhs,
integer, intent(in) nnew,
real(r8), dimension(lbi:,lbj:), intent(in) umask,
real(r8), dimension(lbi:,lbj:), intent(in) vmask,
real(r8), dimension(lbi:,lbj:,:), intent(in) hz,
real(r8), dimension(lbi:,lbj:), intent(in) pmon_u,
real(r8), dimension(lbi:,lbj:), intent(in) pnom_v,
real(r8), dimension(lbi:,lbj:), intent(in) pm,
real(r8), dimension(lbi:,lbj:), intent(in) pn,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(inout) tl_t )
private

Definition at line 70 of file tl_t3drelax.F.

80!***********************************************************************
81!
82 USE mod_param
83 USE mod_scalars
84!
85! Imported variable declarations.
86!
87 integer, intent(in) :: ng, tile
88 integer, intent(in) :: LBi, UBi, LBj, UBj
89 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
90 integer, intent(in) :: nrhs, nnew
91
92# ifdef ASSUMED_SHAPE
93# ifdef MASKING
94 real(r8), intent(in) :: umask(LBi:,LBj:)
95 real(r8), intent(in) :: vmask(LBi:,LBj:)
96# endif
97 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
98 real(r8), intent(in) :: pmon_u(LBi:,LBj:)
99 real(r8), intent(in) :: pnom_v(LBi:,LBj:)
100 real(r8), intent(in) :: pm(LBi:,LBj:)
101 real(r8), intent(in) :: pn(LBi:,LBj:)
102
103 real(r8), intent(inout) :: tl_t(LBi:,LBj:,:,:,:)
104# else
105# ifdef MASKING
106 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
107 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
108# endif
109 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
110 real(r8), intent(in) :: pmon_u(LBi:UBi,LBj:UBj)
111 real(r8), intent(in) :: pnom_v(LBi:UBi,LBj:UBj)
112 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
113 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
114
115 real(r8), intent(inout) :: tl_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
116# endif
117!
118! Local variable declarations.
119!
120 integer :: i, itrc, j, k
121
122 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FE
123 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FX
124
125# include "set_bounds.h"
126!
127!----------------------------------------------------------------------
128! Compute horizontal diffusion relaxation of tracers between current
129! and previous representer tangent linear Picard iteration trajectory
130! (basic state).
131!----------------------------------------------------------------------
132!
133! This is the tangent linear of the relaxation terms that appear in
134! the representer model. Since the representer model is linear, we
135! DO NOT need to include the tangent linear of the Hz terms below
136! because Hz is computed from the background trajectory (AMM).
137!
138 DO itrc=1,nt(ng)
139 IF (tl_tdiff(itrc,ng).gt.0.0_r8) THEN
140 DO k=1,n(ng)
141!
142! Compute XI- and ETA-components of diffusive tracer flux (T m3/s).
143!
144 DO j=jstr,jend
145 DO i=istr,iend+1
146 fx(i,j)=0.5_r8*tl_tdiff(itrc,ng)*pmon_u(i,j)* &
147 & (hz(i,j,k)+hz(i-1,j,k))* &
148 & (tl_t(i ,j,k,nrhs,itrc)- &
149 & tl_t(i-1,j,k,nrhs,itrc))
150# ifdef MASKING
151 fx(i,j)=fx(i,j)*umask(i,j)
152# endif
153 END DO
154 END DO
155 DO j=jstr,jend+1
156 DO i=istr,iend
157 fe(i,j)=0.5_r8*tl_tdiff(itrc,ng)*pnom_v(i,j)* &
158 & (hz(i,j,k)+hz(i,j-1,k))* &
159 & (tl_t(i,j ,k,nrhs,itrc)- &
160 & tl_t(i,j-1,k,nrhs,itrc))
161# ifdef MASKING
162 fe(i,j)=fe(i,j)*vmask(i,j)
163# endif
164 END DO
165 END DO
166!
167! Time-step horizontal diffusion relaxation term (m Tunits).
168!
169 DO j=jstr,jend
170 DO i=istr,iend
171 tl_t(i,j,k,nnew,itrc)=tl_t(i,j,k,nnew,itrc)+ &
172 & dt(ng)*pm(i,j)*pn(i,j)* &
173 & (fx(i+1,j)-fx(i,j)+ &
174 & fe(i,j+1)-fe(i,j))
175 END DO
176 END DO
177 END DO
178 END IF
179 END DO
180!
181 RETURN
integer, dimension(:), allocatable n
Definition mod_param.F:479
integer, dimension(:), allocatable nt
Definition mod_param.F:489
real(r8), dimension(:,:), allocatable tl_tdiff
real(dp), dimension(:), allocatable dt

References mod_scalars::dt, and mod_scalars::tl_tdiff.

Referenced by tl_t3drelax().

Here is the caller graph for this function: