MASTER

OAK RIDGE NATIONAL LABORATORY

 

operated by

UNION CARBIDE CORPORATION B
for the

U.S. ATOMIC ENERGY COMMISSION
ORNL- TM- 1070

 

STABILITY ANALYSIS OF THE MOLTEN-SALT

REACTOR EXPERIMENT

S. J. Ball
T. W. Kerlin

 

NOTICE This document contains information of a preliminary nature
and was prepared primarily for internal use at the Oak Ridge National
Laboratory. it is subject.to revision or correction and therefore does
not represent a final report.
—_———— e

e ieiiiieo o~ LEGA e e m e e
— LEGAL NOTICE 7

This report was prepared as an account of Government sponsored work. MNeither the United States,

nor the Commission, nor any person acting on behalf of the Commission:

A. Makes any worranty or representation, expressed or implied, with respect to the accuracy,
completaness, or usefulness of the information contained in this report, or that the use of

any informotion, apparatus, method, or process disclosed in this report may not infringe

privately owned rights; or
B. Assumes any liabilities with respect to the use of, or for damages resulting from the use of
ony information, apparatus, method, or process disclosed in this report.
As used in the above, ‘‘person acting on beha!f of the Commission’’ includes any employee or
contractor of the Commission, or employee of such contractor, to the extent that such employee
or contractor of the Commission, or employee of such contractor prepares, disseminates, or
provides access to, any information pursuant to his employment or contract with the Commission,

or his employment with such contractor.

 

 

 

 
et e g

wh ’

 

 

 

ORNL-TM~1070

*

Contract Nd. W-7405~eng-26

 

1IN MCERAR SOTENCE ABSTRAINS

L3 o

-

 

STABILITY ANALYSIS OF THE MOLTEN-SALT
REACTOR EXPERIMENT

S. J. Ball o
Instrumentation and Controls Division

T. W. Kerlin
Reactor Division

" LEGAL NOTICE

Thie report was prepared as an account of Government sponsored work, Neither the United
Btates, nor the Commission, nor any person acting on behalf of the Commission:
A. Makes any warranty or Tepresentation, expressed or implied, with respect to the accy-
' racy, gompletenens, or usefulness of the information contained in this report, or that the uge
;" of any information, apparatus, method, or process disclosed In this feport may not infrings
: privately owned rights; or - .
- B. Assumes any Uabilities with respect o the uge of, or for damages resulting from the
use of any tnformation, spparatus, method, or process disclossd in this report,
As used In the above, ‘‘person acting on behalf of the Commission’ tncludes sny sm-~
© ployee or contractor of the Commission, or empioyee of such contractor, fo the extent that
* such employee or contractor of the Commiesion, or employee of guch contractor prepares, -
diaseminates, or provides access to, any information pursuant to his employment or contract
- -with the Commission, or his employment with such contractor, ’

 DECEMBER 1965

' OAK RIDGE NATTONAL LABORATORY
~ Osk Ridge, Tennessee
. ...° operated by . o
' UNION CARBIDE CORPORATION
. for the |
U.S. ATOMIC ENERGY COMMISSION

 

 

 

 
 

B Sy T R e S ARy ST S e ol S L T e e AT A T

e

   

¢

w e crp e
¥ o S
» P P
=

3

%

- o * ek e s

S T N Ay

= 4.
e
.ffi“w

o)

 

oo

ui) l

i A e

.8. COnCluSlonS .Q..,ll..V!....VQ.V.Ill..l..\‘.’Il.. .......... . 0

iii
CONTENTS

AbstraCt 8 & & 0 8 0 & W PSS RS ES RN RSN Nt 2 bR a0 e e

’ l. IntrOduCtion v e . ¢ & 9 8 2 s s 2 ® 8 0 8 2 8 0 5 0 9 B U NS P L e e h
2 . Descripticn Of the LiSRE * 8 2 8 28 0 0N BB O SRS RS R e s RN
3. Review of Studies of MSRE DynamicCsS sceeeeescancasaas cons

4. Description of Thedretical MOGElS veveraneenrsrononeanas

Core Fluid Flow and Heat Transfer .......cevvvnvnninnns
NeutronKineticS 8 0 & F e et B PR RS SERRS T ENSS OSSN
Heat Exchanger ana REAIBEOT wevernnnrsunrnonsaseoensens

Fluid Transport and Heat Transfer in Connecting Piping

XenonBehaVior ® & 8 9 0 " P 4 088 B 0 BB BB SRS S g e ® 9 8 & & & 0 & 09

Delayed POWEr «eveeevcecerscssensossessncans treasne csenes
5. Stability Analysis ResultsS seeessacsnvccecnssnaans ceaan

Transient ReSpPONSE scesecrcccisorsssssrssscncssssescsnnsse

Closed-Loop Frequency ReSPONSE sseesesescessocrnsasensns

Open-Loop Frequency Response ceeeveeetceenesseevannenss

Pole Configuration eeevseescieevnvessnessacesassasennans

6. Interpretation of Results s.iivecieerensvesensarnerocans
Explanation of the Inherent Stability Characteristics .

Interpretation of Early Results ¢.eceeeee Ceesienasan -

7. Perturbations in the Model and the Design Parameters ..

Effects of Model Changes cheerrsaana teeasarascasennne P

Effects of Parameter . Changes et e eeeeererrae s
Effects of Design Uncertainties Fasesasesns trsascsernae

REferenceS '.'...i"l.'-.."C"‘O“._.Q‘"I.'V....V..-...l'...'..rfl'llQl-

~ Appendix A. Model Equatlons ....;;......;.;..;.......}.....
f:-ApPendix B. Coefflclents USed in the Model Equatlons ...;.

* 4080

L B B BN J

s an s

8 000 P

...l.

Appeidix C. General Descrlptlon of MSRE FrequenCy—Response‘

. _code_V'......I...'....I:l_l...l.....'.....l"'.'......'I.'.......

Appendix Dl Stablllty EX'tI‘ema C&lCUl&tlon . o. re e e e s ."o.c. *a

" %00

NN N

14
16
16
16
16
18
18
19
23
25
27
_27
31
34
34

. 35

40
42
45

49

67

71

75

 

 
 

-) ’afi .

'1‘{ ™
|
" b~ —’

o\l o

STABILITY ANALYSIS OF THE MOLTEN-SALT
REACTOR EXPERIMENT

S. J. Ball T. W. Kerlin

Abstract

A detalled analys1s shows that the Molten-Salt Reactor
Experiment is inherently stable. It has sluggish transient
response at low power, but this creates no safety or opera-
tional problems. The study included analysis of the tran-
sient response, frequency response, and pole configuration.
The effects of changes in the mathematical model for the
system and in the characteristic parameters were studied.

A systematic analysis was also made to find the set of
parameters, within the estimated uncertainty range of the
design values, that gives the least stable condition. The
system was found to be inherently stable for this condition,
as well as for the design condition.

The system stability was underestimated in earlier
studies of MSRE transient behavior. This was partly due to
the approximate model previously used. The estimates of
the values for the system parameters in the earlier studies
also led to less stable predictions than current best values.

The stability increases as the power level increases and
is largely determined by the relative reactivity contribu-
tions of the prompt feedback and the delayed feedback. The
large heat cgpacities of system components, low heat transfer

- coefficients, and fuel c1rculat10n cause the delayed reac-
t1v1ty feedback.

1. Introduction

‘Investigations of inherent stability constitute an essential part of

- a reactorrevaluationQ; ThlS 1s partlcularly true for a new type of system,
7'sucfi‘as the MSRE. The flrst con31derat10n in such an ana1y31s is to de-
eetermlne whether the system possesses 1nherent self—destructlon tendencies.

'iiLOther less 1mportant’but 51gn1f1cant con31derat10ns are the influence of
'rlnherent characterlstlcs on control system requlrements and the pos31—

éb111ty of conductlng experiments that require constant condltlons for ex-

tended perlods.

 
 

-

Several spproaches may be used for stability analysis. A complete
studyrof power reactor dynamics wouid take into account the inherent non-
linearity of the reactivity feedback. It is not difficult to calculate
- the transient response of nonlinear systemS'with analog or digital com-
puters. On the other hand, it is not currently possible to study the
stablllty of multlcomponent nonlinear systems in a general fashion. The
usual method is to 11nearlze the feedback terms in the system equatlons
and use the well-developed methods of 11near-feedback theory for stablllty
analysis. ThlS leads to the use of the frequency response (Bode plots or
Nyquist plots) or root locus for stability analysis. This study 1nc1uded
nonlinear tran51ent-response calculations and linearized frequency-response
and root-locus calculations.

The stability of a dynamic system can depend on a delicate balance
of the effects of many components. This balance may be altered by changes
in the mathematical model for the system or by changes in the values of
the parameters that characterize the system. Since neither perfect models
nor exact parameters can be obtained, the effect of changes in each of
these on predicted stability should be determined, as was done'in this
study.

- The transient and frequency responses obtained in a stability analy-
sis are also needed for comparison with results of dynamic tests on the
system. The dynamic tests may indicate that modifications in the theo-
retical model or in the design data are needed. Such modifications can
provide a confirmed model that may be used for interpreting any changes
possibly observed in the MSRE dynamic behavior in subsequent operation
and for predicting, with confidence, the stability of other similar

systems.

2. Description of the MSRE -

The MSRE is a graphite-moderated, circulating-fuel reactor with fluo-
ride ealts of uranium, lithium, beryllium, and zirconium as the fuel..l
The basic flow diagram is shown in Fig. 1. The molten fuel-bearing salt
enters the core matrix at the bottom and passes up through the core in

channels machined out of 2-in. graphite blocks. The 10 Mw of heat

.

"."Jf 4

©§ pPow
 

 

 

 

 

 

 

 

i
.xr‘L >

B
- f
.Y

&

 

 

 

 

 

 

 

PRIMARY
HEAT
EXCHANGER

 

 

 

ORNL-DWG 65-9809

PUMP

—r

 

RADIATOR

 

 

 

 

 

Fig. 1. MSRE Basic Flow Diagram.

 

 
 

 

generated in the fuel and transferred from the graphite raises the fuel
temperature from 1175°F at the inlet to 1225°F at the outlet. When the
system operates at lower power, the flow rate is the Sameras at 10 Mw,
and the temperature rise through the core decreases. The hot fuel salt
travels to the primary heat exchanger, where it transfers heat to a non-
fueled secondary salt before reentering the core. The heated secondary
salt travels to an air-cooled radiator before returning to the primary
heat exchanger. | | ' |
Dynamically, the two most important characteristics of the MSRE are
that the core is heéerogeneous and that the fuel circulates. Since this
combination of important characteristics is uncommon, a detailed study
of stability was required. The fuel circulation acts to reduce the ef-
fective delayed-neutron fraction;rto reduce the rate of fuel temperature
change during a power change, and to introduce delayed fuel-temperature
and neutron-production effects. The heterogeneity introduces a delayed

feedback effect due to graphite temperature changes.

The MSRE also has a large ratio of heat capacity to power production.

This indicates that temperatures will change slowly with power changes.
This also suggests that the effects of the negative temperature coeffi-
cients will appear slowly, and the system will be sluggish. This type
of behavior, which is more pronounced at low power, is evident in the
results of this study.

Another factor that contributes to the sluggish time response is the
heat sink — the air radiator. An approximate time constant for heating

and cooling the entire primary and secondary system was found by consider-

ing all the salt, graphite, and metal as one lumped heat capacity that
dumps heat through a resistance into the air (eink), as indicated in Fig.
2. TFor the reactor operating at 10 Mw with a mean reactor temperature
of about 1200°F and a sink temperature of about 200°F, the effective re-

sistance must be

1200°F — 200°F

o
T5 T = 100°F/Mw .

Thus the overall time constant is

v

6)jf- >
 

 

 

 

 

ORNL-DWG 65-9810
REACTOR HEAT CAPACITY 2% 12 Mw-sec/°F
S RESISTANCE VARIES WITH AIR
‘s FLOW RATE

 

SINK MEAN AIR TEMPERATURE VARIES
TEMPERATURE WITH AIR FLOW RATE

 

 

 

Fig. 2. MSRE Heat Transfer System with Primary and Secondary Sys-
tem Considered as One Lump.

 

&) h,.. v | .

 

 

 
 

12 M.sec o 100

T = ;200 sec

I

= 20 min .

For the reactor operating at 1 Mw, the sink temperature incresses to about
400°F. This is due to a reduction in cooling air flow provided at low
power to keep the fuel temperature at 1225°F at the core exit. In this
case the resistance is |
1200°F — 400°F oo
1M - 800°F/Mw ,

and the overall time constant becomes

Mw.sec P .
12 5 X 800 = 2600 sec |
=23 hr .

This very long time-response behavior would not be as pronounced with a
heat sink such as a steam generator, where the sink temperatures would

be considerably higher.

3. Review of Studies of MSRE Dynamics

Three types of studies of MSRE dynamics were previously made:
(1) transient-behavior analyses of the system during normal operatibnr
with an automatic controller, (2) sbnormal-transient and accident studies,
and (3)ftransient-behavior analyses of the system without external con-
trq}. fThe automatic rod control system operates in either a neutron-flux

control mode, for low-power operation, or in a temperature control mode

at higher powers.2 The predicted response of the reactor under servo éon—

trol for large changes in load demand indicated that the system is both
stable and controllable. The abnormal-transient and accident studies
showed that credible transients are not dangerous.>

The behavior of the reactor without sérvo control was initially in-
vestigated in 1960 and 1961 by Burke.* 7 A subsequent controller study
by Ball® in 1962 indicated that the system had greater inherent stability
 

 

‘ffi

"

o

1

A

than predicted by Burke. Figure 3 shows comparable transient responses
for the two cases. The differences in predicted response are due to dif-
ferences both in the model and in the parameters used and will be dis-
cussed in detail in Section 6.

