EPA-R2-73-259
  MAY 1973               Environmental Protection Technology Series
    Heat and  Water Vapor Exchange
    Between Water Surface
    and  Atmosphere
                               \
                         ^^^^    »•
                                Office of Research and Monitoring
                                U.S. Environmental Protection Agency
                                Washington, D.C. 20460

-------
RESEARCH REPORTING SERIES
Research reports of the Office ot Research and
Monitoring, Environmental Protection Agency, have
been qrouped into five series. These five broad
categories were established to facilitate furtter
development and application of environmental
technology. Elimination of traditional grouping
was consciously planned to foster technology
transfer and a maximum interface in related
fields. The five series are:
1. Environmental Health Effects Research
2. Environmental Protection Technology
3. Ecological Research
4. Environmental Monitoring
5. Socioeconomic Environmental Studies
This report has been assigned to the ENVIRONMENTAL
PROTECTION TECHNOLOGY series. This series
describes research performed to develop and
demonstrate instrumentation, equipment and
methodology to repair or prevent environmental
degradation from point and non-point sources of
pollution. TillS work provides the new or improved
technology required for the control and treatment
of pollution sources to meet environmental quality
standards.

-------
                                                   EPA-R2-73-259
                                                   May  1973
          HEAT AND WATER VAPOR EXCHANGE BETWEEN

               WATER  SURFACE AND ATMOSPHERE
                             By
                     Wilfried Brutsaert
                     Cornell University
                  Ithaca, New York,  14850
                     Project 16130  DIP
                   Program Element  1B1032
                      Project Officer

                   Dr. Bruce A. Tichenor
 Pacific Northwest Environmental  Research Laboratory
          National  Environmental Research Center
                 Corvallis, Oregon 97330
                        Prepared for
            OFFICE OF RESEARCH AND MONITORING
           U.S. ENVIRONMENTAL PROTECTION AGENCY
                  WASHINGTON, D.C. 20460
For sale by the Superintendent of Documents, U.S. Government Printing Office, Washington, D.C. 20402
            Price 90 cents domestic postpaid or 66 cents GPO Bookstore

-------
EPA Review Notice
This report has been reviewed by the Environ-
mental Protection Agency and approved for
publication. Approval does not signify that
the contents necessarily reflect the views
and policies of the Environmental Protection
Agency, nor does mention of trade names or
commerical products constitute endorsement
or recommendation for use.
13.

-------
ABSTRACT
The physical and mathematical aspects of simultaneous
turbulent heat and water vapor exchange between a large
open water body and the surrounding atmosphere were
studied. Thus analytical and numerical solutions were
developed for various conditions of fetch, surface
roughness, atmospheric stability, etc., that are likely
to be of physical importance.
One of the main findings was that in spite of some
theoretical limitations the semi-empirical turbulent
diffusion model provides a method for the prediction
of heat and water vapor transfer, that should be useful
for engineering calculations. Although the interaction
between momentum, sensible heat and water vapor is
considerable, for evaporation or cooling calculation
purposes it is probably permissible to uncouple their
transfer mechanisms provided the water surface tempera-
ture is known and the averaging period under considera-
tion is of the order of, say, a week.
In addition the validity of the semi-empirical assumption
in the near—water surface layer was analyzed by determin-
ing the anisotropy of the eddy diffusivity, the effects of
radiative transfer and of water wave action on the eddy
diffusivity. Finally, a practical method was developed
to determine evapotranspiration from the surrounding land
surface based on the geostrophic drag concept.
This report was submitted in fulfillment of Project
Number 16130 DIP under the (partial) sponsorship of the
Environmental Protection Agency.
1 3.

-------
CONTENTS
Section Page
I Conclusions 1
II Recommendations 3
III Introduction 5
IV Turbulent Diffusion Model 7
V Discussion 11
VI Prediction of Evaporation from a
Large Water Body 21
VII Acknowledgements 25
VIII References 27
IX Publications 29
X Appendix 31
V

-------
SECTION I
CON CLUS IONS
1. The semi-empirical turbulent diffusion model with
Reynolds’ analogy provides a method for the prediction
of heat and water vapor transfer from a water surface
which should be adequate for engineering calculations.
2. For practical evaporation and or sensible cooling
calculations with averaging periods of, say, one week,
or more, it is permissible to uncouple momentum, heat
and vapor transfer. For “average” conditions, namely
a near—neutral atmosphere and a roughness of the water
surface of approximately 0.02 cm, the theoretical
calculations are in agreement with the empirical formula
given by equations (16) and (17).
3. The eddy diffusivity used in the semi-empirical
approach is highly anisotropic. However, for calcula-
tions of the fluxes near the water surface only the
vertical component has to be included. The effect of
the longitudinal, lateral and off—diagonal components
is for all practical purposes negligible.
4. With proper values of the parameters the power law
for the wind profile in the atmospheric surface layer
gives good results in the prediction of turbulent trans-
fer. The calculations are in agreement with solutions
obtainable by means of the more complicated laws such
as the logarithmic, log-linear and the other generally
accepted profile laws.
5. The effect of the difference in surface roughness
of the water and the surrounding land on heat and vapor
transfer near the surface is very small, since the air
flow and the surface shear stress adjust rather quickly
to the roughness of the new surface. Because presently
available methods are not quite adequate to predict flow
in the higher regions of the boundary layer, this points
to the importance of a better understanding of the
transfer mechanisms right at the surface.
1

-------
6. The discontinuity of surface temperature of the
evaporating water body and the surrounding land has
a strong effect on the turbulent fluxes. For example,
a relatively cold and neutral air moving over a
relatively warmer water surface may experience extreme
instabilities.
7. Under conditions normally encountered in the lower
atmosphere the effect of molecular and radiative trans-
port mechanisms and of water vapor stratification on
turbulent transfer is not always negligible. Inclusion
of these effects in the semi—empirical turbulence theory
produces expressions for the wind velocity profile and
other turbulence parameters which are in good agreement
with available experimental data.
8. The turbulent transfer over a water surface disturbed
by waves is affected not only by atmospheric conditions
but also by sea-state. Except under conditions of develop-
ing swell, this effect is a function of c/ui, where c is
the phase velocity of the dominant wave and u* the friction
velocity.
9. Evapotranspiration from large areas can be determined
satisfactorily by means of the geostrophic drag concept.
The method can be applied by using standard meteorological
data obtained at the surface and in the atmospheric
boundary layer.
2

-------
SECTION II
RECOMMENDAT IONS
The recommendations of the project are best formulated
by summarizing areas of research which in the light of
this study appear to require further work. The findings
of the present study allow the prediction of heat and
water vapor exchange for periods larger than one week.
However, for more accurate short term calculations
the following topics should be studied in greater detail.
1. The transfer mechanism in the air right at an air—wet
surface interface. This involves a better definition of
the flow, the nature of the so—called laminar sublayer
and the turbulence at the interface; also, a determination
of the relationship between the flux of heat and water
vapor from the surface into the atmosphere and such factors
as the roughness of the surface, the viscosity of the
air, the molecular diffusivities of the diffusing substances
and the shear stress of the wind.
2. The study of simultaneous heat and water vapor transfer
for an inversion, that is under extremely stable conditions
in the neighborhood of the critical Richardson number. The
nature of this type of phenomenon dictates that it be con-
sidered under unsteady conditions, and that radiative
transfer be included in the general formulation of the
problem.
3. The inclusion of radiation in the formulation of the
transfer mechanisms in the atmospheric surface layer. The
equations for radiative transfer are, however, so compli-
cated that this may necessitate the development of a
simpler semi—empirical theory analogous to that used to
obtain engineerng solutions to turbulence problems.
4. The effect of water waves on the turbulent transfer
near the surface. Specifically, this will encompass a
better characterization or parameterization of the waves,
the definition of the surface roughness in terms of wave
field parameters, and the development of similarity func-
tions to describe the turbulence under conditions of devel-
oping swell preferably on the basis of a suitable wind wave
generating theory.
3

-------
SECTION III
INTRODUCTION
Large inland water bodies, lakes or manmade reservoirs
constitute a valuable resource available to human society.
To preserve this natural resource and to allow its optimum
utilization a better understanding of the physical proper-
ties of these large water bodies and of their interaction
with the environment is required. Of critical importance
in this respect is the heat and mass exchange at the
interface between the water and the surrounding atmosphere.
Large reservoirs, lakes and the continental shelf have
also become increasingly attractive as sites for construc-
tion of power generating plants and related industries
because of the abundance of condenser cooling water. The
advent of nuclear power and the rapid increase of such
plants make the need for cooling water especially severe
and critical. Moreover, the fact that water entering the
evaporation process becomes nonrecoverable and is thus a
“consumptive” use is fundamental to water resource planning
and management. Mass transfer and evaporation must be
predicted as accurately as possible since such data are
indispensable to design the capacity of manmade reservoirs
or to assess the value of natural water bodies for municipal
and industrial water supply, navigation, irrigation of
agricultural lands and recreation.
Therefore the overall objective of this project was to gain
a better understanding of the simultaneous turbulent exchange
of momentum, heat and water vapor between a large water body
and the surrounding atmosphere for given conditions of local
atmospheric advection and of solar and other energy inputs
into the water body. The results of the theoretical analy-
sis should allow a better evaluation of some of the practi-
cal methods presently available to predict heat transfer
and evaporation which are based on simpler but cruder
models. A by—product of this investigation is a clearer
statement of some of the conditions at the water surface
which are essential in the solution of lake circulation
and stratification problems. Finally, it should be possi-
ble to apply some of the findings of the project in the
solution of certain aspects of air pollution problems.
S

