AEDC-TR- 76-85

EQUILIBRIUM COMPOSITIONS AND COMPRESSIBILITY

FACTORS FOR AIR IN THE TEMPERATURE

RANGE 1500 TO 13000 KELVIN

AND FOR DENSITIES FROM 10

6

TO 10

3

TIMES

STANDARD DENSITY

NATIONAL BUREAU OF STANDARDS

HEAT DIVISION

WASHINGTON, D. C. 20234

June 1976

Final Report for Period July 1974 - June 1975

Approved for public release; distribution unlimited.

Prepared for

DIRECTORATE OF TECHNOLOGY (DYR)

ARNOLD ENGINEERING DEVELOPMENT CENTER

AIR FORCE SYSTEMS COMMAND

ARNOLD AIR FORCE STATION, TENNESSEE 37389
NOTICES

When U. S. Government drawings specifications, or other data are used for

any purpose other than a definitely related Government procurement

operation, the Government thereby incurs no responsibility nor any

obligation whatsoever, and the fact that the Government may have

formulated, furnished, or in any way supplied the said drawings,

specifications, or other data, is not to be regarded by implication or

otherwise, or in any manner licensing the holder or any other person or

corporation, or conveying any rights or permission to manufacture, use, or

sell any patented invention that may in any way be related thereto.

Qualified users may obtain copies of this report from the Defense

Documentation Center.

References to named commercial products in this report are not to be

considered in any sense as an endorsement of the product by the United

States Air Force or the Government.

This final report was submitted by the National Bureau of Standards, Heat

Division, Washington, D. C. 20234, under contract F40600-72-M-0002 with the

Arnold Engineering Development Center, Arnold Air ForCe Station, Tennessee

37389. Mr. Elton R. Thompson (DYR) was the AEDC project monitor.

This report has been reviewed by the Information Office (01) and is releasable to the National Technical Information Service (NTIS). At NTIS, it will be

available to the general public, including foreign nations.

APPROVAL STATEMENT

This technical report has been reviewed and is approved for publication.

FOR THE COMMANDER

ELTON R. THOMPSON ROBERT O. DIETZ

Research & Development Director of Technology

Division

Directorate of Technology
UNCLASSIFIED

REPORT DOCUMENTATION PAGE READ INSTRUCTIONS

BEFORE COMPLETING FORM

REPORT NUMBER

AEDC-TR-76-85

2 GOVT ACCESSION NO. 3 RECIPIENT'S CATALOG NUMBER

i. t,t l e r-n d subuti^ EQUILIBRIUM COMPOSITIONS AND

COMPRESSIBILITY FACTORS FOR AIR IN THE

TEMPERATURE RANGE 1500 TO 13000 KELVIN AND

FOR DENSITIES FROM 10" 6 TO 10 3 TIMES

STANDARD DENSITY

5. TYPE OP REPORT a PERIOD COVERED

Final Report - July 1974

- June 1975

6. PERFORMING ORG. REPORT NUMBER

7. AuTHORfsj

Max Klein and Lester Haar

8. CONTRACT OR GRANT NUMBER'S,

F40600-72-M-0002

9. PERFORMING ORGANIZATION NAME AND ADDRESS National Bureau of Standards

Heat Division

Washington, D.C. 20234

10. PROGRAM ELEMENT. PROJECT TASK AREA 4 WORK UNIT NUMBERS

Program Element 61102F

12. REPORT DATE June 1976

II. CONTROLLING OFFICE NAME AND ADDRESS Arnold Engineering Development Center (DY)

Air Force Systems Command | ,3 number of pages

Arnold Air Force Station, Tennessee 37389 75

1«. MONITORING AGENCY NAME & ADDRESS'// dillerent Irom Controlling Ollice) 15. SECURITY CLASS, 'ol this report; UNCLASSIFIED

15a. DECLASSIFICATION DOWNGRADING

SCH ED; N/A

16- DISTRIBUTION STATEMENT 'ol this Report)

Approved for public release; distribution unlimited,

'7. DISTRIBUTION STATEMENT (o! the abstract entered in Block 20, it dillerent Irorr. Report;

18. SUPPLEMENTARY NOTES

Available in DDC

19. KEY WORDS 'Continue on reverse side it necessary and identity by block number) air

gases

mixtures

thermodynamic properties

high density

20. ABSTRACT f Continue on reverse side if necessary and identify by block number) A formalism for computing thermodynamic properties and equilibrium

compositions of reacting gas mixtures at high densities is developed. This formalism makes use of the Haar-Shenker version of the

augmented hard sphere equation of state. The formalism is used for

computing compressibility factors and equilibrium compositions of

air in the temperature range 2,000 K - 13,000 K for densities from

that appropriate to air at standard conditions (0°C, 1 atmosphere)

DD ,j2n

M

73 1473 EDITION OF 1 NOV 65 IS OBSOLETE

UNCLASSIFIED
UNCLASSIFIED

20. ABSTRACT (Continued)

to densities 1000 times as great. These tables are consistent with

the earlier NBS tables which covered the range up to 100 times

that density. Some interesting features of the dependence of composition on density are discussed.

UNCLASSIFIED
AEDC-TR-76-35

PREFACE

This research was conducted by the National Bureau of Standards (NBS),

Heat Division, basically under sponsorship of the Arnold Engineering Development Center (AEDC), Air Force Systems Command (AFSC), under Air Force

Contract F40600-72-M-0002. The AEDC project monitor for this work was

Elton R. Thompson, AEDC DYR. The work was also supported in part by

the Atomic Energy Commission under Interagency Agreement AT(49-16)3003,

Modification No. 2. This work covers essentially the period July 1974 to

•June 197 5. Principal investigators and authors of this report were Max Klein

and Lester Haar of NBS.

Klein, M. , Haar, L. . Equilibrium compositions

arc. c j-cressicility factors for air in one

temperature range 1500,-tG 1^023 kelvia and

for densities from 10~ c to 10^ times

standard density., AETJC-TR-7'6-65 , 77 pages

(Available free: the National Technical

Information -Service, Springfield^ VA 22161,

June 1976)

.

221

1

AEDC-TR-76-85

TABLE OF CONTENTS

I. INTRODUCTION 5 II. THE EQUATION OF STATE OF HAAR AND SHENKER 11

A. Background 11

B. Derivation for a Pure Substance 14

C. Comparison with Experimental Data and Other Theories 19

1. Sensitivity to the selection of the attractive part of

the effective potential 19

2. Relationship to other fundamental equations of state ... 20

3. Comparisons of the HS equation of state with experiment

and with other theories 22

4. Summary » 25

III. THE EXTENSION OF THE HAAR-SHENKER EQUATION OF STATE TO MIXTURES ... 26

IV. THERMODYNAMIC FUNCTIONS IN THE ONE FLUID HS MODEL 30

A. Helmholtz Free Energy, Chemical Potential and Activity

Coefficients 30

B. Entropy, Enthalpy and Gibbs Free Energy 36

V. NUMERICAL METHODS 38

VI. RESULTS 40

A. The Dependence of Species Concentrations on Density 40

B. Density Effects on the Equation of State 45

C. Suggested Possibilities for Improvements in the Tables 47

VII. TABLES

1. Q* Values for the (12,6) Hard Core Potential 51

2. y'^ Values for Several Densities for Representative Species. . . 52

3. Comparisons of the Compressibility Factors, PV/RT for Air for

Four Approximations (Including that of Grabau and Brahinsky) as Functions of Density at Various Temperatures 53

4. Compressibility Factors, PV/RT, for Air for Three Approximations

as Functions of Density at Various Temperatures 54

VIII. REFERENCES 56

IX. FIGURES

1. Reduced Virial Coefficients B^ Vs Reduced Temperature T . The

Maxima Tend to Lower Temperatures with Increasing N, Supporting

the Argument made in the text that, except at low temperatures,

the higher virials are representable by hard-sphere interactions

without attraction 59

3
AEDC-TR-76-85

2. Effective potentials vs molecular separation in arbitrary units

for several values of m and T* = 2.5. It is seen that though the

shape of the "bowl" is quite sensitive to the characteristic

parameter m, the sphere diameter at this temperature is only

weakly sensitive to it 60

3. The sphere diameter in reduced units vs reduced temperature for

effective potentials with m=9 and m=10 and for the ZBH theory.. 61

4. Compressibility factor for argon vs reduced density p(Np^) = p for T = 119. 8K. Curves 1, this work; Curves 2 and 4, ZBH

theory from references 22 and 23, respectively; Curves 3, ZBH

theory based on Monte Carlo calculations of perturbation

terms'^; Curves 5, PVT experimental data^. The "computer

experiments"^ are designated by circules, Monte Carlo; and

squares, molecular dynamics . 62

5. Compressibility factor vs density for Argon for T = 161. 7K. . . 63

6. Compressibility factor vs density for Argon for T = 239. 6K. . . 64

7. Compressibility factor vs density for Argon for T = 328. 25K.. . 65

8. Compressibility factor vs density for Argon for T = 373K. ... 66

9. Compressibility factor vs density for Argon for T = 673K. ... 67

10. Sensitivity of Equation of State for ZBH theory to choice of

potential function. Compressibility factor is plotted vs density in amagat. For Curve 1, the reference potential is

the Lennard-Jones (18,6); for Curve 2, the Lennard-Jones (12,6).

The numerical results are from calculations furnished by

Toxvaerd23 using potential parameters e/k = 119.8, a = 3.405

for the latter and equation (12) for the former 68

11. Sensitivity of present theory to choice of potential functions.

Potential functions are described in caption 10. The small

difference between Curves 1 and 2 is mostly due to the slight

differences in second virials produced by the two potentials

in the vicinity of the isotherm . 69

12. Compressibility factor vs density in amagat for nitrogen

T = 123. 15K. Curves 1, present theory, Curves 2, experimental

PVT measurements. 35 70

13. Compressibility factor vs density in amagat for Nitrogen

T = 173. 15K. . 71

14. Compressibility factor vs density in amagat for Nitrogen

T = 273. 15K 72

15. Compressibility factor vs density in amagat for Nitrogen

T = 673K 73

16. Concentration vs density for several important species for

T = 3000K 74

17. Concentration vs density for several important species for

T = 9000K 75

4
AEDC-TR-76-85

I, INTRODUCTION

Aerodynamic design requirements demand the ability to calculate tables of

thermodynamic properties and equilibrium compositions of air for thermodynamic

states far from normal laboratory conditions. Earlier design requirements led

to the production of tables at very high temperatures ranging in density from

quite low densities to intermediate densities. The NBS tables- (which were published under AEDC support) are examples of a response to this earlier need.

These tables covered the temperature range 1500 K < T < 15,000 K and the

density range 10

0 < c/c < 10^, where p is the density at the standard condi-

- o o tions of T = 273.15 K and P = 1 bar. More recently, there have been indications^ of needs for tables covering considerably higher densities at these

elevated temperatures

.

In this project, we were requested to attempt to develop methods for the

calculation of the properties of air for densities up to 1000 times normal sea level density. If successful, these tables, taken together with the NBS tables

already published, would result in tables being available for the enormous

—6 3 density range 10 < c/c

q

< 10 (i.e. nine orders of magnitude) at temperatures

well outside those associated with laboratory experiments. Such tables would

then represent an extrapolation from ordinary conditions of over two orders of

magnitude in temperature and at the same time, an extrapolation of at least one order of magnitude in density from the earlier tables

.

Combined extrapolations of such magnitudes in temperature and density pose

very difficult problems. They must of necessity be based on fundamental

properties of the systems under study, properties which might not vary over the

range of the extrapolation. These properties must also be used in a framework

of fundamental theoretical methods valid over the entire range of temperatures

and densities of interest. In addition, specific numerical methods must be

developed and computer programs produced for transforming these numerical methods

into the actual calculation of tables. These numerical problems include the

often difficult task of producing methods for the solution of sets of non-linear

3 algebraic equations, with such methods being required to produce solutions for a wide variety of values for the unknown parameters.

5
AEDC-TR-76-85

It should be obvious that, because of the extent of the extrapolations

required, empirical correlation methods are not appropriate here.

Such engineering methods involve the least squares fitting of measured data to

equations containing arbitrary parameters. Such parameters generally do not

have physical significance and can therefore only be used with extreme caution

if at all in extrapolations beyond the range of experimental data on which their

values are based. In fact, using arbitrary parameters obtained by least squares

fitting, in such long extrapolations beyond the range of their data base

generally produces highly erratic behavior. Furthermore, considerable additional difficulty could be anticipated in applying such methods to mixtures. Mixing

rules for non-physical parameters in the context of empirical correlations are generally arbitrary and often not useful even within the range of the data base

and hence their behavior in extrapolations is unpredictable.

The original plan on which this research was based envisioned the use of

the then developing and fundamentally based integral equation theory for the

4 equation of state of fluids. At the time this work began, that theory, though

very promising for one component systems, had not been applied to mixtures nor

had any attempt been made to integrate it into the context of chemical or phase

equilibrium. It was not clear how one might calculate the free energy and

chemical potentials of a multicomponent mixture within this approach. The study

of the possibility of using the integral equation approach was, in fact, to

constitute the major part of the research program supported by this contract

with the probability of success not entirely clear. A rather complex and

purely numerical approach was envisioned which, although it might of necessity

be very complicated in a numerical sense, might nevertheless be expected

ultimately to be made to work. Because of its expected complexity, the basic

numerical parts of the problem were postponed in favor of a close study of the

details of the integral equation approach itself. Initially, this involved

examining integral equation methods as a means for representing the equation of

state of pure substances before looking at them as possible methods for mixtures.

Any statistical mechanical method which could be used in this work