There are two 1mportant aspects of the MSRE's inherent stablllty
characteristics that were observed in the earlier studies. First, the
reactor tends to become less stable at lower powers, and second, the
period of oscillation is very long and increases with lower powers. As
shown in Fig. 3, the period is about 9 min at 1 Mw, so any tendency of
the system to oscillate can be easily controlled. Also, since the system
is self-staebilizing at higher powers, it would not tend to run away, or
as in this case, creep away. The most objectionable aspect of inherent
oscillations would'be'their interference with tests planned for the re-

actor without automatic control.

 

4. Description of Theoretical Models

Several different models have been used in the dynamic studies of
the MSRE. Also, because the best estimates of parameter wvalues were modi-
fied periodically, each study was based on a different set of parameters.
Since the models and parameters are both important factors in the predic-
tion of stability, their influence on predicted behavior was examined in
this study. Tables 1, 2, and 3 summarize the various models and parameter
sets used. Table 1 liSte the parameters for each of the three studies,
Table 2 indicates how each part of the reactor was represented in the

three dlfferent models, and Table 3 indicates which model was used for

~ each study. . The three models are referred to subsequently as the "Re-
duced," "Intermediate," and;"Complete models, as designated in Table 2.
'The'medels'arefdeseribed'ifi;fihis section, and the equations used are given

-1n.Append1x A. The'COeffieientsrfor-each case are lieted'in_Appendix B.

Core Fluld Flow and Heat Transfer ,

 

A typlcal scheme for representing the thermal dynamics of the MSRE
core is shown in Flg. 4. The arYTrows 1nd1cat1ng heat transfer require

additional explanation. It was desired to base the calculation of heat

 
 

 

ORNL-DWG 65-9811

ITIAL STUDIES {1960, 1961)

1962 STUDY

FLUX POWER (Mw)

       

0 200 400 600 800
TIME (sec)

Fig. 3. Respbnse of MSRE Without Controller to Decrease in Load
Demand.

 
 

o

%)

r*

N

ry

Table 1.

 

Sunmary of Parameter Values

Used in MSRE Kinetics Studies

 

 

transit time, sec

Burke Ball Present
Parameter 1961 1962 Study
Fuel reactivity coefficient, °F! =3.3 x 1075 -2.8 x 10”° —4.84 x 10°°
Graphite reactivity coefficient, —6.0 X 107° —6.0 X 107° =3.7 x 1079
o
-
Fuel heat capacity, Mw:sec/°F 4,78 4.78 4.19
Effective core size, ft3 20.3 24.85 22.5
Heat trensfer coefficient from 0.02 0.0135 0.02
fuel to graphite, Mw/°F |
 Fraction of power generated in fuel 0.9% 0.9 0.934
Delayed power fraction (gamma 0.064 0.064 0.0564
heating) -
Delayed power tlme constant, sec 12 12 188
Core transit time, sec 7.63 9.342 8.46
Graphite heat capacity, Mw.sec/°F  3.75 3.528 3.58
Nuclear data |
Prompt-neutron lifetime (sec) 0.0003 0.00038 0.00024
Total delayed-neutron fraction  0.0064 0.0064 0.00666
Effective delayed-neutron 0.0036 0.0041 (0.0036)%
fraction for one-group
approximation
Effective decay constant for 0.0838 0.0838 (0.133)%
.one-§roup approximatiOn{m"
se AR
Fuel transit tlme in- external 14,37 17.03 16.73
primary circuit, sec . o | -
Total secondary 1oop coolant _-.”' 724.2_ 24.2 21.48

 

qgix groups used, see Appendlx B fbr individual delayed-neutron

fractlons (B) and decay constants (K)

 
 

 

10

Table 2. Description of Models Used
in MSRE Kinetics Studies

 

Reduced Intermediate Complete

 

Model Model - Model
Number of core regions‘ 1 9 9
TNumber of delayed-neutron groups 1 1 6
Dynamic circulating precursorsa No No _Yes.
Fluid transport lagéb ' First Fourth-order Pure

order  Padé ~ delay

Fluid-to-pipe heat transfer No Yes o Yes
Number of heat exchanger lumps 1 1 - 10
Nunber of radiator lumps 1 T 10
Xenon reactivity No No | : Yes

 

8Tn the first two models, the reduced effective delayed-neutron
fraction due to fuel circulation was assumed equal to the steady-
state value. In the third model, the transient equations were treated
explicitly (see Appendix A for details).

Prye Laplace transform of a time lag, 7, is e '~. The first
order approximation is 1/(1 + ts). The fourth order Padé aspproxima-
tion is the ratio of two fourth-order polynomials in ts, which gives
s better approximation of e™'® (see Appendix A).

Table 3. Models Used in the Various
MSRE Kinetics Studies

 

Study Model Used

 

Burke 1961 analog (refs. 4-7) Reduced

Ball 1962 analog (ref. 8) Intermediate
1965 frequency response Complete
1965 transient response Intermediate
1965 extrema determination® Reduced
1965 eigenvalue calculations Intermediate
1965 frequency response with Complete

extrema. data

 

®The worst possible conbination of pa-
rameters was used as described in Section 7.

 

€y
o

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

11
ORNL-DWG 65-9812
GRAPHITE
Ts
HEAT TRANSFER
FUEL-TO-GRAPHITE ' DRIVING FORCE
» HEAT TRANSFER :
»
TEts IN =g=d FUEL_LUMP 1 FUEL_LUMP 2 | ——T2. our
-~ ] Try 7 Tr2 —
F1
QPERFECTLY MIXED
SUBREGIONS
Fig. 4. Model of Reactor Core Region; Nuclear Power Produced in
all Three Subregions. !
f !
Loa
f}-l

)

v

M

 
 

12

transfer rate between the graphite and the fuel stream on the difference
of their average temperatures. The outlet of the first lump or "well-
stirred tank" in the fuel stream is taken as the fluid average temperaturé.
Thus a dotted arrow is shown from this point to the graphite to represent
the driving force for heat transfer. Hdwever, all the mass of the fluid
is in the lumps, and the heat transferred is distributed equally between
the lumps. Therefore solid arrows.are shown from the graphite to each
fluid lump to indicate actual transfer of heat. '

This model was developed by E. R. Mann* and has distinct'advantages
over the usual model for representing the fluid by a single lump in which
the following algebraic relationship is used to define the mean fuel tem-
perature:

TF inlet + TF outlet

) TF mean = > .

 

The outlet temperature of the model is given by

TF outlet = ZTF mean ---TF inlet .

Since the mean temperature variable represents a substantial heat
capacity (in liquid systems), it does not respond instantaneously to
changes in inlet temperature. Thus a rapid increase in inlet temperature
would cause a decrease in outlet temperature — clearly a nonphysical re-
sult. With certain limitations on the length of the flow path,9 Mann's
model avoids this difficulty.

The reduced MSRE model used one region to represent the entire core,
and the nuclear average temperatures were taken as the average graphite
temperature (TE) and the temperature of the first fuel lump (Efl)' The
nuclear average temperature is defined as the temperature that will give
the reactivity feedback effect when multiplied by the appropriate tem-
perature coefficient of reactivity.

The intermediate and complete models employ the nine-region core

model shown schematically in Fig. 5. - Each region contains two fuel lumps

 

*¥0ak Ridge National Laboratory; now deceased.

&y

‘\
 

3,

 

13

UNCLASSIFIED
ORNL~DWG 63-7349R

 (TRdour

¢ 4

 

 

 

 

 

s -
fl

 

 

 

 

 

 

 

 

 

 

 

L T
T YFiow

 

- Fig. 5. Schematié-'.'Diagrazn of Nine-Region Core Model.

 
 

 

14

end one graphite lump, as shown in Fig. 4. This gives a total of 18 lumps
(or nodes) to represent the fuel and nine lumps to represent the graphite.
The nuclear power distributiofifiand nuclear importances for all 27 lumps
were calculated with the digital code EQUIPOISE-3A, which solves steady-
state, two-group, neutron-diffusion equations in two dimensions.

Tests were made on the MSRE full-scale core hydraulic mockup'® to
check the validity of the theoretical mpdéls of core fluid transport. A
salt solution was injected'suddeniy'into“the1200—gpfi‘water stream at the
reactor vessel inlet of the mockup, and the response at therreactdr outlet
was measured by a conduétivity probe. The frequency response of the sys-
tem was computed from the time response'by Samlon's method*? for a sam-
pling rate of 2.5/sec. The equivelent mixing characteristics of the
theoretical models are computed from the transfer function of core outlet-
to-inlet temperature by omitting heat transfer to the grgphite and adding
pure delays for the time for fluid transport from the point of salt in-
jection to the core inlet and from the core outlet to the conductivity
probe location. A comparison of the experimental, one—region, and nine-
region transfer functions is shown by frequencyhrespbnse plots in Fig. 6.
Both theoretical curves compare favorably with the experimental curve,
especially in the range of frequencies important in the stability study
(0.01 to 0.1 radians/sec). The relatively large attenuation of the mag-
nitude ratios at frequencies as low as 0.1 to 1.0 radians/sec is due to
a considerable amount of axial mixing, which is to be expected at the
low Reynolds number of the core fluid flow (~1000). This test indicates
that the models used for core fuel circulation in the stability analyses

are adequate.

Neutron Kinetics

 

The standard one-point, nonlinear, neutron kinetics equations with
one gverage delayed-neutron group were used in all the analog and digital
transient response studies. Linearized equations were used for all the
other studies. In the studies of a nine-region core model, weighted
values of nuclear importance for each of the 27 lumps were used to compute

the thermal feedback reactivity. 8Six delayed-neutron groups and the

iy

 
 

i
L
i
i
2
i
i
i
i
i

 

 

- ORNL-DWG 65-9813
1.0

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

0.5
Q o2
o=
<
@
oo
: o o :
: S -REGION
£ ‘CORE” MODEL
3 3 o005 .
NINE-REGION _-
CORE MODEL
= o
- 0.02"
0.0
0
-60 '
q§~'~ i
- T
120 ™
_ N : :
| o o o N\ | —EXPERIMENTAL
; . -180 s - ey ,
' , B , 1] : \\
8 -240, AN\
ul
2 -300 \
% z : 11 -\
: -360 \\
_ . _ A4\ | 1-ONE-REGION
: _ )\A( CORE MODEL
g o o -420 , \
L . | | Nwe-ReGION /] \'\
- ‘ ) o . : . -CORE MODEL | 3
- L S .. -480 :
~ -sq0 L——L L 4 1] : i\ _ _
> 001 002  005. 01 02 05 1 2 5. 10
: e . ,-FREQUEch(ramunssed ' ‘ Coamte L
Fig. 6. Comparlson of Frequency Response of Experlmental and Theo-
ret:n.cal Reactor M:inng Characterlstlcs. _ o e
|
;

 
 

a : 16

 

dynamlc effects of the C1rculat10n of the precursors around the prlmary

loop were included in the complete model. o _ , ' i

- Heat Exchanger and Radlator -

The lumplng scheme used to represent heat transfer in both the heat
‘exchanger and the radiator 1s:shown in Fig. 7. As with the core model,
two 1umps are used to represenf each fluid flow path. The.reduced and
 intermediate models both used one sectlon as shown. The complete model

used ten of these sectlons connected sequentlally.

Fluid Transport and Heat Trsnsfér in Connecting Piping

The reduced model used;singie;well-stirredétank approximations for
-fluid transport in the piping from the core to the heat exchanger; from
 the heat exchanger to the core, from the heat exchanger to the radiator,
and from the radiator back to the heat ekchanger. Since the flow is
highly turbulent (primary system, Re m 240,000; secondary system,

Re w 120, 000), there is relatively little axial mixing, and thus a plug
flow model is probably superior to the well—stlrred—tank model. Fourth-

 

order Padé approx1mations were used in the intermediate model and pure
delays in the complete model (see Appendix A). Heat transfer to the pip-

-ing and vessels was also includedein the complete model.

Xenon Behavior

The transient poisoning effects of xenon in the core were considered

only in the complete model. The equations include iodine decay into xe-

' ’,I’

non, xenon decay and burnup, and xenon absorption into the graphite.

Delayed Power
In all three models, the delayed-gamma portlon of the nuclear power

was approximated by a first-order lag.
 

 

 

 

)

B

Fig.

 

17

ORNL-DWG 65-9814

 

 

 

_
£ IN ]

 

 

 

 

 

 

HEAT TRANSFER

 

TT Ls— TUBE

 

 

 

' HEAT TRANSFER
I DRIVING FORCE

 

 

 

r2’ OUT ——

 

T
2* IN

e e

 

 

 

 

 

 

 

 

WELL-STIRRED
TANK (TYPICAL)