-------
This objective was accomplished by the analysis of some
of the turbulent diffusion mechanisms in the surface
layer of the atmosphere. The basic theoretical frame-
work of this analysis is outlined in the next Section.
6

-------
SECTION IV
TURBULENT DIFFUSION MODEL
A conservative substance admixed in a fluid is transferred
relative to a fixed coordinate system first through con-
vective flow of the fluid itself and second through
molecular motion relative to the convective motion of the
fluid. Under fully turbulent flow conditions the molec-
ular motion may be neglected as compared to the motion
resulting from the fluid convection. Decomposing the
velocity field and the concentration of a scalar admixture
in a time mean and a fluctuation, one can write the
equation of continuity for the admixture as follows
= - I ( i ) + ( ) + (WE)] (1)
where C and c are the mean and fluctuating concentration
of the admixture, respectively. Similarly, U, V and W
are the mean velocity components and u, v and w are the
fluctuating velocity components in the mean wind direction
x, in the lateral y and in the vertical z, respectively.
It is very difficult to solve equations such as (1) because
of the problem of closure. This means that for any problem
of turbulent transfer there are always more unknowns than
equations. It is for this reason that the semi-empirical
theory of turbulence was developed. The fundamental
assumption of this approach is that the turbulent flux
components (üE), (GE) and (WE) can be related linearly to
the mean concentration gradients. For the present problem
of diffusion from a water surface at ground level it is
permissible to assume that the mean lateral velocity is
zero, V = 0, and that the mean vertical velocity W is at
least an order of magnitude smaller than the horizontal
velocity U. Thus equation (1) can be reduced to a
relatively simple partial differential equation in C only,
(e.g. Brutsaert, 1970)
+K
dtax’xxax xzaz
+ — (K — .) + —( (K + x (2)
ay yy 3y z xz ox zz az
7

-------
where (d/dt) = U(a/ax) + W(a/az) + ( / t) and where the
K—terms are the components of the eddy diffusivity tensor
for the admixture under consideration.
In the case of heat and vapor transfer the conservation
equations may be readily obtained merely by replacing C
by c T and q, respectively, where Cp is the specific
heat at constant pressure, T the absolute temperature
and g the specific humidity.
Conservation of momentum is formally describable by the
Reynolds equations. However, in the atmospheric surface
layer the flow is governed mainly by frictional forces
due to the proximity of the earth’s surface. Also, the
mean flow is practically (but not quite) plan-parallel,
i.e. as mentioned U >> W. All this means that for the
present problem it is permissible to consider the mean
horizontal wind as a scalar so that (2) also describes
the conservation of momentum if C is replaced by U.
Thus simultaneous turbulent momentum, heat and water
vapor diffusion near a water surface is governed by
the equation of continuity
au÷aW_ 0 (3)
ax az
and by the equations of conservation of horizontal
momentum, sensible heat and water vapor, namely
dU a( aU
dt x xx x xz az
a m 3U (Km + Km au) (4)
+ ( ) + a z xz ax z z
dT a h T h T
(K — + b—)
a xh T a h aT h T (5)
+ C yy 5 + xz + R
8

-------
a +KV )
dt i xx ax xz 3z
+ (K , ) + (K + K (6)
where the superscripts m, h and v refer to momentum, heat
and water vapor, respectively.
In general, these four equations must be solved simultane-
ously, subject to, among others, the following boundary
conditions
UW0 atz=O,allxandy. (7)
LFvo + FHO +R = (lY)Rs + l LD + Fw (8)
at z = 0 and x and y on the wet surface
in which “e is the latent heat of evaporation, F and
FHO the vertical fluxes of water vapor and sensi Ie heat,
respectively, R and RLfl the up- and downward long wave
radiation fluxes, respecEively, R the down coming short
wave radiation, F the heat flux released from below the
surface and y and l the albedo’s for short and long
wave radiation, respectively.
Clearly, the complete mathematical solution of this problem,
although not presenting any insurmountable difficulties, is
quite involved and handy practical for design calculations
of heat transfer and evaporation. Therefore, the major
effort in this project consisted of determining the relative
importance of the various terms of equation (2) and the
sensitivity of the solutions to different boundary condi-
tions and also of analyzing the degree of interaction among
the transfer of heat, water vapor and momentum. This was
done primarily to develop and evaluate approximate but
simpler methods. In addition, however, because of the
uncertainties of certain assumptions concerning the various
physical mechanisms, a large amount of attention was also
devoted to a better formulation of turbulent and other
transfer phenomena.
9

-------
SECTION V
DISCUSS ION
In this section a description and a discussion is given
of the studies that were carried out with the basic
turbulent diffusion model presented in Section IV.
Relative Im ortance of Individual Terms in Equation (2) .
Under certain conditions, which are not unco on on a long
term basis, it is permissible to assume, first, that
momentum transfer is strictly vertical, so that equation
(4) can be replaced by a one—dimensional wind profile
equation with constant roughness and atmospheric stability
over the land and water surface with U = U(z) and W = 0;
and second, that the temperature of the water surface is
constant, so that the heat and vapor fields are separable
and each can be analyzed by itself. Although this special
situation occurs only rarely at any given instant, it is
of interest because it allows the determination of the
relative importance of the individual terms in the turbu-
lent diffusion equation (2) without the need to consider
the complicated interaction among equations (4), (5) and
(6).
For this purpose several methods of solution were devised
for the turbulent diffusion equation (2) subject to a
steady concentration boundary condition at z = 0.
First a numerical procedure was developed for the case
of a steady velocity profile field under near neutral
conditions, which is homogeneous in the horizontal. In
the finite difference formulation of this type of equation
two difficulties arise caused by the automatic satisfac-
tion of one of the boundary conditions at the surface and
by the infinite size of the solution domain. The numerical
scheme was designed to overcome these difficulties by the
use of appropriate transformations. The results of
numerical experiments showed that the effect of the longi-
tudinal diffusion term and of the off-diagonal term are
both negative but usually negligible for evaporating or
cooling surfaces with an area of at least a few square
meters. It was also found that with suitable parameters
for roughness and stability the power law can be as useful
as the more complicated logarithmic law. The power law
was determined to be
11

-------
U = (5.5/7m) u (Z/Z 0 )m (9)
where m is close to (1/7), u the friction velocity and
z the roughness length. The details of this numerical
study are presented in publications 1 and 2 (Section IX).
A diagonalized form of equation (2) without the K terms
was also solved analytically for the same type of 9 oundary
conditions by means of a regular perturbation method with
a small parameter. This solution confirmed essentially
the conclusions of the numerical study since it showed that
the effect of longitudinal diffusion is negative but
negligible for evaporating surfaces larger than say 100 cm
in the direction of the wind. The effect of the lateral
diffusivity KJ 7 had already earlier been shown to be
negligible. This means then that for the problem of
evaporation or heat transfer from a surface at ground level
equation (2) (and thus equations (4), (5) and (6) as well)
may be simplified to
dC — a (K aC) (10)
dt 3z
in which the subscripts zz have been dropped from the
diffusivity. The details of the perturbation analysis have
been presented in publication 3.
Comparison of Solution of Equation (10) with Empirical Data.
Several empirical mass transfer formulae are available in
the literature for the determination of evaporation from
lakes and from pans at ground level. These were compared
with the theoretical solution of equation (10) for steady
conditions obtained by assuming the validity of Reynolds’
analogy and of the power law of equation (9) for the wind
profile. It was found that although this theoretical
solution, which was first given by Sutton, has its limita-
tions, it is a useful tool for the study of vapor or heat
flux under conditions of lateral advection. Conversely
the empirical formulae, in particular that proposed by
Rarbeck, should be quite adequate for many design problems
since their general functional form is compatible with
this simple theoretical model. For large water surfaces
the wind profile power of equation (9) appears to be
m = 1/8, and the roughness close to z 0 = 0.02 cm. The
details of this analysis are described in publication 4.
12

-------
Interaction between heat and vapor in forced convection .
After it had been determined that equation (10) with the
assumption of equation (9) can provide solutions in agree-
ment with experimental data, by analogy with (10) the
following simplified forms of equations (5) and (6), were
considered
a
U - (K -) (11)
U = .1 (K ) (12)
ax az
in which again equation (9) for U was assumed to replace
(4) and the Reynolds analogy was assumed valid. The link
between these two equations was assumed to be provided by
the boundary condition (7) whereas the other boundary
conditions were taken as those used in the Sutton solution.
The fact that equations (11) and (12) are linked leads
naturally to the necessity of solving them simultaneously.
The solution was obtained analytically by the construction
of Green’s function which, when incorporated in the bound-
ary condtions, produces two integral equations. These in
turn were solved by transformation into two algebraic
equations by means of the Laplace Transformation. The
results of this solution showed, how for a simple case
sensible heat and water vapor transfer and also the water
surface temperature may be expected to depend on the meteo-
rological conditions and on the rate of change of energy
content of the water body and the surrounding land. Due
to advection the water surface temperature and the turbu-
lent fluxes vary in the downwind direction. In contrast
to what is sometimes assumed for a deep lake the surface
temperature is usually quite different from the wet-bulb
temperature. The main conclusion of this study was, that
the error introduced by the use of an experimentally
obtained average temperature, in the calculation of evapo-
ration and turbulent sensible heat transfer, is probably
quite small, and negligible for practical purposes. Of
course, this provides one more indication of the validity
of Sutton’s solution and Harbeck’s formula in situations
where horizontal momentum advection is negligible. The
details of this interaction analysis are given in
publication 5.
13

-------
Investigation of Various Semi-Em iricaI Formulations of
the Turbulence . Because the validity of the theoretiáal
calculations carried out by means of the solution of
equations (4), (5) and (6) is so critically dependent on,
first, the mathematical suitability, and second the
corespondence to real physical mechanisms of the semi-
empirical turbulent diffusion models under various con-
ditions, this aspect of the problem also had to be
investigated.
In the solutions described above intensive use was made
of the power law for the wind profile as given by equation
(9). Actually,this law has been widely used in the
solution of numerous turbulent diffusion problems.
However, because it is admittedly only approximate and
without much of a theoretical model to justify it, such
solutions have not always received the attention they
deserve among atmospheric scientists. Nevertheless, as
pointed out above, evaporation under neutral conditions
can be calculated by means of the power law equally well
as by means of the logarithmic law. Therefore it was
felt desirable, that the power law also be made applica-
ble under non—neutral conditions. This was accomplished
by determining C and m in equation (9) as functions of
the roughness and of the (Monin-) Obukhov stability length
L, through comparison with the log—linear law for stable
conditions, and Obukhov’s formula (KEYPS) for unstable
conditions. As was to be expected, for neutral conditions
the results show that m and C are approximately equal to
1/7 and 6, respectively. C increases (decreases) slightly
whereas m decreases (increases) with decreasing (increas-
ing) stability. The results are tabulated in publication
6.
The next step in this research consisted of evaluating
the effect of certain important factors on the wind profile
that have usually been neglected in the past, namely
radiative and latent heat transfer, the molecular viscosity
of the air and the molecular diffusivities of heat and
water vapor in the air. The basic theoretical approach
was analogous to Obukhov’s. However, the (Monin-) Obukhov
stability length was adjusted to obtain a more nearly
height-independent similarity. Even though the Obukhov
model contains some strong assumptions, the inclusion of
the effects of radiation and of the molecular diffusivities
on the dynamics of the flow yields a profile law, which in
contrast to other previously proposed laws is valid over
the whole range of atmospheric stability conditions over
14

-------
rough as well as over smooth surfaces. For example, the
resulting solution appears to agree well with Obuknov’s
(KEYPS) profile for unstable conditions or neutral con-
ditions over a rough wall or over a smooth wall at some
distance from the wall. Very close to the smooth wall,
where viscosity is taken into account, the resulting
solution reduces to the linear profile of the laminar
sublayer. Hence it does not have a singularity right
at z 0 for a smooth surface with z 0 = 0, and it is in
good agreement with the generally accepted concept of
the “law of the wall”. Under stable conditions the
resulting solution is in good agreement with the experi-
mentally verified log-linear law. This seems to indicate,
that under unstable conditions the effects of radiation
and of the molecular viscosity and diffusivities tend to
cancel each other and that very little error is intro-
duced if both are neglected. However, under stable con-
ditions these effects probably reinforce each other and
must be considered. This new formulation of the eddy
diffusivities and the wind profile should be especially
useful in turbulent transfer calculations and numerical
modeling of the atmospheric surface layer on digital
computer. The initial concept of introducing the molecu-
lar viscosity into the mixing length under neutral con-
ditions is described in publication 7 and 8.
The project dealt primarily with turbulence as an
efficient mechanism of heat and vapor transport. Because
radiative transfer has for some time been known to play
an important role in the maintenance of turbulence under
inversion conditions very close to the critical Richardson
number this matter needed further investigation. Thus an
expression was derived relating the critical flux
Richardson number with the critical (gradient) Richardson
number. In contrast to an earlier analysis by Townsend,
which was restricted to the atmosphere well outside the
earth’s boundary layer, the treatment of the problem was
applied specifically to turbulent transfer near the
surface and it took account of the effect of evaporation
on the stability. The effect of radiation on the rate of
destruction of the mean square of the temperature was
obtained by considering the radiative flux divergence in
a stratified atmosphere and by using a simple functional
relationship to represent empirical emissivity data. It
was found that under conditions normally encountered in
the lower atmosphere, the effect of evaporation and
radiation on turbulence is not negligible. The critical
Richardson number, which is the criterion for the
existence or cessation of turbulence, is sensitive to
15

-------
both factors. There is no definite critical Richardson
number but it falls in a range between 0.25, below which
turbulence is highly probable, and somewhat larger than
0.5, above which turbulence is unlikely, This critical
Richardson number can be expressed in terms of evapora-
tion, radiation and the ratio (aw/u ) of the variance of
the vertical velocity fluctuations and the friction veloc-
ity; this ratio, (a /u ) in turn, also appears not to have
a definite value. Evaporation and radiation cause it to
be larger than unity under neutral conditions. These
results show that the proposed model for turbulence and
radiation together with the assumption Kh = K 1 T I is consis-
tent with, or at least not contradicted by, the available
experimental evidence. However, more research is needed
in this area. The details of this analysis have been
given in publication 9.
Finally, there has been increasing evidence of the fact
that there exists a strong interaction between the turbu—
lent flow of the air and motion of the underlying water
waves. Because this interaction is a highly nonlinear
process with a complicated feedback mechanism it is still
not very well understood. Therefore, first a review was
made of recently published data. It appeared that Volkov’s
measurements taken on the Western Mediterranean were
especially suitable. By assuming a similarity function
suggested by Kitaygorodskiy a relationship was derived
between the actual eddy diffusivity Km and the apparent
eddy diffusivity, K , estimated purely from the wind
velocity profile. ‘rhis was done in terms of (c/ut) where
c is the phase velocity of the dominant wave and u the
friction velocity. It was found that the eddy diffusivity
near the air-water interface, determined in the usual way
from the mean velocity profile, may lead to considerable
error depending on the sea state. Thus measurements of
the mean wind velocity or humidity profiles lose a great
deal of their usefulness unless additional measurements
are made of the sea—state. At any rate, a knowledge of
the similarity function allows a correction to be applied
to the turbulent flux obtained by means of the bulk
aerodynamic method on the basis of profile measurements.
This matter is described in detail in publication 10.
In publication 11, a theoretical derivation is presented
of this similarity function for the wind profile gradient
and for the standard deviations of the velocity fluctua-
tions. The model is mainly based on the assumption that
the relative speed of the air with respect to the celerity
of the waves is the most important factor governing the
16

-------
magnitudes of the wave induced velocity fluctuations.
Additional assumptions are that the mean wind profile is
logarithmic and that the Reynolds stresses are constant
in the vertical. The results were in satisfactory
agreement with available experimental data for the case
of decreasing swell. For developing swell further
research will be required.
The Computation of Evapotranspiration from Large Land Areas .
In the determination of heat and vapor losses from a lake
or reservoir it is desirable that the evapotranspiration
from the surrounding land be known. Most methods, that are
available at present, and that are sufficiently practical
to be applied with normally observed meteorological data,
are almost purely empirical without a firm physical base
and quite limited because they were developed for special
hydrologic conditions. There are equations that are
based on a more rational approach, but they provide only
the vertical vapor flux at a given point where they require
a detailed knowledge of the energetics and of the turbulence
structure of the atmosphere.
Therefore, a practical, yet rational procedure was developed
for computing regional evapotranspiration. The assumptions
underlying the model are essentially those used in the
one—dimensional mass transfer methods, but the shear stress
is determined by means of the equation i P C Va o
which p is the density of the air, C is the ge sttophic
drag coefficient as first proposed b Lettau and Vqo the
surface geostrophic wind. The main advantage of tF e method
is that all data needed in the computation, namely the
surface temperature, the pressure of the air, the specific
humidity at two levels in the atmospheric boundary layer
and the surface geostrophic wind speed can be obtained
from standard meteorological data such as Climatological
Data, Daily Weather Maps and Rawinsonde Data which are
published regularly. Generally, good agreement was
obtained between mean monthly evapotranspiration computed
for a number of stations and the corresponding evaporation
data obtained from Class—A pans. Although this testing
is obviously not conclusive, because of the limitations
inherent in pan evaporation to determine actual regional
evapotranspiration, it does offer some encouraging evidence
of the soundness of the model. The details of this study,
with a step by step procedure, are given in publication 12.
17

-------
Heat and Water Vapor Exchange Under Conditions of Unstable
Atmospheric Stability and Momentum Advection . In the study
of the interaction between heat and water vapor in forced
convection described above (reference 5), it was assumed
that the momentum transfer is strictly vertical. Therefore,
equation (4) could be replaced by a wind velocity profile
law U = U ( z), W = 0, and equations (5) and (6) could be
simplified to (11) and (12). Actually in this earlier
study a phenomenon, in which as a result of the sudden
temperature jump at x = 0, the stability of the air would
definitely vary in the x-direction, was treated as one in
which the stability would affect the flow only in an aver-
age way. On a long term basis, as the agreement with the
experimental data shows, the effect of this temperature
disontinuity on the advection of momentum is probably
negligible. However to determine the effect of this
discontinuity under extreme conditions it was decided to
consider the case in which the air moves from over a
relatively cool land surface onto a warmer water surface.
This situation, which may create severe instability in
the air is not unconmon, for example, over a deep lake in
the fall and early winter months.
In recent years the problem of momentum advection in the
atmospheric surface layer has received considerable atten-
tion. A very thorough literature survey on the problem of
a sudden (step—) change of the surface roughness on the
air flow under neutral conditions revealed that very close
to the surface the flow adjusts rather quickly to the new
surface, but that it is still difficult to predict the
flow in the higher regions of the boundary layer. This
means that a thorough understanding of the transport mecha-
nisms in the iumtediate neighborhood of the surface is of
utmost importance in the determination of the vertical
transfer rates. Although there are several methods for
analyzing the change in roughness problem, the method
consisting of solving equation (4) but with only the ver-
tical diffusion term, namely
(13)
in conjunction with (3), seemed to be preferable. The
reasons for this are mathematical and computational
simplicity and the fact that the mixing length models are
still adequate within the surface layer.
18

-------
Thus it was concluded that the problem of simultaneous heat
and water vapor transfer under conditions of momentum advec—
tion, as a result of a temperature and humidity discontinu-
ity at x = 0, may be described by equations (13), (3) and
simplified forms of (5) and (6), namely
U + W = (Kh ) (14)
u 2. + w = (K” 2.) (15)
3z z
Note that in contrast to equations (11) and (12) the terms
in W are retained because of the effect of momentum advec—
tion. The diffusivities were assumed to be given by empiri-
cal but mathematically convenient functions of the stability
of the Dyer - Businger type.
Numerical solutions of this set of equations showed that the
sudden instability of the flow initially can greatly increase
the efficiency of the turbulent exchange processes at the
surface and that the vapor flux can greatly contribute to
the atmospheric stabiltiy. The solution was shown to be
relatively insensitive to the exact form of the Monin-Obukhov
similarity functions or the exactness of Reynolds’ analogy or
von Karman’s constant. In the formulation of the problem it
is essential to keep the eddy diffusivities bounded at the
value corresponding to the top of the atmospheric surface
layer. Blackadar’s modification of the mixing length pro-
vides a suitable way of accomplishing this. For large
fetches and near neutral conditions the solution becomes
similar to that obtained in Sutton’s problem and with
Harbeck’s empirical formula. The details of this study are
given in Section X, Appendix.
19

-------
SECTION VI
PREDICTION OF EVAPORATION FROM A LARGE WATER BODY
In view of the theoretical findings described above, evapo-
ration from a free water surface can be calculated, with
a degree of accuracy sufficient for most engineering
problems, by means of the following equation
E=Ea+buz(eo_ea) (16)
where E is the average evaporation in cm/sec, e 0 and ea
are the vapor pressure in millibars at the water surface
and in the air unaffected by the water surface (i.e. upwind
from the water body), respectively, U is the wind speed
in cm/sec at an elevation z, b is the so—called mass trans-
fer coefficient and Ea is the evapotranspiratiOn from the
land surrounding the water body. Similarly sensible heat
flux in(cal/cm 2 sec) may be estimated from
H = Ha + (b p c /O.622) U o — Ta) (16k)
where Cp = 0.24 (cal/g °K) is the specific heat at constant
pressure, T and Ta the temperature of the water surface
and the air, respectively, in °C and p the pressure in mb.
From a practical point of view the main problem is the
determination of the coefficient b. In general, it is a
fairly complicated function of the stability of the air,
the direction of the wind, the topography and the nature
of the surrounding land, the size and the geometry of the
evaporating surface, and numerous other factors which can
only partly be included in the presently available
theoretical models. Therefore, as it is very difficult to
ascertain b on a purely theoretical basis, for a particular
water body, lake or even stream, it is probably best
determined experimentally. Harbeck and Meyers (1970)
suggested that this be accomplished by means of the energy
budget method, which albeit more accurate can only be
used for a short intensive measurement period on account
of the generally more expensive equipment and data process-
ing involved.
21

-------
When, however, b cannot be determined experimentally, for
averaging periods of a few weeks to a month Harbeck’s
(1962) coefficient should be quite adequate. Thus equation
(16) can be applied with U at 2m above the surface and
with (Brutsaert and Yeh, 1J70)
b = 5.32 10 A 0 ’° 5 (17)
where A is the water surface area in cm 2 . This equation
was shown to agree with the forced convection model under
near neutral conditions at 20°C over a surface with a
roughness of approximately z 0 = 0.02 cm.
Equation (16) with (17) takes account of the effect of
horizontal advection of the wind, since ea refers to the
humidity upwind from the water body and since A is included
in the formulation. Sometimes, when horizontal advection is
small and under near neutral conditions the local evapora-
tion can be calculated by means of the drag coefficient on
the basis of local variables, namely
B = CE p U (g 0 - (18)
where q is the specific humidity at the water surface and
q 2 at tEe level z. The local heat flux may be estimated
by the analogous equation
H = CE p c U - T ) (18k)
For z = lOm, and z = 0.02 cm one obtains by means of the
logarithmic wind p%file CE 1.36 l0 . This is in
excellent agreement with the value 1.2 i0 recently
observed over the ocean east of Barbados (Pond, et al . ,
1971). For an assumed mean temperature of 20°C, an air
pressure of 1013.2 nib and a density of the air p = 1.20 10
g/cc equation (18) can be written in a form similar to (16),
i.e. approximately
E = 1.0 10 U 10 (er, — e 10 ) (19)
22

-------
where now U 10 and e 10 are the wind and the local vapor
pressure at 10 m. above the surface.
Note finally that the variables used in any of these
transfer equations should preferably be means taken over
periods of a day or less. If, as is often done, monthly
means are used to calculate the monthly mean evaporation,
serious errors may result. It can easily be shown that
the difference between the mean of the products of two
variables and the product of the means equals the covari—
ance between the two variables. For example Jobson (1972)
made an analysis of this type of error by analyzing data
obtained at Lake Hefner, Oklahoma, and found that using
daily means resulted in an error of 0.088% with a variance
of 29%, and monthly means in an error of 5.7% with a
variance of 30%.
23

-------
SECTION VII
ACKNOWLEDGEMENTS
The preparation of this report was arranged through the
School of Civil and Environmental Engineering, Cornell
University. Prof. C. D. Gates, Head of the Department
of Water Resources Engineering was Acting Project
Director from September 1969 to August 1970 during the
absence of W. H. Brutsaert.
As listed in Section VII, Publications, Dr. G. T. Yeh,
Dr. R. N. Weisman, Dr. I. M. Cheng and Mr. J. A. Mawdsley
contributed substantially to the project. Mr. El-Sahragty
performed a literature survey and some initial calculations.
The support of the project by the Environmental Protection
Agency, and the help provided by Messrs. W. A. Cawley,
D. G. Stephan, B. F. LaPlante, D. Oyster, W. J. Lacy and
B. A. Tichenor, the Grant Project Officer, is acknowledged
with 5incere thanks.
25

-------
SECTION VII
REFERENCES
Brutsaert, W. On the anisotropy of the eddy dIffusivity.
Jour. Meteorol. Soc. Japan, Ser. 2, 48, 411-416, 1970.
Brutsaert, W. and G. T. Yeh. Implications of a type of
empirical evaporation formula for lakes and pans. Water
Resources Res. 6, 1202—1208, 1970.
Harbeck, G. Earl, Jr. A practical field technique for
measuring reservoir evaporation utilizing mass-transfer
theory. U. S. Geological Survey Prof. Paper 272-E, 1962.
Harbeck, G. Earl, Jr. and J. Stuart Meyers. Present day
evaporation measurement techniques. Jour. Hydraul. Div.
Proc. ASCE, 96 (HY7), 1381—1390, 1970.
Jobson, Harvey. Effect of using averaged data on the
computed evaporation. Water Resources Res. (A.G.U.) 8,
513—518, 1972.
Pond, S., G. T. Phelps, 3. E. Paquin, G. McBean and
R. W. Stewart. Measurements of the turbulent fluxes of
momentum, moisture and sensible heat over the ocean.
Jour. Atmos. Sci. 28, 901—917, 1971.
27

-------
SECTION IX
PUBLICATIONS
1. Yeh, T. T. and W. Brutsaert, A numerical solution of
the two—dimensional steady state turbulent transfer
equation. Monthly Weather Review 99, 494-500, 1971(a).
2. Yeh, G. T. and W. Brutsaert, Sensitivity of the solution
for heat flux or evaporation to off—diagonal turbulent
diffusivities. Water Res. Research 7, 734—735, 1971(b).
3. Yeh, G. T. and W. Brutsaert, Perturbation solution of
an equation of atmospheric turbulent diffusion. Jour.
Geophys. Res. (Oceans and atmos.) 75, 5173—5178, 1970.
4. Brutsaert, W. and G. T. Yeh, Implications of a type of
empirical evaporation formula for lakes and pans. Water
Res. Research 6, 1202—1208, 1970.
5. Yeh, G. T. and W. Brutsaert, A solution for simultaneous
turbulent heat and vapor transfer between a water sur-
face and the atmosphere. Boundary Layer Meteor. 1,
106—124, 1971(c).
6. Brutsaert, W. and G. T. Yeh, A power wind law for tur-
bulent transfer computations. Water Resources Res. 6,
1387—1391, 1970.
7. Yeh, G. T., A note on a universal velocity profile.
Jour. Hydraul. Engg., Chinese Inst. Hydraul. Engrs.,
Taipei, 15, 1972.
8. Yeh, G. T., Unified formulation of wall turbulence.
Jour. Hydraul. Div., Proc. ASCE, 98, (HY12), 2263-2271,
1972.
9. Brutsaert, W., Radiation, evaporation and the mainte-
nance of turbulence under stable conditions in the
lower atmosphere. Boundary Layer Meteorology 2, 309—
325, 1972.
10. Cheng, 1—Ming and W. Brutsaert, Wave effect and eddy
diffusivity in the air near a water surface. Water
Resources Research, ! 1439-1443, 1972.
29

-------
11. Brutsaert, Wilfried, Similarity functions for turbulence
in neutral air above swell. In preparation, 1973.
12. Mawdsley, John A. and Wilfried Brutsaert Computing
evapotranspiration by geostrophic drag concept. Jour.
Hydraul. Div., Proc. ASCE, 99, (HY1), 99—110, 1973.
30

-------
SECTION X
APPENDIX
Heat and Water Vapor Transfer from a Water
Surface under Unstable Atmospheric Conditions
by
Richard N. Weisman and Wilfried Brutsaert
Formulation of the Model . The problem under consideration
is one of steady turbulent flow of air from over a homo-
geneous land surface but with a different temperature and
specific humidity. The mean motion over the upwind land
surface is plan—parallel and a function of elevation z
only. Furthermore, turbulent diffusion in the direction
lateral to the mean wind is known to be negligible in the
calculation of evaporation from large water surfaces
(Brutsaert and Yeh, 1969). Therefore, the problem may be
reduced to two dimensions like flow over an infinitely
long strip at right angles to the mean wind. As the air
moves over this surface discontinuity in temperature and
humidity, an internal boundary layer develops within
which the mean wind is no longer strictly horizontal as
it is upwind. The horizontal component, u, and the
vertical component, w, of the mean wind may be related
through the continuity equation (3).
Under conditions of fully turbulent flow, molecular
diffusion is negligibly small. Since the effect of tur-
bulent diffusion in the direction of the wind has also
been found to be negligible, sensible heat and water vapor
satisfy the conservation equations (14) and (15), respec-
tively.
Evaporation and heat flux from the water surface are con-
trolled primarily by the motion in the atmospheric surface
layer. Within this layer the effects of the pressure gradi-
ent and of the earth’s rotation, producing the Coriolis ac-
celeration, are neglibible. Also, near the surface the
motion is very close to horizontal, so that w is probably
at least an order magnitude smaller than u. Thus it is
permissible to neglect vertical momentum in accordance with
Prandtl’s boundary layer assumption. This means that
horizontal momentum takes on the properties of a conserva-
tive scalar admixture of the flow and that the horizontal
Reynolds equation reduces to (13) which is similar to
(14) and (15). This also implies that no consideration
31

-------
is given here to convective cells or return currents aloft
associated with heat islands or ineso-scale sea breeze
problems. In orther words, buoyancy effects are not ana-
lyzed by a vertical equation of motion; rather they are
considered only as they affect the vertical profile of
the mean horizontal wind and the eddy diffusivities by
means of the Obukhov stability length.
In the system of equations (3), (13), (14) and (15) use
was made of the semi-empirical theory of turbulence. How-
ever, there are numerous ways to express the required eddy
diffusivities in terms of the mean quantities. Monin and
Obukhov’s (1954) semi-empirical similarity model for the
lower atmosphere suggests that they be of the form
K ” 1 = ii (Al)
Kh , L/H (A2)
Km = U t/M (A3)
where u = It/p is the friction velocity, 2.. the mixing
length or a length scale characteristic of the turbulence
and the •‘s are supposedly universal functions of (z/L).
The parameter L is Obukhov’s (1946) (Businger and Yaglom,
1971) stability length.
To apply equations (Al), (A2), and (A3), suitable expres-
sions must be adopted for the quantities u , 2., L and the
•‘s. Thus, in contrast to the similarity models for turbu-
lent flow over a homogeneous surface in which the friction
velocity refers to the shear stress at the surface, herein
u = u (x,z) refers to the local shear stress. On account
of (A3), this local friction velocity can be obtained from
= 2.1 ’M u/ z (A4)
Close to the surface the mixing length is usually defined
as 2. = kz or also as £ = k(z + z ) where k is von Karman’s
constant and z 0 the roughness leRgth. As will be shown
below, the upper boundary condition at any x consists of
the attainment of the upwind conditions aloft at some
undetermined elevation outside the developing internal
boundary layer. Therefore to prevent the increase of 2.,
and thus the eddy diffusivities, without bounds,
32

-------
Blackadar’s (1962) modification was adopted, viz.
2. = k(z + z 0 )/ [ 1 + k(z + z)/X] (A5)
where A is some length scale, which constitutes the upper
limit for 2.. As originally proposed by Blackadar (1962),
the expression for A is
V
A = .0027 g
fz
0
where Vg is the geostrophic or free stream wind above the
boundary layer, f is the Coriolis parameter which is a
function only of the earth’s rate of rotation and the
latitude. Taylor (1969) proposed an expression for the
scaling height identical to Blackadar’s except with the
constant equal to .0004. If f = 10 4 /sec, z 0 = .02 cm,
and Vg = 500 cm/sec, the dimensionless scaling height has
a value of around 7 x l0 . Values of 5 x l0 and 7 x lO 4
were used for most calculations and they worked well as
will be discussed further.
The evaluation of the effect of the density stratification
due to water vapor was one of the main objectives of the
study. The Obukhov stability length concept was extended
to also include water vapor concentration gradient as
worked out by Zilitinkevich (1966), namely
3
-u p
kg ( (H/T c ) + 0.61 E ]
where T is some reference temperature, taken as the upwind
surface temperature.
Although the exact form of the Monin-Obukhov similarity
functions is still open to much conjecture (e.g. Monin and
Yaglom, 1971), it was decided to adopt a recent formulation
which has a good experimental basis in the absence of
advection and which has the advantage of mathematical
simplicity. The functions were suggested by Businger
(1966), Dyer (1967) and then Dyer and Hicks (1970) as
follows
33

-------
= (1 - B ) 1 ’ 14 (A7)
= = (1 - B ) ” 2 (A8)
where B, a constant, has been measured to be approximately
16. It should be noted that as the form of these 4)-f unc-
tions is the subject of continued controversy, it is not
the purpose of this study to give credibility to one
formulation rather than to another. Nevertheless, it is
one of the objectives to test the sensitivity of the
results to various values of the parameter B and also of
the power in (M). Therefore the following alternative
formulation is also considered.
= (1 - B ) 1 ’ 12 (A9)
in which the power is in closer agreement with -0.46
obtained experimentally by Morgan, et al. (1971), and
which is in accordance with Reynold’s analogy •M •H
Equations (3), (13), (14) and (15) are subj ect to boundary
conditions that may be specified from the following consid-
erations:
(i) the incoming air is neutral and it has a humidity
profile, g (z), which is known and which results from equi-
librium conditions dicated by the heat flux and the evapo-
transpiration from the windward land surface. The neutral
velocity profile upwind from the water body is given by
__ Z+Z 0 kz
u = u(z) = k (in ) + (AlO)
where the subscript a refers to upwind conditions. Equa-
tion (AlO) results from equation (A3) for a constant fric-
tion velocity and for Blackadar’s (1962) mixing length of
equation (A5).
(ii) High above the water surface, that is at an elevation
well above the height zh of the internal boundary layer
the conditions are the same as upwind and unaffected by
the water body.
34

-------
(iii) at the water surface the temperature is known and
uniform; the specific humidity is saturated and a unique
function of this temperature.
These boundary conditions can be written as
x = 0; z > q =
T T (z)
a
u = u(z) (All)
x > 0; z > Zh; q = q (z)
T T(Z)
u U(Z) (A12)
x> 0; z0; q=q
T=T
u=w=0 (A12)
where the subscript w refers to conditions over the water
and the superscript o to surface conditions at z 0.
Solution of the System . In order to give the equations
and their solution a more universal meaning, they are made
dimensionless with the following non-dimensional quantities:
= U/Ut, % = W/U*, =
= x/z 0 , = z/z 0 , A = A/z 0 , £ = Liz 0 , =
T- T _____
0 0
Tw Pa
H E
PCpU*a (T - T ) - P l*a - q )
Also, the vertical axis is transforented logarthmica].ly to
= in ( + 1) as suggested by Taylor and Delage (1971).
35

