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

Functions/Subroutines

subroutine, public ini_hmixcoef (ng, tile, model)
 
subroutine ini_hmixcoef_tile (ng, tile, model, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, grdscl, diff2, diff4, visc2_p, visc2_r, visc4_p, visc4_r, diffusion2, diffusion4, viscosity2, viscosity4)
 

Function/Subroutine Documentation

◆ ini_hmixcoef()

subroutine, public ini_hmixcoef_mod::ini_hmixcoef ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) model )

Definition at line 29 of file ini_hmixcoef.F.

30!***********************************************************************
31!
32 USE mod_param
33 USE mod_grid
34 USE mod_mixing
35 USE mod_ncparam
36 USE mod_scalars
37!
38! Imported variable declarations.
39!
40 integer, intent(in) :: ng, tile, model
41!
42! Local variable declarations.
43!
44#ifdef SOLVE3D
45 real(r8) :: diffusion2(MT), diffusion4(MT)
46#endif
47 real(r8) :: viscosity2, viscosity4
48!
49#include "tile.h"
50
51 CALL ini_hmixcoef_tile (ng, tile, model, &
52 & lbi, ubi, lbj, ubj, &
53 & imins, imaxs, jmins, jmaxs, &
54 & grid(ng) % grdscl, &
55#ifdef SOLVE3D
56# ifdef TS_DIF2
57 & mixing(ng) % diff2, &
58# endif
59# ifdef TS_DIF4
60 & mixing(ng) % diff4, &
61# endif
62#endif
63#ifdef UV_VIS2
64 & mixing(ng) % visc2_p, &
65 & mixing(ng) % visc2_r, &
66#endif
67#ifdef UV_VIS4
68 & mixing(ng) % visc4_p, &
69 & mixing(ng) % visc4_r, &
70#endif
71#ifdef SOLVE3D
72 & diffusion2, diffusion4, &
73#endif
74 & viscosity2, viscosity4)
75
76 RETURN
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
type(t_mixing), dimension(:), allocatable mixing
Definition mod_mixing.F:399

References mod_grid::grid, ini_hmixcoef_tile(), and mod_mixing::mixing.

Referenced by ad_initial(), roms_kernel_mod::adm_initial(), initial(), roms_kernel_mod::nlm_initial(), rp_initial(), tl_initial(), and roms_kernel_mod::tlm_initial().

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

◆ ini_hmixcoef_tile()

subroutine ini_hmixcoef_mod::ini_hmixcoef_tile ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) model,
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,
real(r8), dimension(lbi:,lbj:), intent(in) grdscl,
real(r8), dimension(lbi:,lbj:,:), intent(inout) diff2,
real(r8), dimension(lbi:,lbj:,:), intent(inout) diff4,
real(r8), dimension(lbi:,lbj:), intent(inout) visc2_p,
real(r8), dimension(lbi:,lbj:), intent(inout) visc2_r,
real(r8), dimension(lbi:,lbj:), intent(inout) visc4_p,
real(r8), dimension(lbi:,lbj:), intent(inout) visc4_r,
real(r8), dimension(mt), intent(out) diffusion2,
real(r8), dimension(mt), intent(out) diffusion4,
real(r8), intent(out) viscosity2,
real(r8), intent(out) viscosity4 )
private

Definition at line 80 of file ini_hmixcoef.F.