(HOLDUP TIME = _’25')

Fig. 7. Model of Heat Exchanger and Radiator. (See discussion of
4 for explanation of arrows.)

 
 

 

 

18

5. Stability Analysis Results

 

Data were obtained with the best available estimates of the system
parameters for analysis by the transient-response, closed-loop frequency-
response, open-loop frequency-response, and pole-configuration (root
locus) methods. The advantages in using various analytical methods are
that (1) comparison of the results provides a means of checking for com-
putational errors; (2) some methods are more useful than others for spe-
cific purposes; for example, thé pole-configuration analysis'gives a
clear answer to the question of stability, but frequency—responSE methods
are needed to determine the physical causes of the calculated behavior;
and (3) certain methods ‘are more meaningful to a reader than others, de-
pending on his technical badkground. The differences between the results
and earlier results*’ are discussed in Section 6, and the effects of '
changes in the mathematical model and the system‘parafieters are exafiined
in Section 7.

The results show that the MSRE has satisfactory inherent stability
characteristics. Its inherent response to a perturbation at low power
is characterized by a slow return to steady state after a series of low-
frequency oscillations. This undesirable but certainly safe behavior at

low power can easily be smoothed out by the control system.2

Transient Response

 

The time response of a system to a perturbation is a useful and
easily interpreted measure of dynamic performance. It is not as fiseful
in showing the reasons behind the observed behavior as some of the other
methods discussed below, but it has the advantage of being a physically
observable (and therefore familiar) process.

The time response of the reactor power to a step change in keff'was
calculated. The IBM-7090 code MATEXP'? was used. MATEXP solves general,
nonlinear, ordinary differential equations of the form

ax

a'é'"‘ Ax + MA(x) x + £(t) , | (1)

 
 

 

n"l

19

where
X = the solution vector,
A = system matrix (constant square matrix with real coefficients),
AA(X) = a matrix whose elements are deviations from the values in A
[thus M (x) x includes all nonlinear effects],
f(t) = foreing funetion vector.

The A matrix was developed for the intermediate model and resulted in a
59 x 59 matrix. | |

The transient response of the neutron level to a step change in k off
of +O Ol% is shown in Fig. 8 for initial power levels of 10 and 1 Mw.
The slow response is evident. Figure 8 also clearly shows that the re-

actor takes longer to return to steady state in a 1-Mw transient than in

. a 10-Mw transient. It is also clear that the power does not diverge

(i.e., the system is stable).

It should be emphasized that this transient analysis included the
nonlinearities inherent in the neutron kinetics equations. The fact that
the results of the other analyses presented below, which are based on
linear models, agree in substance with these results verifies the adequacy

of the linear analyses for small perturbations.

 

Closed-Loop Frequency Response

The closed-loop transfer function is defined as the Laplace trans-
form of a selected output variable divided by the Laplace transform of
a selected input'variablé. If the system is stable, it is possible to
replace the Laplace transform varlable, s, with jw, where j =+/~1 and
w is the angular frequencyrof,a sinusoidal input. A transfer function,
G(w), evaluated at a7particular w is a complex quantity;' The afi@litude
of G(w) physically represents the gain, or the ratio of the ‘amplitude of
the output sinusoid to the amplltude of the input sinusoid. The phase

:~ang1e of G(w) represents the phase difference between the input and out-
"_put sinusoids. A 1ogar1thmlc plot of amplltude ratlo and phase angle as

a fUnctlon of w is called a Bode plot or frequencyhresponse plot.
The,relatlonshlp between the frequency response and_the time response

due to a sinusoidal input is useful conceptually and experimentally. How-

ever, it may be shown that the Bode plot for a linear system also provides

 
 

 

 

 

 

 

8N, CHANGE IN NEUTRON LEVEL (Mw)

 

 

 

 

 

 

 

 

20
ORNL-DWG 65-9815
- 0.4
03 \L ' —
02 - _
h"”a = {Mw ' . - .

01 \\\\\\\Eh : | o - -

0 —— ~ e '
=01

0 100 200 300 400
TIME (sec)

Fig. 8. MSRE Transient Response to a +0.01% &k Step Réaétifiity In-
put when Operating at 1 and 10 Mw.

9p

 
 

wl

 

21

qualitative stability information that is not restricted to any particu-
lar form in input.®? This information is largely contained in the peaks
in the amplitude ratio curve. High narrow pesks indicate that the system
is less stable than is indicated by lower broader peaks.

The closed-loop freguency response was oalculated'for N (neutron
level fluctuations in megawatts) as a function of dk (change in input
k

eff)'
response calculations; see Appendix C) and the complete model were used

‘The MSFR code (a special-purpose code for the MSRE frequency-

for this calculation. 'Thetresults for several power levels are shown in
Fig. 9. Fewer phase-angle curves than magnitude curves are shown in order
to avoid cluttering the plot. | |

Several observations .can be based on the information of Fig. 9.
First, the peaks of the magnitude curves get taller and sharper with lower
power. This indicates that the system is more oscillatory at low power.
Also the peaks in the low#power curves rise above the no-feedback curve.
This indicates that the feedback is regenerative at these power levels.
Also, since the frequency at which the magnitude ratio has a peak approxi-
mately correSponds to the frequency at which the system will naturally
oscillate in response to a disturbance, the low-power oscillations are
much lower in frequency than the 10-Mw oscillations. The periocds of os-
cillation range from 22 min at 0.1 Mw to 1.3 min at 10 Mw.

Figure 9 shows that the peak of the 10-Mw magnitude ratio curve is
very broad and indicates that any oscillation would be small and quickly
damped out. The dip in this curve at 0,25 radians/sec is due to the 25-

sec fuel c1rculat10n time in the primary loop [i.e., (2r radians/cycle)/

| d(25 sec/cycle) 0. 25 radlans/sec] Here a fuel—temperature perturbatlon

'1n the core 1s relnfOrced by the perturbation generated one cycle earller

that traveled around the loop.o Because of the negatlve fuel ~-temperature

'coeff101ent of react1v1ty, the power perturbatlon is attenuated.-

The relatlvely low galns shown at low frequenc1es can ‘be attributed

- to the large change in- steady—state ‘core temperatures that would result

'from a relatlvely small change iobnuclear power with therradlator air

flow rate remalning constant.” ThiS'meane that only & small change in
power is required to brlng about cancellation of an input 8k perturbation

by a change in the nuclear average temperatures.

 
 

22

ORNL-DWG 65-9816

 

iaN
Mfi*

 

000t 2 5 001 2 5 O 2 5 1 2 5 10 2 5 ) 100
' FREQUENCY (rodians/sec) '

90
80

70
PHASE OF -2
£ OF %,8%

60
50
40
30
20

10

PHASE {deg)
o

NO FEEDBACK

 

=90
00001 2 5 oo 2 5 001 2 5 01 2 5 i 2 5 10

_. FREQUENCY (radions/sec)

Fig. 9. Complete Model Analysis of MSRE Frequency Response fbr/Seva
eral Power Levels.

"

 

)
 

)

C

)

 

23

Open-Loop Frequency Response

A simplified block diagram representation of a reactor as a closed-
loop feedback system is shown in Fig. 10. The forward 100p, G, represents
the nuclear kinetics transfer function with no temperature effects, and
the feedback loop, H, represents the temperature (and reactivity)'changes
due to nuclear power changes. _ ,

The closed-loop equation is found by solving for N in terms of dk:
= G ok, - GHN
in

or

N G o
8k, 1 +GH ' (2)
in

The quantity:GH is called-the open-loop transfer function, and represents
the response at point b of Fig. 10 if the loop is broken at point b.
Equation (2) shows that the denominator of the closed-loop transfer func-
tion vanishes if GH = —l. Also, according to the Nyquist stability cri-
terion, the system is unstable if the phase of GH is more negative than
-180° when the magnitude ratio is unity. Thus it is clear that the open-
loop frequency response contains information about system dynasmics that
are important in stability analyses.

A useful measure of system.stability is %he ‘phase margin. It is de-
fined as the dlfference between 180° and the open-loop phase angle when -

the gain is 1.0. A dlscu351on of the phase margln and its uses may be

-found in sultable references on servomechanlsm,theory. 13’ For this appli-
- catlon, it suffices to’ note that smaller phase marglns 1nd1cate reduced
1gstdb111ty. A general rule of thumb in automatlc control practice is that

:a phase margln of at least 30° is de51rable Phase marglns less than 20°
':1nd1cate llghtly damped osczllatlon and poor control._ Zero degrees in-
f.dlcates an osc1llat1ng system, and negatlve phase marglns 1nd1cate 1nsta—

f*blllty.

= The phase margln as a functlon of reactor power level is shown in

"Flg. 11. . The phase margin decreases ‘as the power decreases and goes below

30° at about 0.5 Mw. However, the phase margin is still positive (12°)
at 0.1 Mw. These small phase margins at low power suggest slowly damped

 
 

 

 

 

 

 

 

 

 

24
ORNL-DWG 65-9817
o ~__FORWARD LOOP -
4 3K NET s
| - ' NUCLEAR POWER, Mw
BREAK LOOP HERE
FOR OPEN-=LOOP
ANALYSIS——a.),
H

 

 

 

3k FEEDBACK
- N FEFDBACK LOOP

Fig. 10. Reactor as a Closed~Loop Feedback System.

i ORNL-DWG 65-9818
; 100

PHASE MARGIN

20

10D OF

 

005 01 02 o5 { T2 5 10 20
N, NUCLEAR POWER (Mw)

Fig. 11. Period of Oscillation and Phase Margin Versus Power Level
for MSFR Calculation with Complete Model. |

ol

 
 

 

a)

w)

 

25

oscillations, as has been observed in the transient-response and closed-
loop frequency-respohse results. The period of oscillation as a function
of power level is also shown in Fig. 11.

Figure 12 shows Nyquist plots for the open-loop transfer function,
GH, at 1.0 Mw and 10.0 Mw. It is clear that the unstable condition of
(IGH‘ 1 and phase angle —=180°) is avoided. In order for the phase mar-
gin to be a reliasble 1nd1cation of stability, the Nyquist plot must be
"well-mannered" inside the unit circle; that is, it should not approach
the critical (~1.0 +jO) point. Although the curves shown in Fig. 12 be-
have peculiarly in approachlng the origin, they do not get close to the

critical point.

Pole Configuration

The denominator (l-+ GH) of the closed-loop transfer function of a
lumped-parameter system is a polynomial in the_Laplace transform variable,
s. The roots of this polynomial (called the characteristic polynomial)
are the poles of the sYStem_trahsfer function. The poles are equal to
the eigenvalues of the system matrix A in Eq. (1). A necessary and suf-
ficient condition for linear'stability is that the poles all have nega-
tive real parts. Thus, it is useful to know the location of the poles in
the complex plane and the dependence of their location on power level.

The poles were calculated for the intermediate model of the MSRE
(see Section 4). The matrix used was the linearized version of the
59 X 59 matrix used in the:frahsienf analysis. The calculations were

performed with a mbdificatibn’bf'the general matrix eigenvalfie code QRl4

~ for the IBM=7090. The results are shown in Flg. 13 for several of the
”,fmajor poles. A1l the other poles lie far to the left of the ones shown
- It is clear fromeFlg 13 that the system is. stable at all power levels.

The set of complex poles that goes to Z€ro as the power decreases is the

' ’jfset prlmarlly respons1b1e for the calculated dynamlc performance- The

imaginary part of thls set approx1mately represents the- natural frequency

of osc111at10n of the system follow1ng a dlsturbance.. The frequency of

0501llation decreases as. the power decreases, as observed prev1ously.

 
ORNL-DWG 65-9820
IMAGINARY

0 02 04 06 REAL
o2

 
 
  
 
 
 
  

+-02

.y

- -0.6

 

Ng = 1OMw +-16

IMAGINARY
0 REAL

w* 0.00013—ayg