(particularly one consistent with the earlier tables) including the integral

equation approach, requires the use of intermolecular potential functions. For

the temperatures of interest in aerodynamic problems, such intermolecular
AEDC-TR-76-85

potential functions would certainly be used beyond the range in which data fits were carried out to obtain their paraneters. This means that methods used to obtain these parameters would need to be understood very thoroughly. Furthermore, because of the implicit way in which the potential functions appear in

statistical mechanical theories, the effect of uncertainties in their parameters

on the accuracy of calculations in which they are used is not always clear. For

these reasons, much of the time in this research program had to be spent in the

study of the means by which intermolecular potential functions are obtained

from experimental data of various kinds and in the study of the effect on calculated tables of extrapolating incorrect intermolecular functions . Although

this extrapolation of the potential parameters beyond the range of the data

could be expected to introduce errors, these errors could be expected to be far

less than those obtained from extrapolation of least squared fits to data such

as are common in the usual empirical approaches. Furthermore, one might expect

the general behavior of the functions calculated to be reasonable. This follows

from the fact that the essential behavior of the potential function would be

correct and from the fact that the "actual" potential functions would not vary

with temperature so that the relationship between the potentials chosen and the

"real" ones would not contain any hidden surprises on extrapolation to higher

temperatures

.

Although the present work emphasizes the high density region, it had to be

designed to retain the earlier tables. Since these earlier calculations already

covered densities up to 100 times normal sea level density, the extension from

100 to 1000 times normal density needed to be done using methods which included

the earlier approach. At the highest densities, the earlier tables included a correction for the second virial coefficient based on the intermolecular forces

between the molecules. Hence, almost any approach is consistent with the

earlier tables if it makes use of a statistical mechanical theory based on these

same intermolecular forces provided that this theory has the same (and correct)

linear term in its low density expansion.

As the work progressed, alternative approaches were examined in an attempt

to avoid the numerical complexity expected in the integral equation approach.

These other methods were rejected either because their complexity did not offer

any advantage over the integral equation approach or, more usually, because they
AEDC-TR-76-85

required the extrapolation of parameters which had no physical basis. One

alternative which appeared to hold promise, resulted from research by one of us (Lester Haar) . An equation of state for single component fluids was developed

which seemed to us capable of being used in this work in a computational

procedure far simpler than that projected for the integral equation method.

Also, it was based on fundamental ideas and was thought to provide a natural

extension of the previous NBS tables. The approach is based on an equation of

state due to Kaar and Shenker , (HS) , which was developed for the extrapolation

of low density experimental data of one component system to high densities along

isotherms at ordinary temperatures. According to the assumptions on which it was

based, the HS equation was expected to work well at temperatures above the

critical temperature and to improve with increasing temperature and this was indeed found to be the case by them. They also found indications that the method

could be used to produce engineering calculations of reasonable accuracy even below the critical temperature into the liquid range.

The HS equation uses second virial coefficient data at each temperature as

a basis for extrapolating to high densities. In our studies of the relationship

between second virial coefficients and intermolecular potentials, we saw that

the determination of intermolecular forces from second virial coefficients, if

done properly, could be used to produce a temperature extrapolation of the

second virial coeffcient to higher temperatures even as high as required in

aerodynamic calculations. Thus, given an intermolecular potential function for

a particular substance, a complete high temperature PVT surface could, in

principle, be produced for that substance, with the potential function being

used to calculate second virial coefficients at all temperatures of interest

and these virial coefficients, in turn, serving as the basis for the HS equation

at all densities.

Before it could be used in this work, the HS equation (which was developed

for a pure substance) had to be extended to mixtures in the context of a reacting gas. This we have done and the details of the formalism will be

presented in this report. The resulting formalism was used by us to produce

the attached set of tables. A difficulty is expected to arise in this formalism

whenever there is a need for intermolecular potential functions between

8
AEDC-TR-76-85

fragments which are present in the mixture at high temperatures but which do not

exist at ordinary temperatures. Potential functions between such molecular

fragments can, however, often be characterized approximately by means of the

potential function which describes the vibrational spectra of the combined

fragments or, failing such spectra, from estimates of molecular sizes. Since

estimates of the intermolecular forces between the various species had to be

made for the earlier tables (including those species not present at ordinary

temperatures), the basic data needed for this approach were actually already

available at the start of this work. In order to ensure consistency in the

initial calculations, we decided to use these same estimates of these forces in

the calculations on which we report here.

Although the formal framework of the approach used by us is a valid one, the

values contained in the attached tables should, for a number of reasons, only be

considered an interim set of values :-

1. The estimates of the intermolecular potentials used are based on an old analysis by Wbolley in terms of the (12,6) potential function. As a result of our earlier work under AEDC support, better potentials are now known

and should be used for the major species and these same potentials should also

be used in place of the (12,6) for the new estimates of those associated with

the minor species. As already mentioned, we used these older estimates so as to ensure consistency between these tables and the earlier tables.

2. In several cases, certain unknown intermolecular potential functions were arbitrarily taken as equivalent to others that were known. The effect of

this could be important at the very highest densities.

3. As we shall describe below, we have used an ad hoc approach in the

determination of hard sphere diameters at those temperatures for which the

slope of the second virial coefficient is negative.

4. In this work only the compressibility factor and equilibrium compositions have been calculated.

Each of these compromises has been made in order to produce results to test the

approach, leaving the production of more extensive tables for possible later work.

9
AEDC-TR-76-85

g Grabau and Brahinsky , (GB) , attempted an extrapolation of the earlier NBS

3 tables to 10 times normal density for air and nitrogen. Their method was essentially semi-empirical. Their method was based on an equation having, in

part, a fundamental basis and which contained a dependence on ideas partially

related to the use of intermolecular forces. Their method could not be applied

to a gas of varying composition (i.e. a reacting gas), for which reason they

restricted their air calculations to temperatues below 6000K. Their method

would also be expected to lose accuracy very rapidly at any temperature as the

density increased due to a need for successive subtractions and because of their

neglect of higher order terms. Since their method was based mainly on a graphical extrapolation, even qualitative estimates of its adequacy could not be made

by them.

In this report, we compare our results with the extrapolation of Grabau

and Brahinsky. These comparisons are interesting and will be discussed below.

The agreement between the two extrapolations at intermediate density at the

temperatures which they have in common is expected since the GB extrapolations

were based on the earlier NBS tables with which our results are entirely

consistent. The two results show considerable disagreement at the highest

densities for reasons discussed below.

Our report contains plots at two representative temperatures of species

concentrations for some important species. These plots show some interesting

and possibly unexpected density behavior as discussed below. Information on composition behavior is not available in the approach of Grabau and Brahinsky

since their method is applied to the thermodynamic properties only. Of particular interest in our plots are the dependences of concentrations on density as produced in three widely used approximations - the ideal gas, the second virial

coefficient gas and the HS gas. Regardless of the ultimate accuracy of our own model, the differences among these three approximations can be expected to be

indicative of the magnitude of the effect of the analogous three approximations

in any other (and possibly more accurate) theoretical model.

Many unexpected technical numerical and computer programming problems were encountered in adding on the extreme high density end of the calculation.

These had to do with numerical difficulties in obtaining solutions and will be

discussed only in passing in this report.

10
AEDC-TR-76-85

II. THE EQUATION OF STATE OF HAAR AND SHENKER

A. BACKGROUND

The first step in the calculation of the thermodynamic properties of a re- acting mixture involves the computation of the composition corresponding to

reaction equilibrium. Within the framework of the formalism previously used by

us, this corresponds to the solution of the equations associated with the law of mass action, as modified to take account of any non-ideal effects, subject to

the conservation laws for nuclear types. The mass action relations can be

written

C. = K.(p/p o r

Wi n V. WJl C, Vik (1)

1 1

I

1 k

k

where the C, are the concentrations of the reference species, v., the

K- ~ co^ ik

stoichiometric coefficients for the reaction, K. = K.(T/T ) with K. the

' 1 l o 1 equilibrium constant, -co^ = Ev^-l the net production of particles across the

reaction, p the density at standard conditions (for P = one atmosphere

° (O pressure and for T = T = 273.15 kelvins) and where y* is the effective

° th th activity coefficient for the I non-ideal effect for the i reaction. A

detailed discussion of this equation and its derivation are contained in

reference (3) and will not be repeated here. The y^ associated with the DebyeHuckel theory and that associated with the second virial coefficient are contained in Appendix B of reference (3) . The derivation of the y^ for the HS

equation constitutes a major part of this report. For the moment, it will be

enough to state that, in principle , a y!^ can be obtained for any equation of

state but, in practice, the procedure is very complicated and not always clear.

The understanding that formal expressions for the y! for the HS equation can be

obtained constitutes the motivation for the following extended discussion of

that equation of state, of its development and of tests of its validity.

11
AEDC-TR-76-85

Haar and Shenker (HS) , developed an extremely simple equation of state

based on the virial expansion and on the behavior of the virial coefficients at

high temperatures. The HS equation requires only a knowledge of the second

virial coefficient and its first derivative at each temperature and from this

the behavior at all densities is obtained. As already mentioned, the HS

equation can be used, in a very simple and straightforward manner, to extrapolate PVT data in both the temperature and density directions. It is also

expected to improve with increasing temperatures since the validity of the

assumptions on which it is based improves with increasing temperature. It is

therefore particularly well suited for these calculations.

The development of the equation of state was motivated by earlier work by

9 Haar and Levelt Sengers who showed that only two parameters are required to

correlate thermodynamic properties for a number of simple non-polar gases along

individual isotherms . The two parameters of Haar and Levelt Sengers had to be

different for each temperature and were obtained by fits to experimental data on each isotherm. These two parameters could, in principle, be determined from

experimental data at each temperature in a number of different ways.

The equation of state of Haar and Shenker also has two parameters for each

isotherm with the values of the parameters being determined from values of the

experimental second virial coefficient. Since there are two parameters at each

temperature, two properties of the second virial coefficient are required at

each temperature to determine their values. Haar and Shenker chose for these

the second virial coefficient and its first derivative. The parameters so obtained can be regarded as being an effective temperature dependent molecular

size, which sets the scale of density on the isotherm, and an effective temperature dependent molecular well depth which sets the temperature scale on which

the isotherms are assigned. Each of these parameters is associated with the

experimental system. In the following derivation these two quantities will be

considered to be only slowly varying functions of the temperature which greatly

simplifies their calculation. This assumption is certainly valid in the range

of temperatures of interest in this work. In this way Haar and Shenker

developed an equation of state for the correlation and prediction of high

density data at temperatures above critical using two temperature dependent

parameters. The parameters chosen have a fundamental basis and a simple method

12
AEDC-TR-76-85

for their determination was devised. The fact that the parameters so chosen

could also be associated with notions of corresponding states promised general

success

.

The repulsive energy between two molecules in a gas can generally be

characterized by a very steep function at small intermolecular separations. To a first approximation this repulsion can be described by the interaction between

hard spheres. Since repulsive effects are known to dominate at high temper10

atures , it is reasonable to take the equation of state of a gas of hard spheres

(for which there are only repulsive effects) as a starting point in the development of an equation of soaoe for any gas at hit -

, temperatures. To facilitate

this, Haar and Shenker express the actual equation of state as the sum of a hard sphere contribution (to be calculated by a method as yet unspecified) plus

the difference between the actual equation of state and this hard sphere

contribution. This involves no approximations since the two parts add up identically to the actual equation of state regardless of how the hard sphere

contribution is handled. The first approximating assumption consists in taking,

for the hard sphere contribution, the result obtained from Percus-Yevick theory^

using the compressibility equation of state. It has been shown that this

representation differs only slightly from exact hard sphere theory up to

12 densities approaching 2/3 the close packing density . A second (and more

serious) approximation involves the choice of method for obtaining the hard

sphere diameters needed in the Percus-Yevick theory.

The derivation of the HS equation of state starts from the assumption that

the N body potential of the fluid can be represented by a sum of pair-potentials

13 and that the Ursell-Mayer virial expansion in the density is valid for all

potential functions of interest. The virial series for any potential is then

transformed into a rapidly convergent expansion about the hard sphere series.

Finally, the equation's parameters are fixed by imposing as boundary conditions

the requirement that the first correction to the ideal gas be valid. This last

follows automatically when the second virial coefficient is used to determine

the parameters and is the basis of the consistency between our model and the

earlier NBS tables.

13
AEDC-TR-76-85

B. DERIVATION FOR A PURE SUBSTANCE

The Ursell-Mayer expansion for the equation of state is written,

n=l

PS where Z is the compressibility factor, defined by Z = — ; p the reciprocal

volume, P the pressure, 8 = 1/kT, and N the number of molecules in the system.

The density expansion (2) for the compressibility factor is certainly valid in

the gas phase and can be considered to be an exact representation of the com- pressibility factor Z. The B in (2) are the virial coefficients and are well14 n known integrals, obtained from statistical mechanics, involving the intermolecular interactions among 2,3,4, etc. molecules respectively. We now

formally develop each of these virials about that for a hard sphere of some

Cas yet arbitrary) diameter via the identity

B = B

h - S ' + B - Bh - S>

, (3) n n n n

where B^*

S * is the virial coefficient for the hard sphere. Using (3), we

can rewrite eq. (2)

i " -, n n n=l n=l

where, consistent with the formulation, B^" = B^ = 1. It should be noted that

(4) is an identity and so does not involve any new assumptions. (4) is therefore still an exact representation for Z. The first sum on the right-hand side

of Eq. (4) is the equation of state that would be obtained for a gas made up of

identical hard sphere molecules. The hard-sphere gas has been studied extensively in computer "experiments" via molecular dynamics"'""' and Monte Carlo"^ cal14
AEDC-TR-76-85

dilations and theoretically via the Percus-Yevick approximation . The hardsphere summation in (4) can therefore be considered to be known for any particular sphere size. The task inherent in the evaluation of (4) is to obtain a simple representation for the perturbation terms, i.e. the second sum on the

right of Eq. (4), and to obtain a proper hard sphere diameter to use in the

first term.

We now present an argument which shows that the perturbation terms, i.e.

the second sum on the right-hand side of (4) , converge rapidly above the gasliquid critical temperature so that at such temperatures only the term linear in

density needs to be considered. The prospect that this convergence might

17 persist to somewhat lower temperatures is implied in work by Woolley

To illustrate our argument, in Fig. 1 we plot a few of the lower virial

coefficients for the Lennard-Jones 0-2,6.) pair-potential. The reduced virials

