ROMS
Loading...
Searching...
No Matches
ana_humid.h
Go to the documentation of this file.
1!!
2 SUBROUTINE ana_humid (ng, tile, model)
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! This routine sets surface air humidity (moisture) using an !
12! analytical expression. There three types of humidity: !
13! !
14! 1) Absolute humidity: density of water vapor. !
15! 2) Specific humidity: ratio of the mass of water vapor to !
16! the mass of moist air cointaining the vapor (g/kg) !
17! 3) Relative humidity: ratio of the actual mixing ratio to !
18! saturation mixing ratio of the air at given temperature !
19! and pressure (percentage). !
20! !
21!=======================================================================
22!
23 USE mod_param
24 USE mod_forces
25 USE mod_ncparam
26!
27! Imported variable declarations.
28!
29 integer, intent(in) :: ng, tile, model
30!
31! Local variable declarations.
32!
33 character (len=*), parameter :: MyFile = &
34 & __FILE__
35!
36#include "tile.h"
37!
38 CALL ana_humid_tile (ng, tile, model, &
39 & lbi, ubi, lbj, ubj, &
40 & imins, imaxs, jmins, jmaxs, &
41 & forces(ng) % Hair)
42!
43! Set analytical header file name used.
44!
45#ifdef DISTRIBUTE
46 IF (lanafile) THEN
47#else
48 IF (lanafile.and.(tile.eq.0)) THEN
49#endif
50 ananame( 9)=myfile
51 END IF
52!
53 RETURN
54 END SUBROUTINE ana_humid
55!
56!***********************************************************************
57 SUBROUTINE ana_humid_tile (ng, tile, model, &
58 & LBi, UBi, LBj, UBj, &
59 & IminS, ImaxS, JminS, JmaxS, &
60 & Hair)
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
74 integer, intent(in) :: LBi, UBi, LBj, UBj
75 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
76!
77#ifdef ASSUMED_SHAPE
78 real(r8), intent(out) :: Hair(LBi:,LBj:)
79#else
80 real(r8), intent(out) :: Hair(LBi:UBi,LBj:UBj)
81#endif
82!
83! Local variable declarations.
84!
85 integer :: i, j
86
87#include "set_bounds.h"
88!
89!-----------------------------------------------------------------------
90! Set analytical surface air humidity.
91!-----------------------------------------------------------------------
92!
93#if defined BENCHMARK
94 DO j=jstrt,jendt
95 DO i=istrt,iendt
96 hair(i,j)=0.8_r8
97 END DO
98 END DO
99#elif defined BL_TEST
100 DO j=jstrt,jendt
101 DO i=istrt,iendt
102 hair(i,j)=0.776_r8
103 END DO
104 END DO
105#else
106 ana_humidity.h: no values provided for hair.
107#endif
108!
109! Exchange boundary data.
110!
111 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
112 CALL exchange_r2d_tile (ng, tile, &
113 & lbi, ubi, lbj, ubj, &
114 & hair)
115 END IF
116
117#ifdef DISTRIBUTE
118 CALL mp_exchange2d (ng, tile, model, 1, &
119 & lbi, ubi, lbj, ubj, &
120 & nghostpoints, &
121 & ewperiodic(ng), nsperiodic(ng), &
122 & hair)
123#endif
124!
125 RETURN
126 END SUBROUTINE ana_humid_tile
subroutine ana_humid_tile(ng, tile, model, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, hair)
Definition ana_humid.h:61
subroutine ana_humid(ng, tile, model)
Definition ana_humid.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
subroutine mp_exchange2d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, nghost, ew_periodic, ns_periodic, a, b, c, d)