+-10000
Ny 10 Mw/{ [

~ Fig. 12 .' Nyquist Diagrains for Complete Model
No = 1 and 10 Mw.

 

IMAGINARY PART OF POLE {sec!)

0.07
0.06
0.05
0.04
0.03
0.02 .
10 Mw

0.0 !0 NEGLIGIBLE POWER
ol 01 DE PENDENCE

-0.01
-0.02
-0.03
-0.04
-0.05
-0.06

-0.07

Fig. 13. Major Poles for the MSRE.

*-—0
-0.07 -0.06 -005 ~004 -0.03 -0.02 -001

 

92

 
 

 

 

 

 

27

6. Interpretation of Results

Explanation of the Inherent Stability Characteristics

A physical explanation of the predicted stability characteristics
is presented in this section, and an'attem@t is made to explain the rea-
sons for the changes in inherent stability with power level. The reasons
for the behavior are not intuitively obvious. Typically a feedback system
will become more stable when the open-loop gain is reduced. The MSRE,
however, becomes less stable at lower powers. In the discussion of open-
loop frequency response.(Seet. 5) it was noted that the forward-loop
transfer function G represefits the nuclear kinetics (N/8k) with no
temperature-feedback effects, and from the equations (Appendix A), the
gain of G is directly propertional to power level. Thus the MSRE is not
"typical" but has the characteristics of what is known as a "conditionally
stable" system.l’

The MSRE analysis isicomplicated fuither by the complexity of the
feedback loop H, which represents the reactivity effects due to fuel and
graphite temperature changes resulting from changes in nuclear power. A
more detailed breakdown of the components of H is given in Fig. 14. This
core thermal model has two inputs, the nuclear power N and core inlet tem-
perature T 3’ and three outputs, nuclear average fuel and graphite tem-

peratures T? and T¥, and the core outlet temperature Tco' The block "Ex-

2
ternal Loops" reprgsents the primary loop external to the core, the heat
exchanger, the seCOndary?lbop,tand the radiator. All the parameters are
treated as perturbatlon quantltles or dev1at10ns from their steady-state
values. Also the radlator alr flow rate 1s adgusted S0 that with a given
steady-state power level, the core cutlet temperature is 1225 F. This
means that the feedback loop transfer function H also varies W1th power
level., | o _t , | |
If we 1ook at the effect of perturbatlons in power, N, on the core
1nlet temperature, Tc ,rwe can see that the effects of dlfferent air flcw
rates are only apparent at 10w frequenc1es, as in Fig. 15 which Shows
the Bode plots for T /N at Nb 1 and 10 Mw. It is important to note
that at low power and at low frequency, the magnitude of the temperature

change is large, and it lags the input N considerably. For example,

 
 

 

ORNL-DWG 65-9822

 

 

 

 

 

 

 

 

 

 

 

 

 

OPEN LOOP N
Nrgk NUCLEAR POWER, Mw
N a = TEMPERATURE
- COEFFICIENT OF
REACTIVITY
L_Ts, CORE F = FUEL
CORE OUTLET THERMAL | CORE INLET
TEMPERATURE, oo MODEL. TEMPERATURE, T, G = GRAPHITE
- x = NUCLEAR AVERAGE
EXTERNAL
== LOOPS -

 

 

 

 

 

 

 

Fig. 14. MSRE as'Closed-Loop Control System.

A a

¥ &

1000

Tei/N, RATIO OF CORE INLET TEMPERATURE
TO NUCLEAR POWER (°F/Mw)

 

08 2 5 404 2 5 403 2 5 g2 2 5 o

8¢

PHASE (deg)
®
o

=100
=120

-140

 

-160
105 2 5 o4 2 5 103 2 5 102 2 5 4o

| FREQUENCY (radians/sec}
Fig. 15. Closed-Loop Frequency-Response

Diagram of Core Inlet Temperature as a Function
of Nuclear Power.

 
 

 

)

 

29

at Np = 1 Mv and w = 0.0005 radians/sec, the magnitude ratio is 170°F/Mw
and the phase lag is 80°. The block diagram of Fig. 14 indicates that

the nuclear average temperature perturbations in T; and Té can be con-
sidered to be caused by the two separate inputs N and Tci' For example,
the open-loop transfer function Tf/N (withTci constant) is 5.0°F/MW at
steady state, and there is little attenuation and phase lag up to rela-
tively high frequencies, as in Fig. 16, which shows the open-loop transfer
functions of the.nuclear average temperatures as functions of N and Tci'
Returning to the example case of Ng = 1 Mv and w = 0.0005 radians/sec,

we note that the prompt feedback effect of 5 F/MW from TF/N is very small
compared with temperature changes of the entire system represented by a
Tci/l\! of 170°F/Mw at —80°. (uote that TF/TCi = 1.0 at 0.0005 radians/sec.)
The important point here is that for low power levels over a wide range

of low frequencies, the large gain of the frequency response of overall
system temperature relative to power dominates the feedback loop H, and
its phase angle approaches —90°. | '

The no-feedback curve;in Fig. 9 shows that at frequencies below about
0.005 radians/sec, the_forward-loop transfer function N/8k (open loop)
also has a phase approaching —90° and a gain curve with a —1 slope (i.e.,
one-decade attenuation per”decade increase in freQuency). With both G
and H hauing phase angles approaching —90°, the phase of the product GH
will approach —180°. If the magnitude ratio of G were such that [GH| = 1.0
under these conditions, the system would approach instability. From the
Bode plot of Fig. 9, it can be seen that at a power of 0.1 Mw, IGH\ 1.0

at- 0. 0045 radlans/sec (22 mln/cycle), since the peak in the closed-loop
’.occcurs there. At lower powers and. consequently lower galns G, lgH | ap-

aproaches 1.0 at even Lower frequenc1es, where the phase of GH is closer
. to —180°. ThlS accounts for the 1ess stable conditions at the lower

L*pcwers and lower frequencles.,

At the higher powers, IGHI approaches 1. O at hlgher frequencles

where the effect of the prompt thermal feedback is s1gn1f1cant For ex-

“ample, the peak in. the 10-Mw closed-loop Bode plot of N/fik Fig. 9, occurs
‘at O.CW8.rad1ans/sec. At thls frequency, IT /NI has & value of 2.0

(Fig. 15) compared with a T%/N_of 4ot at —15° (Fig. 16). Consequently,
the prompt fuel temperature feedback term has a dominant stabilizing effect.

 
 

30

ORNL-DWG 65-9824

 

MAGNITUDE RATIO

 

 

 

 

 

 

 

 

 

PHASE (deq)

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

-105 \W v
-120 S
-135 '
000t 0002 0005 OOf 002 - 005 O 02 . 08 - 1.

FREQUfiNCY _(rodiuns/sec)

Fig. 16. Open-Loop Frequency-Résponse Diagrams of Nuclear Average
Temperatures as Functions of Nuclear Power and Core Inlet Temperature.

-

 
 

 

 

31

The relative importances of the various components of feedback re-
activity are shown in the directed-line diagrams of Fig. 17 for power
levels of 1 and 10 Mw and at the frequencies at which the oscillations
occur. The vectors labeleqz—fikF and —6kG represent the products of the
nuclear temperature components and the reactivity coefficients that re-
sult from a unit vector inpfit 8k, -+ The vectors labeled "prompt" are
the effects due to the nuclear power input based on no change in core
inlet temperature. Those labeled "loop" are caused by changes in core

inlet temperature alone. For example:

 

 

BkF prompt N T*
= £ 8k/°F)
5K = B N Oy
in in | closed loop open loop
and
- T*
5kG loop N Tci G /o
ok, Bk B T %, (8k/°F)
in in | closed loop closed loop “ci |open loop

The net 8k vector is the sum of the input and feedback vectors. For the
1-Mw case, ok net is greater than Bkin; this indicates regenerative feed-
back and shows up on the closed-loop Bode plot (Fig. 9) as a peak with a
greater magnitude ratio than that of the no-feedback curve.

The increased steblllzlng effect of the prompt fuel temperature term
in going from 1 to 10 Mw is also evident. These plots clearly show the
diminished effect of the graphite at the higher frequency.

_ In both cases, too, the plots show that a more negatlve graphlte

temperature coefflclent would tend to increase the net Bk vector and

‘Thence destablllze the system.

'Interpretatlon of Early Results

_ The prev1ously published results of 8, dynamic study performed in
1961 predlcted that the MSRE would be less stable than 1s predicted in

thls study This is partly'because of dlfferences 1n design parameters

-~ and partly ‘because of dlfferences in the models used. These differences

were reviewed in Section 4 of this report.

 
32

ORNL-DWG €5-9825

Py = 1OMw ,
" w= 00778 RADIANS/SEC

T g8k NET

\ '-Skin

/\’_-—w | |
~Bk. TOTAL '
B - -8k, PROMPT

~8k LOOP

I}Sk-NET

Fp= tMw
w = 0.0145 RADIANS/SEC

-8k PROMPT

 

} o 8“'in‘

~84g LOOP ~Bk; PROMPT

Fig. 17. Directed-Line Diagrams of Reactivity Feedback Effects

1 and 10 Mw for Complete Model.

at

o

 
 

 

33

The most significant parsmeter changes from 1961 to 1965 were the
values for the fuel and graphlte temperature coefficients of react1v1ty
and the changes in the fuel heat capacity. Table 1 (Sect 4) shows that
the new fuel coeff1c1ent is more negative, the new graphlte coefficient
is less negatlve, and the fuel heat capacity is smaller than was thought
to be the case in 1961 ‘All these changes contribute to the more stable

'behav1or calculated w1th the current data.. (The destablllzlng effect of

a more negative graphlte temperature coefficient is discussed in Section
7.) - |
The most important change, however, is the use of a multiregion core

model and the calculation of the nuclear average temperature. In the

| single-region'core, T¥ is equal to the temperature of the first of the

_ F Pe:
two fuel lumps or the average core temperature (in steady state). In

the nine-region core,.T* is_computed by multiplying'each of the 18 fuel-

| F | |
lump temperatures by a nuclear importance factor, I. In the single-

region model, the steady-state value of TF/N (with T 5 constant) is only
2.8°F/Mv compered w1th 5,0 F/Mw for the nine-region core model. This
dlfference occurs because in the nine-region model, the downstream fuel-
lump temperatures are affected not only by the power 1nput to that lump
but also by the change 1n the lump's inlet temperature due to heating of
the upstream lumps. This point may be illustrated by noting the differ-
ence between two'single-region models, where in one case the nuclear im-
portance of the-first lumptis 1.0 and in'the other case the importance

of each 1ump'is-015 As an example, say the core outlet temperature in-

_creases 5 F/Mw The change 1n T* for a l- Mw 1nput would be .

AT%," Il ATJ_ + Ig ATZ o

':Intheifirst'case L = 1.0 and.Aflh 2.5°F, 50 AT% = 2.5°F. In the sec-
’ond'case;”Ii = Ip = 0. 5 A$1 = 2. 5 F, and AT '5°F 80 AT; = 3, 75°F, or

‘a 50% greater change than in case ‘one. For many more lumps, thls effect

1s even greater.-

' As was shown above, the prompt fuel react1V1ty feedback effect was
the dominant stablllzlng mechanlsm at both 1 and 10 Mw, so the original

single-region core model_would give pessimistic results.

 
 

34

- 7. Perturbations in the Model and the Design Parameters

Every mathematical analysis of a physical system is subject to some
uncertainty because of two quéStiOns:' How good isrthe mathematical model,
~ and how’accurate_are the values oftheaparameters.used in the model? The
- influence of changes in the assumed model were therefore inrestigatea,
and the sen51t1v1ty to parameter varlations of the results based on both
the reduced and conmplete models was determlned An ana1y31s was also
performed to determine the worst expected stablllty performance W1th1n

the estimated range of uncertainty of the system parameters.

Effects of Model Changes

 

The erfects of changlng the mathematlcal model of the system.were
determlned by comparlng the phase marglns with the reference case as
each part of the model was changed in turn. The follow1ng changes were
made: o | ' |
| 1. Core Representation. A single-region core model was used in-
stead of the nine-region core used in the complete model. |

2. Delayed-Neutron Groups. A single delayed-neutron group was
used instead of the six-group representation in the complete model.

3. Fixed Effective B's. The usual constant-delay-fraction delayed-

neutron equations were used with an effective delay fractibn, B, for each |

precursor. The effective B was obtained byycalculating the delayed-
neutron contribution that is reduced due to fluid flow in'the steady-
state case. This is in contrast to the explicit %reatment of dynamic
circulation effects in the complete model (see Appendix A). |

4, First-Order Transport Lags. The Laplace transform of a pure
delay, e '°, was used in the complete model. The first-order well-
stirred-tank approximation, 37%3;53 was used in the modified model.

5. Single-Section Heat Exchanger and Rediator. A single section
was used to represent the heat exchanger and radiator rather than the
ten-section representation in the complete model. | |

- 6. Xenon. The xenon equations were omitted in contrast to the ex-

plicit xenon treatment in the complete model.

 
 

e

 

»n

 

35

The results are shown in Table 4. The only significant effect is
that due to a change in the core model. The results for the one-region
core model indicate considerably less stability than for the nine-region
core modei. This difference is due primarily to the different way in
which the nuclear average temperature of the fuel is calculated, as was

discussed in detail in Section 6.

Table 4. Effects of Model Changes on Stability

 

Phase Margin Change in Phase Margin

 

 

(aeg) from Reference Case®
Model Chenge 7 g | (deg)
At 1 Mw At 10 Mw At 1 Mw At 10 Mw
Reference case (complete' 41 99
~ model) [
One core region o 27.5 0 56.5 -13.5 ~2.5
One delayed-neutron 38 () -3
group :
Fixed effective B's 40 o8 -1 -1
First~-order transport 41 100 0 +1
lags _
Single-section heat 4] 98.5 0 —0.5
exchanger and radiator '
Xenon omitted - 41 100.5 0 +1.5

 

Reference case is complete model with current data. An increase

in phase margln 1nd1cates greater stablllty.

bNyqulst plot does not cross unlt circle near frequency of oscil-

'~lat10n.;

'Effects-of'Paramfiter,changesé

B Frequency-response sen31t1v1t1es and pole sen51t1v1t1es were cal-

.f ulated Frequencyhresponse sens1t1V1t1es are deflned as fractlonal

'-;e:changes in magnltude or" phase per fractlonal change in a system parameter.

The_magnltude frequencyeresponse sen31t1V1t1es were calculated_for Sev-
eral important parameters with the MSFR code (see Appendix C) for the

complete model. The senéitivitieS'were obtained by differences between

 

 
 

 

36

the results of calculations at the design point and those of calculations
with a single parameter changed slightly. The results of these calcu-
- lations are shown in Fig. 18. Calculations were also performed on the
system with the reduced model with a new computer code calied SFR (Sensi-
tivity of the Frequency Response). This code calculates magnitude and
phase sensitivities for a general system by matrix methods. This calcu-
lation was restricted to the reduced syétem.representation because of
the large cost of calculations for very large matrices. The results of
this calculation are given in Fig. 19. In Figs. 18 and 19, a positive
sensitivity indicatés that-a dééreasé_in the sysfiem parameter will de- -
crease the magnitude of the ffequency response. The situation is re-
versed for negative‘sensitivitiés.r |

'The sensitivities shown in Figs. 18 and 19 can be used to estimate
'the'effects 6f posSible futfiié updating of the MSRE design parametefs
used in this study. 1In addifiion, they support other general observations
obtained by other means. For instance, Fig. 19 shows that the sensitivi-
ties to loop effects, such as primary and secondary loop transit times,
are important relative to core effects. This indicates that the external
loops strongly influence the system dynamics, as was concluded in Sec-
tion 6.

Similar information may be obtained from pole sensitivities (or
eigenvalue sensitivities). These are defined as the fractional change
-of a system pole due to a fractional change in a system parameter. The
sensitivity of the dominant pole (the pole whose position in the complex
plane determines the main characteristics of the dynamic béhavior) is
usually the only one of interest.

The dominant pole sensitivities for a nunber of system parameters
are shown in Table 5 for power levels of 1 and 10 Mw. These results may
be used to estimate the effect of future updating of the MSRE design
~ parameters, and they also furnish some insight as to the causes of the
calculated dynamic behavior. For instance, it is noted that the sensi-
tivity to changes in the graphite temperature coefficient is only about
'one-fourteenth the sensitivity to changes ifi the fuel temperature coef-
ficient at 10 Mw. At 1 Mw, the graphite effect is about one-third as (:}
large as the fuel effect. This indicates that a decrease in power level

e e e —T— e

 
 

»

 

Wy

 

RATIO OF FRACTIONAL CHANGE IN AMPLITUDE TO
FRACTIONAL CHANGE IN PARAMETER

 

 

37

ORNL-DWG 65-9826
POWER =10

FUEL TEMPERATURE COEFFICIENT
OF REACTIVITY

GRAPHITE TEMPERATURE CQEFFICIENT
OF REACTIVITY

FUEL HEAT CAPACITY

4
Q.00  0.002 005 oo 002 0.05 ol 02 05 1 2 5 10

FREQUENCY (radians/sec)

Fig, 18. Frequency-Response Sensitivities of Complete Model.

 
 

38

ORNL-DWG 65-9827

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