B„, B„, and B. with

2 3 4

* ..2 „ 3,n-l

B =B / [-r ttNo J , n n 3

are plotted against the reduced temperature T = (3e) , where a is the LennardJones length parameter and e the well depth. At low temperatures the contribution from the attractive part of the potential (regions of negative energy) is

important and the lower virials tend to large negative values. At the higher

temperatures, the repulsive part of the potential (regions of positive energy)

tends to yield the dominant contribution and the virials become positive. Thus,

in Fig. 1 we see that B^ monotonically increases, with decreasing slope, from

large negative values to positive ones, finally passing through a maximum. It

is apparent that as the temperature increases, the contribution of attraction

decreases relative to that of repulsion, even though the contribution of

repulsion is itself slowly decreasing with increasing temperature. We suggest

that at some temperature, not much beyond the temperature where B achieves its

*

* maximum, the contribution of attraction to B is relatively small. This occurs

* r at a reduced temperature somewhat above T =30. Thus if the hard-sphere

contribution is chosen appropriately, we would expect that the perturbation

15
AEDC-TR-76-85

terms linear in density would be quite small at temperatures somewhat above

T =30. Likewise we expect that for the third, fourth, etc. virials the

contribution of attraction becomes relatively small at temperatures starting

somewhat above the temperatures for which these virial coefficients achieve

their (initial) maxima. The maximum for the third virial coefficient occurs at

about T =1.25 corresponding to a temperature near the liquid-vapor critical

point for simple substances. The temperature at which the fourth virial

coefficient would achieve its maximum is also in this neighborhood. (This is

18 also true for the fifth virial coefficient not shown in the figure ).We invoke

a corresponding states argument and assume that the disappearance of attractive

contributions to the third and higher virial coefficients above the critical

temperature should be a general property of any simple gas. Thus, we suggest

that at temperatures somewhat above the critical temperature the second virial

coefficient includes all of the major effects due to molecular attractions, and

that the higher virial coefficients are primarily determined by the repulsive

interaction. For each substance, at any given temperature , we therefore represent the repulsion between the molecules by that in a gas made up of identical

hard spheres whose diameter is somehow chosen so as to be appropriate to that

substance at that temperature. We then use this hard sphere gas to represent

the total contribution of each of the virials above the second at that temperature. This certainly should be a good model at the temperatures of interest in

Based on the preceding arguments, for temperatures of interest here, the

second sum in the equation of state, Eq. (4), can be truncated after the term

linear in the attractive contribution. The hard sphere part, on the other hand,

contains contributions to all orders of the density. To represent this hard

sphere part, we employ the results of the Percus-Yevick theory"^ using the

compressibility form for the equation of state. This is a good approximation

at low densities, and is in error by, at most, 5% at densities approaching 2/3

close packing of hard spheres, a density well beyond those of interest here.

The equation of state is written, therefore,

*

our work.

2 + 4y (^ - 1) , (5)

16
AEDC-TR-76-85

where 4y = bp, and b is the hard sphere second virial coefficient,

tt N a

3 (6)

a being the temperature dependent hard sphere diameter. For a given temperature

Eq. (5) is a two constant equation of state, these constants being the hard

sphere diameter a and the well depth associated with the representation of . An important feature to be used below is that the equation (5) is easily integrable in closed form to yield a free energy.

As already mentioned, to evaluate the two equation of state parameters for

a particular gas at each temperature we shall employ the numerical values, for

that gas, of the second virial coefficient and its first temperature derivative.

To bring out the connection between these two parameters and the molecular

diameter and the intermolecular well depth, we introduce an effective intermolecular interaction that has the general features of a typical pair potential,

except that it is specifically characterized by a hard sphere cut-off at some diameter a(T) at which point it is joined to an attractive bowl of well depth

e(T) CNote that we have explicitly indicated the temperature dependence of these

quantities.) The purpose of this effective function is to provide a means for

transforming the repulsive and attractive parts of the "actual" potential of

the gas into an explicit hard sphere diameter and a well depth. Typical

functions of this kind are shown in Fig. 2.

A simple numerical method has been developed for the determination of the

parameters a and e at each temperature. The method starts with tables of

the effective potential function (i.e. the function with the hard sphere cutoff). At any given temperature, these reduced quantities are required to yield

the values associated with the experimental system being described. Thus the

conditions on these reduced quantities are,

reduced second virial coefficients B and their first derivative T —, for

dT

B exp

(T) = b B (kT/e) (7)

17
AEDC-TR-76-85

and dBexD^ * dB* (kT/e)

T —^2- = b T 2£_ (8)

dT

where T is the temperature of interest. On the right-hand side, we have

indicated the dependence on T by kT/e to emphasize the fact that e/k is an unknown quantity. The right-hand sides in (7) and (8) are the reduced quantities as calculated for the effective potential function (i.e. the potential with

the hard core). It should be noted that the right hand side of equation (8) contains

da

the implicit assumptions that — and de/dT can be neglected in the calculation of

a and e. On dividing Eq. (8) by Eq. (7), there results

T ^exp_ T* dB*. _ 0

* B dT * * (9)

exp B dT

* It is a simple matter to produce a table of values for the quantity Q as a

* * function of T for the effective potential function. By way of illustration Q

values are listed in Table 1 for a particular effective potential. The same quantity (i.e. the left hand side of (9)) is then calculated from experimental

data as a function of T. (9) is then solved for e/k at a given experimental

* temperature T. This is done by starting with the experimental value of Q at

* that value of T, and, by interpolation, finding that value of Q and the value

of T associated with it in the table of values calculated for the effective

potential. e/k is then calculated at that temperature from the relation

e/k = T/T*

This procedure guarantees that Eq. (9) is automatically satisfied at the chosen

experimental value of T. This e/k value is then used in (7) to obtain b and

3 therefore a . By carrying out this procedure at each temperature, a(T) and e(T) are determined for all experimental points. The equation of state is then calculated at each temperature T using b(T) and B (T ) in equation (5). (The use

of B^ (T ) is equivalent to the use of e(T))

18
AEDC-TR-76-85

One of the difficulties associated with the above procedure involves the

calculation of dB/dT from experimental data which are generally not smooth and

not presented at convenient temperature intervals. A reasonable way to do this

is first to fit the experimental B(T) data to a realistic (i.e. one which contains

a dependence on r for small r) potential function (as opposed to the effective

potential which has a hard core cut-off for small r) and to calculate smooth tables

exp

of B and —— values using that function,

exp dT 6

The procedure outlined for solving Eq. (9) does not work for those temperature at which dB /dT < 0. Because the effective potential has a hard core

exp r repulsion, its second virial coefficient does not have a region of negative

slope and solution of (9) becomes impossible since the negative value of Q

associated with the experimental system is being sought in the table for the

effective potential which contains only positive values. Under such conditions, we have proceeded by neglecting the attractive contribution by setting B equal

to unity in (7). This leads to a negligible discontinuity in b(T) but in a non- negligible one in its temperature derivative.

C. COMPARISONS WITH EXPERIMENTAL DATA AND OTHER THEORIES

1. Sensitivity to the selection of the attractive part of the effective

potential.

It should be obvious that the effective potential plays no fundamental role

in the HS theory but is used only for computational convenience, being used to extract an effective hard sphere diameter from the second virial coefficient.

It clearly should not be allowed to introduce any of its own character into the

calculation. For this reason, we shall precede a detailed comparison with

experiment, by an examination of the sensitivity of our method to any particular

choice for the shape of the bowl used in the effective potential. In Fig. 3 a

k plot of the sphere diameter a(T) vs. T is shown for the two effective potential

functions corresponding to m=9 and m=10 . As expected from the behavior of a typical intermolecular potential function at small separation, the effective

sphere diameter is a monotonically decreasing function of the temperature. The

scale of the abscissa is normalized so that a(T) = 1 for the effective potential

function given by m=12. The parameters a and z were calculated from equations

(7) and (9) for each of these effective potentials using values for B„ and dB/dT

19
AEDC-TR-76-85

calculated with a Lennard-Jones (12,6) potential in place of experimental data.

The choice of a particular effective potential function, that is, the choice for

the shape of the potential bowl in Figure 2, is far from unique, for, in

addition to the dependence on m, the bowl shape could be affected by adjusting

the exponent in the attractive term.

At a reduced temperature T =3.0, a(T), and hence the equation of state, was found to be relatively insensitive to m, that is, to the shape of the bowl

appended to the hard sphere core. At higher temperatures, the sphere diameter

becomes completely independent of m. In such a case, the procedure used for the

evaluation of the sphere size can be further simplified, as will be discussed

below. On the other hand, as the temperature is reduced, the hard-sphere size

tends to become increasingly sensitive to m. This results in a useful procedure

for determining an optimum value for m for a given substance. This consists of

comparing experimental PVT data with those predicted by the HS equation at a low temperature for effective potentials characterized by several values of m

until a best fit is obtained. Because of the insensitivity to m already

described, the value of m chosen can obviously be used at higher temperatures.

The method has thus been modified to produce a good low temperature fit. Haar

and Shenker used an isotherm near 3/2 times the critical temperature for the low

temperature fitting.

We include, in the next section, a comparison of the predictions of the HS

equation with other theories and with experimental • data for argon and nitrogen.

It is easily determined that m=9 and m=10, respectively, are reasonable values

of m for these fluids. Since the sensitivity of the equation of state

properties to m is weak except at low temperatures, m=9 can also be taken for

nitrogen when the temperatures of interest do not extend much below twice the

critical.

2 . Relationship to other fundamental equations of state.

It is useful to compare the HS equation with equations of state which have

the appearance of being more fundamentally based. The latter are invariably

much more complicated than is the HS equation so would have to produce far

superior results to it to warrant their choice over it. In this comparison,

20
AEDC-TR-76-85

we shall include an example of an integral equation (the Percus-Yevick) as well

as the perturbation theory approach of Barker and Henderson. Although the latter

is developed totally within the language of statistical mechanics, it is not

unrelated to the Haar-Shenker approach which has been described with an emphasis

on phenomenological language.

A major objective of statistical mechanics is to predict the properties for

real fluids at high densities from the known properties of the dilute gas. To

accomplish this, it is usual practice to reduce exact theories (such as that

associated with Eq. (2)) to theories in which interactions among several

particles are pairwise additive. The properties of the dense fluid can then be

formulated in terms of the detailed structure of the potential functions which

describe the forces between pairs of particles. Considerable progress has been

realized using this approach. Relevant to our work are the expansions in

density based on Percus-Yevick (PY) theory and expansions in reciprocal temper19 20

ature using the Zwanzig theory as modified by Barker and Henderson (ZBH)

.

These have been tried for several potential models including the hard and soft

spheres, the square well, and the Lennard-Jones potential. The ZBH temperature

perturbation theory appears to be the more successful when compared to results

of computer experiments, particularly at liquid temperatures but also for gases

at high densities and at high temperatures.

Though the ZBH theory is a physically satisfying approach and does compare

well (except at low liquid temperatures) with results of "computer experiments",

we note several practical limitations. The most serious of these is associated

with the application of the theory to real fluids. In the ZBH theory the

perturbation terms are obtained as an expansion about the hard sphere. Barker

20

and Henderson have shown that the theory can yield quantitative results when

the sphere diameter is expressed in terms of a kind of Boltzmann average of the

molecular separation, where the average is taken over the positive energy region

of the pair potential. Sphere sizes calculated in this way are also contained

in Figure 3. It is a fact, however, that the shape of the pair potential even

21 for simple systems is quite uncertain. In this connection it has been proven

that for realistic (non-monotone) potentials it is not possible to obtain an unambiguous representation of the pair potential from second virial coef f icients,,

Also at high densities the sphere volume tends to affect the equation of state

21
AEDC-TR-76-85

properties somewhat like an excluded volume, so that ambiguities in the sphere

diameter are amplified in their effect on the equation of state at high

densities. We shall present an example of this below. Thus, though the ZBH

theory is theoretically satisfying, its application would seem to be limited to

situations where the pair potential is known, as is the case for "computer

experiments". The fact that the use of this equation of state to produce thermodynamic tables requires a complicated numerical integration poses a second

limitation to the ZBH approach. Because of this, results have so far been

obtained as a theoretical end in themselves and are therefore of somewhat

limited utility to the engineer or scientist who desires a simple analytic

representation of the equation of state as a predictive tool.

It should be appreciated that the hard-sphere diameters are only conve- nient artifices. Physical interpretation is meaningful only in the context of

the particular overall theory. However, since the ZBH and the present theory

are quite sensitive to the hard-sphere diameter, agreement between them at

least to within several percent for a(T) would be necessary for the two

approaches to yield comparable equation of state properties at high densities.

3. Comparisons of the HS equation of state with experiment and with other

theories .

