Rational Bureau of Standard

Library, E- in - BldgWAR 2 4 1971

NBS TECHNICAL NOTE 396

Data Analysis

for Isoperibol Laser Caiorimetry
NBS TECHNICAL PUBLICATIONS

PERIODICALS NONPERIODICALS

JOURNAL OF RESEARCH reports National

Bureau of Standards research and development in physics, mathematics, chemistry, and engineering.

Comprehensive scientific papers give complete details of the work, including laboratory data, experimental

procedures, and theoretical and mathematical analyses. Illustrated with photographs, drawings, and

charts.

Published in three sections, available separately: • Physics and Chemistry

Papers of interest primarily to scientists working in these fields. This section covers a broad range of

physical and chemical research, with major emphasis on standards of physical measurement, fundamental

constants, and properties of matter. Issued six times

a year. Annual subscription: Domestic, $9.50; for- eign, $11.75*. • Mathematical Sciences

Studies and compilations designed mainly for the mathematician and theoretical physicist. Topics in mathematical statistics, theory of experiment design,

numerical analysis, theoretical physics and chemistry, logical design and programming of computers

and computer systems. Short numerical tables. Issued quarterly. Annual subscription: Domestic,

$5.00; foreign, $6.25*. • Engineering and Instrumentation

Reporting results of interest chiefly to the engineer

and the applied scientist. This section includes many

of the new developments in instrumentation resulting

from the Bureau's work in physical measurement,

data processing, and development of test methods.

It will also cover some of the work in acoustics,

applied mechanics, building research, and cryogenic

engineering. Issued quarterly. Annual subscription:

Domestic, $5.00; foreign, $6.25*.

TECHNICAL NEWS BULLETIN

The best single source of information concerning the Bureau's research, developmental, cooperative and

publication activities, this monthly publication is designed for the industry-oriented individual whose

daily work involves intimate contact with science and

technology

for engineers, chemists, physicists, re- search managers, product-development managers, and

company executives. Annual subscription: Domestic,

$3.00; foreign, $4.00*.

• Difference in price is due to extra cost of foreign mailing.

Applied Mathematics Series. Mathematical tables, manuals, and studies. Building Science Series. Research results, test methods, and performance criteria of building materials, components, systems, and structures. Handbooks. Recommended codes of engineering

and industrial practice (including safety codes) de- veloped in cooperation with interested industries,

professional organizations, and regulatory bodies. Special Publications. Proceedings of NBS conferences, bibliographies, annual reports, wall charts,

pamphlets, etc. Monographs. Major contributions to the technical

literature on various subjects related to the Bureau's

scientific and technical activities. National Standard Reference Data Series.

NSRDS provides quantitive data on the physical

and chemical properties of materials, compiled from

the world's literature and critically evaluated.

Product Standards. Provide requirements for sizes, types, quality and methods for testing various industrial products. These standards are developed cooperatively with interested Government and industry groups

and provide the basis for common understanding of product characteristics for both buyers and sellers. Their use is voluntary.

Technical Notes. This series consists of communications and reports (covering both other agency and

NBS-sponsored work) of limited or transitory interest. Federal Information Processing Standards Pub

lications. This series is the official publication within

the Federal Government for information on standards

adopted and promulgated under the Public Law

89-306, and Bureau of the Budget Circular A-86

entitled, Standardization of Data Elements and Codes

in Data Systems.

CLEARINGHOUSE

The Clearinghouse for Federal Scientific and

Technical Information, operated by NBS, supplies

unclassified information related to Government-gen

erated science and technology in defense, space, atomic energy, and other national programs. For

further information on Clearinghouse services, write: Clearinghouse

U.S. Department of Commerce

Springfield, Virginia 22151

Order NBS publications from: Superintendent of Documents

Government Printing Office

Washington, D.C. 20402
UNITED STATES DEPARTMENT OF COMMERCE

Maurice H. Stans, Secretary

NATIONAL BUREAU OF STANDARDS • Lewis M. Branscomb, Director

NBS TECHNICAL NOTE 396

ISSUED FEBRUARY 1971