-------
Equations (3), (13), (14), (15) and (A5) to (A8) become
A A A
- + e = 0 (A14)
ax
A A A A A
U + we = e (AlS)
A A A A
e — - (A16)
ax ac
A A A A A
A39 A_C9 _ aH
u A+we e — (A17)
A A
ax ac
A
A A
A te au
(A18)
U*$ 4
A
C A
A Lute
E= — (A19)
A A
A
A Lute ae
____ — (A20)
*v aZ
A
A keC
(A21)
1 + keC/A
A
(1 B Z)_1/4 (A22)
IJ
A
$ = $ = (1 — B ) 112 (A23)
V H
L
36

-------
A (A24)
A*H+B E
The boundary conditions, equations (All), (A12), and (A13)
are now
A A A A
x = 0; 1 > 0; Q ( ) = 0
A A
o R) = 0
A A
u R) = + (A25)
>° > h’ 0 (i;) =0
A A
o ( ) = 0
A
u (ç) = + (A26)
A
A A A
x>o; 0; 0=1
6%
6=1
A
u=w=0 (A27)
The iat, A, denoting non—dimensional quantities is now
dropped for convenience.
Two ØJ.mensionless parameters appear in the stability length,
L, a a result of the non—dimensionalization process, namely
0
T -T kgz
= — a ) 2 (A28)
T
kgz
B = -.61 - q ) (A29)
• U
*a
37