In this section the HS equation of state is used to calculate compressibility factors for the real fluids argon and molecular nitrogen. The results

are compared with PVT experimental data, and, in the case of argon, with

results of "computer experiments" and with ZBH theory. The experimental second

virial coefficients used as input data to obtain the required parameters are smoothed values calculated from model pair potentials as determined from the

experimental second virials.

The Figs. 4-8 contain isotherms calculated for argon using the HS

equation of state and the ZBH theory. The compressibility factor is plotted

versus the reduced density, p = Na 3 p, for isotherms from 119. 8K (=.8 critical)

to 673. 15K. Curves labeled #1 represent the compressibility factors predicted

by the HS equation of state. The curves #2, #3 and #4 represent results of

second order temperature perturbation ZBH theory, for a Lennard-Jones gas with

22
AEDC-TR-76-85

parameters e/k = 119. 8K and o = 3.405 A as reported in references 22, 23 and 24,

respectively: the curves #2 and #4 are obtained by numerical integration of

approximate expressions for the perturbation terms and are essentially equivalent treatments; the curves #3 refer to an "exact" calculation of the second

order perturbation terms via Monte Carlo techniques. The curves #5 represent

experimental compressibility factors measured for argon. The data points on Figs. 4, 5 and 6 are results of "computer experiments" for the Lennard-Jones

gas with the above parameters, the circled points referring to Monte Carlo

results and the boxed points to molecular dynamics results. The "computer

experiments" could be uncertain by 5 to 10%.

Since use of the Lennard-Jones potential with the above parameters does not

25

produce second virial coefficients which fit the data for argon below 200K , Figs. 4-6 are interesting primarily as comparisons of theory with the results of

the "computer experiments." These figures show that the present theory, as well

as the ZBH theory, are only qualitative at these temperatures. However, the

ZBH curves #3 yield a somewhat closer approximation to the results of "computer

experiments" for liquid densities. As previously stated, the Haar-Shenker

equation of state tends to degrade at low temperatures. But it is apparent from

Fig. 4 that the theory is still at least qualitatively good at temperatures

even as low as .8 of the critical temperature. In fact it is only for curves #3

(which involve extensive numerical calculations to evaluate the second order

perturbation terms) that the temperature perturbation results are a significant

improvement over the HS equation of state.

An explanation for the fact that the "wrong" potential gives the correct

results in the ZBH theory for argon has been offered by Barker, Henderson and

26 Smith. They argue, that the (12,6) pair-potential with the above parameters

happens to be an "effective potential" that, to first order for argon, accounts

for high density non-additive effects.

The Figs. 6-8 include the temperature region for which the pair potential

used produces a good fit to the experimental second virial data. Comparisons

of our results with experimental PVT measurements for these temperatures are therefore more meaningful. The HS equation of state (curves #1) tends to follow

the PVT experimental data (curves #5) fairly closely at the lower densities, up

23
AEDC-TR-76-85

to about p = 0.6 (which is approximately twice the critical density). At

higher densities it tends to yield values that are low, but in most cases only

by less than 5%. The results are roughly comparable to those for the more complicated ZBH theory. The ZBH theory (curves #2, //3, #4,), however, tends to underestimate at densities near the critical and to overestimate at high

densities. The results of the "computer experiments" seem to scatter among the

various theories. Thus, in Fig. 6. at p* = .75, the Monte Carlo results tend

k to favor the HS theory, at p = .9, they favor the ZBH; while at p* = 0.55 they

fall between the results of the two theories, in fact, almost on the PVT

experimental curve #5.

We have stated above that uncertainties in the pair potential for real

fluids limit the utility of the ZBH theory as a tool for predicting the equation

of state properties for such fluids. To illustrate this we compare the equation

of state of argon for a particular isotherm with the ZBH and HS theories, in

which different inverse power representations for the pair potential are used,

21 each of which produces an equivalent fit to the experimental second virials. In Figs. 10 and 11 equation of state results are presented for two representations

25 for the pair potential for argon: the (18,6) with the parameters

e/k = 160.87 K

a - 3.261 A .

as curves #1; the (12,6) with parameters given earlier, curves #2. Fig. 10

includes the results for the ZBH theory; Fig. 11 those of the HS equation. Both

figures refer to the isotherm 239. 8K. In these figures the compressibility

factor is plotted against the density in amagats. All the numerical data for

23

Fig. 10 were furnished by Toxvaerd. The results for Fig. 10 show the two curves near coincidence up to a density about 2/3 the critical (critical density for

argon = 300 amagats) but sharply divergent at higher densities. By comparison,

the results for the HS theory in Fig. 11 are relatively potential independent.

The reason for this independence is obvious. The HS equation depends only on

the experimental second virial coefficient and its first derivative. Since the same values for the second virial coefficient are obtained for either potential

for a range of temperatures near this one (T = 2.0), the calculated second

24
AEDC-TR-76-85

virials and their first derivatives are equivalent. Since the details of the

repulsive branches of the potentials are different, however, the calculation of

the hard sphere diameter according to the prescription of Barker and Henderson

yields different results for the two potentials. This difference, in turn,

produces a difference in the equation of state predicted.

Finally, in Figs. 12-15 we present results of using the HS theory for calculating the equation of state of molecular nitrogen. The compressibility factor is

plotted versus the fluid density in amagat units. Here, the experimental

second virial coefficients were represented by smoothed tables calculated using

the Lennard-Jones pair potential with parameters e/k = 95.781K, and o = 3.712 A. The curves #1 refer to the present theory; the curves #2 to experimental PVT

measurements. The latter extend to 10,000 atm in Fig. 15. The sphere diameters

are obtained from Fig. 3 with m=10. The comparison of theory with experiment

for nitrogen is quite similar to that for argon. As with argon, the HS theory

is in good agreement with experiment at low densities, but at the higher

densities the experimental PVT isotherms tend to be slightly steeper.

4. Summary .

Based on the virial expansion and on the behavior of the virial coefficients

at high temperature, Haar and Shenker derived a quantitative yet simple equation

of state which is valid for real fluids over a density range from the dilute

gas to densities approaching that of the solid at temperatures above twice

critical, and which requires only a knowledge of the second virial coefficient

and its first derivative at each temperature. This equation of state is much

simpler than the ZBH temperature perturbation theory and furthermore does not

require reference to the precise details of the pair potentials. The equation

proved to be quite successful in comparisons with experimental data. These

comparisons were naturally carried out for ordinary temperatures since experimental data exist only for such temperatures. Figures 2 and 3 contain such

comparisons for argon and nitrogen.

Through the use of the intermolecular potential function, the HS equation

of state can be used, in a very simple and straightforward manner, to extrapolate PVT data in both the temperature and density directions. Thus, the inter25
AEDC-TR-76-85

molecular potential is used to extrapolate the second virial coefficient as a function of temperature and this, in turn, is used with the HS equation to cover

all fluid densities. This approach is therefore particularly well suited for

aerodynamic calculations since, as already mentioned, these require temperature

extrapolations of up to a factor of ten for densities up to those which, at low

temperatures, correspond to the liquid. Because of the decrease in the attractive contribution to the second virial coefficient with temperature, this

equation of state should improve with increasing temperature. As a result, the

comparisons made at ordinary temperatures should easily be sufficient for

estimating the expected adequacy of the theory at high temperatures. Since this

equation of state depends only on the second virial coefficient which, in turn,

is determined once the intermolecular potential is known, it becomes possible

to develop an entire PVT surface given this intermolecular potential function,

or equivalently , given a sufficient (generally small) amount of low density PVT

data for the substance at ordinary temperatures, from which data the intermolecular potential can be obtained.

The derivation of the equation of state is based on an expansion in

density, where the reference state is a gas of hard spheres. We have presented

a plausibility argument which indicates that, when the sphere diameter is chosen

appropriately, the terms that account for the differences between the properties

of the actual fluid and those calculated for a fluid of hard spheres are sharply

attenuated at temperatures above the liquid-vapor critical temperature, for

densities up to that of the solid.

It has been shown that the temperature dependence of the equation of state

at temperatures above the critical is determined by two parameters which depend

on temperature and that these parameters can be obtained from the second virial

coefficient and its first temperature derivative.

III. THE EXTENSION OF THE HAAR-SHENKER EQUATION OF STATE TO MIXTURES

Haar and Shenker developed their equation of state for use with pure substances. Our needs are, of course, for a theory applicable to mixtures since

at aerodynamic temperatures even pure nitrogen becomes a mixture as a result of

dissociation and ionization. There are several ways in which a theory for pure

26
AEDC-TR-76-85

fluids can be extended to mixtures. Perhaps the simplest are what have been

called the one and two fluid van der Waals' theories. In the one fluid model

the mixture behaves as if it were a single component fluid (which can be called

the equivalent fluid) with any parameters used to describe the fluid being

averaged over the composition. This averaging is carried out by taking the

parameters associated with the individual constituents and suitably weighting

them to the extent that the constituents are present in the fluid. For example,

a fluid made up of a mixture of hard spheres would be described as a one fluid van der Waals' model by the equations associated with a single component hard

sphere fluid but with the single relevant parameter (that associated with the

molecular diameter) averaged over the composition. A natural way of doing this

is to take for the hard sphere volume

where a is the diameter of the molecules of the equivalent fluid and a.. eq

th

ij

the diameter for the interaction between a molecule of the i species and one

of the species in the actual fluid. For hard spheres a.. — = (a. + a.).

ij 2 1 3

The equation of state for this mixture is easily derived. From the virial

theorem, it is possible to derive a general equation of state for a fluid in

28 terms of the distribution of pairs of particles in the fluid. Because of the

abruptness of the hard sphere interaction, this equation, for a one component

hard sphere fluid, reduces to

where g(j) is the probability that a pair of molecules will be found a distance

a apart. For a mixture of such spheres, the equation of state becomes

o

3 = T X, X

eq i

pa

3 g(a) (10)

X. X. a

3

, i 1 ij (11)

27
AEDC-TR-76-85

ttl where g. . (a..) is the probability that a pair of molecules, one of the i

1J 1J th

species and one of the j species, are separated by a distance a . In the one fluid model, then, (11) is replaced by (10) in which

§ij (o..) = g (a)

and ( 12 ) o

3 = y x. x. o?.

This model is consistent with the results of the density expansion in the

Ursell-Mayer virial equation of state.

A number of models other than (12) can be devised for using (10) in place

of (11) for a mixture. Each of these, though reasonable, does not lead from

(11) to (10) in a natural way. One might, for instance, average the diameter

rather than the volume so that

a = y^ X. X. a. . = y X. a. . since a.. = - (a . + a ) i j xj x xi ij 2

v xx 23

i. j i

One might also average the volumes over like species only, i.e. a

3 = E X_^ CT |^-

We shall consider (12) as the only reasonable model, particularly since it is

the only one that leads to results which are satisfactory.

A two fluid theory can be obtained from the approximation

, ( 0 ) = i U . (a. .) + g. . (a.,) I 'i: iy 2 L 11 11 j J J

If, at the same time, one uses the fact that, for hard spheres a . . = r (o.+ a.)

ij 2 x j (11) becomes, for this model,

28
AEDC TR-76-85

RT ~ 1 + ~ P L X

i

°- 8±i

(o. . ) (13)

i

In order for this to be of the pure fluid form, it must be assumed that g is

the distribution function for a pure fluid of hard spheres of volume a. given

by

2

Thus, if this were a binary system, it would appear to be made up of two pure

fluids of different diameters, namely a

3 = (Tjj + X

2 a

3 2 and a

3 = X

1 ar^ +

X

2 °22' Hence the name two fluid theory.

A fluid of hard spheres is a highly idealized model for an actual fluid so might be thought to be quite useless for testing theories of fluids. There are, however, properties for fluids of hard spheres as calculated by computer

simulation methods. If such results obtained for hard sphere mixtures are taken

as "experimental" data, it is then reasonable to compare them with (12) and (13)

to see which is the better approximation to a mixture of hard spheres. This

29 has been done by Henderson and Leonard who found the one fluid theory to be,

by far, the superior.

In a later paper, these same authors carried out a similar comparison"^ for

a fluid whose molecules interact in accordance with a Lennard-Jones (12,6)

potential function. Such a potential contains both attraction and repulsion,

and is often characterized by a molecular diameter and by a potential well depth.

The one fluid model now follows from ea

3 = E X. X. e.. a3 ... Comparisons were ±a i J ij ij

made for an equimolar mixture and the one fluid model was found to be superior,

particularly in the prediction of the excess free energy and heat of mixing.

The intercomparison for the excess volume of mixing was somewhat ambiguous,

however.

29
AEDC-TR-76-85

It would thus seem clear that the most reasonable way in which a theory for

a pure fluid can be applied to that for a mixture is via a one fluid model.

This we have chosen to do for the HS equation of state. This simplifies our

task since we need only write down the pure fluid form of the equation of state

which we have already derived and interpret the parameters in the form of the one fluid theory.

By (5), the HS equation of state for a pure fluid is written

pv = i±v+£ + 4/5 _ A

y (14)

(l-y)

RT „ ,3 \b

where y = . Then Eq . (14) becomes an equation of state for a mixture

according to the one fluid model, if

b = Y\ X. X. b. . (15)

i J

1 2 13

b = t x. x.

4 A 1 J B. (16)

It should be remembered that Eq. (9) must be solved for b„ for each of the

interactions included.

IV. THERMODYNAMIC FUNCTIONS IN THE ONE FLUID HS MODEL

In this section of the report we shall be concerned only with the real gas

contributions to the various thermodynamic properties. The ideal gas contributions and certain precautions required in their calculation are contained in

reference (3)

.

A. HELMHOLTZ FREE ENERGY, CHEMICAL POTENTIAL, AND ACTIVITY COEFFICIENTS.

The combination of the equation of state (14) and the one fluid model, (15)

30
AEDC-TR-76-85

and (16), is sufficient for the calculation of the equation of state of a real

gas mixture at high densities. With such an equation of state it then becomes

possible to calculate the real gas part of the Helmholtz free energy from

from which the chemical potential can be calculated using

u

i \dnj T, p, n_. when n. refers to all n., j ^ i.

J J

When the ideal and real parts of the chemical potential are separated it

becomes possible to identify the activity coefficient for each species. By

properly combining the activity coefficients for each species taking part in a chemical reaction, one can define an effective activity coefficient for that

chemical reaction. This effective activity coefficient directly modifies the

equilibrium constant to produce the effect of non-ideality on the chemical

reaction. To see how this goes, let

A „.

where is the ideal gas part of the chemical potential and Ay^ the real

part. It follows, then, that the activity coefficient for this species is

given by

Ay.

w = lny

±

The effective activity coefficient for species i in a chemical reaction with

stoichiometric coefficients v . , is then given by

31
AEDC-TR-76-85

which is equivalent to

y] = exp [(E V. . Au. - Ay ) /RT]

j

3 (17)

We shall now proceed to the detailed calculation of the free energy and

chemical potential for the HS equation of state.

If at each volume the reference state for the free energy is taken as that

of the ideal gas at that volume then

V

A(V) = A(«) -f PdV

A

(0)

(V) = A

(o)

(~) - / P

<•> dV

But, since A (oo) = (oo) it follows that

A(V) - A

(o)

(V) = - / (P - P

(0)

) dV

nZRT j ^(o) nRT , „ . ^, r ^ , Since P = —— and P = ——, where Z is the compressibility factor and n the

total number of moles, it follows that

A(V) - A

(o)

(V) = - nRT

OOf {1-1) ^

The real part of the chemical potential is then given by

32
AEDC-TR-76-85

RT

y .-

i

(o)

i RT

The evaluation of the common, integral is obviously necessary to progress in the

derivation of both the free energy and the chemical potential.

According to the HS theory,

Cl-y) 3 b

>^2 - 2

>T+4) + (| - 1) 4y

Cl-y) 3 b