05
0 1T i
' "]
_o5 |SECONDARY LOOP TRANSIT TIME
05 -
_ . i
FUEL-SALT HEAT EXCHANGER
HEAT TRANSFER COEFFICIENT
o -05 S 2
e 05 .
ul -
W : B : SECONDARY-SALT RADIATOR ‘
d'“_‘ 'GlflAIlDHET:I'l;EMPERATURE COEFFICIENT HEAT TR\I(\NSSFER COEFF(I)CIENT ' .
£5 05 L1l L L LTE -05 ' :
ag 15 0
2__’& . T ITY -‘“'\\
gg .0 -05 "\L,
Zw CORE RESIDENCE TIME FOR FUEL SALT
9 1.0 o el : C s,
2z 1.0
°%
_lo //
W 05 |- paa
53 0 .
Qe N
Q9 s GRAPHITE HEAT CAPACITY o5 L_PRIMARY EXTERNAL LOOP TRANSIT TIME
wE 05 T 0.5
e
2 FRACTION OF POWER \ [ L—T SECONDARY-SALT HEAT RS ey
o« 0.5 L GENERATED IN THE FUEL N1 05 L EXCHANGER HEAT TRANSFER COEFFICIENT
0.5 FGRAPHITE-TO-FUEL HEAT 05
o |_TRANSFER COEFFICIENT o T .
T , MASS FLOW RATE TIMES| 1T
SPECIFIC HEAT OF SECONDARY SALT
-0-5 _0.5 1 i L bbb 31 1 1 i1l
0001 2 5 oot 2 5 o1 2 0.001 2 5 oo01 2 5 ot 2
»
FREQUENCY (radians/sec) FREQUENCY (radians/sec)
Fig. 19. Frequency-Response Sensitivities of Reduced Model at a ,

Power of 10 Mw Based on Current Data.

 

 
 

 

 

39

Table 5. Pole Sensitivities

 

x| oA

 

 

| | R
Parameter = x (see footnote a)
At 10Me At 1 Mw
Fraction of power that is prompt —0.944 -2.515
Neutron lifetime 0.00944 0.0129
Fuel temperature coeffieient of reactivity 0.858 1.701
Graphlte temperature eoefflclent of ~0.0627 —0.493
react1v1ty - :
Fraction of power generated in the fuel -0.7328 0.561
Graphite-to-fuel heat transfer coefficient 0.0434 0.177
Fuel heat capacity 1.024 1.315
Moderator heat capacity ~0.616 -0.359
‘Fuel-salt heat exchanger heat transfer 0.0157 0.254
coefficient o _
Fuel transit time in core -0.606 —0.787
Fuel transit time in external primary 0.659 0.804
circuit
Secondary-salt heat exchanger heat transfer 0.00708 0.449
coeff1c1ent
Secondary-salt loop transit tlme —0.305 —0.622
Secondary-salt radlator heat transfer —0.0155 —-0.0754
- coefficient = - S o o . S
' Heat exchanger heat capa01ty | :HO;OO745 b ~0.009%69
"Effeetlve precursor decay constant I e—0;304 ;,-,,—O 726
_ ”:fleme constant for delayed-gamma em1s51on ;—0.00858,j1 0. 0536
S ,;Total delayed-neutron fraction , '0.0103_” 0.159
' Q;EEffectlve delayed-neutron fractlon 1f0;788.-'_}50.221

 

U o8 is the real part;of the dominant pole.-
,_f0.01865_for;lO]Mwfiand;é03001818_for 1 Mw. -

The values are

 
 

40

causes modifications in the dynamic behavior that accentuate the relative
effect of graphite temperature fgg@bgck.fl It is also noted that the vari-
ous heat transfer coefficients have a much léréér relative effect at 1 Mw
than at 10 Mw; this indicates that the coupling between system components

has a larger influence at low power than at high power.

Effects of Design Uncertainties

 

A new methodl® for automatically finding the least_stafile condition
in the range of uncertainty in the design parameters was used. A com-
putef code for the IBM-7090 was used for the calculations. The method
is described in some detail in Appendix D. Thé least steble condition
is found by using the steepest-ascent (or gradlent-progectlon) method of
nonlinear prograrming. . _ ,

The quantity that is maximized is the real part of the dominant pole
of the system transfer function (or equivalently, the dominant eigenvalue
of the system matrix). Less stable conditions are accompanied by less
negative values for the real part of the dominant pole, and instability
is accompanied by a pole with a positive real part. The maximization
involves a step-wise determination of the particular combination of sys-
tem parameters within the uncertainty range that causes the réai'part of
the dominant pole to have its least negative value. If the maximized -
pole has a negative real part, instability is not possible in the uncer-
tainty range. If the maximized pole has a positive real part, instability
is possible in the uncertainty range if all the system parameters deviate
from the design point in a particular way. |

A key factor in the stability extrema analysis is the availability
of the appropriate ranges to assign to the system parameters. The ranges
appropriate for the MSRE were furnished by Engel.17 It was decided to
use a‘wide range on the important nuclear parameters (neutron lifetime,
fuel temperature coefficient of reactivity, and graphite fempérature co=
efficient of reactivity). These parameters were sllowed to range between
two-thirds and three-halves of the design values. All other ranges were
assigned by cons1der1ng the method of evaluating them and the probable
effects of aging in-the reactor environment. The ranges of the 16 system

parameters chosen for this study are given in T&ble76.

 
 

‘i

 

41

 

 

 

fEffectxve delayed—neutrone fe'

fractlon

Table 6. Ranges on System Parameters for Extrems Calculations
: Ranges
. Parameter ‘ -
Low Design Value High
Neutron lifetime, sec. 1.6 X 1074 2.4 x 1074 3.6 x 1074
Fuel temperature coefficient, -7.26 X 10°° —4.84 x 10°°. =3.23 x 10°°
8k/k- °F o |
Grephite temperature coeffa,- —5.55 x 1077 ~3.70 x 107° =2.47 x 10°°
cient, 8k/k-°F - | |
Fraction of power generated - 0.92 0.9335 0.95
in fuel o '
Grephite-to-fuel heat trens-v 0.013 0.02 0.03
fer coefficient, Mw/°F :
Fuel heat capacity; Mw-sec[fF 1.13 1.50 1.910
Graphite heat capac1ty, 3.4 3.58 3.76
Mw-sec/°F | |
Fuel-salt heat exchenger heat 0.1613 0.2424 0.3636
transfer coefficient, Mw/°F
Fuel transit time in core, = 6.96 . 8.46 10.25
sec : B
Fuel transit time in external 15.75 17.03 18.60
prlmary circuit, sec |
Secondary~salt heat exchanger 0.1001 0.129 0.1686
heat transfer coefficient, - | |
Mw/F |
._Secondary 1oop transit time,ire 2l.7 24 .2 32.7
. sec , o s - -
--Heat exchanger heat capac1ty,i”10,0738  ; ) 0.0738 -+ 0.4216
 Mw-sec/°F s St e T e
' Effective precursor decay o 0d11 0.133 . 0.15
~constant, sec” -1 _W__,,+-5 DR | . :
Time constant for delayed 120 188 270
gamma, emission, sec T S e . S '
0.0032 0.0036 =  0.0040

 

 
42

The reduced model was used with current parameters for locating the
‘least stable condition in the uncertainty range. This gave results at
a much lower cost than with a more complete nmodel. This was considered
adequate because other calculations showed that the reduced modei_predicts
lower stability than the complete model. The reasons for this were ex-
ploréd in Section 6. Experience with other calculations also shdwed:that
changes in system parameters gave qualitatively the same type of changes
in the system performance with either model. | | |

The set of parameters for the least stable condition is listed in B
Table 7. The least negative value of the dominant eigenvalue calculated
with the reduced model changes from.(u0.0187 + 0.0474 j) sec! at the de-
sign point to (—~0.00460 + 0.0330 j) sec™! at the worst condition for 10
Mw. For 1 Mw, the change is from (—0.00182 # 0.0153 j) sec™® to
(+0.000574 + 0.0134 j) sec. fThis indicates that instability is im-
possible in the uncertainty range at 10 Mw but that the reduced model |
predicts an instability at 1 Mw for a combination of parameters within
the uncertainty range. This condition gives a transient with a doubling
time of about 1/2 hr and a period of oscillation of about 8 min per cycle.

It is evident that the calculated instability at the extreme case
for 1 Mw is due to the inherent pessimism of the reduced model"(sge_Sect.
6). This was verified by using the MSFR code for the complete model with
the pafameters describiné the extreme case. It was found thét the phase
margin for 10 Mw was 75° for the extreme condition (vefSus'99° for the
design condition), and the phase mfirgin for 1 Mw was 21° for the extreme
condition (Yersus 41° for the design condition). Thus, it is concluded
that the best available methods indicate that the MSRE will be stable
throughout the expected range of system parameters. |

8. Conclusions

This study indicates that the MSRE will be inherently stsble for all
operating conditions. Low-power transients without control will persist
for a long time, but they will eventually die out because of inherent

feedback. Other studies have shown that this sluggish response at low

 
 

 

43

 

‘Table 7. Values of System Parameters at the
Least Stable Condition
Parameter At 1 Mw

At 10 Mw

 

Neutron lifetime, sec 3.6 X 1074 (H)® 3.6 x 107% (H)

 

 

Fuel temperature coefflcient
8k/k°F ,

Graphite temperature coefficient,
8k/k- °F

~3.23 x 10~% ()

~5.55 x 1073 (L)

- =3.23 x 1077 (H)

~5.55 x 10> (L)

Fraction of power generated in fuel 0.95 (H) 0.95 (H)
Graphite-to-fuel heat transfer 0.03 () 0.02535
coefficient, Mw/°F
Fuel heat capacity, Mw.sec/°F 1.91 (H) 1.91 (H)
Graphite heat capacity, Mw.sec/°F 3.40 (L) 3.40 (L)
Fuel-salt heat exchanger heat - 0.3636 (H) 0.3636 (H)
transfer coefficient, Mw/°F |
Fuel transit time in core, sec 6.9 (L) 6.9 (L.)
Fuel transit time in external 10.25 (H) 10.25 (H)
primary circuit, sec
Secondary-salt heat exchanger heat 0.1686 (H) 0.1686 (H)
transfer coefficient, Mw/°F
Secondary loop transit time, sec 2L.7 (1) 21.7 (L)
"Heat exchanger heat capacity, 0.4216 (H) 0.0738 (L)
Mw-sec/°F
Effective precursor decay constant, 0.11 - (L) 0.15 (H)
sec™ , : : ‘
. Time constant for delayed gamma -~ 120.0 (L)  120.0 (L)
. emission, sec . - S |
o Effective delayed-neutron fractiOnf - 0.0040 | (H)_t770.0040' '(H)

 

®Letters in parentheses indlcate whether parameters are at high

values (H) or low values (L)

 
44

power can be eliminated by the control system, which suppfesses tran-
sients rapidly. | | |

The theoretical treatment used in this study included all known

effects that were considered to be capable of significfintly,influencing.

the system dynemics. BEven so, for safety and also for obtaining basic

reactor information, system stability will be checked experimentally in

dynamic tests, which will begin with zero power and which will continue
- through full-power operation. ‘

1]

 
 

 

10.

11,

12.

 

45

References

R. C.VRobertson, MSREJDesign and Operations Report, Part 1: Descrifi—
tion of Reactor besigfi,fiSAEC-Report ORNL-TM-728, Oak Ridge National
Laboratory, January 1965.

Oak Ridge NationalnL&boratory, MSRE Semiann. Progr. Rept. July 31,
1964, pp. 135-140, USAEC Report ORNL-3708.

P. N. Haubenreich et al., MSRE Design and Operations Report, Part 3:
Nuclear Analysis, p. 156, USAEC Report ORNL-TM-730, Oak Ridge Na-
tional Laboratory, Feb; 3, 1963.

. O. W. Burke, MSRE — Preliminary Analog Computer Study: Flow Accident

in Primary System, Oak Ridge National Laboratory unpublished internal
document, June 1960,

0. W. Burke, MSRE — Analog Computer Simulation of a JT.oss of Flow
Accident in the Secondary System and a Simulation of a Controller
Used to Hold the Reactor Power Constant at Low Power Levels, Oak
Ridge National Laboratory unpublished internal document, November
1960. S

0. W. Burke, MSRE — An Analog Computer Simulation of the System for
Various Conditions — Progress Report No. 1, Oak Ridge National Lab-
oratory unpublished internal document, March 196l.

0. W. Burke, MSRE ~ Analog Computer Simulation of the System with a
Servo Controller, Oak Ridge National Laboratory unpublished internal
document, December 1961. . |

S. J. Ball, MSRE Rod Controller Similation, Oak Rldge National Lab-

oratory qnpubllshed'lnternal document, September 1962.

S. J. Ball, Simulatibn'df Plug-Flow Systems, Instruments and Control

_S stems, 36(2): 133-140 (Fébruary 1963).

Oak Ridge Natlonal Laboratory; MSRE Semiann. Progr. Rept. July 31,
1964, pp. 167*174 USAEC Report ORNL-3708.

H. A. Samulon, Spectrum.Analy51s of Tran31ent ResPonse Curves,

Proct. I.R.E., 39: 175-186 (February 1951). |
‘S. J. Ball and R.'K;,Adams, MATEXP — A General Purpose_Digital Coms

puter Code for Solving:Noniinear, Ordinary Differential Equations

by the Matrix Exponential Method, report in preparation.

 
 

13.

14.

15.

16.

17.

18.

19.

20.

46

J. E. Gibson, Nonlinear Automatic Control, Chapt. 1, McGraw-Hill,
N.Y., 1963. ”

F. P. Imad and J. E. Van Ness, Eigenvalues by the GR Transform,
SHARE-1578, Share Distribution Agency, IBM Program Dzstribution,-
White Plains, N.Y., October 1963. L .

H. Chestnut and R. W. Mayer, Servomechanlsms and Regulating System .

 

Design, Vol. I, p. 346, Wiley, N.Y., 1951 , ,
T. W. Kerlin and J. L. Lucius, Stabillty Extrema in Nuclear Power

 

Systems with Design Uncertainties, Trans. Am. Micl. Soc.,,7(2).
221-222 (Novenmber 1964). - - » |

J. R. Engel, Oak Ridge National L&boratory, personal communlcatlon

to T. W. Kerlin, Oak Ridge National Leboratory, Jan. 7, 1965.

S. J. Ball, Approximate Models for DlstributediParameter Heat—TTansfer
Systems, ISA Trens. ., 3: 3847 (January 1964)

B. Parlett, Applicétions of Laguerre's Method to the Matrix Eigen-
value Problem, Armed Forces Services Technical Information Agency
Report AD-276-159, Stanford University, May 17, 1962.