Nat. Bur. Stand. (U.S.), Tech. Note 396, 34 pages (Feb. 197i; CODEN: NBTNA

Data Analysis for Isoperibol Laser Calorimetry

E. D. West

Quantum Electronics Division

Institute for Basic Standards

National Bureau of Standards

Boulder, Colorado 80302

V

°»e au of

NBS Technical Notes are designed to supplement the Bureau's regular publications program. They provide

a means for making available scientific data that are

of transient or limited interest. Technical Notes may

be listed or referred to in the open literature.

This work was partly supported by the Night Vision Laboratory, U.S. Army

Electronics Command, monitored by S. E. Smathers, and by the Advanced

Research Project Agency of the Department of Defense, monitored by Lt. Col. John M. MacCallum under ARPA Order No. 891.

For sale by the Superintendent of Documents, U.S. Government Printing Office, Washington, D.C. 20402

(Order by SD Catalog No. C 13.46:396), Price 40 cents.

Data Analysis for Isoperibol Laser Calorimetry

E. D. West

Isoperibol calorimeters (those operating in a constant-temperature environment) are used to measure

the power and energy in laser beams relative to electrical

standards. The derivation of the basic formula is reviewed.

Two methods are presented for analyzing the data taken at

equal time intervals: (1) An approximate manual method

with criteria for avoiding significant errors of approximation

and (2) A least squares method for use with automatic digital

computers.

Key Words: Calorimetry; laser; laser calorimetry;

laser energy; laser power.

Laser calorimetry affords a method of comparing the power and

energy in laser beams directly to electrical power and energy. Laser

calorimetry thus provides an accurate and convenient means of calibrating devices for laser power and energy measurements in terms of

the absolute joule or watt. The electrical standards required -- standard resistors and standard cells -- are far more accurate than

is now required for laser power and energy measurement; they are

easy to use and available in most laboratories or obtainable at reasonable cost.

Calorimeters used at the National Bureau of Standards for these

1, 2 . . measurements operate m a constant-temperature environment and

are termed isoperibol calorimeters. The simplified theory of opera3, 4, 5, 6

tion of these calorimeters is given in various books. The same
measurement method can be derived from a slightly less simplistic

7 model and more rigorously from detailed considerations of the heat

flow problem in the calorimeter and its relationship to the first law of

o thermodynamics. The detailed theory shows that the natural response

of the calorimeter to a sudden energy input is an infinite series of

exponential terms decaying with time. This theory shows how heat

exchange due to these higher order terms is taken into account by

proper attention to calorimeter design and operation.

A temperature-time curve for an isoperibol calorimeter is shown

in Fig. 1. This curve is typical for calorimeters operating at a temperature above the average temperature of the surroundings. In the

first part of the experiment, the temperature is decreasing as a single

exponential function of time. Then an electrical or laser input raises

the temperature. Some time after this input the temperature again

decreases as a single exponential function of time. Although this

curve is certainly the usual one for laser calorimetry, the data analysis in this Note applies equally well to any method of operating the

calorimeter -- including processes which absorb heat.

1. THE WORKING EQUATION

The first law of thermodynamics states that the work W done on

a calorimeter by a laser beam or in an electrical calibration is given

by the increase in the energy AU stored in the calorimeter minus the

heat Q transferred from the constant-temperature surroundings:

W = AU - Q. (1)

The experimental problem is to determine the heat exchange and the

stored energy. Of course, it must be shown that one is applying the
o

pQ• H

rH(LPh

CO • H

09MC

N >

(U

m en O rO

en ° —< CD

i> +->

a) 5c

rH a;

<D • H g

4->U

fH +J (1) uo

Oh

—1 X —

a

(U 'J c

>. >

•i-h

rH

n l+H

rt

u

C

—

rH •*>> ILI ai +j to

tti rt 3 P M tti

rHC

<+H

rt XI

tu

rH3

rt- o

Oh cd

(D pQ • rHa

CO

ex) cxJ <p

Oh ^ o

rH "£ ° rrt O r-1

T3

r-H (Tj C! Ph U (Tj

3UfUVd3dlAI31
first law properly for each calorimeter. For example, tests must be

made to show that various combinations of AU and Q give the same

value for the work W and that AU = Q during rating periods (W = during rating periods; that is, there is no electrical heating or laser

input. ) Such tests are made to demonstrate the adequacy of a particular

design and are essential to the estimation of the accuracy of a particular calorimeter. Discussions of design and testing calorimeters are

scattered through the literature, but some of these have been collec9-H

ted. Restrictions which the theory imposes on the experimental

measurements are derived from the detailed consideration of the heat

o flow problems. This paper will not treat the problems of design and

testing, but will instead be concerned with the analysis of data from

proven isoperibol calorimeters.

To use Eq. (1) for measurement, the quantities AU and Q must

be determinable in terms of observed quantities. In the linear theory

o of isoperibol calorimeters, the change in internal energy of the calorimeter is shown to be proportional to the change in temperature

AU = C(T - Tj), (2)

where T - T represents the change in the observed temperature.

(T may be looked on as a reference temperature at which U has the

value of U . ) The constant of proportionality C will be different for

each particular calorimeter. This constant is a sum of heat capacities

Q

weighted according to the temperature gradient in the calorimeter.

The procedure described below takes this weighting into account.

The determination of the heat transferred from the environment

is such a complicated problem that the development given here is

more illustrative than rigorous. We begin with the observation that
sometime after an electrical or laser input to the calorimeter, the

observed temperature T changes according to the relation

dT/dt = - e (T - TJ (3)

where T is the observed temperature for which dT/dt = and e is

oo termed the cooling constant. Equation (3) is a representation of

Newton's law of cooling. Multiplication by the constant C gives the

relation

C dT/dt = - (Ce) (T - TJ. (4)

It is apparent from Eq. (Z) that the left-hand side of Eq. (4) is the rate

of change of the internal energy. Since the electrical or laser work

term must be zero when Eq. (3) holds, the right-hand side of Eq. (4)

must represent heat transfer as required by Eq. (1). For a calorimeter behaving as in Fig. 1, the internal energy is decreasing and

heat is flowing from the calorimeter because T > T . The constant ° OO Ce has the dimensions of a heat transfer coefficient and is frequently

so called; in this instance, T is merely an observed temperature and

not the temperature of the calorimeter surface from which heat is lost.

When the existence of a temperature gradient in the calorimeter is

o acknowledged, Ce loses this physical meaning. It has been shown

that e is the smallest eigenvalue in the heat flow problem.

When Eq. (3) holds, the calorimeter is said to be in a rating

period. Another criterion for the rating period is obtained by integration of Eq. (3) between limits and placing it in the exponential form

(T - T ) e~

£lt ~ ^ I oo T - T = (T - T \ *-e (t ~ t l) (5)

oo
where the constant of integration is obtained from knowledge of the

temperature T at the time t . Eq. (5) is reached in practical calorimeters because the other eigenvalues are appreciably greater than e and the corresponding exponentials usually become negligible in a reasonable time. The heat exchange due to these higher order terms

can be accounted for if the calorimeter is properly calibrated electrically and care is taken that temperature observations extend from an

initial rating period, in which Eq. (3) [or (5)] holds through the "main"

period, in which the laser or electrical energy is added, into a final

rating period in which Eq. (3) [or (5) with new integration constants T

and t ] again holds. If the data are analyzed by methods outlined F

below, the measurement does not require that dT/dt is small. These

methods therefore speed up the measurement by eliminating the

"infinite" time required for T to equal T in Eq. (5). A principal

source of systematic error is that the heat transfer due to the higher

order terms depends on the location of the source. Since the laser

energy obviously cannot be absorbed exactly where the electrical

energy is generated, the calorimeter must be designed so that the two

sources appear to be equivalent as far as accounting for the heat transfer by the use of the observed temperature. This problem is discussed

in reference 8. With these restrictions the small heat transfer q due

to the transients is shown to be proportional to the electrical work on

the calorimeter heater or to the electromagnetic work of the laser beam

q = kW. (6)

The balance of the heat exchange, and probably the greater part of it,

is obtained by integrating the right-hand side of Eq. (4). The quantity

Q in Eq. (1) is then the sum
Q = - C£ /(T - T ) dt + kW. (7)

The heat quantity Q is usually negative in laser work because T is

usually greater than T ; that is, the calorimeter is usually heated

above its equilibrium temperature. Substituting from Eqs. (2) and (7)

into Eq. (1) and rearranging, we obtain the working equation for isoperibol calorimetry

t W = C

[TF -TI+e/ (T-Tjdt], 1 - k

L F I c J

t

v oo' "M' (8)

where T and T are the observed temperatures at times t and t respectively. The quantity in brackets in Eq. (8) is usually called the

corrected temperature rise. It is a constant when Eq. (3) holds, as is

apparent from integrating Eq. (3) to obtain

TF- TF= e J

< T F- T»

)dt '

*F

which shows that the decrease in T is just equal to the increase in the

integral term. One must be on guard at this point that the data really

fit Eq. (3), that is, that deviations from the smooth curve are really

random. To check on this point, the computer program described

below prints out the deviations, which can be examined to see if the

positive and negative deviations are grouped or in patterns. The calculations can also be done manually, but manual calculations require

considerable time.

Temperature is not observed directly, but is calculated from

some other observed quantity. The temperature is frequently a linear
function of an observed quantity x, such as a thermopile voltage. In

such a case, that is, when

T = a + bx, (9)

it is not necessary to convert to temperature. This linear relationship is frequently an adequate representation of the temperature and

is used in our work to relate the temperature to the output of a difference thermocouple or thermopile when the temperature difference is

small. Substitution for T from Eq. (9) into Eq. (8) will show that the

quantity a cancels and the quantity b can be factored and placed outside

the bracket. Putting -—— = E, which is customarily called the energy

equivalent of the calorimeter, the working form of Eq. (6) becomes

W = E [x F- Xj + , J* ( X . XJ dt]

(10)

The quantities x, x , x , t, t , t , may be observed directly; the

F I F I quantities c and x are calculated from observations of x and t during

oo ° rating periods using Eq. (5). Although x may be observed directly,

it is commonly not done.

Equation (10) is identical to Eq. (8) so that x is merely a new

temperature scale with a different zero and a different sized degree.

The quantity E in Eq. (8) or (10) is to be determined from a known

work quantity, an electrical calibration, for example. It may accordingly be regarded as the calibration factor. Note that it is not a sum

7 8

of the heat capacities of the "parts" of the calorimeter.

The quantities E, x , and e are constants which are properties

of the calorimetric system. It is prudent to maintain control charts

on each of these quantities. The calibration factor E should remain
constant within the required accuracy, although a predictable drift

with time is almost equally useful. The quantity x is affected somewhat by changes in the room because the calorimeter exchanges some

heat with the room through the opening for the laser beam. This

quantity might therefore be expected to show small variations with

time, perhaps related to room temperature. The cooling constant e

is markedly affected by the gas pressure in the calorimeter and affords

an operating check on the vacuum system if the calorimeter is evacuated. The calibration factor E should not be appreciably affected if x and e are different for different experiments, but it is incumbent on

the user to show that this is true.

The direct computation of Eq. (1 0) can be carried out with a manually-operated desk calculator. Historically, much attention was

given to simplifying the problem of evaluating the integral, but, when

readings of x or T are spaced equally in time, the direct manual computation may be less time-consuming than establishing the validity of

a shortcut method. The straightforward calculation may avoid some

of the pitfalls of the shortcut methods. Actual data lack the smoothness of the analytical equations and various schemes are employed to

reduce the effects on the calculated quantity W of the "noise" or "scatter" in the observations. The computation outlined in the next

section uses averaging of a number of points to smooth the data. If

the temperature-time data are recorded on a chart, the smoothing

and the integration can be done graphically. Data analysis by the

method of least squares with a digital computer provides the best

smoothing and the most information about the experiment. The least

squares procedure used at NBS is discussed in later sections of this

Note.

The question of whether to use a manual or an automatic digital
method should be resolved, as usual, on the basis of how much work is

to be done and what accuracy is required. The extra work of programming the computer can be offset by its greater speed and accuracy, if enough work is to be done. The computer should probably be used for

the most accuracy work, in order to obtain a more thorough analysis

of the results, even though the difference in computed values may be

slight.

2. DIRECT MANUAL EVALUATION

OF THE TEMPERATURE-TIME DATA

The method of data analysis presented in this section is suggested

for use when automatic digital computers are not available and when

the user does not wish to develop his own methods. The basis of the

method can be readily understood from a simple illustration of how the

constants required in Eq. (8) or (10) are extracted from discrete and

noisy temperature observations spaced equally in time.

Temperature-time data plotted in Fig. 1 indicate the typical

course of an experiment with an isoperibol calorimeter. The smooth

curve is sketched to emphasize that the individual data points do not

slavishly follow the theoretical single exponential. Points to be used

in our illustration are given numbers from to N + 6. Ordinarily,

seven points in each rating period will be too few, but they are sufficient to illustrate the method. From these points we must obtain

T , T , e and T^ for Eq. (8) and the integral from t to t . Equation (3) is the basic equation for determining the four constants. The derivative dT/dt can be approximated by AT/ At using any

pair of points in a rating period, but a derivative based on adjacent

points will likely be a poor approximation because of the scatter in the

data. Spacing the points farther apart will reduce the effects of

scatter but may get into the problem of how well the slope of a chord

10
approximates the derivative of an exponential curve. This problem

will be considered at the end of this section.

To apply Eq. (3) to the data in Fig. 1, the slope of the chord is

taken as an approximation to the derivative at the time corresponding

to the midpoint of the. chord. The following approximation to Eq. (3) is

obtained from the points, spaced at equal time intervals z, and numbered to 4 and point 2, which is midway between

T

4

' T

, , (ID —-

:

= - e (T - T . 4 z 2 °°

The time interval is just four times the fixed interval z between

observations. From the first seven points, two more equations can be

written:

T

5

- T

~ii = - S(T

3

" T»>

T

6

- T

2

~4i = " £ (T

4

" T»»-

The average of these equations is

1 /

T

4

+T

5

+T

6

T +T l + T

2\

iz\ 3 3 ;

-T , + T^ + T, T + T, + T \ T^ + T + T.

. ,_? I i+ eT . d2)

3

oo

If the initial rating period is divided into a first section preceding the

midpoint and a second section following the midpoint (see Fig. 1), the

preceding equation represents the average temperature of the first

section minus the average temperature of the second section divided

by the time interval and equal to the cooling constant e times the

11
average temperature of the interval less the quantity € T . Seven points are usually not enough to get the noise out of the

data. We therefore need to generalize Eq. (12) to deal with a larger

number of points. The number of points is not arbitrary but is restricted by Eq. (11) which requires that there be an observed temperature (e. g. , T ) midway between two observed temperatures (e. g. , T

and T ). This method of computation works if we take 2K equal intervals in each rating period where K is an odd integer. Then the more

general form of Eq. (12) for the initial rating period can be written

2K K-l

E V E T

i i.K+1 i=0 ^ (3K - ! ) (14)

= -€ y T. +e K T . i=J(K+l)

For the final rating period starting at point number N and extending for

2K time intervals, as in Fig. 1, the corresponding equation can be

written

N+2K N+K-l

i=N+ K+ l

X i=N N+ i(3K-l)

= - e V T. + 6KT . (K+l) z i=N+|(K+l)

l °o

Equations (14) and (15) are solved simultaneously to obtain the following

formula for calculating the cooling constant e

12
N+2K N+K-l 2K K-l

Z T i" Z T

i" E Ti + E Ti _1 i=N+K+l i=N i= K+ 1 i=

(K+l) z i(3K-l) N+i(3K-l) (16)

E v E

i=J(K+l) i=N+J(K+l)

T.

i

This value of the cooling constant is inserted in Eq. (14) or (15), which

is then solved for T , giving two of the quantities required by Eq. (8).

The scheme of computation indicated in Fig,, 1 requires the temperature T.. and T at the midpoints of the rating periods and integration

1 F

between those points. We use the average temperatures obtained bydividing the sums from the right hand side of Eqs. (14) and (15) by K,

the number of temperatures in the sums

i(3K-l)

Tl4 E T,

i=i(K+l)

N+i(3K-l) (1?)

k E

i=N+i(K+l)

T

F

: K 2^ T.

The integration required in Eqs. (8) and (10) is carried out

between the midpoints of the rating periods; that is, the time t = t and t = t,, . are limits on the integration in Eqs. (8) and (10). The F N+K s integration is performed by the trapezoidal rule, which approximates

the curve between adjacent points by a straight line. By this rule, the

integral I is equal to the product of the length of the time interval z

13
and the sum of one -half the first and last temperatures plus all the

temperatures in between:

N+K -1

(« ^

K

+ X! T. + \ T ) *-* i 2 N+K/

i=K+l

The analytical formula of Eq. (8) or (10) can now be represented in

terms of the temperature observations

W

N+i(3K-l)

i=N+i(K+l)

i(3K-l)

i+|(K+l)

T. |+ el - e N z T

OO (18)

There is a problem, as mentioned earlier, of possible errors

from taking the time interval z too long or from taking the rating

periods too long. The problem arises because the theoretical exponential curve has been approximated by a straight line to obtain the derivatives in Eqs. (14) and (15) and the average temperatures T and T

F I in Eqs. (17) and (18). Neglecting the scatter in the data, we will consider just the error in these approximations.

The temperature in the initial rating period is given by Eq. (5).

To remove any reference to a particular rating period or part of a

rating period, we write the equation in the form

T = Ae oo et

(19)

where A is the value of T - T at the time chosen for t - 0. A shift of

oo the time scale (from initial to final rating period, for example) will

merely require a change in the value of A and will not affect our

argument on the relative magnitude of the errors. We first consider

14
the error in the derivative in Eq. (11). For chosen value of K an

approximate value of the derivative is computed by the left-hand side

of Eq. (11).

\ T - T dT\ K+l

)

dt / app (K+ 1 ) z

Substitution into this equation from Eq. (19) gives

dT \

dt j

A

. - e (K+ 1 ) z A(e - 1) (20) aPP (K+l) z From Eq. (19) the actual derivative at the midpoint of this interval is

dT _ -£e(K+l) z (21)

To evaluate the error, we expand the exponentials in a series

and subtract Eq. (20) from Eq. (21) obtaining

Error = e A t^llfjL . ^(K+nV +

\ . \ 24 48 * ')

Dividing by Eq. (21), we obtain the fractional error

^«W.2 -^.W.3

t...

Fractional Error = . . -e(K+l) z/2

e Evidently the fractional error increases as larger values are chosen

for K and z. Proper choice of K and z may make the error as small

as required. As an example of a practical case, for K = 9, z = 5 sec. , and the calorimeter cooling constant = . 002 sec , the fractional error

in the derivative is only . 0004, which is acceptable for most work.

The error in the derivatives causes an error in Eq. (16).

15
Usually the derivative in the initial rating period is small, so that the

fractional error in the cooling constant is about the same as the fractional error in the derivative. The corresponding error in Eq. (18) can

be calculated when the integral I and the temperature T are known.

The average temperatures also appear in Eq. (18) and these also

have small errors introduced by the method of computation. In Eq. (17)

we average K successive points. The temperature of the n

tJl point T

is obtained according to Eq. (19):

T = T + Ae" enZ

n °o

Numbering these points from to K-l and writing the exponential in

series form we obtain the average temperature T in the following

form

K-l K-l

Tav4£vT% 4J:[,.,z+ ^...l

(22) n= n= o

The temperature T at the midpoint, when t = ]rz(K-l), is also obtained

from Eq. (19)

A r (K-l) €z , (K-l)

2

e 2 z

2

"I T = Too +A

l!- —T— + 5 •••_!• U3)

From Eqs. (22) and (23) we find the error

K-l

T - T = A

(K-l) e z i V^ 2 2 2

av ft " 2 Z-f n e z + • n=0

r

2

or, after simplifying and using 2_\ n = r(r+l) (2r + l)/6,

n=0

16
2 2 T - T = A 6 z av

[(K-l) 2

/ 8 - (K-l) (2K-1)/12 + . ..]

Eq. (18) requires the temperature difference T - T . The error in

F I this difference is obtained from the preceding equation:

(A F

-A

I

) e Z z

2

[(K-l)

2

/ 8 - (K-l) (2K-D/12+ . .]

Fractional Error = 1

.__„ ,

. (A -AJ e-J^'K" 1 ' F 1

For a practical case in which K = 9, z = 5 sec, and e = . 002 sec , the

-4

fractional error is about 3X10 . The preceding discussion indicates

that there is a compromise between taking a large number of points to

avoid systematic errors in the computation. That the manual method

can be quite adequate is shown by comparison with results of the least

squares method described in the next section. For a calorimeter used

at the National Bureau of Standards having a cooling constant of . 0024

_1

sec and with K = 9, z = 2 sec, the manual method gives 4287 for the

bracket in Eq. (18), compared to 4285 for the least squares method.

Using twice as many points to obtain the average temperatures T and

F

T , the manual method gives 4288. The precision of the calorimeter

is a few tenths of a percent, so these results would probably be acceptable. The main disadvantage of the manual method, aside from the

time involved in the computation, is that one must be continually on

guard against these systematic errors. For similar repetitive measurements the systematic errors need be examined only once and the

manual method should prove quite satisfactory.

3. DERIVATION OF THE LEAST SQUARES SOLUTIONS.

General discussions of the method of least squares for the

treatment of experimental data can be found in many books on

17
12

experimental physics or statistics. Natrella " deals carefully with

the problem and gives many worked examples to illustrate the applications. This treatment is recommended reading for those who would

13

develop their own methods of data analysis. Young gives a readable

introduction for those who would like a better understanding of the

power and elegance of the method and the basis of the development

which follows.

For isoperibol calorimetry, the relationships required by the

least squares data analysis are derived from Eq. (5). Using subscripts

I and F to designate temperatures and time in the initial or "fore"

rating period and the final or "after" rating period, we have the two

equations

(24)

T - T~= (TF

- TJ e"

€(t

" V Final

These equations contain the quantities T , T , T , and e which are

required in the calculation of Eq. (10). These are the quantities for

which we want to make the best possible evaluation from the data

available.

The actual data points will not lie exactly on the smooth curves

of Eq. (24), so we define residuals r for each observation as the difference between the observed temperature at the time of the observation

and the corresponding value calculated from Eqs. (24). In order to

simplify the work, the substitutions A = T - TT and B = T - T„ are

made to give the following expressions for the residuals -- r. for the

i point in the fore period and r. for the j ttl point in the after period.

18
r = T. - T + A e"

€(ti" *l) , ... , ,_ c

. l l oo Initial (25a-)

VVT«> + B ."<*-tF>

Fina! (25b)

Note that the difference T - T required by Eq. (8) is equal to the F 1 difference A-B. The least squares calculation is set up to minimize

the sum S of the squares of the residuals. If there are S. time intervals

in the initial rating period and mtime intervals in the final rating

period,

1 2

m

2

i=o

l j=0 J (26)

Although the A and B contain T , which is also to be determined from

the data, it can be shown that substitution of these parameters does not

alter the values obtained for TT , 17 , e, and T by the method of least

I f oo ' squares. Ordinarily one should plan the experiment so that i and m

are approximately equal, say within a factor of two of one another. If

I is very different from m, the estimates of T and e will depend primarily on the larger set of data, as is apparent from Eq. (13). In most

of the work in this laboratory, the uncertainties in the temperature

(in microvolts) and time (in seconds) are the same for both periods, so

that we usually take i - m. However, the development which follows

will allow for i ^ m, to allow for flexibility. (The computer program

in the next section requires that both I and m be even numbers. ) For the sum S to be a minimum the partial derivatives of S with

respect c to the four constants -- A, B, T , and e -- shall be zero.

OO Substituting from (25) into (26), taking the derivatives and putting

Z. = t. - t„ and Z = t -t„ we obtain the equations

1 i I J J F

19
^=2 V

i=0

dA

:' €Zi (T - T + A e' €Zi

)

i oo (27a)

m

££ := 2 V e"

£Z

J (T. - T + B e"

eZ

J) (27 b)

as

3TOO

- 2

"i

i=0

p

pl=-2

m

T + A e"

eZi

) + V (T.

oo *-^ 1 j=0 J T + B e J) OO

i=0

A Z. e' €Zi (T. - T + A e~ eZi

) 1 l

°° m

B Z. e"

€Z

J (T. - T + B e~ €Z

J)

jr ^ j °°

(27c)

(27d)

The least squares problem is now reduced to finding values of

A, B, T , and e which make all four partial derivatives simultaneouslyequal to zero. This task is complicated by the appearance of € in the

argument of an exponential, but the solution can be found by successive

approximations. We will now develop the equations needed for this

procedure. To reduce the number of symbols which must be written

in the subsequent development, we define the following symbols in terms

m

s = 2^ e J L

mE-eZ- V

of the ob served quantities.

I

3 o

l

s 5

= i:(e- €Zi

) 2

6

20
i m

.,-Et, s 8 =Et.

I m

i m

i m

s

Setting Eqs. (27) equal to zero, which is the condition for a minimum

value of S, substituting the above definitions, and rearranging, we

obtain the following four equations:

s

l Too

" s 3 (28a)

A =

S 5

^

s 2

T<*> " S 4 (28b)

s 6

3

7

+s

8

- Sl s 3

/s

5

-s

2

s 4

/ s 6 Too = 2

1 (28c)

i+m+2-s,/s_- 15 s / s , 2 6

A S

9

" A Too S 1X + A

2

s I3 + B s 1Q

- B Tm

. + B

Z

. = (28d)

21
These equations can be solved for A, B, T , and e. The procedure

used is to make a first estimate of e by the method outlined in the preceding section and to use this value to estimate T in Eq. (28c). Both

quantities are then used in Eqs. (28a) and (28b) to estimate A and B.

Newton's method (discussed in books on calculus) is then used to estimate a new value for e from Eq. (28d), based on the values just obtained

for T , A and B. If we call the left-hand side of Eq. (28d) F(e), an

(i+l)th estimate of € is obtained from the relationship

6. - e. "i+l^ €

i + F. - F

' F

i>

i i1

where the derivative is approximated by Ae/AF. The value of e is

arbitrarily taken to be 1. 001 e and is used to calculate F . After the

new value of e is obtained, new values of T . A and B are estimated

and the procedure is repeated until successive estimates of e differ by

a negligible amount.

4. THE COMPUTER PROGRAM

The computer program which carries out the isoperibol calculation is written in a version of BASIC language for use on a time-sharing,

remote access computer. This version of BASIC does not require the

use of LET in a replacement statement (e. g. , line 145). The lines and

their functions are as follows:

80 DATA 2, 1>40>40> 15

90 READ Z>N2,L*M,N5

The data supplied in 80 and identified with variables in 90 are: Z - constant length of time interval (normally seconds, but may be arbitrary if constant); N2 = the number of successive blocks of data to be

used from the after period to give N2 separate solutions; L = number of

22
equal time intervals in the fore period; M = number of time intervals

in the after period; N5 = the number of points beyond the maximum

temperature where the final rating period is to be started. The number

N2 will ordinarily be one. Putting N2 = 2 will cause the problem to be

done twice and will require new values of L, M, and N5, which will be

read in line 715. This procedure checks whether different combinations

of AU and Q give the same answer.

100 OPEN /C/* INPUT

Instructs the computer to open the data file / C/ in which the temperature data have been stored. The data must be in the order in which

they were taken. The temperature may be in any scale or on an arbitrary X scale (See Eq. (9)). The last entry in the data file is a fake

100, 000 used for identification. For this program, the data must be

taken at equal time intervals and no temperature observation may be

omitted from the file /C/.

110 DIM TC200)* I (200)* SC 14)

A dimension statement reserving storage for 200 temperatures, 200 integrals, and 14 summations.

120 FOR 1=0 TO 200

130 INPUT FILE T

140 IF T=100000 THEN 170

145 N1 = I 150 T(I)=T

160 NEXT I Reads the data into storage and assigns each temperature a subscript.

Line 140 ends the read-in when it sees the false 100, 000. Line 145

counts the number of data points and line 175 prints the total.

23
170 I(0)=0

175 PRINT"TOTAL POINTS = "N1 + 1 180 FOR 1=1 TO Ml

190 1(1 > = I(I-l>+Z*T(I-l>/2+Z*TCI>/2

200 NEXT I

Sets the integral to zero and computes by the trapezoidal approximation and stores the integral from the first point to each point used in

the calculation.

210 T9=TC0>

220 FOR I=L TO Nl

230 IF T(I)<T9 THEN 250

240 T9=T<I>

245 N9=I

250 NEXT I 252 PRINT "MAX TEMP IS T ("N9")"

253 PRINT

254 PRINT"TIME CORR RISE COOLG CONST CONV TEMP INIT T FINAL T' 255 N=N9+N5

Selects the maximum temperature by a series of comparisons, assigns

it the number N9, and prints it out. In 255 the point N, which is the

beginning of the final rating period, is taken N5 points beyond the

maximum. This number will depend on the particular calorimeter. If

point N consistently shows a large deviation from the smooth curve

(deviations are calculated beginning with Line 660) the after period

should be started farther beyond the maximum. The maximum temperature is stored in 240 and printed out in 252.

256 T1 = 257 T2=0

258 T3=0

259 T4=0

260 FOR 1=0 TO (L/2-1)

2 62 Tl=Tl+4*CTCL-I>-T<I>3/CLt2+2*L>/Z

24
264 NEXT I 266 FOR 1=0 TO L

268 T2=T2 + T(D/(L+1)

270 NEXT I 272 FOR 1=0 TO (M/2-1)

27 4 T3=T3+4*CT(N+M-I )-T(N+I ) 3 /(M t 2+2*M) /Z

276 NEXT I 278 FOR 1=0 TO M

280 T4=T4+T(N+I )/(M+l)

282 NEXT I 300 EC 1) = (T1-T3)/CT4-T2)

Following the procedure very similar to that outlined in Section 2,

these lines make a first estimate of e from the data, using Eq. (3).

The average derivative is Tl in the initial rating period and T3 in the

final rating period. The average temperature is T2 in the initial rating

period and T4 in the final rating period.

3 10 E(0)=1 .001*E< 1 )

Computes a fake cooling constant for later use in Newton's method of

approximation.

315 N4=0

Sets the index N4 to zero. N4 counts the number of complete least

squares solutions using different blocks of data.

320 FOR 1=0 TO 9

Starts the successive approximations and limits the number of passes

to ten. The computation for I = is not based on the data, but it provides values for the next approximation.

330 FOR J= 1 TO 14

340 S(J>=0

3 50 NEXT J

Sets all sums to zero. The subscripts correspond to those in Eqs. (28).

25
360 FOR J=0 TO L

370 X=EXPC-Z*ECI >*J)

380 SC1>=SC1>+X

390 SC3>=SC3>+TCJ)*X

400 S(5)=S(5)+Xt2

410 S(7)=SC7)+T( J)

420 SC9>=SC9>+Z*J*TCJ)*X

430 S( 1 1) = SC 1 1>+Z*J*X

440 S(13)=S( 13>+Z*J*Xt2

442 NEXT J

These lines compute the summations for the value of E(I). These

summations are based on the data in the initial rating period.

444 FOR J=0 TO M

446 X=EXPC-Z*E(I)*J]

448 S(2)=S(2)+X

450 SC4>=SC4>+TCN+J)*X

455 S(6)=S(6)+Xt2

460 SC8>=SC8>+TCN+J>

470 SC 10)=S( 10)+Z*J*TCN+J)*X

474 S(12)=S( 12)+Z*J*X

476 SC 14>=S( 14>+Z*J*Xt2

480 NEXT J

These lines compute the corresponding summations from the data in

the final rating period.

490 T8 = S(7)+S(8>-SC3)*S( 1 ) /SC 5) - SC 4> *SC 2) /SC 6.) 49 5 T7 = L + M + 2-SO> t2/SC5)-SC2) T2/SC6)

500 T8=T8/T7

The value of T is computed according to Eq. (28c).

oo

510 A=CT8*SC1)-SC3)]/SC5)

520 B=CT8*SC2)-SC4)]/SC6)

The values of A and B are computed according to Eqs. (28a) and (28b).

26
530 F(I >=A*S(9)+B*S< 10)-A*T8*S< 1 1 ) -B*T8*SC 12> + At2*SC 1 3> +B t2*S< 1 4)

This line computes the value of the left-hand side of Eq. (28d), which

we seek to make zero. The quantity F is given a subscript because two

successive approximations are used in Line 580.

550 IF I<1 THEN 590

560 IF F(I)=F(I-1) THEN 800

570 IF 1=9 THEN 592

580 E(I+1)=E<I)-CE(I-1)-E(I)3*F(I)/CF(I-1)-F(I)]

582 IF ABSCE(I+1>-E(I) 3<lE-7 THEN 740

590 NEXT I 592 PRINT "RAN TEN TIMES"

These lines compute the next approximation to the cooling constant

E(I +1). Line 550 avoids the computation for 1=0 because it cannot

be computed until two successive values of E and F are available. Line

560 prevents the possibility of division by zero in 580. Line 570 stops

the number of approximations at nine, but completes the calculation

using E(9). Line 580 computes the next approximation E(I + 1) by

Newton's method, using the chord AE/AF to approximate the derivative.

Line 582 terminates the approximations provided that 3 successive

-7

values differ by less than 10 . This criterion may have to be altered

to suit a particular calorimeter or time scale. It is selected for a

typical cooling constant of . 002 sec ; it would probably not be suitable

for very much smaller cooling constants. Line 592 provides printout

to indicate a possible fault.

605 T5=A-B+EU >*CI CN)-T8*Z*N]

610 PRINT Z*N*T5*1E-6*INT<1E6*ECI>>ME-3*INT<1E3*T8>*T8-A;T8-B

Line 605 calculates the corrected temperature rise, which is the

quantity in brackets in Eq. (8).

27
664 U=0

668 FOR J=0 TO L

670 C=T8-A*EXPC-E(I >*Z*J3

674 U=U+CC-T(J>] t2

678 PRINT T(J>*C; 1 E- 4*1 NTC < TC J> -C> * 1 E4)

680 NEXT J

684 FOR J=N TO N+M

686 Z1=Z*CJ-N>

688 C=T8-B*EXPC-E«I)*Z1]

690 U=U+CC-T< J)3 »2

694 PRINT T(J)5C; 1E-4*INT( (T(J>-C)*1E4)

696 NEXT J

698 PRINT "RESIDUAL STD DEV= "SQRC U/CL+M-2)

)

These lines compute the smooth exponential curves for the initial and

final rating periods and the deviations of the experimental points from

these curves. The number of points L + M + Z is decreased by four in

698 because four constants have been determined from the data.

700 IF N4=N2-1 THEN 730

7 10 N4=N4+1

7 15 READ L*M,N5

720 GO TO 320

7 30 STOP

Line 700 determines whether another calculation has been asked for.

If so, line 710 increases the index N4 and line 715 reads the parameters

which determine the data to be used in the next computation.

740 IF ABSCE(I)-E(I-l)]<lE-7 THEN 605

750 GO TO 590

These lines complete the criterion for terminating the successive

approximations.

800 PRINT "FCI >=F(I-l)"j"F(I)= M F(I)

810 GO TO 605

Line 800 serves as an alert to a possible fault.

28
5. SAMPLE CALCULATIONS WITH THE COMPUTER PROGRAM.

As an illustration of its operation and the interpretation of

results, two computations have been made on the data in Table I, which

are from a typical experiment with a small calorimeter -- one of a group

designated the C-series and having energy equivalents corresponding to

1 or 2 joules per kelvin and cooling constants between . 00Z and . 003

-1

sec

The print-out for the first computation is shown in Table II. This computation took the first 40 points for the initial rating period

and 40 points for the final rating period starting 15 points past the

maximum temperature. The first two lines of the print-out give information to judge whether the data are sufficient. In this case, the computation used 55 points beyond the maximum, which occurred at the

49 point. There were thus 104 points used, comfortably less than

the total of 1 1 9 available.

The next line shows the length of time between the initial and

final temperatures, the corrected rise, which is the desired quantity

in the bracket in Eq. (8), the cooling constant s, the convergence temperature T , and the initial and final temperatures. Note that T is

well below any observed temperature, which is consistent with the

decrease in temperature during the initial rating period.

The balance of the print-out consists of a column of the input

data, a column of the calculated temperatures for the times corresponding to the input data, and a column of the differences. The differences

indicate how well the computation fits the data for the rating periods.

It is extremely useful in discovering points which have been entered

incorrectly, or an improper fit of the data. The last entry is the

standard deviation of the observations from the smooth curves. Points

with deviations three or more times this quantity should be checked.

29
A long run of deviations of the same sign indicates a poor fit

of the data. An example is given in Table III, which shows the deviations for the same set of data, but with the final rating period starting

only one point past the maximum. Note the long series of points which

lie above the smooth curves and followed by the long series of points

which lie below the smooth curves. The standard deviation also indicates a poor fit, but it is a less sensitive indicator.

30
ooin
>ccvjoor"-vO',Tis*ri/-

) vor-oocomooo — toocvjpCVJCVJCVJCVJOOnOTC)— OCOnO ssssnnnnnwftiw i » i i + + + + + + +

u

<D

(-> CDH

•-Iu

o

I—

I

U

.—

i

o

CDo

CO

1—1d

oo

-,r<str>cv)r-\0'T-*voo<s Nor^coooios — co o o o <>

i i » i + + + ++ + +

t-COOs^tSvOU") — CO CV) T ^

- >or-r-cocoojcoif>t*- —> in ts>

I I ! I + + ++ + +

wcv]<\)OJr-r-u">r>.-iscor~

i i i i + + + ++ + + +

NOr-t*»cotSoor-co — <r co co tststsscoconcoconcvjcv]

i i i i + + + ++ + + m «- r- co <m co — \oo — «— o

vOr-l*-CO'<TCOOtS3<\jNOtS)'<T wcv]Cvicv]or-in<TOJtsori i i i + + + + + +

W

O

P

WtMOjajoooNOTcvjcsorStStS»tS»-'C0C*5(*)<DC0OJC\J

I I I I + + + »+ + + +

"TSsoojvot-tscvj ^•^woo

OJWWCVJ'TOOO'^'CVJSONt^ ssisissnfooncocMW

I I I I + + + ++* + +

\ot-f-cocjco<rN0co — <To

WCVJWCVJCVJCDvO^CV) — or-

<§<StSQ<SCOC0C0C0C0CVJOJ

I » • I !++ + + COOlD — vO^DOO o r» ^" c* s

CVJWWOJWCOsOTCVJ — ooos

S)StS(S<SC*)C0f0C0C0CVJCV3®

I I I I I+ + ++ + + +®

31
Table II

Print-out from Data Analysis Program

TOTAL POINTS = 1 19

MAX TEMP AT PT NO 49

TIME CORR RISE COOLG CONST CONV TEMP I NIT T FINAL T

128 4285.322 .2374E-02 -397.253 -262.815 3590.65C

-263 -262.815 -.185 3591 3590.653 .3465

-263 -263.4521 .4521 3572 3571 .757 .2425

-264 -264.086 .086 3553 3552.95 •49E-01

-265 -264.7171 -.283 3535 3534.233 .7664

-265 -265.345 .345 3516 3515.604 .3951

-266 -265.9701 -.3E-01 3497 3497.064 - .0644

-267 -266.5921 - .40 79 3479 3478.61

1

.3881

-267 -267.21 13 .21 12 3 460 3460.247 -.2467

-268 -267.8274 -.1726 3442 3441 .968 •314E-01

-268 -268.4407 .4407 3424 3423.777 .223

-269 -269.051

1

.51E-01 3406 3405.672 .3283

-270 -269.6585 -.3415 3387 3387.652 -.6522

-270 -270.2631 .263 3369 3369.718 -.718

-271 -270.8648 -.1352 3351 3351 .869 -.8688

-271 -271 .4636 .4636 3334 3334.104 -.1042

-272 -272.0597 .597E-01 3316 3316.424 -.4238

-273 -272.6529 -.3471 3299 3298.827 .1728

-273 -273.2433 .2433 3281 3281 .313 - . 3 1 39

-274 -273.8309 -.1691 3264 3263.884 .1 164

-275 -274.4157 -.5843 3246 3246.536 -.5359

-275 -274.9978 -.23E-02 3229 3229.27 -.2704

-276 -275.577 - .423 3212 3212.087 -.0867

-276 -276.1535 .1535 3195 3194.984 .156E-01

-277 -276.7274 -.2727 3178 3177.963 •368E-01

-277 -277.2985 .2984 3161 3161 .023 -.226E-01

-278 -277.8669 - . 1332 3144 3144. 162 -.1623

-278 -278.4325 .4325 3127 3127.382 -.3819

-279 -278.9955 -.45E-02 31 11 31 10.681 .319

-280 -279.5558 - .4442 3094 3094.059 -.592E-01

-280 -280.1 136 .1135 3078 3077.516 .4838

-281 -280.6686 -.3314 3061 3061 .051 -.516E-01

-281 -281 .221 .221 3045 3044.664 .335

-282 -281 .7708 -.2292 3028 3028.356 -.356

-282 -282.318 .318 3012 3012.124 -.1243

-283 -282.8626 -.1374 2996 2995.969 •305E-01

-283 -283.4046 .4046 2980 2979.891 .1087

-284 -283.944 -.56E-01 2964 2963.889 .1 107

-284 -284.4809 .4809 2948 2947.963 .37E-01

-285 -285.0153 •153E-01 2932 2932.1 12 -. 1 123

-286 -285.5471 -.4529 29 17 2916.337 .6634

-286 -286.0764 .0764 2901 2900.636 .36 43

RESIDUAL STD DEV= .3331032

32
Table III Print-out for Rating Period Started Too Soon.

3844 3860.785 -16.7851

3832 3840.762 -8.7619

3817 3820.833 -3.8328

3801 3800.997 •26E-02

3783 3781 .255 1 .7448

3764 3761 .606 2.39 42

3745 3742.049 2.9512

3726 3722.583 3.4163

3707 3703.21 3.7899

3688 3683.927 4.0725

3668 3664.735 3.2644

3649 3645.634 3.3662

3630 3626.621 3.3781

361 1 3607.699 3.3007

359 1 3588.865 2.1343

3572 3570.12 1 .8795

3553 3551 .463 1 .5365

3535 3532.894 2. 1058

3516 3514.412 1 .5879

349 7 3 49 6.017 .9831

3479 3477.708 1 .2918

3460 3459.485 .5145

3442 3441 .348 .6515

3424 3423.297 .7033

3406 3405.33 .6702

3387 3387 .447 -.4473

3369 3369.649 -.6489

3351 3351 .934 -.9341

3334 3334.303 -.3026

3316 3316.753 -.754

3299 3299.288 - .2878

3281 3281 .903 -.9 38

3264 3264.601 - .6014

3246 3247.38 -1 .3804

3229 3230.24 -1 .2403

3212 3213. 18 -1 .1808

3195 3196 .201 -1 .2014

3178 3179 .301 -1 .3019

3161 3162.48 1 -1 .4818

3144 3145.741 -1 .7407

3127 3129 .078 -2.0784

RESIDUAL STD DEV = 2.591687

33
REFERENCES

1. D. A. Jennings, IEEE Transactions on Instr. and Meas. IM 15, 161

(1966).

2. D. A. Jennings, E. D. West, K. M. Evenson, A. L. Rasmussen,

and W. R. Simmons, NBS Technical Note 382 (1969).

3. J. Coops, R. S. Jessup, and K. Van Nes in Reference 9.

4. J. M. Sturtevant in Technique of Organic Chemistry, A. Weissberger,

ed. , Interscience Publishers, Inc. , New York (1959).

5. W. P. White, The Modern Calorimeter, Chemical Catalog Co. , New York (1928).

6. W. A. Roth and F. Becker, Kalorimetrische Methoden zur

Bestimmung Chemischer Reaktionswarmen, F. Vieweg and Son,

Braunschweig, Germany (1956).

7. E. D. West and K. L. Churney, J. Appl. Physics 39, 4206 (1968).

8. E. D. West and K. L. Churney, J. Appl. Phys. 41_, 2705 (1970.

9. Experimental Thermochemistry, F. D. Rossini, ed. , Interscience

Publishers, Inc., New York (1 956).

10. Experimental Thermochemistry, Vol. II, H. A. Skinner, ed. , Interscience Publishers, Inc., New York (1962).

11. Experimental Thermodynamics, Vol. I, J. P. McCullough and

D. W. Scott, ed. , Plenum Press, New York (1 967 ).

12. M. G. Natrella, Experimental Statistics, Nat'l. Bur. Stds.

Handbook 91, U. S. Gov't. Printing Office, Washington, D. C. (1963).

13. H. D. Young, Statistical Treatment of Experimental Data, McGrawHill Book Co. , New York (1962).

34

GPO 833 - 296

U.S. DEPARTMENT OF COMMERCE

WASHINGTON, D.C. 20230

OFFICIAL BUSINESS

PENALTY FOR PRIVATE USE. $300

POSTAGE AND FEES PAID

U.S. DEPARTMENT OF COMMERC