Since y» = 7-=:, — = - & . Thus

4v' V y

Cl-y) 3 „ b 1 '

In 1^, let 1-y = x so that dy = -dx and

h

.j-^'a-x) 2 -2(i-x)+4 dx

1 x1 x

Since, obviously,

33
AEDC-TR-76-85

I

2 = - 4y(E

-1)

then

/

CZ-l) f- = -

| [Ojr

y

) 2 -1] + ta (l-y) - 4y (5 -1) Q8)

so that

A(V) - A

(Q)

(V)

RT % [ (IT>

2 -U -n to (l-y) + 4ny(? - 1)

The non-ideal contribution to the chemical potential is then given by

Ay

RT

The first term has already been evaluated. The second is obtained by differentiating (18). Thus

3 (-2)(-l) 1 B .J 9y . 3 .B, 2 " ,3 -lTy

-

X)

an"." ^ to.^ (l-y) 7

) l i

At this point use must be made of the one fluid model. Thus

1 3b 1 3n. 4V 3n. 4V

l l

2 E n. b..

(E n.)

2

2 En. n. b., .

(E n.)

1_ 2 4V n

En. b. . - b

34
AEDC-TR-76-85

Also

9 1 3B B_ 3b

"ST CB/b) " b 3n

± b

2 3n

±

3B 9b Since, obviously, -— will have the same form as — but with B replacing b, dn^ dn^ we can immediately write down

dn. <V 2 n b

I n. B. . - r E n b

I

] lJ b . j ij

Combining all results obtained to here yields

^i _ 3 RT 2

£n Ci-y) +4 C- - Dy

b

3 1-y b

d-y)

22

b

I n . b . . ....I jJ

, 2 + 4y 5

"En. B.

.

*2 En. b. . B j ij

For simplicity we write this as

En. b. . Ay. j u

RT"

= y o + y

l nb

E n. B. . + y

2 ^-^b—

where

35
AEDC-TR-76-85

y o =

2 -1

Xi-y)'

- In (1-y) -

^i- 2 ?\— 3 +

1=7 " 4

1 ! [Cl-y)

J 7

j

P

2 = 8y

Substitution in (17) then yields

= exp - cj u + —— E n./Ev..b, n -b.„)

i o nb

% %\. ij j£ i^/

+ ^ n„ / E v.. (B.„- B. n )\ /RT

nb

t I L xj

C19)

where -oj. = E v.. -1 is the net decrease in the number of particles in the i

reaction defined earlier. With this expression for y

' it becomes possible to

calculate the effect of the HS theory on the equilibrium composition of the

mixture

.

B. ENTROPY, ENTHALPY AND GIBBS FREE ENERGY

It is now possible to calculate all of the thermodynamic properties predicted for the HS model. The entropy follows from

so that the real gas part is given by

36
AE DC-TR-76-85

AS

R

VIA

3T

AA T _3 /1A

RT ST \RT

But

3^ +4 M

+ 4ny ty^).

Now

and

3T Cb

} " b 3T b

Z BT

so that -r= TTT = 1A

T \ RT

l(l-y)

+ J- -4 1 » |B

3 1-y | b ST b 3T

Thus

_S_ R

3n

2 c-U 2 -1

1-y

-r. in CL-y) + 4 n y G -D

(Cl-y)

3 ^ i b 3T

4yT _3B b 3T

3b 3B

Trrr and -r= are expected to be quite small because of the high temperatures of

dl dl

interest. It is therefore possible that the last two terms might be negligible

with respect to the first three. This would need to be examined, however.

This expression for the entropy and the earlier one for the equation of

state can be used to calculate the enthalpy and Gibbs free energy. Specific

37
AEDC-TR-76-85

heats can be calculated from these functions either by numerical differentiation

or by means of algebraic expressions which can be derived by differentiation of

3 the above expressions .

V. NUMERICAL METHODS

The primary step in calculating the thermodynamic properties consists of

the determination of the equilibrium composition of the mixture. This requires

the solution of equations (1) subject to the constraint that the total mass of

each nuclear type be conserved. The effective activity coefficients, , are given by equation (19) for the contribution of the HS model to the equilibrium

constant and in appendix B of reference (3) for that of the Debye-Huckel theory

of ionic solutions. There are two lower level approximations to the HS model

which we have already considered in earlier work. We have already used the

ideal gas approximation, in which all y

' are taken equal to unity. We have

also made use of the second virial coefficient approximation for which the y^*

of the HS model are replaced by the term linear in the density in their density

expansions (which term appears in appendix B of reference (3)).

The method used for the solution of the equations for the concentrations in

both the ideal gas and second virial gas has been described in some detail else3 where. The procedure used here was required to be applicable to all three

approximations, i.e. the ideal gas, second virial and HS. Although based on these earlier methods , there was a considerable modification of certain details

of the original approach to allow for greater flexibility in obtaining solutions.

In order to permit reference species to be chosen at will simply through

modification of the input data (as described in reference (3)), the calculation

of the equilibrium constants was made part of the computer program. This led to

the discovery of an error in the equilibrium constant used for 0^ in the

previous NBS tables. Our correction of this error caused the search method used

for finding the electron concentration to be unstable at low temperatures. A

not insignificant amount of time was spend in isolating this problem and in

correcting it. The problem was associated with orders of magnitude increase in

the concentration of 0^ which occurred after the correction was inserted. As a result, the electron concentration now became orders of magnitude smaller than

the concentration of 0_ and, in fact, was now calculated as a small difference

38
AEDC-TR-76-85

between the concentrations of NO and CL. This often led to negative interim

values for the electron concentration during the non-linear search procedure.

Since the electron concentration has only a minor effect on the values obtained

for the other reference species, the problem was solved by instituting a grid

search for the electron concentrations, over positive values only , whenever a negative guess value was obtained for that concentration, leaving all other

concentrations fixed at their current values.

At low temperatures, it was also necessary to produce a very strong

"damping" of the electron concentration by means of a q value smaller than

(but in the neighborhood of) unity. This reduced the "natural" excursions in

the electron concentration from iteration to iteration.

Initially, precisely the same search procedure was used for the HS model

as was used earlier for the ideal and second virial gases except that the

for the HS model now had to be computed at each iteration. Problems arose,

however, because these y.' became quite large for certain species at the highest

densities. In fact, y.' values approaching 10 were encountered. This produced

instabilities in the search procedure used. This problem was solved by what

might be called a dual level search procedure. In this procedure an initial set

of y

1 were computed based on the initial guesses. These y.1 were maintained

constant until a set of concentrations was obtained which satisfied the mass balance equations for the gas with these initial y

' values. With these

concentrations, a new set of y^' was computed. These new y.' were now held

constant and a second set of concentrations obtained which satisfied the mass balance equations for this second set of y.' values. The procedure was carried

out repeatedly until the largest change in the y^' when recomputed was less than

a given tolerance. To reduce the possibility of producing numerical instabilities, only a fraction of the change calculated for each y.

1 was used in any new

iteration and, in any case, the y

' were not allowed to change by more than some arbitrary factor. In order to reduce the computer time required, a coarse tolerance was placed on the values of the concentrations accepted as correct for

the first few sets of values (generally incorrect) obtained for the y

1 in this

procedure.

The present computer program has, to a large extent, retained the ability

39
AEDC-TR-76-85

of the earlier program to obtain solutions essentially independent of the

initial guess values for the concentrations of the reference species. The

electron concentration is perhaps an exception to this but then mainly only

when the computation on an isotherm begins at quite high densities.

Since the few problems which we have encountered with obtaining solutions

must depend on the specific search procedure used by us, it would be useful to study the process of solution for other search procedures. Because of its

simplicity, we would place particular emphasis on a study of the direct search

31 method of Hooke and Jeeves

Once the species concentrations are obtained in the above manner (with the mass balance equations satisfied to the tolerances specified for each reference

species) , the thermodynamic properties of the mixture can be calculated in a straightforward manner.

VI . RESULTS

In this section we discuss a number of special features of the results

obtained, reserving our discussion to those results which illustrate the effects

of density. We compare results obtained by us for the ideal gas, second virial

coefficient and HS models with each other and, where appropriate, with the

results of the GB extrapolation. We shall also compare the concentrations

predicted for various species among the three approximations used by us. The

8 Grabau-Brahinsky model does not include the calculation of the species concentrations so cannot be included in that part of the discussion.

Certain interesting features of the predicted dependence of the concentrations on density are also discussed. Of particular interest is the demonstration

32 of the "strength" of Le Chatelier's principle

A. THE DEPENDENCE OF SPECIES CONCENTRATIONS ON DENSITY.

The formalism which we have developed has the capability of predicting

density effects on species concentrations. This ability is important for two reasons. It is obviously important when the concentrations of the species them40
AEDC-TR-76-85

•

selves are needed, as might be the case when specific species have particularly

interesting properties. Examples might be the electrically charged species

(which obviously affect the electrical conductivity very strongly) and any

species which radiate in a useful part of the spectrum. Being able to predict

density effects on concentrations is also important in the calculation of the

density effects on the thermodynamic properties. The expressions for the thermodynamic properties of the reacting mixture are sums over concentrations times

properties of individual species plus sums of products of concentrations times

properties of pairs of species. In our model, the latter enter into the calculation of the average (i.e. one fluid) molecular diameter. The thermodynamic

properties of the mixture therefore depend on temperature and pressure through

the dependence of the concentrations on these state parameters.

Density effects on the thermodynamic properties of air have generally been

calculated by others on the assumption that the species concentrations either

do not change with density or simply obey the ideal gas mass action law. This

was essentially the assumption of Grabau and Brahinsky. As we shall show, this

kind of assumption breaks down, particularly at the highest densities where the

identities of the major species in the mixture change, when density effects are specifically taken into account. This can change the entire character of the

gas being studied.

Care must be taken here not to place too much emphasis on the actual

numerical values obtained by us for the concentrations. The species concentrations at high densities are expected to be sensitive functions of the assumptions

made in any calculation of this kind. In our model, they can be expected to depend very strongly on the intermolecular forces used (especially at the highest densities) as well as on the method of calculating the hard sphere diameter

at temperatures for which the second virial coefficient has a negative slope.

This problem is associated with the fact that the effective hard core potential

produces a second virial coefficient whose slope is always positive. The

relationship of the actual numerical values obtained by us to actual air might

also be expected to depend strongly on our omission of such species as ^0^, C^ t 0^ and 3ome of the molecular ions.

Because of Le Chatelier's principle, errors introduced into the calculation

41
AEDC-TR-76-85

of the concentrations for any one of these reasons should be less than that

which might be estimated from a direct calculation of the apparent effect of the

omission. This principle is essentially a statement of the competition among

the various chemical reactions, because of which species concentrations "resist"

change. Thus, suppose a correction to the equilibrium constant in a certain

reaction produces an increase in the concentration of a particular species.

This increase will cause a competing reaction containing that same species to move in such a direction as to reduce the concentration of this particular

species. The overall effect of the original correction to the equilibrium

constant is thereby reduced. An example of this is described below. Because

of this and because any errors resulting from our various assumptions can be

expected to affect the different species perhaps randomly according to sign, the

ttl overall effects of our zero order assumptions on the properties computed should

actually be much smaller than might be expected from an examination of the

effects of the various separate approximating assumptions. This error reduction

might be considered to be a decided advantage in favor of the use of a microscopic molecular model such as ours.