-------
These two parameters must be specified to solve the system
of equations. Both parameters contain a term resembling an
inverse Froude number, kgz 0 /u providing some indication
of the relative importance of gravity and inertia effects,
in this case, free and forced convection. A contains a
fractional change in surface temperature while B contains
the change in surface specific humidity at the leading edge.
A and B are referred to as the temperature change stabil-
ity parameter and the specific humidity change stability
parameter, respectively. Both represent a measure of the
discontinuity at the leading edge.
As an illustration, consider a surface with a roughness
length of = .02 cm. If the wind speed at 2 m is 100
cm/sec, the logarithmic velocity profile yields a friction
velocity of 5 cm/sec. An average temperature change at a
leading edge might be taken as .05 and, if the upwind con-
dition is dry, the surface specific humidity jump might
be .03 kg/kg. Thus, putting the above values into equations
(A28) and (A29) yields the average values: A -.015 and
= —.005.
In the numerical solution, A was given values equal to 0,
—0.001, —0.01, —0.05, and —0.1 and B was given values equal
to 0, —0.001, -0.003, and —0.01. The higher values of A
and B correspond to larger surface temperature and specific
humidity changes at the leading edge, to a lower upwind
velocity, or to a higher roughness value. These ranges of
values for the stability parameters are probably larger
than would be met in nature and still be applicable to the
problem considered here.
The equations (A14) to (A24) with boundary conditions (A25),
(A26), and (A27) were solved using an explicit Euler scheme.
More sophisticated and accurate methods would incur an
increase in computer time and probably are not warranted for
the complex problem considered here. For a complete descrip-
tion of the numerical procedure, the reader is referred to
Weisman (1973).
Central differences were used for vertical gradients in all
the equations except in the continuity equation (A14) which
was formulated using backward differences as suggested
first by Peterson (1969). A difficulty arises in using
central differences because it thereby becomes impossible
to define the gradients and, therefore, the fluxes at the
surface.
38