104!***********************************************************************
105!
106 USE mod_param
107 USE mod_mixing
108 USE mod_scalars
109!
111#ifdef DISTRIBUTE
113# ifdef SOLVE3D
115# endif
116#endif
117!
118! Imported variable declarations.
119!
120 integer, intent(in) :: ng, tile, model
121 integer, intent(in) :: LBi, UBi, LBj, UBj
122 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
123!
124#ifdef SOLVE3D
125 real(r8), intent(out) :: diffusion2(MT), diffusion4(MT)
126#endif
127 real(r8), intent(out) :: viscosity2, viscosity4
128!
129#ifdef ASSUMED_SHAPE
130 real(r8), intent(in) :: grdscl(LBi:,LBj:)
131# ifdef SOLVE3D
132# ifdef TS_DIF2
133 real(r8), intent(inout) :: diff2(LBi:,LBj:,:)
134# endif
135# ifdef TS_DIF4
136 real(r8), intent(inout) :: diff4(LBi:,LBj:,:)
137# endif
138# endif
139# if defined UV_VIS2
140 real(r8), intent(inout) :: visc2_p(LBi:,LBj:)
141 real(r8), intent(inout) :: visc2_r(LBi:,LBj:)
142# endif
143# ifdef UV_VIS4
144 real(r8), intent(inout) :: visc4_p(LBi:,LBj:)
145 real(r8), intent(inout) :: visc4_r(LBi:,LBj:)
146# endif
147#else
148 real(r8), intent(in) :: grdscl(LBi:UBi,LBj:UBj)
149# ifdef SOLVE3D
150# ifdef TS_DIF2
151 real(r8), intent(inout) :: diff2(LBi:UBi,LBj:UBj,NT(ng))
152# endif
153# ifdef TS_DIF4
154 real(r8), intent(inout) :: diff4(LBi:UBi,LBj:UBj,NT(ng))
155# endif
156# endif
157# if defined UV_VIS2
158 real(r8), intent(inout) :: visc2_p(LBi:UBi,LBj:UBj)
159 real(r8), intent(inout) :: visc2_r(LBi:UBi,LBj:UBj)
160# endif
161# ifdef UV_VIS4
162 real(r8), intent(inout) :: visc4_p(LBi:UBi,LBj:UBj)
163 real(r8), intent(inout) :: visc4_r(LBi:UBi,LBj:UBj)
164# endif
165#endif
166!
167! Local variable declarations.
168!
169 integer :: Imin, Imax, Jmin, Jmax
170 integer :: i, j
171#ifdef SOLVE3D
172 integer :: itrc
173#endif
174 real(r8) :: cff
175
176#include "set_bounds.h"
177!
178!-----------------------------------------------------------------------
179! Set horizontal, constant, mixing coefficient according to model flag.
180!-----------------------------------------------------------------------
181!
182 IF (model.eq.inlm) THEN
183 viscosity2=nl_visc2(ng)
184 viscosity4=nl_visc4(ng)
185#ifdef SOLVE3D
186 DO itrc=1,nt(ng)
187 diffusion2(itrc)=nl_tnu2(itrc,ng)
188 diffusion4(itrc)=nl_tnu4(itrc,ng)
189 END DO
190#endif
191#if defined TANGENT || defined TL_IOMS
192 ELSE IF ((model.eq.itlm).or.(model.eq.irpm)) THEN
193 viscosity2=tl_visc2(ng)
194 viscosity4=tl_visc4(ng)
195# ifdef SOLVE3D
196 DO itrc=1,nt(ng)
197 diffusion2(itrc)=tl_tnu2(itrc,ng)
198 diffusion4(itrc)=tl_tnu4(itrc,ng)
199 END DO
200# endif
201#endif
202#ifdef ADJOINT
203 ELSE IF (model.eq.iadm) THEN
204 viscosity2=ad_visc2(ng)
205 viscosity4=ad_visc4(ng)
206# ifdef SOLVE3D
207 DO itrc=1,nt(ng)
208 diffusion2(itrc)=ad_tnu2(itrc,ng)
209 diffusion4(itrc)=ad_tnu4(itrc,ng)
210 END DO
211# endif
212#endif
213 END IF
214!
215! Update generic values.
216!
217 IF (domain(ng)%SouthWest_Test(tile)) THEN
218 visc2(ng)=viscosity2
219 visc4(ng)=viscosity4
220#ifdef SOLVE3D
221 DO itrc=1,nt(ng)
222 tnu2(itrc,ng)=diffusion2(itrc)
223 tnu4(itrc,ng)=diffusion4(itrc)
224 END DO
225#endif
226 END IF
227
228#if defined UV_VIS2 || defined UV_VIS4 || \
229 ((defined ts_dif2 || defined ts_dif4) && defined solve3d)
230!
231!-----------------------------------------------------------------------
232! Initialize horizontal mixing arrays to constant mixing coefficient.
233!-----------------------------------------------------------------------
234!
235# ifdef _OPENMP
236 IF (domain(ng)%Western_Edge(tile)) THEN
237 imin=bounds(ng)%LBi(tile)
238 ELSE
239 imin=istr
240 END IF
241 IF (domain(ng)%Eastern_Edge(tile)) THEN
242 imax=bounds(ng)%UBi(tile)
243 ELSE
244 imax=iend
245 END IF
246 IF (domain(ng)%Southern_Edge(tile)) THEN
247 jmin=bounds(ng)%LBj(tile)
248 ELSE
249 jmin=jstr
250 END IF
251 IF (domain(ng)%Northern_Edge(tile)) THEN
252 jmax=bounds(ng)%UBj(tile)
253 ELSE
254 jmax=jend
255 END IF
256# else
257 imin=bounds(ng)%LBi(tile)
258 imax=bounds(ng)%UBi(tile)
259 jmin=bounds(ng)%LBj(tile)
260 jmax=bounds(ng)%UBj(tile)
261# endif
262!
263# if defined UV_VIS2
264 DO j=jmin,jmax
265 DO i=imin,imax
266 visc2_p(i,j)=viscosity2
267 visc2_r(i,j)=viscosity2
268 END DO
269 END DO
270# endif
271# ifdef UV_VIS4
272 DO j=jmin,jmax
273 DO i=imin,imax
274 visc4_p(i,j)=viscosity4
275 visc4_r(i,j)=viscosity4
276 END DO
277 END DO
278# endif
279# ifdef SOLVE3D
280# ifdef TS_DIF2
281 DO itrc=1,nt(ng)
282 DO j=jmin,jmax
283 DO i=imin,imax
284 diff2(i,j,itrc)=diffusion2(itrc)
285 END DO
286 END DO
287 END DO
288# endif
289# ifdef TS_DIF4
290 DO itrc=1,nt(ng)
291 DO j=jmin,jmax
292 DO i=imin,imax
293 diff4(i,j,itrc)=diffusion4(itrc)
294 END DO
295 END DO
296 END DO
297# endif
298# endif
299#endif
300
301#ifdef VISC_GRID
302!
303!-----------------------------------------------------------------------
304! Scale horizontal viscosity according to the grid size. The values
305! used during initialization above are over-written.
306!-----------------------------------------------------------------------
307!
308# ifdef UV_VIS2
309 cff=viscosity2/grdmax(ng)
310 DO j=jstrt,jendt
311 DO i=istrt,iendt
312 visc2_r(i,j)=cff*grdscl(i,j)
313 END DO
314 END DO
315 cff=0.25_r8*cff
316 DO j=jstrp,jendt
317 DO i=istrp,iendt
318 visc2_p(i,j)=cff*(grdscl(i-1,j-1)+grdscl(i,j-1)+ &
319 & grdscl(i-1,j )+grdscl(i,j ))
320 END DO
321 END DO
322# endif
323# ifdef UV_VIS4
324 cff=viscosity4/(grdmax(ng)**3)
325 DO j=jstrt,jendt
326 DO i=istrt,iendt
327 visc4_r(i,j)=cff*grdscl(i,j)**3
328 END DO
329 END DO
330 cff=0.25_r8*cff
331 DO j=jstrp,jendt
332 DO i=istrp,iendt
333 visc4_p(i,j)=cff*(grdscl(i,j )**3+grdscl(i-1,j )**3+ &
334 & grdscl(i,j-1)**3+grdscl(i-1,j-1)**3)
335 END DO
336 END DO
337# endif
338#endif
339
340#if defined DIFF_GRID && defined SOLVE3D
341!
342!-----------------------------------------------------------------------
343! Scale horizontal diffusion according to the grid size.
344!-----------------------------------------------------------------------
345!
346# ifdef TS_DIF2
347 DO itrc=1,nt(ng)
348 cff=diffusion2(itrc)/grdmax(ng)
349 DO j=jstrt,jendt
350 DO i=istrt,iendt
351 diff2(i,j,itrc)=cff*grdscl(i,j)
352 END DO
353 END DO
354 END DO
355# endif
356# ifdef TS_DIF4
357 DO itrc=1,nt(ng)
358 cff=diffusion4(itrc)/(grdmax(ng)**3)
359 DO j=jstrt,jendt
360 DO i=istrt,iendt
361 diff4(i,j,itrc)=cff*grdscl(i,j)**3
362 END DO
363 END DO
364 END DO
365# endif
366#endif
367
368#if !defined ANA_SPONGE && \
369 ( defined uv_vis2 || defined uv_vis4 || \
370 ((defined ts_dif2 || defined ts_dif4) && defined solve3d))
371!
372!-----------------------------------------------------------------------
373! Increase horizontal mixing coefficients in the sponge areas using
374! the nondimentional factors read from application Grid NetCDF file.
375!-----------------------------------------------------------------------
376!
377 IF (luvsponge(ng)) THEN
378# ifdef UV_VIS2
379 DO i=istrt,iendt
380 DO j=jstrt,jendt
381 visc2_r(i,j)=abs(mixing(ng)%visc_factor(i,j))* &
382 & visc2_r(i,j)
383 END DO
384 END DO
385 DO i=istrp,iendt
386 DO j=jstrp,jendt
387 visc2_p(i,j)=0.25_r8* &
388 & abs(mixing(ng)%visc_factor(i-1,j-1)+ &
389 & mixing(ng)%visc_factor(i ,j-1)+ &
390 & mixing(ng)%visc_factor(i-1,j )+ &
391 & mixing(ng)%visc_factor(i ,j ))* &
392 & visc2_p(i,j)
393 END DO
394 END DO
395# endif
396# ifdef UV_VIS4
397 DO i=istrt,iendt
398 DO j=jstrt,jendt
399 visc4_r(i,j)=abs(mixing(ng)%visc_factor(i,j))* &
400 & visc4_r(i,j)
401 END DO
402 END DO
403 DO i=istrp,iendt
404 DO j=jstrp,jendt
405 visc4_p(i,j)=0.25_r8* &
406 & abs(mixing(ng)%visc_factor(i-1,j-1)+ &
407 & mixing(ng)%visc_factor(i ,j-1)+ &
408 & mixing(ng)%visc_factor(i-1,j )+ &
409 & mixing(ng)%visc_factor(i ,j ))* &
410 & visc4_p(i,j)
411 END DO
412 END DO
413# endif
414 END IF
415
416# ifdef SOLVE3D
417# ifdef TS_DIF2
418 DO itrc=1,nt(ng)
419 IF (ltracersponge(itrc,ng)) THEN
420 DO j=jstrt,jendt
421 DO i=istrt,iendt
422 diff2(i,j,itrc)=abs(mixing(ng)%diff_factor(i,j))* &
423 & diff2(i,j,itrc)
424 END DO
425 END DO
426 END IF
427 END DO
428# endif
429# ifdef TS_DIF4
430 DO itrc=1,nt(ng)
431 IF (ltracersponge(itrc,ng)) THEN
432 DO j=jstrt,jendt
433 DO i=istrt,iendt
434 diff4(i,j,itrc)=abs(mixing(ng)%diff_factor(i,j))* &
435 & diff4(i,j,itrc)
436 END DO
437 END DO
438 END IF
439 END DO
440# endif
441# endif
442#endif
443!
444!-----------------------------------------------------------------------
445! Exchange boundary data.
446!-----------------------------------------------------------------------
447!
448#ifdef UV_VIS2
449 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
450 CALL exchange_r2d_tile (ng, tile, &
451 & lbi, ubi, lbj, ubj, &
452 & visc2_r)
453 CALL exchange_p2d_tile (ng, tile, &
454 & lbi, ubi, lbj, ubj, &
455 & visc2_p)
456 END IF
457#endif
458#ifdef UV_VIS4
459 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
460 CALL exchange_r2d_tile (ng, tile, &
461 & lbi, ubi, lbj, ubj, &
462 & visc4_r)
463 CALL exchange_p2d_tile (ng, tile, &
464 & lbi, ubi, lbj, ubj, &
465 & visc4_p)
466 END IF
467#endif
468
469#ifdef SOLVE3D
470# ifdef TS_DIF2
471 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
472 DO itrc=1,nt(ng)
473 CALL exchange_r2d_tile (ng, tile, &
474 & lbi, ubi, lbj, ubj, &
475 & diff2(:,:,itrc))
476 END DO
477 END IF
478# endif
479# ifdef TS_DIF4
480 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
481 DO itrc=1,nt(ng)
482 CALL exchange_r2d_tile (ng, tile, &
483 & lbi, ubi, lbj, ubj, &
484 & diff4(:,:,itrc))
485 END DO
486 END IF
487# endif
488#endif
489
490#ifdef DISTRIBUTE
491# ifdef UV_VIS2
492 CALL mp_exchange2d (ng, tile, model, 2, &
493 & lbi, ubi, lbj, ubj, &
494 & nghostpoints, &
495 & ewperiodic(ng), nsperiodic(ng), &
496 & visc2_r, visc2_p)
497# endif
498# ifdef UV_VIS4
499 CALL mp_exchange2d (ng, tile, model, 2, &
500 & lbi, ubi, lbj, ubj, &
501 & nghostpoints, &
502 & ewperiodic(ng), nsperiodic(ng), &
503 & visc4_r, visc4_p)
504# endif
505# ifdef SOLVE3D
506# ifdef TS_DIF2
507 CALL mp_exchange3d (ng, tile, model, 1, &
508 & lbi, ubi, lbj, ubj, 1, nt(ng), &
509 & nghostpoints, &
510 & ewperiodic(ng), nsperiodic(ng), &
511 & diff2)
512# endif
513# ifdef TS_DIF4
514 CALL mp_exchange3d (ng, tile, model, 1, &
515 & lbi, ubi, lbj, ubj, 1, nt(ng), &
516 & nghostpoints, &
517 & ewperiodic(ng), nsperiodic(ng), &
518 & diff4)
519# endif
520# endif
521#endif
522
523 RETURN
subroutine exchange_r2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
subroutine exchange_p2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
Definition exchange_2d.F:66
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 nghostpoints
Definition mod_param.F:710
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
integer, dimension(:), allocatable nt
Definition mod_param.F:489
real(r8), dimension(:,:), allocatable tl_tnu2
real(r8), dimension(:,:), allocatable tnu2
logical, dimension(:), allocatable luvsponge
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic
real(r8), dimension(:,:), allocatable nl_tnu2
logical, dimension(:,:), allocatable ltracersponge
real(r8), dimension(:,:), allocatable tl_tnu4
real(r8), dimension(:), allocatable visc2
real(r8), dimension(:), allocatable ad_visc4
real(r8), dimension(:), allocatable visc4
real(r8), dimension(:), allocatable tl_visc4
real(dp), dimension(:), allocatable grdmax
real(r8), dimension(:,:), allocatable ad_tnu2
real(r8), dimension(:), allocatable nl_visc4
real(r8), dimension(:), allocatable ad_visc2
real(r8), dimension(:,:), allocatable ad_tnu4
real(r8), dimension(:), allocatable tl_visc2
real(r8), dimension(:,:), allocatable nl_tnu4
real(r8), dimension(:), allocatable nl_visc2
real(r8), dimension(:,:), allocatable tnu4
subroutine mp_exchange2d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine mp_exchange3d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, lbk, ubk, nghost, ew_periodic, ns_periodic, a, b, c, d)

References mod_scalars::ad_tnu2, mod_scalars::ad_tnu4, mod_scalars::ad_visc2, mod_scalars::ad_visc4, mod_param::bounds, mod_param::domain, mod_scalars::ewperiodic, exchange_2d_mod::exchange_p2d_tile(), exchange_2d_mod::exchange_r2d_tile(), mod_scalars::grdmax, mod_param::iadm, mod_param::inlm, mod_param::irpm, mod_param::itlm, mod_scalars::ltracersponge, mod_scalars::luvsponge, mod_mixing::mixing, mp_exchange_mod::mp_exchange2d(), mp_exchange_mod::mp_exchange3d(), mod_param::nghostpoints, mod_scalars::nl_tnu2, mod_scalars::nl_tnu4, mod_scalars::nl_visc2, mod_scalars::nl_visc4, mod_scalars::nsperiodic, mod_scalars::tl_tnu2, mod_scalars::tl_tnu4, mod_scalars::tl_visc2, mod_scalars::tl_visc4, mod_scalars::tnu2, mod_scalars::tnu4, mod_scalars::visc2, and mod_scalars::visc4.

Referenced by ini_hmixcoef().

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