The operation of Le Chatelier's principle can be seen by comparing the

change in various concentrations actually obtained in the calculation when the

density is changed at constant temperature with that which might have been

-0)

.

expected from the change in the density factor (p/p

o

) 1 in equation (1).

Consider, for example, going from log p/p

q = 2.0 to 3.0 at T = 3000K. For the

species N„0, a concentration enhancement by a factor of 1.6x10^ could have been

4 expected based on the value y. = 4.3x10 at log p/p = 3.0 and to. = 0.5, whereas

the actual enhancement obtained

X involved only a factorox

of 1.7x10^. There

was, therefore, a reduction by a factor of 100.

Figures 16 and 17 contain plots of concentrations against density for a number of species for the temperatures 3000K and 9000K. The concentrations are those which were calculated using the full density effect with the HS equation,

and those based on the ideal gas.

The dependence on density of the concentration of oxygen and the effect of

this dependence on the concentration of the other species are particularly

interesting. The rapid decrease in the concentration of molecular oxygen at

42
AEDC TR-76-85

the highest densities is very dramatic for the HS gas, especially when compared

to its behavior in the other two approximations. As a result, the gas at 3000K

and 1000 times normal density consists mainly of ^ and NO rather than of ^ and

O2 as predicted in the ideal gas approximation. Secondary effects from this

reduction in 0^ concentration can be seen in the marked reduction in the

concentrations of 0~, 0* and 0 predicted by the HS theory as compared to the

predictions of the ideal gas approximation.

The reduction in 0£ concentration "results" from the enhancement in the

production of NO, NO£ and ^0 and the associated requirement for the production

of atomic oxygen for the formation of these species. The cause of the enhancement of these species can be seen in Table 2 which contains the y.' values for

all species at T = 3000K for several densities. In the concentration units used

by us, it is necessary to multiply the y by (p/p ) 1 to obtain the full

density effect along an isotherm. As is pointed out below, in some instances

this produces enhancement factors of 10^ for the equilibrium constant.

An interesting effect, which might be called a second order effect of the

reduction in 0£ concentration, was seen among the carbon containing compounds.

According to the y^' values of Table 2

} the concentration of C should drop

drastically with increasing density. This is especially true when the factor

-1 -3

Cp/p ) is added (bringing in an additional factor of 10 at the highest

o density). Instead, the concentration of C increases slightly with density.

This comes about because the reduction in 0^ concentration causes a decrease in

the production of CO2. This "frees" carbon atoms which become available for

the enhancement of the CO and C concentrations.

Another second order effect is that of the slight increase in electron

concentration with increasing density, quite the opposite of the ideal gas

behavior. This results from the decrease in the concentrations of O2 and 0

which results in additional free electrons. The increase in the electron

concentration, in turn, results in a drop in the concentration of N0+

, the main

electron producer, causing the increase in electron concentration to be somewhat

reduced from that expected purely on the basis of reduction in O2, another

example of Le Chatelier's principle.

43
AEDC-TR-76-85

Those species for which y^' has the value of unity in Table 1 are either

reference species, for which no equilibrium constant is needed, or are ions for

which no estimates of virial coefficients were available. It should be noted

that the y^

1 values for the species ^0 and N0^ begin to approach 10^ at the

highest densities. This effect is further enhanced by the factor 33 for each

-1/2 arising from the factor (p/p Q

) . The overall density effect for these

species at p/p

q = 1000 as compared to p/P

Q = 1«0 is then such as to multipy the

equilibrium constant by a factor of over 10^. It is therefore obvious why there

is such a strong dependence of concentration on density.

The competition at the highest densities among NO, NC^ and ^0 (and

especially between the last two) is particularly interesting. As the density

increases, all three species' concentrations increase. At 3000K this occurs mainly at the expense of molecular oxygen. As the density increases, however,

the enhancement of N^O begins to proceed at such a pace as to "require" oxygen

atoms from other reactions so that, ultimately, the increase in N^O concentration takes place at the expense of the concentrations of NC^, and NO. At

9000K, the initial enhancement occurs at the expense of the oxygen atom con- centration but, at intermediate densities, produces a reduction in molecular

oxygen concentration. Eventually, the increased N£0 concentration occurs as a result of a reduction in NO and NO^ concentrations.

The absence of estimates of the virials for interactions involving ionic

species obviously leads to errors at the highest temperatures where charged

species become non-negligible. Although these charged species do not dominate

in our approximation, it is conceivable that they might become major constituents through enhancements caused by large y

1 values which might be obtained.

This effect might be expected to be smaller for the ionic species than it was

for the neutrals NO, N0„ etc. for two reasons. First of all, even at high temperatures, repulsion between ions of like sign is reduced by the effects of

attraction between those of opposite sign. Secondly, our inclusion of the DebyeHuckel limiting law (see reference (3)) already Includes part of the interaction

between charged particles.

44
AEDC-TR-76-85

B. DENSITY EFFECTS ON THE EQUATION OF STATE.

Table 3 contains the compressibility factor as a function of density at

selected temperatures for all three approximations as well as for that of Grabau

and Brahinsky. The ideal gas approximation contains no correction to the

equation of state for the effect of density. Since there is a "density effect"

in the law of mass action as written in the units used by us there is a variation in the calculated compressibility factor with density. By our

definition, the compressibility factor for the ideal gas is simply the total

number of moles of the mixture. According to the equations of chemical equilibrium for the ideal gas as written in our units, the effective equilibrium

constant for a particular reaction in the ideal gas approximation is given by

eff . nij -co.. K. = K. (T/T.) 1 (p/p ) 1 and so increases with density for those reactions

i l (J o in which a net decrease in the number of particles results in the production,

from the reference species, of the species associated with the reaction. Such

a reaction, for example, is

N

2 0 = N2 +

| 02

where two molecules are produced for every three which react. On the other

hand, a reaction for which there is a net increase in particles has an equilibrium constant which decreases with density. An example is

0 = ~2 °2

where two oxygen atoms are produced for each molecule of molecular oxygen. The

net result of this combination of enhanced effective equilibrium constant with

density for reactions in which the number of particles is decreased and decreased effective equilibrium constant with density for those in which it is increased is to produce a decrease in the total number of particles with density - a result which is clearly visible in Table 3.

With very few exceptions (e.g. some of the ionic species) all virial

45
AEDC-TR-76-85

coefficients are expected to be positive at the temperatures of interest (and,

in fact, even at much lower temperatures). It is for this reason that the HS

model is expected to be particularly appropriate here. Because of this positive

behavior, at these temperatures density corrections to the equation of state

must be positive at all densities and are required to increase with density.

Any approximation which has this character must therefore produce a density

correction to the ideal gas which has the correct sign, if not the proper

magnitude. The simplest such correction is associated with taking only the

second virial coefficient. According to Table 3, such an approximation does

produce an increase in compressibility factor with density for the intermolecular

potential functions used here since these produce positive second virial

coefficients. The decrease in the magnitude of this correction with increasing

temperature at constant density is caused by a reduction with increasing temperature in the magnitude of the second virial coefficient for many of the interactions used since the temperatures of interest are above those at which the

second virial coefficients for these interactions exhibit maxima.

The GB approximation is strongly dependent on our second virial coefficients, since those authors made use of our earlier results for the second

virial gas to tie down their extrapolations at 5000 and 6000 kelvins. For this

reason, their predictions should be in close agreement with the results for the

second virial gas up to the densities at which that approximation is expected

to be valid or to densities somewhat below 100 amagats. This is essentially

the behavior exhibited in Table 3. Since the GB model is based on the results

of our second virial approximation, this agreement is not a test of the GB model

but rather serves as a test of the computer program developed by those authors

as well as of our own and as a test of the GB input data as obtained from our virials. Since their approximation includes an estimate for the effect of third

virial coefficients and since the HS approximation does also, and since both of

these are of the same sign, the predictions of the GB calculation should be

expected to agree with results for the HS model to slightly higher densities

than the second virial gas (and this is also exhibited in Table 3).

Our present results are entirely compatible with the predictions of the

second virial gas. This is no more than expected since both calculations were

based on the same second virial coefficients and since our model reduces

46
AEDC-TR-76-85

exactly to the second virial gas in the limit of low density when all but

linear density terms are neglected.

As the density increases beyond approximately 100 amagats, our results for

the HS model begin to deviate very rapidly from the second virial results,

becoming over twice as large at a density of 1000 times normal. It should be

noted that there is a similar relationship between the results for the GB

approximation and those for the second virial gas except that since those two

approaches involve very similar approximations, differences between their

results are considerably smaller than between the HS and second virial

approaches, especially at the highest densities.

C. SUGGESTED POSSIBILITIES FOR IMPROVEMENTS IN THE TABLES.

As mentioned above, a number of approximations were made in order to

expedite the completion of these calculations. These approximations were

considered to be sufficiently minor so as not to affect a study of the effect

of density corrections on the concentrations and of the feasibility of carrying

out such calculations. Because of the nature of the calculation and because of

the applicability of Le Chatelier's principle, most of the approximations should

not be expected to affect the accuracy of the compressibility factors drastically. In this section we shall describe ways in which these approximations

might be relaxed in order to produce more accurate tables. Please note that

the order in which the approximations are discussed is not necessarily related

to the order of their importance.

In this calculation we include only the compressibility factor from among

the thermodynamic properties. Thermodynamic properties can be obtained either

through numerical operations on our tables of compressibility factors or through direct calculation of the properties from the equations given in the

text. For some of the properties, such direct calculations require a knowledge

of the temperature derivatives of the temperature dependent parameters associated with the effective potential function, i.e. b(T) and e(T), particularly the

first. This, in turn, requires an improvement in the method used for the calculation of the hard sphere diameters for temperatures above that at which the

second virial coefficient attains its maximum. While the method used by us

47
AEDC-TR-76-85

produces a relatively trivial discontinuity in the value of b(T) at the changeover temperature, it does produce a discontinuity in its temperature derivative

which is much larger. It is clearly possible that these discontinuities in the

db. . values of jj will also be reflected in dZ/dT values obtained in a numerical

differentiation. Hence, even though it might be possible to neglect the temperature derivative of b(T) as being small, it is probably very necessary to improve on the method of calculating b(T) at high temperatures to ensure that

spurious discontinuities in slope are not introduced.

There are a number of methods which could be used for the calculation of

b(T) for temperatures above that at which the second virial attains a maximum

and these should be investigated. Since the effect of attraction can be

totally neglected at these temperatures, the most promising method might be one

in which the problem of finding the two parameters b(T) and e(T) at each temperature is replaced by that of obtaining b(T) only. The present method

essentially does this but in a very arbitrary manner, and must be modified so

as to produce a smooth table of values for db/dT.

An obvious improvement in the tables will also result when the intermolecular potentials used by us are replaced by improved ones. The most

important of these have already been determined by us in earlier work under

this contract. In an earlier report, we estimated the possible effect of this

on calculated tables. Although the effect was shown to be considerably smaller

than were the substantial differences reported here between our results and those

of the GB model, they were nevertheless found to be not negligible. Such a study needs to be made within the context of the HS theory, it being otherwise

impossible to place meaningful estimates of precision on our results. For the

second virial coefficient alone, the ratio of the value predicted for the (18,6)

potential to that predicted for the (12,6) at a temperature of 5000° when the

value at 500° is correctly predicted by both, is given approximately by

dT

/"5000 V

V500 /

2 18 12

= (10)

1

18

48
AEDC-TR-76-85

for a difference of approximately 14%.

Those properties which depend on db/dT would probably only be modified

slightly by changing the potentials since db/dT is small for almost all

27

potentials at the temperatures of interest. In our earlier report this was

shown to be true for the first density corrections to such properties.

Improvement of the tables through the improvement of certain of the intermolecular potentials used (e.g. pair interactions involving NC^) would take

considerable additional effort. This would require a literature search for

experimental data for second virial coefficients and viscosity data for the

relevant species and the determination of parameters for intermolecular forces

using such data.

Related to this but somewhat broader in scope is. >the need for a detailed

study of mixing rules by means of which potential functions which describe the

interactions between unlike species are inferred from those which describe the

interactions between like species. This is particularly important for such

pairs as NO which are major constituents under the conditions of interest.

One of the unexpected problems which we met in this work had to do with

the need for having enough interaction virials for the description of the interactions between a given major species and other important species. Because the

net effect of these interactions on the equilibrium constant generally appears

as a smaller difference between larger quantities, such quantities cannot be

arbitrarily neglected. We solved this problem partly by a shift to other

reference species and partly by arbitrary approximation of the unknown interactions. A study needs to be made to establish a criterion for determining

when such interactions can be neglected. The importance of this can be seen from the sheer number of possible interactions which can be needed in a calculation of this kind. Thus, in a mixture of n constituents there are n(n+l)/2 pair interactions. A mixture of 30 species therefore has 465 possible

pair interactions!! The determination of these would constitute a tremendous

job made particularly difficult in our case by the facts that data are not

available for the appropriate binary mixtures (since many of these species can- not be handled at ordinary temperatures while other species are not available

49
AEDC-TR-76-85

in more than trivial concentrations at equilibrium at such temperatures). This

is, for example, true of the atomic species. Another major problem comes from

the fact that of the 465 pairs, 435 involve interactions between unlike species.

The possibility of there being data on the 435 binary mixtures needed at a sufficient number of temperatures and for a sufficient number of relative

concentrations from which to infer potentials of interaction is extremely small.

Thus, many of the interactions need to be estimated by whatever means is

available. Clearly, any reduction in the number of pair interactions needed

and in the accuracy with which the remaining ones are needed produces a comparable direct reduction in the amount of work involved in calculations of

this kind. Thus, criteria need to be established, within the framework of our model, by means of which it can be determined when a particular pair interaction can be neglected and when a pair interaction contributes sufficiently

little so that it can be approximated in a rather cavalier fashion.

50
AEDC-TR-76-85

Table 1. Q* Values fo

nit 1

n an — u . ouu /

u

J.J

n q ^ n ^qq ^

i

"5 A J • O

i nn n ^ q 7n

—U. J5 /Uj 7 J . / Xl « z9nu _n— U . J^77^7 / / J /

•3J ft

• O

i ^n _n "ifi^in j.j/ 1 Q Xi An4U —_nU . J^AAA7 40 0 / 4A • un

Xi JsnU —n c;9« c —U ;i . JiOJl A4 . A4

i fin ^n—U . jUO" ^nscn

J

A A x . / u _n ^ftftl

8

A ft 4. O

i sn _—nU. Afifisfi 4UUJU J^ . nu

iX . Qn-7 U _n AAlftQ JR 9 • Z

?z • nnuu —U.Pi A9nfi9 4Z U O Z J . 4

Z? • xuin u « jyo

t

j

A J . D

Z9 • 9niU - n 1 79 A u . J / z jo ^ fl J • o Z.9 JU^n a n

Z9 . An4U - n "39 9 9 A U . JZ Z Z O

A 9 0 . Z

9z • ju^n n 9QAAA A / 2.60 -0.27075 6.6

2.70 -0.24456 6.8

2.80 -0.21812 7.0

2.90 -0.19147 7.2

3.00 -0.16461 7.4

3.10 -0.13758 7.6

3.20 -0. 11038 7.8

3.30 -0.08304 8.0

3.40 -0.05556 8.2

the (12,6) Hard Core Potential

n* T* n-k

— n n97Q c — U i . UZ / Z7 J ft A X1 • J^8^88 O J OO

_n— U . nnn9A U UUZ 4 fto • 8o i sm 1

1

X • J uxxx

n n97SQ q-7 . nU X1 JSSQ7fi J J 1 o nU . nsssn U J J JU q-7 « Z9 X1 . filftAA U1044

nU . ns^si U O J J X 9 A X1 . fi771 U / / X J5 n i i

i

u . X XXJsqJ 9-? fi • O 1 7':

!

1 X • / J Ji88oo

U * Z_ Z. *+ U '-T 9y • 8o X1 . 7QAfi^ / 74UJ

n 98i s

U . Z OXJ J 1 8S1A1

nU • JJOUi ^ftfii 1 1 0 Z9 . 141 A7/ JS JT nu * j^9Sft7 y _j (j / XX1 1 • J5 9 9QA74

nU « AS^98 t J -> L. U 12 .

0

9 AA909

nU s Jsinfts J.UU J X19z_ • J5 A9 « JU? S89S9J7

n sfiftsn Ji ^i • nU 9 7^fi79 / JU; .7 nu . fi9fi9Q u z o Z

" X J . J Z9 . 88A9fi OOHZU

nu • fiftAi UU41;7 14i a . nU J . U J X / O nU • 7A91A / HZXH 141 A . JS J • X 7 -7 J 4

n ftnm q iJ s1 . nu J? • J19fiQA Z u 7 H

nuft ^8^9

o

jo jz 1J1 ^ R . J J • 4 / 4 J O 0.91651 16.0 3.62225

0.97476 16.5 3.76994

1.03307 17.0 3.91767

1.09143 17.5 4.06542

1.14984 18.0 4.21319

1.20829 18.5 4.36098

1.20829 19.0 4.50879

1.32532 19.5 4.65661

20.0 4.80446

51
H oo o o O

r

O CMO o o iO

O CM 00 o o o -d- CO o o lo <t o o o 00 o o f—

O o

lo oo

<r rH rH

CO CM <r CM

LO <r oo o o o rH

rH o CO o o o LO

10

O CO o o o v£> CO

ON

O •o iO o o u-i rH o lo

CM

O o O

H

H sO rH o

<u •HV

<Do.

(1) >

•HU

crj o

J-l Q_ C

CU Q.

If) 0) 00O