-------
By making certain gross assumptions, equations (A15), (A16),
and (A17) can be simplified to heat-type equations for
which an approximate stability analysis can be made. The
maximum allowable downwind step size derived from the
stability criteria is a function of the velocity at the
first vertical grid point and the distance to that grid
point from the surface. For reasonable vertical grid spac-
ings, the required horizontal step size is quite small which
necessitates an unreasonable amount of computer time to
move a physically relevant distance downwind. This and the
aforementioned difficulty were resolved with the so-called
constant flux wall layer.
In this scheme it is assumed that a thin layer of air exists
next to the surface in which the vapor viux, Ew, heat flux,
H. , and friction velocity, u , are constant with height.
This constant flux layer was used by Peterson (1969) and
Taylor (1970); Taylor and Delage (1971) have discussed this
method.
The first difficulty, namely that of defining the fluxes at
the surface, is resolved since the flux equations (A18),
(A19), and (A20) can be integrated analytically to obtain
the profiles of humidity, temperature, and velocity within
the constant flux layer; the fluxes may then be calculated
by knowing the values of humidity, temperature, and velocity
at the top of the layer. These values were obtained from
the Euler method.
If now the height of the constant flux layer is allowed to
grow as the internal boundary layer grows, the second diffi-
culty, that of small forward step sizes, is also resolved.
Thus, the forward step sizes can be increased considerably
allowing a gradually faster numerical solution as it moves
downwind. Note also that though there is little physical
justification for a constant layer under advecting condi-
tions, it can be shown from equations (A18), (A19) and
(A20) that the fluxes become constant as z approaches zero.
Moreover, as the air mass moves downwind, equilbrium condi-
tions become established in which horizontal gradients tend
to zero and the vertical fluxes are constant with height.
Results . Averaging the local vapor flux at the surface,
E ,(x), over the fetch of the water, x 0 , one obtains for
the total average evaporation
39