F. E. Hohn, Elementary Matrix Algebra, MacMillen, N.Y., 1958.

 

%

 
APPENDICES

 

 

 

 
 

 

Ean
 

 

49

Appendix A. Model Equations

Core Thermel Dynamics Equations

The differential thermal dynamics equations for a single-core region
are given below (see Fig. 5).
First fuel lump:

dT 1 K hA
g = £ Gl = =
== (T, .. =T, ) +eiP 4 — (T, —-T_) . (A.1)
dat Try Fl,in F1 MCpl T KC_l + Kéz MCpl G F1
Second fuel lump:
4T 1 K hA
F2 = = % G2 - =
— == (T —~T,) + 55— P + e oe— (T, ~ T ) . (A.2)
dat Tgo F1 F2 MCPZ T KGl + KG2 MCP2 G Fl
Graphite lump:
::G = MEA. (Tfl "EE) + KG;C+ o2 Fp - (a.3)
PG PG
In these equations,
EFl = mean fuel temperature in first well-stirred tank, or
lump, °F,
t = time, sec,
Tg = transit (or holdup) time for fuel in first lump, sec,
Tr 4y = inlet fuel temperature to first lump, °F,
) _
' Kl = fraction of. total power generated in first fuel lump,
qul-= heat capac1ty of flrst lump, MW-sec/°F,
P& = total power, M,
Kfil = fractlon of total power generated in graphite adjacent to

- first fuel lump,

K =_fract10n:ofitotalwpowér generated in graphite adjacent to
secdndffiel-lg@p,

mean heat trensfer coefficient times area for fuel-to-

graphite heat transfer, Mw/°F,

mean graphite temperature in section, °F,

Qral
i

 

 
 

50

i&z = mean fuel temperature in second lump, °F,
Tpp = transit time for fuel in second lump, sec,

Ké = fraction of total power generated in second lump,
MQPZ = heat capacity of second lump, Mwesec/°F,

MC o = heat capacity of graphite, Mw-sec/°F.

Nuclear importances:

3k = | -
s @4
5=\akzfi ~ . (A.5)
2" 2 or, e 0

| . 3k = . ‘
Bk, = IGB-TEATG, | (A.6)

where
Skl 2.a°= changes in effective reactor multiplicetion due to tem-
272

perature change in fuel lumps 1 and 2 and the graphite,
respectively,

IFl,F2,G = importance factors for fuel lumps 1 and 2 and the
graphite, respectively; note that |

nine sections IFl F2 :

> (I,) = 1.0,

nine sections

g-]%- E o = total fuel temperature coefficient of reactivity,

F
F 5k/k- °F,

%%}-Ecx = total graphite temperature coefficient of reactivity,

G G
8k/k- °F. |
In the nine-region core model, the individual regions are combined

as shown in Fig. 5. The nuclear average fuel and graphite temperatures,

u.

 
 

 

51

reactivity feedback, and core outlet temperatures are computed as func-
tions of nuclear power and core inlet temperature for both the analog and
frequency~-response models. -

The core transfer function equations solved in the MSRE frequency-'
response code are as follows.

Single Core Region. The equations are obtained by substituting the
Leplace transform variable, s, for %t in Egqs. (A.1) through (A.6)7and
solving for TFl’ FQ’ and dk in terms of the inputs TFl, and PT.

Tt is assumed (without loss of generallty) that the wvariables are written
as deviations from steady state. Thus the Laplace transformed egquations

that follow do not contain initial conditions:

 

 

 

 

‘ o * Ke2
1 K, Ky mA T oW |,
T TLin T TR R e oA |
F1 Lo pL  far T e Moy 8 T
G
T z (a.7)
= - . 3 Ao
F1 | | hA
1 nA
LA Kil Ry MC
TFL Ke2 SM;’ hA
G
T Ji(s) TFl,in + J,(s) I (A.7a)

 

 

8 + —

 

’ (A"a)

 

 

To & I5(s) Ty s+ 9,(8) By (A.8a)

 
 

 

 

 

 

 

 

52
o hA
11 tE KEZK MEA MCPGhA = 1) 3y(s)
_ F2 \Gl G2 p2 \s + -fib—-—
. :k G 2
T, = ‘ -
F2 < +_;}_ Fl,in | |
F2 (A.9)
*1 : K('.z hA Mgw _ K, M(ip 5 Mcp2
T 7K MC —i — 1] 98 v
2 Gl G2 P2 s + MC : p2 g + 'fl'c':""
| ¢ | G —,f)
+ = . _ " — »
s +-—i— T
Tp2
| TF2 g H (s) ¥1,in + Ba(8) Bp s (A.9)

 

= [IFl J,(s) + I, Hl(S)] %1%; + I J5(s) %I-.;; »  (A.10)

?Fl,in
-f-’i [ (s) + I }{2(3)] + 1,9 (s) (A.11)
i;T Ty 92 F2 T &
Bk = Hy(s) %Fl,l + H,(s) %T . (a12)

Multiregion Core

The overall transfer functions for an axial section of core consist-
ing of several regions in series are. complicated by the fact that the in-
puts to the upper (or downstreem) regions are affected by the response
of the lower regions. A block diagram illustrating the coupling in terms
of the transfer functions HI_A(S) is shown in Fig. A.1l.

The general forms of the coupled equatlons of n regions in series

!:—5)

>

ci J=l

e

 
 

 

ORNL-DWG 65-9828

 

  
  

- CORE INLET
TEMPERATURE, 7 Teo,

———crc.

—= CORE OUTLET
TEMPERATURE, 7,

 

3k

Fig. A.l1. Series Connection of Single-Core Regions.

 

 

 
 

54

5 | | |
co _ |
—=H (H]_,2 H o Hl,n) + 1, 5 (51’3 H 4o Hl,n) + .ee, (A.14)
Py | .
sk
'1_\'--:= HB,l + }I].,l {H3’2 + Hl’z [H3,3 + H.]_’B (H3,4 + eve 3 (A.lS)
ci |
A n
Bk ' |
75; " :j; H4"j ¥ H2}l H3’2 * (Hz;l H1)2 + H2;2) H3,3

+ (Hz,l' H1,2 31’3 + H2’2 1—[1,3.+ H2,3) %,4 +vee . (A.16)

The mean value of the core outlet temperature, E%o for m axial sections

in parallel is
T = jgl (FF,) Ty s s (A.17)

where FFj is the. fraction of the total flow in the jth axial section.
The total Bk is simply the sum of all the individual contributions.

The calculation of nuclear aversge temperature transfer functions
wag added to the MSRE frequency-response code as an afterthought; conse-
quently there is some repetition in the calculations. The transfer func-

tions of nuclear average temperature contributions from each core region

 

 

are

T .

Fl,in
z £ = I, I5(s) & Hy(s) , | (A.19)
Fl,in

T

5" Iy Ja(s) + I, Hy(s) = Hy(s) , (4.20)

T |

UQ

 
 

 

u!

'l-3>| >
k=i %

 

55
Ag' | |
~. = IG sz(-s) E HB(S) ’ _ | (A.21)
PT

where the asterisk indicéééS’a'hublear:average temperature.

The equations for the total nuclear average temperatures of the nine-
region core model“are'defived'the,same way as the general equations for
8k, Eqs. (A.15) and (A.16). The sécond subscripts refer to the nine core

regions, as designated in Fig. 5. The equations are

. 9
L= Z: +H, , H
3 A Y7,5 T V2,2 5,3

(H, - Hi 3t 8 5) Hy ) +H, 5 Hy o

(s ot iy ) By o + 1y g Hy o, (A.22)

=l

9
= H + B
géé 8,5 * 2,2 H,3
+ (Hiz”2 Hl,3 3) + H2 5 6 ¢ ;,,'

+ (3215,3i:6.+_H2;6)_Eé,7 +_H2;8_H5:9’ _(A.23)

CH W
.a'---1,2z,:5,sfi5’j ' }11'2-:'?’3 :
1,2 Hi 3 5 4 + H1 5 Hs 6

o+ Hi 5 Hi 6 H 7 + Hl g Hs 9 , ‘(A324)

 
56

= L He,5 * P12 He,3

+H12H13 64+H15 6,6

+H15H16 6'7+HIB 6, (“_"25)

Neutron Kinetics

 

Nonlinear Equations

 

Neutron balance:

2 =2 k(1 -y 1] +§, NG e (ae)

Precursor balance:

ac knBi‘ : _ . '
—_— - A, C. . . o (A;27)7*

In these expressions

neutron population,

n
t
1* = prompt-neutron lifetime, sec,

time, sec,

k = reactor multiplication,
'BT = total delayed-neutron fraction,

B. = effective delayed-neutron fractlon for ith precursor group,r

i
N = decay constant for ith precursor, ~
Ci = ith precursor population.

~ For the one=-group approximation, the effectlve B in the precursor |
balance equation was simply the sum of the B's for the six grOups. The
average decay constant N was calculated from Eq. (A.28):

 
{rpopulation in

 

Mo
o
'—lo

J=te
il
I.-..l

. | (A.28)

 

>
!

=t

i o
>1+FD

Linéérized:EquatiOnsiwith Circulating Precursor Dynamics. The dif-

ferential equation for the precursor population in the core is

-\ T
. { % c.tt) C.(t -1 € 1L
dcl(t) =k81 n( ) - 7\5_ Ci(t) - i( ) + l( TL) ’ (A.29)

* _
dt l ° Tc

 

 

formation decay

rate of change rate of ] [rate of]
of precursor nl =
in core in core

core

e rate of reentry by precursors
;j-ia o l;-that left core 1y, seconds ago
caving and decayed to a fraction ?

Lcore T
1'L of their previous value

where

Bl total delayed-neutron fractlon fbr 1th precursor group,
Tc = core holdup time, sec, -

flézL = loop holduP tlme, sec.r'”_7°;"

For thlS treatment, we assumed that the core is a- well-stlrred tank

~ and that the precursor transport around the loop is a pure delay. We
;QObtained the llnearlzea neutron kinetlcs eqpatlon used in the frequency-
fiforesponse calculatlon from Eqs. (A 26) and (A 29)

 
 

58

6

 

 

 

. AB!
1=ty [ T oeny
| ' i=1 ' ~T. (s+A\
. s+7\j'_+%—[l-e L i] |
l n ) | k .
—§= — °» — » (A.30)
n : : .
0 1 - (1-8) ky + SI* ~k/ >, 13
| i=1 1 -1 (8+)\,)
8+ N +=—|l—e
where

ny = steadyhstate nuclear power, Mw,

ko = steady-state multlpllcation constant

and the 01rcumf1ex (~) indicates 8 perturbation quantlty) that is, o

The wvalue of the critlcal reactOr multiplication factor ko is comp
puted by setting dn/dt and dci/dt equal to zero in Eqs. (A.26) and

(A.29):

 

 