r-l

0) aiU

fo

CO

.1) •H

i-i •H

COc

CD Q

rH

CTJ

r-l CU>

CU COo

4-1

CDcu3

crj >- CO

>-

CO

oo

LO

ooCO

CO r-~ LO o o CO CO o LO

00 CT\ rH o o O o CM o CM o o sD LO o rH LO o sD LO o o rH O o r-~ O o o o rH rH CO 00 rH O CM

CO CO

CO CM o LO O o <r CO o CM LO o \D O O CM o o rH CM

LO CO o O en rH o OS o o o <r CM o vO vO o O O rH rH C7\ vO rH o LO

rH rH

LO CI rH o o CO

O o CO rH

vO uO rH o o CM CO o O o

On CO rH o o rH LO o oo o

<r oo oo o o -3" CO o CM o o o rH rH rH rH rH o rH

CTN 00 CO o o LO LO o LO

\D <£> o o o LO CM CO CM rH o o LO <r o as CM OS o o o o o CO o o o o rH rH rH rH H o rH

ro CM o O v£> LO o vO

CO CM rH o O rH o O 00 ON

CTn

O o O o O CO as

as CT\ o o o o O CT\ as O O o rH rH rH rH i—1 O o

CO CM rH o o rH o O as as

l-> as o o o o O CO a-, CTN OS as o o o o o a\ as

<3s OA as o o o o o as as O o o rH rH rH rH rH o o o2

CM CM O O

E2

O

O

a 5 2
CO ,

,__)

Q o H rH m ^0 ON

i-t

H CN CN <r Ln

PQo rH i—i

rH rH rH

H CN CO

~? CO

~~j rTi CJJ CO

LM CN * <r —J O ?l rH rH CO CO CN

nr1 | i—i CN

_

QOO

H CN

_J Lm o *T Ul 00 00 CO O CO

*^ C7< CO

!-T\ 55

2j

CO

u i | en

rH , 1 rH 1 . { rH OJ

4£

H CO 0

<4H ro ~?

m

—.

_^ J~| _^ IS

fc-l to oo CO oo

PS

M 03 nttU CN CN

CQ (U qj rH

QQ o Q

O

Q 1 1 5 ' t 1 1 1 1 rH 1 ) 0

It) rrt CO

!ii 3 vD

O ca

4-1 x> ts o cfl GJ

a) o 03 o rH CN cr^ "T

4J o O

O CO CO

•H CO 'vD rH

i—1

4J o r—1

CO i—!

-H CO rri TO o o rH CN CO .O J2 Cm

•H 4-1 iH p—1

rH rH rH <H CN

CO

U

CO 00 CtJ OJ c u •H o o CN <T o r> CO vo CN CO CO £ 3 •H CO CO 0^ rH IT) o 0 iH o —i

LT> lO UO CO o u c C/j o O rH CN| CO CO

ti CD * <u M PS rH rH 1—i

rH i—1 1—1

CN LO

4=

4J

CO

O

<4-l c

w>NH

M o i—1 o Ul

O 0 *2

v«0 LM ^2

CN CO •H c o co vO 00 pH CN CO

CO 4J o o o i—ICN o co

C cd •H o rH CN rr) lO o> 0 E *J m c

CO •H

u CM iH rH rH rH rH rH rH CM

•H

X

a o « M ptj H \© vO \D rH

m rH

m

ex. iH rH CN <r rH $ w i—1

CMO o> c^ <0^ CO 0 03 O O o C^ a> CT\ CTn o 0) O o o as a» o> M TJ 3 t—l

iH rH rH o o o o o • o

CO

iH

O

O c o r~ r» o o O o o CO CO O o

Q.

O o o -3- rH i-H o

iH o o 00 i—l CO o o

i-H o LO in ON

m o 1—1 '—1 m orH

O

Q.

Q.

O o o O o O o

60

O o o CM <r CO o O

iH o CM CM CN CN cn

53
AEDC-TR-76-85

C/3

CE

no on o O r^. \o t—I co r-» in co co in OONrH-NTCMONONONOrHininr^

OOOOi—IN-*0ONHr--Nifi

CMcMrHcocMrONOOcor^ONr-~NO

<f «) H -J M O i-1 ifl in* H n o\