-------
1
W= — / E ,dx
00
Since the extreme case is the most interesting one, the
results presented are for upstream conditions neutral
and dry. For these conditions, the total average heat
flux is defined by
x
1 °
HF=— I Hdx
xo 0
Similarly, dimensionless temperature, 0, and specific
humidity, Q, are numerically equal for neutral, dry upwind
conditions.
Figure 1 is a curve of total average vapor flux, W, versus
the fetch of the water surface, x 0 , for various values of
the A stability parameter while B is held at a constant
average value. Figure 2 is the same as figure 1 except
that A is kept constant at an average value and B is the
parameter. Figures 3 and 4 plot the above results for A
equal to zero and B varying, and equal to zero with A
varying.
As can be seen from the figures 1 through 4, the total
average evaporative flux decreases as the fetch increases.
Also, evaporation is greater for the larger values of A
and B ; that is, when the instability created by the
temperature and/or specific humidity change at the leading
edge is large, the evaporation is large. The effect of A ,
the surface temperature stability parameter, is quite sig-
nificant. The effect of the surface specific humidity
change parameter, B , is smaller than that of A , but there
is still a twenty percent increase in evaporation at
= i0 7 from B = -.001 to B = -.01 as seen in figure 2.
From the effect of B , the specific humidity stability
parameter, on the evaporation, as seen in figures 2 and 4,
it can be stated that the concentration buoyancy term in
the modified Obukhov length, equation (A6), has a signifi-
cant impact on the momentum exchange. Hence, even when no
temperature discontinuity occurs between the land and
water surfaces, the air flow does not remain neutral; rather
the flow becomes unstable and it does affect the evaporative
flux significantly.
40

-------
Figure 1: Dimensionless vapor flux, W, versus dimensionless fetch,
x 0 , for various values of stability parameter A , and B constant.
100
—1
10
I - ’
io2
io2
lOLl
106
io8
xo

-------
i i I T 11F T 1 11
A -O.O1
vv
B; : O.O1
B -O•OO 1
-2 I 11 ________ — I I JI 1 1 Ii 1 t 1 Il U
10 10 10
Figure 2: Dimensionless vapor flux, W, versus dimensionless fetch, X 0
x 0 , for various values of stability parameter B , and A constant.

-------
‘1 I I
B O
-1 _______ -___________ A -O.O1
10 /
/ ,
-Q.QQ1
-21_IIJ I LH I IHI ____ IH
io2 io6 io 8
Figure 3: Dimensionless vapor flux, W, versus dimensionless fetch,
x 0 , for various values of stability parameter B , and A zero.

-------
100
10’
108
w
Ax O
-0.01
B 0
1
io2
-0.001
io 4
Figure 4: Dimensionless vapor flux, W, versus dimensionless fetch,
x 0 , for various values of stability parameter B , and A zero.
io6
xo

-------
Both figures 3 and 4 show that, by allowing momentum to be
advected, the solution is greatly affected. When A and B
are both zero, the modified Obukhov length, L, goes to
infinity and all the 4’s become unity. Momentum is no
longer advected and the velocity profile remains unchanged
downwind. Note that the case A = B = 0 is the Sutton-
type problem. The main difference with the Sutton solution
is that the present case is not solved for the power law
of the wind profile but for Blackadar’s log-linear profile
equation (AlO). The related solution with the logarithmic
profile has been treated by Yeh and Brutsaert (197la).
Figure 5 is a plot of total average evaporative flux, W,
versus fetch, x 0 , for smaller evaporating surfaces. One
curve represents the calculation performed with
PM = (1 - B z/L)l/ 2 , equation (A9), instead of
= (1 - B z/L)]-/ 4 , equation (A7). As can be seen in
figure 5, the solutions tend to merge at values of the
fetch, x 0 , greater than about lOs; for larger values of x 0 ,
no difference in the curves can be seen between solutions
using the —1/2 exponent as compared with the -1/4 exponent
in the •M function. Although only a few runs were made
using the —1/2 exponent and compared with the solution
using the —1/4 exponent in ctM, it is probably safe to infer
that no difference would occur in the solutions for any A ,
B combination above x 0 = io .
A reason for the above result is that, beyond the initial
high instabilities at the leading edge, the stability con-
dition downwind approaches neutral quickly and the velocity
profile changes only slowly. This can be seen by expand-
ing both equations (A7) and (A9) in a Taylor series. Thus,
for small z/L, the two expressions for M are about equal.
If z/L is large (at the leading edge), the two expressions
for 4 M are different and the results are significantly
affected by the form of the dimensionless shear stress.
Figure 6 is a graph of surface shear stress, Tw, versus
downwind distance, x. As can be seen, the shear stress
increases downwind, greatly for the more unstable condition
of A* = —.1 and less so for smaller A . The shear stress
should not increase indefinitely. Taylor (1970a), using
Prandtl’s mixing length, 2 .. = k (z + z 0 ) , and giving no
consideration to the problem of vapor and its effect on
stability, solved a similar problem to the one considered
here and found that the shear stress did increase indef i-
nitely. By using Blac]cadar’s mixing length, equation (A5),
45