k, critical = e~ (A3)
| 1=B,+ ) —2
T 1 -hiTL
=l A\, +=—[{1-e
1 'I.’c

Heat Exchanger snd Radiator Equations

The coefficients for the heat exchanger and radiator gquatiohs are
given in terms of time constants and dimensiohiesé parsmeters., "A"de-
tailed discussion of this model is given in ref. 18;; The equatibns,
based on the model shown in Fig. A.2, are | '

at ?*l? (Tl,in - Tl)- +'€'f (T_T ~T) ST (_A'Bz)
aT 2 | - .

1,out - o= T : .

at 51 (Ty = T, out) * T* (Tp =) » (8.33)

e @ -T) e -T), (e

e A = -

ifi;-

 
 

 

59

‘o

ORNL-DWG 65-8829

~/ASSIERRRERESY

_SHELL s °

 
 

  

A NSRRI

O=(

A

 

AN

        

.. TTUBET _
Heat Exchanger and Radiastor Tube Model.

Fig. A.2.

 
 

60

- \&'"
dT, 2

0 1 o
- D e - — - -
T (T, 3o~ T) + 5 (T, T,) + F:# (T, - T,) , (A.35)

2

45 out
-—-4211— £ e (
dt tg

e =y .%o = *
T, - T2,out) + -5-2; (TT - T2) + -1-:-%- (Ts - Tz) ,. (A.36)

eT Y | |
® T (.T2_ - TS.) > (A.37)

where the nomenclature of Fig. 7 gpplies to the temperatures, T, and s

T = mean shell temperature, 6F3 -»

t* = transport time, sec, | | -

T = heat transfer time constantMQb/hA, sec, .

n = section length = hA/WCP, dimensionless. ~ =

The subscripts have the following meanings:

1 = fluid 1 side,
2 = fluid 2 side,
T = tube side,

s = shell side.

Also,

h = heat transfer coefficient, Btu/sec.ft?.°F,

A = heat transfer area, ft?, : l 5
M = mass of tube or shell, 1b, | |

c:p = specific heat, Btu/lb-°F, | | .
W = mass flow rate, 1b/sec.

Since it is the radiator air flow rate that is varied to change the
loed demand, the radiator shell-side coefficients will vary with power
level. The coefficients listed in Appendix B are for a 197-1b/sec air
flow rate, corresponding to 10-Mw power removal at design temperatures.
In all the studies, h

air

was varied as the 0.6 power of the flow rate,
and WI. as the first power. '
eir

 
 

@
‘>

The solutions of the Laplace-transformed transfer function equations

are

>

(2-n)) +nB

1out _ Lr, (A.38)

>

.tls+f‘<".v

. Lin

>

n
: £
2,0ut _ , s L, (a.39)

 

>

t*¥ s + 2
2,in 2

| . DE |n, + ¢ (2-n )]
¢ 1,out [ 1 1 7
| = = — - , - (A.40)
t¥s + 2

>

>

' n
~ BFin, + (2 = n -n)A+----————S A
2,out== [2 . 2 8 Ts s +1 ]

 

, (A.41)
' * ,

>

where L

n S

S
*  —
t28+n2+2+nsw

‘7T s + 1 -
s

1
T St ._'Ll. —===A
) = X T 0 , g
i

 

  

 

 

 
 

62"

- - — .
#ls+nl+2 nlB

‘To compute the transfer functions of -an arbitrary nunber of éqnal lumps
connected in series, we considered first the transfer fUnctions for two
lumps in series, as in Flg. A 3, where for each lump,

o - lyout

5
L Tl,in'

%1n

The transfer functions for the two combined lumps are

 

 

T ec
Alzou.t - 11 , T (A42)
1 -GG/ 1
1,in lcomb. 34
T GG'~ -
.A_ZJQ."E ._......_.g._ , (A.43)
1- GG/
2,in lcomb. 34
A2zout .Gl + 2 3 - © (A.44)
T 3 1-%%
1,in fcomb.
T C G/ ', | I
T o 1—%4 I
2,in Comb. :

To solve for more lumps in serieé,'we set the primed functions eqnal‘to
the respective combined transfer functions and repeated the'computation.

 
 

 

 

> ‘ i

in

Fig. A.3.
Series. -~

 

63

ORNL-DWG 65-9830

 

 

Diagram for Computing Transfer Functions for Two Lumps

 

 
 

 

64
o
Equations for Piping Lags |

- The first-order (or well—stifred-tank) epproximation used in the
reduced model, is given by

dT -
where
T = meen (and outlet) fiuid temperature, '°'F,

T,, = fluid inlet temperature, °F,
7 = fluid holdup time, sec, | ’
and heat transfer to the piping is neglected. o ; »

~ In the intermediate model, fourth-order Padé epproximations were
used. They are series expansions of the Laplace transfbrm.expressions
for a pure delay: '

-

=18 1072 — 53615 + 1201°s2 — 13.5513s3 + ¢%a% | _
e f 1 , - - (A-47)
1072 + 5361s + 1207282 + 13.5513g3 + %%

The heat transfer to the piping end the reactor vessel was approxi-
mated by lead-lag networks. The method for obtaining the coefficients
Iy and L, is described in detail in ref. 9. The general form of the
equation is

~

 

 

out st 1 7 oo
~ = *

Tin Los + 1 | : : _
In the complete model, the heat transfer to the piping and the trans-

port lag were represented by the exact solution to the plug-flow equa-
tions:? |

T - n/{1l+1,8) _
"_\_C_J'llb_ =e T 7l ¢ P ’ ‘ (A.49)

in

 
 

s e e e e A e e A T AT S SR AL 5 oo areetan s ot

"\

where

B
n

Tp

 

65

section length = hA/WCP, dimensionless,
transport lag, sec, -
time constant for heat transfer to pipe = MCP/hA, sec.

Equations for Xenon Behavior

Xenon

was considered only in the present frequeney-response model.

The following differential equations were used:

where

H
i

Ok =

ax

T = K3X/ — (K2 + K3Pg) X (A.50)

ax/ , |

o = KT + KX, = st + KyP , (A.51)
2% = —KglI + KoP , (A.52)
k = =K30X/ — K11X (A.53)

xenon concentration in graphite, atoms/cc,

= xenon concentration in fuel salt, atoms/cc,

iodine concentration in fuel salt, atoms/cec,

nuclear power, Mw,

fchange in reactor multipllcatlon factor due to change in

xenon concentratlons,

K- =

constants,

' Del&yed Power Equatlons e

The equatlon fbr total thermal power, PT’ includes a flrst-order

J

'  1ag approximation of the delayed nuclear power due -to gamma heatlng.

 

*' "dn.'(n"PT) A
dt (l Kd)_'a‘f"' Tg K (A.54)

 
 

 

K. = the fraction of flux power delayed,
n = flux power, Mw, _
P, = total thermal power, Mv, | Lo e

?rg = effective time constant for the delayed power, sec.

The frequency response of the thermal pcwer.in"terms of nis |
(1-K)ts+1 o
= £ 0 - (A.55)
T8 +1 ' ' |
e

 

.‘3)[6 >

 
 

 

 

'fl

67

 

Appendix B. Coefficients Used in the Model Equations

Core Thermal Equations Data

 

 

0.003161

Study Region 1o, (sec) Tpo (sec)r K, K, K.y Keo

Burke 1961 1 3.815000 - 3,815000  0.470000 ° 0.470000 0.030000 0.030000
analog _ '

Ball 1962 1 1.533 1.607 0.01493  (0.01721 0.000946 0.001081
analog 2 2.302 1.574 0.02736 0.0455 0.001685 0.00306

3 1.259 1.259 0.04504 0.04656 0.003029 0.003131
4 1.574 3,064 0.05126 0.04261 0.003477 0.002395
5 2.303 1.574 0.03601 0.06069 0.002216 0.004081
6 1.259 1.259 0.06014 0.06218 0.004044  0.004182
7 1.574 3.066 0.06845 0.05664 0.004603 0.003184
8 2.621 1.525 0.06179 - 0.07707 0.00392 0.00583

9 1.779 2.983 0.09333 0.07311 0.006277 0.004305

1965 frequency- 1 1.386000  1.454000  0.014930 0.017210 0.0009%46 0,001081
response and 2 2.083000 1.424000 0.027360 0.045500 C,001685 0.003060
eigenvalue 3 1.139000 1.139000 . 0.045040 0.04656Q0 0.003029 0.003131
calculations 4 1.424000 = 2,772000 0.051260 0.042610 0.003447 0.002395

5 2.084000 1,424000 ~ 0,036010 0.060690 0.002216 0.004081
6 1.139000 " 1.,139000 - 0.060140 0.062180 0.004044 0.004182
7 1.424000 2.774000 0.068450 0.056640 0.004603 0.003184
8 2.,371000 1.380000 0.061790 0.077070 0.003920 0.005183
9 1.610000 2.700000 0.093330 0.073110 0.006277 0.004305

1965 extrema 1 4.230 4.230 0.470 0.470 0.030 0.030
determination . ‘

1965 frequency- 1 1.14068 1.19664 0.01518  0.01750  0.000695  0.000794
response with 2 1.71431 1.17195 0.02783 0.04627 0.001237 0.002247
extrema data 3 0.93740 0.93740 ~  0.04581 0.04735 0.002224  0.002299

4 1.17195 2.28136 0.05213 0.04333 0.002531L 0.001758
5 1.17513 1.17195 0.03662 0.06172 0.001627 0.0029%9
6 0.93740 0.93740 0.06116 0406324 0.002969 0.003071
7 1.17195 | 2.28300 - 0.06901 0,05760 0.003380 0.002338
8 1.95133 - 1.13574  0.06284 0.07838 0.002878 ©  0.003806
9. 1 2.22210 ~0.,07435 0.004609

1.32503

0.09492

 
 

 

, Mcpl Mcp2 MC e . o
Study Region - hA : \
(Mi-sec/°F)  (Mw.sec/°F) * (Mwesec/°F) (Mw/°F) I Tpo Is

Burke 1961 1 0.763000 0.763000 3.750000 0.020000 1.000000 O. 1.000000
analog ‘

Ball 1962 1 0.0188 0.0197 0.070 0.265 X 10~  0.02168 0.02678  0,04443
analog 2 0.0638 0.0436 Q.2114 0.814 x 10~ 0.02197 0.06519 0.08835

3 0.0349 0.0349 0.1601 0.609 x 10~?  0.07897 0.08438 0.16671
4 0.0436 0.0850 0.2056 0.79% x 10~  0.08249 0.04124  0.12077
5 0.1080 0.0738 0.3576 1.338 x 10~3 0.02254 0.06801L 0.09181L
6 0.0590 0.0590 0.2718 1.031 x 10~3 0.08255 0.08823 0.17429
7 0.0738 0.1437 0.3478 1.343 x 10~3  0.08623 0.04290 0.12612
8 0.2970 0.1726 0.9612 3.685 x 10~2 0.02745 0.05521 0.08408
9 0.2014 0.3380 0.9421 3,624 x 1073 0.06936 0.03473 0.10343

1965 frequency- 1 - 0.015100 0. 015800 0.070000 0.000392 0.021680 0.026780 0.044430
response and 2 0.051200 0.034900 0.211400 0.001204 0.021970 0.065190 0.088350
eigenvalue 3 0.028000 0. 028000 0.160600 0. 000900 0.078970 0.084380 0.166710
calculations 4 0.035000 0.068200 0.205600 0.001174 0.082490 0.041240 0,120770

‘ ‘ 5 0.086600 0.059200 0.357600 0.001977 . 0.022540 0.068010 0.091810
6 0.047300 0.047300 0.271800 0.001525 0.082550 0.088230 0.174290
7 0.059200 0.115200 0.347800 0.001985 0.086230 0.042900 0.126120
8 0.238000 0.138400 0.961200 0.005445 0.027450 0.055210 0,084080
9 0.161500 0,271.000 0.942100 0.005360 0.069360 0.034730 0,103430

1965 extrema 1 0.750 0.750 - 3,58 0.020 1.0 0.0 1.0
determination : '

1965 frequency- 1 0.01581 0.01654 0.070 0.588 x 10~2 0.02168 0.02678 0.04443
response with 2 0.05360 0.03654 0.2114 1.806 x 10- 0.02197 0.06519 0.08835
extrems data 3 0.02931 0.02931 0.1606 1.35 x 10-3 0.07897 0.08438 0.16671

4 0.03664 0.07140 0.2056 1.761 x 102  0.08249 0.04124 0.12077
5 0.09066 0.06197 0.3576 2,965 x 10> 0.06801 0.06801L 0.09181
6 0.04952 0.04952 0.2718 2.287 x 103 0.08255 0.08823 0.17429
7 0.06197 0.12060 0.3478 2,977 X 103  0.08623 0.04290 0.12612
8 0.24915 0.14489 0.9612 g.167 x 10~3  0,02745 0.05521  0.08408
9 0.16907 0.28370 0.9421 8.04 X 1073  0.06936 0.03473 0.10343

C

. . .

)

89

 

 

 

 
 

 

 

69

Flow fractions FF in the four sections (nine-region core only)‘for

all studies were

FF (1) FF (2) FF (3) FF (4)
0.0617 0.1383 0.234 0.566 °

Neutron Kinetics Equations, Date*
1965 frequency-response, eigenvalue, and extrema data for decay

constant %i(sec'l):

N N M A, As s
3.01 1.14 0.301 0.111 0.0305 0.0124 °

 

195 frequency-response data for total delayed-neutron fraction for

ith precursor group, B;:

L B2 B3 B4 Bs Bé
0.00028 0.000766 0.002628 0.001307 0.001457 0.000223 °

1965 eigenvalue and extrema data for effective delayed-neutron frac-

 

tion for ith précursor groups |
B2 B2 B3 B4 Bs Bs
0.000277 . 0.000718 0.001698 0.000499 0.000373 0.000052 °

 

Heat ExchangerAand Radiator Data

Radiator air-side data given for 10-Mw conditions: heat exchanger
fluid 1 = coolant, fluid 2 = fuel; radiator fluid 1 = coolant, filuid 2

= air.

*

 

o om ng % 2T TR s
, , (sec) (sec) (sec) (sec) (sec)
Burke 1961 and  Heat ex- 0.980 0.906 0  1.95 2.24 1.75  1.165
-+« Ball 1962 - changer - - i o ' S
- analogs ‘Radiator - 0.0882 0,260 O 7.14  0.0L © 2.35  19.7

1965 frequency- Heat ex-" 1.10 = 1.366 0.1363 2.01 2,29 - 0.569  0.304 1l.14
response, eigen- changer = - o : -

value, and ex- Radistor 0.8803 0.2591 0O 6.52 0.01L 2.35 19.7 N
trema determinae T : _ S ' o
tions

195 frequency-  Heat ex- 1.60  1.611  0.1363 1.80 2.29 2.5 1.16 1.14

. response with

extrema data

changer

Radiator

 

0.983 0.2591 O

5.84

0.01

2.35

*¥See Table 1 of Section 4 for additional information.

19.7

 
 

70

 

 

 

 

 

Piping Lag Date
1965 Freguency ' 1965 Frequency
Burke — Ball Response and 91265 ?xtz:mas  Response with
1961 1962 Eigenvalue elermination Extrems Data
Core to heat n : a -
exchanger T 3.09 5.72 5.77 5.77 / .6.30
, | . , :
P .
Heat ex- n 0.155 0.155 o 0.155
changer to T 9.04 9.02 8.67 8.67 | 9.47
core . T 15.15 15.15 , s 15.15
_ : ?
"~ changer to T 5.2 5.2 471 4,70 4.22
radiator p 6.67 - 6.67 | - 6.67 &
Radistor n 0.40 0.40 - 0.40
to heat t 10.11 10.11 8.24 - 8.24 B 7.38_
exchanger Tp 6.67 6.67 L | 6.67
The coefficients of the lead-leg approximation for fluid heat tréns-
fer to the piping apply to the Ball 1962 anslog study:
Heat Exchanger Heat Exchanger Radiator to
to Core to Radiator Heat Exchanger
La 14.35 5.84 5.58
Lo 16.67 7.69 8.33
Xenon Equation Data -
Data for 1965 frequency response and frequency response with extrema ‘ N
data: ' : | -
Ky = 1.587 x 1076 Ky = 2.885 x 107%
Ky = 2.2575 x 10°° Kg = 2.9 X 10°°
K3 = 1.654 x 1076 Ko = 9.47 x.10-%
Ks = 1.0714 x 10™° Ky1 = 1.03 x 10~4
K¢ = 1.059 x 1073
Delayed Power Equation Data o —

 

See Table 1 of Section 4.

 
 

 

 

‘o

fia

 

[

 

71

-Appendix- C. General Description of MSRE
Frequency-Response Code

The MSRE frequency-response code (MSFR)'is written in FORTRAN IV
language for the IBM 7090 computers at the Oak Ridge Central Data Pro-
cessing Facility. This language has builtin capebilities for handling.
complex algebra that result in considerable savinés of programming effort.

MSFR uses transfer function techniques (rather than matrix methods)
to compute frequency response. It exploits the fact that a reactor sys-
tem is made up of separate components, each having a certain nunber of
inputs and outputs, which tie in with adjacent components. The subrou-
tines written for each subsystem were useful in other reactor and process
dynamics calculations. - The MAIN program of MSFR performs input, output,
and supervisory choree, and calls the subroutines. A subroutine called
CLOSED must be written to compute the desired closed~loop transfer func-
tions from.the component transfer functions.

The transfer function approach has several edvanteges over the ma-
trix methods:

1. Input parameters are the physical coefficients of the subsys-

tems, rather than sums and differences. This not only makes generating

'input data easy,_it allows the computer to carry out the sum-difference-

type arithmetic internally. Several matrix type computations for which
the matrix coeffiCiehts'were‘generafied "carefully" with long slide rules
resulted in large errors in the frequency response.

2. The frequency response of dlstrlbuted-parameter models can be
computed exactly with MSFR, whlle most matrlx calculatlons are llmlted

to lumped-parameter models.rigc

3. MSFR calculatlons are much faster. The 7090 can put out between

: ,1000 and 2000 frequency-response pOlnts per minute for the complete

' -model Iyplcal runnlng tlmes for current matrlx calculatlons are much

longer.
_ The matrix technlque has the advantages that spe01a1 programmlng

~is not required for each dlfferent problem, and no algebralc manipulatlons

-;joflthe-equations are requrred. HAlso, matrlx_manlpulatlons can be used

 
72

for optimization calculations, eigenvalue calculations, time-response
calculations, and possibly many'others, all with the same input data.

The advantages of both methods were ekploited in this study.

The following subroutines of MSFR have potential as generally use-
ful packages:
| 1. PWR, which calculates the frequency response (N/ak) of the
" nuclear kinetics equations for up to six groups of precursors, with an
option for including c1rculat1ng precursor dynamics.

- 2. CLMP, which computes the frequency response of & "typlcal" core
region (as noted in Appendix A). Inputs are power and inlet temperature,
and outputs are outlet temperature, nuclear average temperatures, and Bk.

3. COR9, which calculates the overall frequencyaresponses_of the
MSRE nine-lump core model using CLMP outputs. , -

4. LHEX, which calculates the transfer functions of a 1ufiped-
parameter heat exchanger (as in Appendix A), with an input optiofi for
solving for up to 99 typical sections in series in a counterflow con-
figuration. ,

5. PLAG, which computes the frequency response of piping lagé for
an erbitrary number of first-order series lags, & fourth-order Padé ap-
proximstion, or a pure delay, or combinations of these, with heat trans-
fer to the piping.

Figure C.l shows the block diagram used as a guide to compute the
closed-loop transfer functions. Typical outputs of the subroutine CLOSED
are fi/&fi (closed loop), Nyquist stability information, nuclear average
temperatures %;/Sfi and %E/Sfi, and %co/afi and %Ei/Si'

Several commonly used transfer functions are

~

 

 

ci - G1.0G11G13614Ga5 |
~ g GS = GryGglGo + 1, (¢.1)
Teolopen primary 1 = G12G13G14G15
loop
and
N Gy |
~ C.2
sk = G4GsGS ( )

1+ GJ.Gz(-]—_-:E-;G—S' + G3 = GX)

o

 
 

 

73

‘o

ORNL-DWG 65-9819

/f(NO)

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

8k | wrsk | M Pr/N
G1 G2
Bkye/Pr
Gx y THERMAL POWER, Pr
RODS j 1 '
-~ 3k + - 5k /Py . Te/Pr
: 63 [ 1 ~—T7F
¥- . I G16
i - Bk/fei _ CORE Teo/Pr
G4 G5
: ' TC?/PT _rg
ru | TF/Tei (g 3 T COREIN |Teof7ei | T CORE OUT 617
- UF Gig [ 66 : t
A
r*, .
VESSEL |- "
PRIMARY LOOP
| 619 .
4 rpo/rpi THE (PRDIN ¥
- 69
. = Y
rpo/rsi - HEAT rso/rpi
61t |[EXCHANGER| 6t0
| N —
Tug (SECYINY | Too/Tsi Tue (SEC) OUT
: _ Y
%
. " PhE|  sEconparv Loop  |™.PifE|
. . ?:5 ] RADIATOR | 614
| . . o . _ :=' rso/rsi v
Fig. Cl MSFR Reference Block Diagram. :
-

'\

 
 

 

 

T4 | | e

Fach of these closed-loop equations can be written as a single
FORTRAN IV statement, so it is a simple inatter to generate different
functions. An option is also available in MSFR to firint out all the
internal or component transfer functions. FORTRAN;_IV listings, decks,
and input information may be obtained from S. J. Ball.

 
 

)

at

 

75

Appendix D. -Stability Extrema Calculation

For a linear description of a reactor system,_the eigenvalues of
the system matrix must all have negative real parts for stability. A
technique16 was developed that_éystematically seeks out the combination
of system.parameters'that causes the.least stableicondition in the
feasible range (causes the dominant eigenvalue to become as positive as
possible). This technique utilizes a form of the gradient-projection
method to explore the hypersurface that defines the stability index (the
negativeness of the real part of the dominant eigenvalue) as a function
of the system parameters. The upper and lower limits on the expected
ranges of system parameters constitute constraint surfaces that limit
the ares of search on the performance hypersurface.

The real part of the dominant eigenvalue is labeled B. The change
in B due to smell changes in the system parameters, X5 is given by

4R = YB.dX = |VB| 'dk'l cos 6 , | (D.1)
where
d8 = incrementael change in B,
3B S "
vB-el&;-+GZS%2-+ cse

dJ?=eJ_Xm+eng2+...,

a unit vector,

e,

; = angle between the_yectors.

‘Thus the maximum change in B occurs when 6 = O; that is, the changes
infthesystem,parameteré“afe in-théflSame vector direction as the gradient
vector. It is therefcre_eipected that the greatest change in B will
occur when the system.pafdmeters change in proportion to.their corre~
sponding elements in'the.gfadientlvector:

ok, es

l

 
 

76

where ¥ is a real positive coefficient whose magnitude is chosen to
insure 'that constraints are satisfied.

| Tt is clear that the calculated velues of the componeni:s of VB are
the key quantities in implementation of this method. The method for
finding VB can be developed from the characteristic equation for the

systém given in determinant form:

D=IA-‘SII=0, - | (DB)

where

A = the systém matrix,
8 = an eigenvalue of A,
I = the unit diagonal matrix.

We now write D with some arbitrary eigenvalue, 8y factored out:

D= (s —-sk) F(s) = 0, (D.4)

where F(s) is & nonzero determinant if s, is a simple eigenvalue. We
differentiate Eq. (D.4) with respect to an element, 8447 of the matrix,

A, and with respect to s:

aD OF(s) Os

S— = (s-—sk)g————F(s)é—k—, a - (.5)
*1J "13 ®13 -
%182 = (s — sk) BFSS + P(s) . (D.6)

We then evaluate Egs. (D.5) and (D.6) for s = s, and take their ratio to
get Bq. (D.7):

oD
ds, E;'_,j ’s:sk _
= - . (D'7)
a. 13 oD
Js ls:sk

TS e i

 
 

4

i

a"

 

 

77

The deri%ative; BD/B& 4+ is just the cofactor of 8 5 in [A = 81],

1]

“and aD/Bs is the negative of the sum of the cofactors of dlagonal ele-

ments of [A — sI]. Thus Eq. (D.5) may be written

355; = -trig*-'- | (D.8)

ij z: Cff
- =1 -
If we choose sk to be the least negative 31genvalue, the real part of
) is Just B. Slnce alj is real, we may wrlte
® [ c "
o = Re _E-_Q— . (D-g)
i : . _
o Y Cer

f=1

The derivative with respect to a system pareameter, xz,,is easily ob-
tained from Eq. (D.9), since the following relation holds:

. - | (p.10)

The usefulheSS df:Eq.:(D;iO) rests on the ability to calculate the

_rsystem eigenvalue, k? to glve IA - s I] This mabeé readily accom=-
_'plished by using one of the standard eigenvalue com@utation methods,

guch as Parlett's methodlg or the QR method. 14

The cofactors. in Eq. (D 10) could be calculated dlrectly'w1th a

method such as Gaussian eliminstion. However, this tedious procedure may
be circumvented by application of a useful theorem from matrix algebra.

 
 

 

78

It is known that the cofactors of parallel lines in a matrix with order

'n and rank n -~ 1 are prqportional.2° Since [A — SkIJ is & matrix With
these properties, the cofactor calculation may be simplified. For in-
stance, if the cofactors of the first row and first column of [A — skI]

are calculated, all other cofactors are given by

_ % Oy

i3 vl | (D.11)
Use of Eq. (D. 11) to find the cofactors shown in Eq. (D 10) glves a prec-
tical method for finding the derivatives afi/ax needed to carry out the
gradient-projection step shown in Eq. (D.2). -

Gradient methods ere useful for finding local extrema fbr noniinear
prdblems. However, it nmy'be possible for the surface of B versus sys-

- tem parameters to have many peaks. The only tedhnique currently suiteble
for handling this provlem is to use multiple starts. The computer code
developed to implement this method is set up to use maltiple sterts au~
tomatically.

The procedure for carrying out the maximization from a given base
point is to recalculate the eigenvalues for several new parameter sets
specified by steps out the gradient vector. The point that gives the
system with the largest value of 8 is then used as s new starting point.
This is repeated until a maximum within the constrained set of system

parameters is found.

.

 
 

 

£

‘\a

<O W

8.

10-34.
35,
36.
37.
38.

39.

40.

41.

42.

43.

b,

45.

46.

47,

48.

49.

50'

51.

52.

53.

54.

55.

56.

57"

58.
59.
60.
61.
62,

V' 63.

65.
66.

67.
68.

69.
70.
71.

R.

R.

.‘L.

A.
R.

Je
M.

.S.

S.

E .

F.
R.
E.
c.
G.
0.

'So

Ww.
G.
E.
Je

W--

L.
G.

F.

S.

¥.
N.

.J.
E.
. K. Ergen

D.
A.

H-

c.
R.

R'

 

79

Internal Distribution

 

K. Adams

. M. Adamson

G. Affel
G. Alexander
H. Anderson
F. Apple

. I« Baes

M. Baker
A. Baker
J. ‘Ball

E. Beall
S. Bettis

F. Blankensfiip”

Blumberg

G. Bohlmann
J. Borkowski
H. Burger
W. Burke
Cantor

L. Carter
I. Cathers
L. Compere
A. Conlin
H. Cook
T. Corbin
A. Cristy

. L. Crowley

L. Culler

. G. Davis

J. Ditto

. G. Donnelly

A. Doss
E. Dunwoody
R. Engel

P. Epler

E.'Ferguson  -
P. Fraas -

. N. Fray .
A. Friedman

. N. Fry '
.J.

H Frye, Jr.
"H. Gabbard

B. Gallaher‘ '
. R. Grimes
. G. Grindell

H. Guymon

72,
73.
74,.
75.
76.
77.
78.
79,
80.
8l.
82.
83.
84~108.
109.
110.
111.
112.
113.
114.
115,
116,
117,
118,
119.
120,
121.
122,
123,
124,
125,

126.

127.
128,
129.
130.

131.

132.
133.
- 134,
. 135,
136.
137.
138.

- 140.
141.
142,

- 139.

P.

C.

p.
G.
p.
V.
P.
A.

R.
M.

T

S.

D.
'A. I. Krekovisk

J.

R.

Je

C.

C.
H.

W.
H.
C.
L.

W.
R.
H.

H.
A.
H.

B.

J.
M.

J.
H.

H.

ORNL-TM-1070

H. Harley

S. Harrill

N. Haubenreich
M. Herbert

G. Herndon

D. Holt

P. Holz
Houtzeel

. L. Hudson
. R. Kasten

J. Kedl

J. Kelly
W. Kerlin
S. Kirslis
J. Knowles

W. Krewson

. R. Lawrence
C.

G. Lawson
B. Lindauer

« I. Lundin
. N. Lyon

H. Marable
D. Martin
E. Mathews

‘G. MacPherson
. C. McCurdy

B. McDonald

F. McDuffie

K. McGlothlan .
E. McNeese

J. Miller

R. Mixon

L. Moore

G. O'Brien

‘Patriarca
‘R. Payne

M. Perry .
B. Piper

E. Prince
L. Redford
Richardson

. C. Robertson

C. Robinson
C. Roller

* PO ROU.X

C. Savage

 
143.
144.
145,
146.
147.

148.

149.
150.
151.
3 152 .
153.
154.
155,
156.
157.
158.
159.
160.
161.

186~-187.
188.
189.
190.
191.
192.
193.
194.
195.
196210,

 

H.
A.
D.
H.
J.
M.
A.
P.
W.
I.
R.
H.
R.
A.

5.

R.
G.
D.
We.

D.
R.
Se.
D.
H.
R.
We
Se
M.

 

W. Savage

W. Savolainen
Scott

E. Seagren

H. Shaffer

J. Skinner

N. Smith

G. Smith

F. Spencer
Spiewak
C. Steffy
H. Stone
S. Stone
Teboada,

R. Tallackson
E. Thoma
M. Tolson
B. Trauger
C. Ulrich

80

162.
163.

164.

165.
166.
167 .
168.
169.
- 170.
171.
172.

173.

174~175.
176177,
178-179.
180~184.

185.

B. H. Webster

A. M. Weinberg

K. W. West

M. E. Whatley

J. C. White

G. D. Whitman

R. E. Whitt

H. D. Wills

J. V. Wilson

L. V. Wilson

M. M. YarOSh

MSRP Director's Offlce,
Rm. 219, 9204-1

Reactor Division lerary
Central Research Library
Document Reference Section

Laboratory Records Departiment

Laboratory Records, RC

External Distribution

 

FO COPE, AEC, ORO

W. Garrison, AEC, Washington, D.C.
H. Hansuer, University of Tennessee, Knoxville
C. Keeton, Tennessee Valley Authority, EGCR Site

M. Roth, AEC, ORO

F. Saxe, North Caroline State University, Raleigh, N.C.

L. Smalley, AEC, ORO

E. Stephenson, University of Arkansas, Fayetteville, Ark.
J. Whitman, AEC, Washington, D.C.
Division of Technical Information Extension

»