i—I <} i-IOCOH)N>C(M(N-JvOP1 i^cMONOO-d-cooinNOcMOcocM NNrlHNni/lCOvOOmcMH

HHHHHHHt-(rMn(<l<t^

COc

oH

4J

CJ5

1^

to

COc

oH

*J

cfl

o

FX

o

o

o

CU

ooNOsomoLri'-iOiHOtH-*

OOOOi—INPIinoOOHM-J

t—If—It—Ir-Ht—ti—ii—I i—I H CNJ (N N N

<7\COCONOvWvDNrH<fCOC<1r^ mnmcotOHNin-jHco^n CMCOCMCMNOmcOCMrHrHOOO NOoocMrHOOOOOOOOO ooooooooooooo

oo

o

mi—i no o co co rH t—i on o r~- -nt on vtCOWOfOiONtCl<NinNrOts H-JOl^<fO^(NvCNHa\N r-~CNCTNr^coOr—!OOOcsi<fCT\ N N rl rI CM CO <T NO ON O r-I CM CO

i—I H rIt-It—li—lr-It—I H N (N (M (M

r^-r^voc0ONin\Dt—it-loovDin-*

o\ on m n r-^^nnocovo^

\Ot—iNcioNNvom-j-jfnntn CMCMrHtHOOOOOOOOO

o co c ro CN <f ~NT O r-» ro co 00 CM On r-» ON H rH CO rH CM o

no uo r- on CM oo CC m o co co rH nO o m CM co ro o CO ON ro o

00 30 o oo <r on i—

l

LO O CN o ON ON rH nO CM in CM < o o no NO o O eH -nT in in m CO CM ro vD ON -nt rO rH o CO <f CO CO CO ON

re on o o O o rH CM <t co CM oo nO CO rH 1—

1

i—

I

H CM ro <r co nO O m CM CM HHHHHHH(N|nn<tm

u 1H O -X) o rH r-~ o ON 'ON r-» m in o •H NO in n nD CO CM in m 00 ro oo

ItH § > CO co NO nO -NT 00 rH CM oo CM CM co o o rH <r nt cm nO o rH -NT oo CO U H o o o o o rH CM co m on o rH CM

•H o C

<J US o

CO

CM rH rH rH rH rH rH rH rH rH CM CM CM CM n o

o •H 1 UH rJ nO O O CO NO ON rH CO rH CM -» -NT m

03 H rH O .—1 o CM ON r-. -d- rH o ON 00 r~ > rH r~ CM rH C> ON ON on ON ON oo co co H rrj o o O c O ON ON ON ON ON ON ON ON

at de o o o o O ON ON ON ON ON ON ON ON > rH rH rH rH rH rH O o o o o o o o

en

o

o

o

NO

I—It—IrHrHrHi—ICMCOCO*d"in

cooocMr^cocMcor--coooO'<rin COrHinNOOOOOrHOOOOOCOOO ONONrH<t00ONONCOr^rHCOONNO ON*NtC10MDCOCOONOH>fO> H H H i—IrHCMCOinoOOrHCMCO

H H H i—I H H iIiHiHCMCMCMCM

HinNrNvB-JrlinHiCNNri

cm in r*^ t—I on i—I in t—l no cm on no r~-*NOoor~NOooo<rcMrHONOO ONNOCMCOin<rcOCOCMCMCMtHrH iHrHrHOOOOOOOOOO

CO CO u c O CU Q

o r-» <r CO oo CM nD oo nD in rH m NO co ro rH O ro O rH ON CM NO co 00

Cfl C|H -NT m <r NO <r CM a <T -T m o o NO CO ON <r CM rH CO CO ON m

fa o rH -cr «c* nD co oc rH ON o NO o m rH nO rH rH o in ON ro CO rH rH 00 O O rH -NT in nO NO 00 oo oo O NO <r <T rH CO ON 'CC CO NO CO rxi CM CM 00

:>N CO o o O o rH CM -NT co 'CO CO rH rH H rH O o H CM <r 00 rH nO

.lit X rH rH rH rH rH rH rH rH CM cn <r UN NO rH rH H rH rH rH rH rH CM CO CO in

ib: u CO ON CM cn CM CN CM CO m ON NO NO ON NO CM CM rH 00 n ro CO sr in

co •H <r in CN CM rH ON r- CM oo rH ON co NO NO NO CO CM ON CO O icn O ON n NO

ca > rH <r <r in CO NO ON O <n <r r» NO CM W m rH NO O CO -NT o CM <r rH o in NO

QJ o o rH -NT <T CM in r-N o rH ro CO -NT rH CO ON NO -NT r— r~ ON O CNJ in O

u o -a o o o O rH CM co in ON o rH CM -NT O rH rH o o rH CN ro m 00 O rH CM -NT o c B

1 o CM rH rH rH rH rH rH rH rH rH CM CM CM CM O rH rH rH rH rH rH rH rH rH CM CM CM CM

Co CM in

II II ro ON m ON o -NT CO ON CO m rH 00 -T o rH co ro -NT CM O O ON r-- r-» ON -NT H O ON ON 'CO 'CO NO in <r <T -NT ro H rH o o O in rH nO co in co H <T CO

-nt rH o ON ON ON CT- ON ON ON ON ON On ON ON -NT r> CM in NO rH nO CM ON co CO r-~ NO

n) Q ON ON 'ON ON ON ON ON ON 'ON ON ON ON -T o r> nT CJ CN rH rH o o o O o

QJ CU o ON ON ON ON ON ON ON ON ON ON ON ON rH rH o o o O O O o o o O o

rH X) XI rH rH o O o o o O O O o c O o rH rH rH rH rH rH rH rH rH rH rH rH rH H

OcMOcoooNONr^r-~NooOrHO ONOocMOcoooomNcfcMino O t—I O NO O "-CT rH rH ON ON CO CM O

rH CO O rH OOOrHCOOr^-NTrHO oininoNcoooNONO

rHrHCMCONOr-^r-~COO

OCMOcooONONr^r-^NOCOrHO ONOocMOooooom-NTCMino OtHONOON-frHtHCTNONCOcMO

rH CO O rH O00rH00Or-^-<TrHO OininoN<nooNONO iHrHcMc0NOr-~r^00O

ooo oinoinocM~rfNOooooaNONO oooooooomomo

OOrHrHCMCMCMCMCMCMCMCMCO

omomocM^tNococooNONO oooooooomomo

OOrHrHCMCMCMCMCMCMCMCMCO

54
AEDC-TR-76-85

u

c

3

'—

*3 — o *0 — *3" CM cn CC 5

- a X CN in - 00

CO - - - 3P ~* £ X = — — — — —

— — - CM CM ~~

^a-

~3" — *3j CO - - >

~ —

-~ -~

CN

-3- CM X -. CM

EN 5 X -- CM CM

-~ CM — — —

C M C i ^ J

m * -j <r c ^ x

X X X CM -O 3» c — — <• x

— O O X vD £> <0 X

cn —: —i — — — —

r-» O C; i/1 —i co

-~ C w — -3 ^

x x x <r o m o

C 3 -~ CN O cn c- c- W O 3* X- W W W V

— X c— r- cn _~ N N I— N 1^ ^ CO < c- m cm o cm O cm O co W <r O CC> X -3 r~ c~ CM

< -3- -3" a CN -T CM in m CM vO

T. M3 CM m — »

X X CM CM v£> X X vO

sf CN CM CO co r» X

CO CM v£J X CT- X in -3- X

CM CM CM CM co CO •3- 3

- :— < 50

> -—

X — X X x <

N < X C

<T C r" ? x c cn o s- X X

X CN — c- 3~ >C X

la i N X N i."M !" o co Is- m m m \C r-^

r">c^<rxrsx<rx<rxcrj5 ^^^OXNX--Jr-ri it, - s> ifl n c it. s jv « n m

-"^"Orsr-~rsir-~mc^xX'-~NC

— X CM X

- -3- (-» — <r

\£) — CO <• CM CO CO co

>» — cn X m X r~ ~ -

4-1 CM X CO i

- -3- X

CO cn CM CS m CO ^3" c

- — CM CM

1-iCO

CO CM

0) CM CN h > CM CM CM

r* CM r~ -. CN CM = CM CM <? vO CM u CM CM CM CM CM

1 <* CM CM X CM CM

(!) CM — r

CM — —

" ? N C CM N

<r X cn o cn O

X co <• vc X *7 N ^ CM X un -3- ~3" LO

o — x < — c —

cr c~ cn X O r~ o — <r O -3- co o cn ^ X C CO -3" •C 3^ O -3- O X

•3- -3

C-XC:XX<-C^—

!

co X co m <r o —

O x<-crr-<-XOCN

C xr~~-n<3-<--3-mr~-

vot^>*o\HHro>fNins>^<r 2->J5'CiflXi/liCSiNCO x<-r^Hssirir-icOin-Jf-Tc?fO XI~-i.1-INNNi-IHHHWH

CN CO CM CM X CO •3- CM m CM T CO CN X CO I

— — CO COON prj ~ 5 — CN J— X CM r~ X

z

Q.Ol in m m ~~ CM X X r in CM •3" DC X 5

CN CM CM CM CM CN CN CM CO CN CM CM CN CM CM CM CN CO

55
AEDC TR-76-85 REFERENCES

1. Hilsenrath, J., and Klein, M. AEDC-TDR-63-161. Arnold Eng. and Develop.

Center, Tullahoma, Tennessee. AD 339389. (1963).

Hilsenrath, J., and Klein, M. AEDC-TDR-63-162 . Arnold Eng. and Develop.

Center, Tullahoma, Tennessee. AD 432210. (1964).

Hilsenrath, J., and Klein, M. AEDC-TR-65-58. Arnold Eng. and Develop.

Center, Tullahoma, Tennessee. AD 612301. (1965).

Hilsenrath, J., and Klein, M. AEDC-TR-66-65 . Arnold Eng. and Develop.

Center, Tullahoma, Tennessee. AD 630461. (1966).

Hilsenrath, J., Messina, C.G. , and Klein, M. AEDC-TR-66-248. Arnold Eng.

and Develop. Center, Tullahoma, Tennessee. AD 644081. (1966).

2. Grose, W.L., and Nealy, J.E., AIAA Jour. 1164 (1970).

3. Klein, Max, Chapter 7 of "Physical Chemistry" - Vol. I - Thermodynamics , edited by H. Eyring, W. Jost and D. Henderson, Academic Press, New York,

N.Y. (1971).

4. Klein, Max and Green, M.S., J. Chem. Phys. 39, 1367 (1963).

5. Klein, M. , AEDC-TR-67-67 . Arnold Engineering and Development Center,

Tullahoma, Tennessee. AD 649463 (1967).

6. Haar, L. and Shenker, S.H., J. Chem. Phys. _55, 4951 (1971).

7. Woolley, H.W. , Report AFSWC-TDR-62-21 , Air Force Special Weapons Command,

Kirtland Air Force Base, New Mexico (March 1962). (AD 276922).

8. Grabau, M. and Brahinsky, H.S., AEDC-TR-66-247 , Arnold Eng. and Develop.

Center, Tullahoma, Tennessee. AD 646172 (1967).

9. Haar, L. and Sengers, J.M.H.L., J. Chem. Phys. 46, 4181 (1967).

10. Weeks, J.D., Chandler, D. and Anderson, H.C., J. Chem. Phys. 54, 5237 (1971).

11. Wertheim, M.S., Phys. Rev. Letters 10, 321 (1963); J. Math. Phys. _5» 6 ^3

(1964); Thiele, E. , J. Chem. Phys. 39, 474 (1963).

12. Klein, Max, Phys. Fl. 7_, 391 (1964); Ree, F.H. and Hoover, W.G., J. Chem.

Phys. 46, 4181 (1967).

13. Mayer, J.E. and Mayer, M. , "Statistical Mechanics," John Wiley and Sons,

New York (1940)

.

14. Hirsctifelder, J.O., Curtiss, C.F. and Bird, R.B., "Molecular Theory of Gases

and Liquids," John Wiley and Sons, New York, N.Y. (1954).

15. Alder, B.J. and Wainwright, T.F., J. Chem. Phys. 33., 1439 (1960).

16. Wood, W.W. , Parker, F.R. and Jacobson, J.D., II Nuovo Cimento, Suppl., IX

50
AEDC-TR-76-85

Series X, 133 (1958).

17. see reference 7 page 154.

18. Rowlinson, J.R. in "Physics of Simple Liquids" edited by Temperley, H.N.V.,

Rowlinson, J.S. and Rushbrooke, G.S., Wiley and Sons, New York, N.Y. (1968).

19. Zwanzig, R.W. , J. Chem. Phys. 22, 1420 (1954).

20. Barker, J. A. and Henderson, D. , J. Chem. Phys. 47, 4714 (1967); Toxvaerd, S.,

J. Chem. Phys. 53, 2389 (1970).

21. Keller, J.B. and Zumino, B., J. Chem. Phys. 30, 135 (1959).

22. Mansoori, G.A., Provine, J. A. and Canfield, F.B., J. Chem. Phys. 51, 5295

(1969). —

23. Toxvaerd, S. , private communication.

24. Henderson, D. and Baker, J. A., in "Physical Chemistry" VIIIA, edited by

Eyring, H. , Henderson, D. and Jost, W. , Academic Press, New York. N.Y.,

(1971), Chapter 6.

25. Sengers, J.M.H.L., Klein, M. and Gallagher, J.S., "Pressure-Volume-Temperatur. Relationships for Gases," Amer. Inst, of Physics Handbook, Third Edition,

McGraw-Hill, New York, N.Y. (1972), p. 4-204.

26. Barker, J. A., Henderson, D. and Smith, W.R., Phys. Rev. Lett. 21, 134 (1968).

27. Reference 5, p. 22.

28. Hill, T.L., "Statistical Mechanics," (McGraw-Hill Book Co., Inc., New York,

N.Y. 1960).

29. Henderson, D. and Leonard, P.J., Proc. Nat. Acad. Sci. USA 6_8, 632 (1971).

30. Henderson, D. and Leonard, P.J., Proc. Nat. Acad. Sci. USA, 68, 2354 (1971).

31. Hooke, R. and Jeeves, T.A. , Comm. A. CM. _8, 212 (1961).

32. Partington, J.R. "An Advanced Treatise on Physical Chemistry," Longmans,

Green and Co., New York, N.Y. (1949), p. 180.

33. Levelt, J.M.H., Physica _26, 361 (1960); van Itterbeek, A., Verbeke, 0. and

Staes, K. , Physica 29, 742 (1963); Robertson, S.L., Babb, S.E., Jr!, and

Scot, G.S., J. Chem. Phys. 50, 2 1 & 0 (1969).

34. Verlet, L., Phys. Rev. 15_9, 98 (1967); Levesque, D. and Verlet, L., Phys.

Rev. 182 , 307 (1969).

35. Michels, A., Lunbeck, R.J. and Wolkers, G.J., Physica jJ, 801 (1951);

Robertson, S.L. and Babb, S.E., Jr., J. Chem. Phys. 50, 4560 (1969).

57

AEDC-TR-76-85

59
AEDC-TR-76-85

Fig, 2, Effective potentials vs molecular separation in arbitrary units for several

VALUES OF M AND FOR T* = 2,5, It IS SEEN THAT THOUGH THE SHAPE OF THE "BOWL"

IS QUITE SENSITIVE TO THE CHARACTERISTIC PARAMETER M, THE SPHERE DIAMETER AT

THIS TEMPERATURE IS ONLY WEAKLY SENSITIVE TO IT,

60
AEDC-TR-76-85

O

61
AEDC-TR-76-85

.2 .4 .6 .8

p*

Fig, 4. Compressibility factor for argon vs reduced density pGV) = p* for

T = 119. SIC. Curves L this work; Curves 2 and l\, ZBH theory from

references 22 and 23, respectively; Curves 3, ZBH theory based on Monte

(Precalculations of perturbation terms Curves 5, FVT experimental

datajj

. The "computer experiments"^ are designated by circles, Monte

Carlo; and squares, molecular dynamics.

62
AEDC-TR-76-85
AEDC-TR-76-85
AEDC-TR-76-85

5.0 ARGON

I i ,

4.0

Z 7/3

3.0

II /o

2.0

1.0

i i 1 1 I I I I L

0 .2 .4 .6 .8

P*

Fig, 7, Compressibility factor vs density for

ARGON FOR T = 328 , 25Kj SEE CAPTION/ FlG . 4,

65
AEDC-TR-76-85
AEDC-TR-76-85

5.0

ARGON

T= 673K

1.0

.4 .6 .8

Fig, 9, Compressibility factor vs density for

argon for T = 673K; see caption, Fig, 4,

67
AEDC-TR-76-85

5.0

4.0

3.0

2.0

0

0

ARGON

T=239.6K

1

200 400 600

p (Amagat)

800

Fig, 10. Sensitivity of equation of state for 3sH theory to choice of potential

function. Compressibility factor is plotted vs density in amasat. For

Curve 1, the reference potential is the Lennard-Jones (18,6); for Curve 2,

the Lennard-Jones (12,0 , The numerical results are from calculations

FURNISHED BY TOXVAERIT USING POTENTIAL PARAMETERS e/K = 119.8, ° = 3.VB

FOR THE LATTEF AND EQUATION (12) FOR THE FORMER.

68
AEDC-TR-76-85

0 200 400 600 800

p (Amagat)

Fig. 11, Sensitivity of present theory to choice of potential functions. Potential

functions are as described in caption for flg. 10. the small difference

between Curves 1 and 2 is mostly due to the slight differences in second

VI RIALS PRODUCED BY THE TWO POTENTIALS IN THE VICINITY OF THE ISOTHERM.

69
AEDC-TR-76-85

z

yo(Amagat)

Fig, 12, Compressibility factor vs density in

amagat for nitrogen t * 123, 15k, curves 1/

present theorxi curves 2, experimental pvt

measurements.

70
AEDC-TR-76-85

z

0 200 400 600

yo(Amagat)

Fig. 13, Compressibility factor vs density in

amagat for nitrogen t = 173,1510 see caption/

Fig. 12.

71
AEDC-TR-76-85

5.0

NITROGEN

T=273.I5K

1.0

0 200 400 600

/o(Amagat)

Fie. 14, Compressibility factor vs density in

amagat for nitrogen t = 273. 1510 see caption,

Fig. 12.

800

72
AEDC-TR-76-85

Fig, 15, Compressibility factor vs density in ^magat for nitrogen t = 675k; see caption,

Z IG, 12,

73
AEDC-TR-76-85

1.0

N.

o

I—u

<cc

10

IDEAL —— REAL

T= 3000 K

I

1.5 1.75 2.0 2.25 2.5

Log (P/P0 )

2.75 3.0

Fig, 16, Concentration versus density for

several important species for t = 3000k,

74
AEDC-TR-76-85

1.5 1.75 2.0 2.25 2.5 2.75 3.0

Log [P/P0)

Fig, 17, Concentration versus density for

several important species for t = 900ok,

75