-------
w ___ __
a’
ic 10
Figure 5: Dimensionless vapor flux, W, versus dimensionless
fetch, x 0 , for given stability parameters, A and B , and for
two ratios of
101
1 T
A -O.O1
pv
B -O.O
(1-B
)
m
101
io2
x

-------
I I F
io 6
Figure 6: Dimensionless surface shear stress, r ,, versus dimension-
less downwind distance, x, for various values of stability parameter
A , and B constant.
x
1.50
too
io 2
108

-------
it appears from figure 6 that this difficulty is partially
resolved since the shear stress curves have an inflection
point and might reach an asympotic value far downwind.
Nevertheless, also in the present solution the shear stress,
does not seem to tend to unity for large x whereas the
w nd again assumes its upwind neutral profile.
In figure 7, the negative inverse of the surface Obukhov
length, that is an Obukhov length with the surface fluxes,
is plotted against downwind distance for several stability
conditions. This curve would indicate that the air mass
is returning to neutrality as it moves downwind.
Figures 8, 9, and 10 show the effects of changing various
physical and numerical constants on figure 1, i.e. total
average vapor flux, W, versus the fetch, x 0 . Recently,
Businger et al. (1971) have concluded from experimental
wind and temperature profiles under a broad range of
stability conditions that von Karman’s k has a value of
0.35. As can be seen in figure 8 the values of vapor flux,
W, using k = 0.35 are only slightly lower than the results
using k = 0.4. On the same curve, points are shown repre-
senting calculations made with various initial vertical
step sizes, h, in the numerical solution. No change in
the solution is visible.
Figure 9 is a curve of total average vapor flux, W, versus
the fetch, x 0 , for specific average stability parameters
in order to test the sensitivity of the solution to the
value of B, the Businger—Dyer constant in equations (A7),
(A8), and (A9). The solution is quite insensitive to the
value of B; as can be seen from figure 9, a jump from
B = 12 to B = 20 results in a twenty percent increase in
evaporation for a wet surface with a fetch of x 0 i 07.
However, it is doubtful that B ever has a value as low as
12 or as high as 20. In the literature, values of 15 or
16 have been reported and used. The conclusion can be
drawn that any realistic value of the Businger-Dyer
constant or the constant in the O’KEYPS equation (Businger
and Yagloin, 1971) has little effect on the calculation of
evaporation.
In figure 10, total average evaporative flux, W, is plotted
against fetch, x 0 , with average stability parameters for
various values of the scaling height, A. Within a band of
48

-------
io2
Figure 7: Dimensionless negative inverse Obukhov length, -l/L,
versus dimensionless distance downwind, x, for various values of
A , and B constant.
1
L
-4
10
io6
io8
x

-------
100 rir t i i s i i t i i i
A -0.01
B : -0.003
k 0 .35 h 0.1
k 0 .40 h 0 .1
0 k 0 .40 : h0.2
10
U’
o
-2
10 — —‘ 1 1 i J _ I liii 1 1 11I _ I i I ;I I _ i It
10 10 10 10/ io8
Figure 8: Dimensionless vapor flux, W, versus dimensionless fetch, X 0
x 0 , for given stability parameters, A and B , and for Karman’s
k = .35, and two values of initial vertical step size, h.

-------
1o I ill III
A* 0i01 -
B -0.001
-i B 20
10
/
Bz12
I J I I I I I I I I I U
io2 ic io 6 io8
Figure 9: Dimensionless vapor flux, W, versus dimensionless fetch, X
x 0 , for given stability parameters, A and B , and for Businger— 0
Dyer cOnstant, B, equal to 12 and 20.

-------
Figure 10: Diu ensionless vapor flux, W, versus dimensionless
fetch, x 0 , for given stability parameters, A and B , and for
various values of scaling height, A, in the Blackadar mixing length.
1 I —
A -O .O1
100
w
151
io- 2
fl.C’vj3
>\ (1)
5
7r — -

/
/5 i0
/
X 2 .5x1O
io2
io 6
io 7
io 8
xo

-------
A-values, 5 x 1O 4 to 1 x 10 , the solution tends to stay
along the same curve for large values of x 0 . As mentioned
previously, Taylor (1970 a,b) used Prandtl’s mixing length
(A ÷ °) for the temperature jump problem only, and his
solution began to increase around x = 1O . Similarly, as
shown in figure 10 the use of finite A-values greater than
1 x causes the solution to rise but this rise takes
place further downwind and less steeply. For values of X
below 5 x i0 4 , the solution remains on the general curve
but then drops off quite suddenly. Again, the smaller the
value of A, the sooner and the more steeply the solution
begins to drop downwind. It is not immediately obvious
why only a range of scaling height values allows a solution
further downwind. However, the values of the scaling
height are well within the range reported in the literature
for reasonable values of the pertinent physical quantities
involved in the calculation of A.
Discussion . For evaporating fetches, x , greater than about
lO in figures 1, 2, 3, and 4, the total average vapor flux,
W, seems to become a straight line function of x 0 on log-log
paper. Therefore, W can be expressed as a simple power
relationship
w = a x 0 (A30)
where a and n are dependent only on a set of stability para-
meters A , B . The relationship of equation (A30) was
determined graphically from the curves in figures 1, 2, 3,
and 4 and for several other A , B combinations. The coef-
ficient, a, is given in Table 1 and the exponent, n, is
given in Table 2 for the various A , B combinations.
The solutions of this study as expressed by equation (A30)
can indirectly be compared to the simpler problem considered
by Sutton (1934). The Sutton solution for a non-changing
velocity profile expressed as a power law, U = C 15
W = b (A31)
where
1
1-2r r l-2r
r(r) r (l-r)m c
53

-------
Table 1 : Values of the coeffficient a in W = a x 0
-
—B
.1
.05
.01
.001
0
.01
.121
.112
.122
.003
.150
.140
.120
.135
.166
.001
.112
.152
.167
0
Table 2: Values of
the
exponent n
.120
.
in
.167 .210
-n
W = a x 0
-
—B
.01
.003
.001
0
.1
.025
.05
.034
.01
.036
.046
.042
.045
.001
.042 -
.081
.082
.093
0
.045
.089
.093
.112
Table 3: Values of
-
—B .1
m =
n/(1-2n)
.05
.01
‘
.001
0
.01
.037
.048
.003
.027
.037
.051
.096
.109
.001
.046
.098
.114
0
.048
.114
.144
54

-------
and
r = 1 +2in (A32)
The r denotes the gamma function, and c and in are the coef-
ficient and exponent in the power law.
Under neutral atmospheric conditions, it has been found
that in in the velocity profile power law equals approxi-
mately 1/7 or 0.143. For superadiabatic conditions, in
decreases. If n in equation (A30) is set equal to r in
equation (A31), equation (A32) can be solved for in. Values
of in are tabulated in Table 3. For A = B = 0.0, m has
a value of 0.144 which corresponds almost exactly to the
neutral value of 1/7. Thus the case = B = 0 is in
agreement with the Sutton solution. Also, Table 3 shows
values of m decreasing for increasingly negative values
of A and B , i.e. increasingly unstable conditions.
This is not unexpected since the parameter m does indeed
decrease for more unstable velocity profiles (e.g.
Brutsaert and Yeh, 1970b).
An empirical formula for lake evaporation, obtained by
Harbeck (1962) from field measurements, can also be
expressed in the form of equation (A30). The Harbeck
formula is
E = 5.32 x 10 (e — ea) u (A ) • (A33)
where A is the area of the lake in cm 2 , u the wind speed in
cm/sec at two meters, e. and ea the vapor pressures at the
water surface and in the upwind air, respectively, in milli-
bars, and E the mean evaporative flux in cm/sec. The vapor
pressure is related to specific humidity for a given
atmospheric pressure, taken here as 1013 nib. Thus,
E = 7.23 x l0 U 2 x 0 1 (A34)
In order to construct the non—dimensional vapor flux on the
left hand side of the above equation, the velocity at two
55

-------
meters, u2, must be converted to friction velocity, u .
If, as a first approximation, it is assumed that the
velocity profile is logarithmic (neutral conditions), the
substitution of the velocity at two meters into equation
(A34) yields for total average vapor flux:
W = [ 1.8 x io2 z 0 in ( 00)] x 0 (A35)
where both W and x are dimensionless but not z . If the
roughness length, ‘ is known, the coefficient 0 in the
brackets can be calculated. On the other hand, from tables
1 and 2 it can be seen that the exponent in the Harbeck
formula corresponds to A and B both between 0 and -0.001.
Therefore, the coefficient should have a value between .167
and .210. Solving for z 0 in the coefficient in equation
(A35) above yields z 0 = .14 cm for the coefficient, a, equal
to .167 and z 0 = .05 cm for a = .210. Hence, an exponent
of —0.1 seems to imply a roughness length approximately
equal to 0.1 cm. However, this value is almost certainly
too large, since the comparision is made assuming a log-
arithmic velocity profile. For an unstable velocity pro-
file, which the non zero values of A and B require, the
coefficient would be larger and the roughness length
derived therefrom would be smaller. This is in agreement
with the earlier findings of Brutsaert and Yeh (1970), who
obtained a roughness length of z 0 = .0213 cm by equating
coefficients in the Harbeck and Sutton solutions.
The present results can readily be applied for engineering
purposes. If the surface temperature and specific humidity
are known upwind and downwind of a discontinuity and the
upwind shear stress or wind profile were known, A* and B
can be determined. Then, from Tables 1 and 2 values of a
and n may be determined or extrapolated. Equation (A30)
should then yield an estimate for the average vapor flux
over a given fetch, x 0 , greater than, say, iO . For longer
periods, such as one week or more, when the atmosphere may
be assumed to te near neutral on the average, it is probably
permissible to calculate lake evaporation by means of
Sutton’s solution or its empirical equivalent, the Harbeck
formula.
56

-------
Literature Cited in Appendix
Blackadar, A. K,, The vertical distribution of wind and
turbulent exchange in a neutral atmosphere, 3. of Geophys.
Res., 67, 3095—3102, 1962.
Brutsaert, W., and G. T. Yeh, Evaporation form an extremely
narrow wet strip at ground level, J. Geophys. Res., Ocean
and Atm., 74, 3431—3433, 1969.
Brutsaert, W., and G. T. Yeh, Implications of a type of
empirical evaporation formula for lakes and pans, Water
Resources Res., 6, No. 4, 1202—1208, 1970a.
Businger, 3. A., Transfer of momentum and heat in the
planetary boundary layer, Proc. Arctic Heat Budget and
Atmospheric Circulation, Rand Corp., RM-5233-NSF, 305,
1966.
Businger, J. A., J. C. Wyngaard, Y. Izumi, and E. F.
Bradley, Flux-profile relationships in the atmospheric
surface layer, J. of Atm. Sci., 28, 181—189, 1971.
Businger, J. A. and A. M. Yaglom, Introduction to Obukhov’s
paper on “Turbulence in an atmosphere with a non—uniform
temperature”, Boundary Layer Meteor. 2, 3—6, 1971.
Dyer, A. J., The turbulent transport of heat and water
vapor in an unstable atmosphere, Quart. 3. Roy. Meteor.
Soc., 93, 501—508, 1967.
Dyer, A. 3., and B. B. Hicks, Flux—gradient relationships
in the constant flux layer, Quart. 3. Roy. Meteor. Soc.
96, 715—721, 1970.
Harbeck, C. E., Jr., A practical field technique for measur-
ing reservoir evaporation utilizing mass—transfer theory,
U.S. Dept. mt., Geol. Survey Prof. Paper, 272-E, 1962.
Laikhtman, D. L., Physics of the boundary layer of the
atmosphere. (Translated form Russian) Israel Program for
Scientific Translations, Jerusalem, 1964.
Monin, A. S. and A. M. Obukhov, Basic laws of turbulent
mixing in the ground layer of the atmosphere, Trudy Geofiz.
Instit. A N — SSR, No. 24, 151 , 163—187, 1954. (German
Translation: Goering, H., Ed., Sammelband zur Statistischen
Theorie der Turbulenz, Akademie Verlag, Berlin, 1958).
57

-------
Monin, A. S. and A. M. Yaglom, Statistical fluid mechanics:
Mechanics of turbulence. (Translated from Russian). The
MIT Press, Cambridge, Mass., 1971.
Morgan, D. L., W. 0. Pruitt, and F. J. Lourence, Analyses of
energy momentum, and mass transfers above vegetative sur-
faces, ECOM 68—GlO—F, 1971.
Obukhov, A. M., Turbulence in an atmosphere with a nonuni-
form temperature, Trudy Instit. Teoret. Geofiz. AN-SSR,
No. 1, 1946. (English Translation: Boundary Layer Meteor.,
2, 7—29, 1971).
Peterson, E. W., A numerical model of the mean wind and
turbulent energy downstream of a change in surface roughness,
Center for Air Environment Studies, Publication No. 102—69,
April, 1969.
Sutton, 0. G., Wind structure and evaporation in a turbulent
atmosphere, Proc. Roy. Soc. A, 146 , 701-722, 1934.
Taylor, P. A., On planetary boundary layer flow under condi-
tions of neutral thermal stability, J. of Atm. Soc., 26,
427—431, 1969.
Taylor, P. A., A model of airflow above changes in surface
heat flux, temperature, and roughness for neutral and
unstable conditions, Bound. Layer Meteor., 1, 18-39, 1970a.
Taylor, P. A., Numerical models of airflow over Lake
Ontario, IFYGL report, To appear as a Canadian Meteorolog-
ical Memoir, Can. Meteor. Service, Toronto, 1970b.
Taylor, P. A., and Y. Delage, A note on finite difference
schemes for the surface and planetary boundary layers,
Bound. Layer Meteor., 2, 108-121, 1971.
Weisman, R. N., A problem in turbulent diffusion: evapora-
tion and cooling of a lake under unstable atmospheric
conditions, Thesis presented Cornell University, Ithaca,
N.Y. in partial fulfillment of the requirements for the
degree of Doctor of Philosophy, 1973.
Yeh, G. T. and W. Brutsaert, Perturbation solution of an
equation of atmospheric turbulent diffusion. Jour. Geophys.
Res. (Oceans and Atmos.) 75, 5173—5178, 1970.
58

-------
Yeh, G. T. and W. Brutsaert, A numerical solution of the
two—dimensional steady state turbulent transfer equation.
Monthly Weather Review, 99, 494—500, 1971a .
Yeh, G. T. and W. Brutsaert, Sensitivity of the solution
for heat flux or evaporation to of fdiagonal turbulent
diffusivities. Water Resources Res. 7, 734—735, 197lb.
Yeh, G. T. and W. Brutsaert, A solution for simultaneous
turbulent heat and vapor transfer between a water surface
and the atmosphere. Boundary Layer Meteorol. 2, 64-82,
1971c.
Zilitinkevich, S. S., Effect of humidity stratification
on hydrostatic stability, Izv. Acad. Sci. Atmos. Oceanic
Physics, 2, 1089—1094, 1966 (English Trans. A.G.U.
pp. 655—658.)
59

-------
SELECTED WATER 1. Report No. 2. 3. Accession No.
RESOURCES ABSTRACTS w
INPUT TRANSACTION FORM
4. Title 5 Report Date
Heat and Water Vapor Exchange between
Water Surface and Atmosphere 6.
_________________________________________________________________________________ 8. Performing Organization
7. Author(s) Brutsaert, Wilfrjed H. ReportNo.
10. Project No.
9. Organization Cornell University, Ithaca, New York, ii. Contract/Grant No.
School of Civil and Environmental Engineering 16130 DIP
23. Type of Report and
Period Covered
12. Sponsoring Organization Environmental Protection Agency Final Report
15. Supplementary Notes
Environmental Protection Agency report number EPA-R2-73-259, May 1973
16. Abstract The physical and mathematical aspects of simultaneous turbulent
heat and water vapor exchange between a large open water body and the
surrounding atmosphere were studied. Thus analytical and numerical solu-
tions were developed for various conditions of fetch, surface roughness,
atmospheric stability, etc., that are likely to be of physical importance.
One of the main findings was that in spite of some theoretical limitations
the semi-empirical turbulent diffusion model provides a method for the
prediction of heat and water vapor transfer, that should be use u1 for
engineering calculations. although the interaction between momentum,
sensible heat and water vapor is considerable, for evaporation or cooling
calculation purposes it is probably permissible to uncouple their transfer
mechanisms provided the water surface temperature is known and the averag-
ing period under consideration is of the order of, say, a week. In addi-
tion the validity of the semi-empirical assumption in the near-water sur-
face layer was analyzed by determining the anisotropy of the eddy diffu-
sivity, the effects of radiative transfer and of water wave action on the
eddy diffusivity. Finally, a practical method was developed to determine
evapotranspiration from the surrounding land surface based on the geo-
strophic drag concept.
17a. Descriptors
* aporatjon, *Heat Transfer, *Ajr. ..Water Interfaces, Estimating Equations,
Mass Transfer, Water Cooling.
I 7b. identifiers
Turbulent Diffusion
17c. CO WRR Field & Group 02D ___________________________________
18. Availability 19. Security Class. 21. No. of Send To:
(Report) Pages
WATER RESOURCES SCIENTIFIC INFORMATION CENTER
20. Secunty Class. 22. rIce U S DEPARTM tNT OF THE I NTERIOR
W. H. Brutsaert (Page) WASHINGTON. D.C 20240
Abstractor Institution Cornell University
wRs,c ioe(ptv JUNE 1971) GPO 9 t3.2 t - /

-------