NUREG/CR-3239

.UCLA-ENG-8257

ýRGI

COMPBRN- A Computer Code for

Modeling Compartment Fires

Manuscript Completed: August 1982

Date Published: May 1983

Prepared by

N. 0. Siu

School of Engineering and Applied Science

'University of California

Los Angeles, CA 90024

Prepared for

Division of Risk Analysis

Office of Nuclear Regulatory Research

U.S. Nuclear Regulatory Commission

Washington, D.C.r 20655

NRC FIN B0493

ABSTRACT

The computer code COMPBRN-deterministically models the behavior of fire

in a compartment.;_ This manual presents information necessary to run

.COMPBRN, including,:descriptions, of required input and resulting output.

Also included are a sample problem and a listing of the'code (written in

FORTRAN for an IBM 3033 computer): This manual is to be used in conjunction with NUREG/CR-2269, "Probabilistic Models for the Behavior of

Compartment Fires," August 1981, which describes the fire models employed by COMPBRN.

- iii-

CONTENTS

ABSTRACT. ".....

ACKNOWLEDGMENTS .. .

Chapter

1. INTRODUCTION . . ..

;2. METHODOLOGY .......

•iii

.vii

page

31

. 3

Theory ......... ... ............ . ...........3

Additional Modeling During Application 5.................5

3. 'EXPLANATION PROGRAM FLOW........... ... . . . 8

4. OUTPUTS .............................. 11

5. INPUTS . . ................... ........ 13

Overall Job Parameters (&STRT)

Title Card ............

Individual Job Parameters '(&SIZE)

Fuel Bed Definition (&FUELB)

Pilot Fuel Parameters (&PILOT),

Physical Fuel Parameters (&FUELT)

Miscellaneous Data (&MISC)

Communication Data (&SEE) ..

Non-communication Data (&NSEE).

Room Data (&ROOM) ........

Initialization Data (&GINIT) . .

Model Variation Factors (&MODVAR)

Output Definition (&OUTF) . ..

6. A SAMPLE PROBLEM ..........

7. REFERENCES .. . . . ... . . . . .

* 13

14

•14

15

* 16

* 18

* 19

19

* 20

21

.22

* 23

24

26 

. . . . . 51

-V -
Appendix page

A. SOME DETAILS ON THE MODELS FOR WALLS AND THERMAL BARRIERS, .... 52

Thermal Barrier Model .... ..................... .......... 52

Wall Model ............... ................ . . .. . 54

B. CODE LISTING . . . . ........... . ...... . .. ...... . . . . 56

LIST OF TABLES

Table page

1. COMPBRN Sample Problem Input . .; ......... . . ... . 31

2. COMPBRN Sample Problem Output . .. . . .............. ... .... 33

LIST OF FIGURES

Figure page

1. Flow Chart for Computational Model ......... .. 4

2. COMPBRN Subroutine Flow . ........... . .... . . . . .... 9

3. Wood Crib ..... ......... ................... 27

4. Compartment Housing Crib ... ........ ...... ............. . .. 28

5. Model of Wood Crib ... ..... ... ....... ...... . . . . 29

6. Crib Mass Burning Rate .. . .. .. ............. . . . .... . . .. . 48

7. Heat Flux from Wall to Detector .. . . •............ . . . 49

8. Hot Gas Layer Temperature ......... . . ............ . 50

"- vi -
ACKNOWLEDGMENTS

The author gratefully acknowledges the aid and guidance provided by

Professor G. Apostolakis in this work. Also acknowledged are Dr. M. Kazarians for his support, Pickard, Lowe, and Garrick, Inc. for providing

the opportunity to apply the COMPBRN program in nuclear power plant risk

analyses, and Ms. S. Rao and Mr. E. Avidor for their helpful discussions, comments, and criticisms.

- vii -

Chapter 1

INTRODUCTION

COMPBRN is a deterministic- fire hazards, computer program designed to be

used in a probabilistic analysis of fire growth. in a particular room.

As described in References (1) through (3), 'this probabilistic analysis

is required when assessing the .risk 'associated with fires in nuclear

power plants. Possible outputs of the program include, the total heat'

release rate of the fire', the temperature and <thickness 'of Ithe hot gas

layer near the compartment ceiling, the mass burning" rate for individual

,fuel elements,.and the thermal'heat'[flux at-user-specified. locations..

The various sub-models used :to construct COMPBRN are obtained from

the .fire research literature. Efforts, are made to keep the model complexity to a low level, -due to the large uncertainties 'inherent in the

modeling' of ýaccidental fires.

Model simplicity is also desirable from an economic point of view,

since computing costs can be greatly magnified in a probabilistic analysis. However, simplifications in modeling naturally lead to some loss

of accuracy for certain types of scenarios'. The assumptions made in

COMPBRN are geared towards the modeling of relatively small fires in

large rooms. Thus, the code's predictions are expected to be most reasonable for fire scenarios involving small fuel loads and for the-early

portions of fire scenarios involving larger fuel loads..

- .1 -
Verification of the code's predictions for conditions typical of nuclear power plant fires is currently in progress.

This manual documents an updated edition of COMPBRN; the earlier program was documented in a draft version of UCLA-ENG-8113. The prJmary

changes are in the computation of the average temperature and extent of

the hot gas layer lying underneath the compartment: ceiling. In the earlier form, these quantities were derived using a Newton-Raphson iteration scheme. Numerical convergence problems were subsequently experienced in a number of casesý; a slower but more stable, scheme is employed

in this version.

The current version of COMPBRN employed at UCLA is compiled with the

FORTRAN.H Extended'compiler, using 650K bytes of storage space on an. BM

3033 machine. The load module obtained after link-editting requires approximately 200K bytes of storage,., depending upon the amount. of array

spade called for by the MAIN subroutine. .-Running times are relatively

short (on the.order of seconds) for most problems; they are directly related to the number of time steps and. the total number of fuel cells

used in the scenario of interest.

- 2 -
Chapter 2

METHODOLOGY

2.1 THEORY

The overall model underlying COMPBRN is described, in Reference (1); d6-

scriptions of the model for a particular fire scenario and of the-probabilistic risk analysis in which the.model was used are given in References (4) and (5). Briefly, the burning rate of a fuel element is used

to determine the heat output rate of that element. The burning rate depends on the physical properties, of the fuel, and on the compartment

,'ventilation rate. LUsing standard shape factor analysis and idealizing

the flame as a cylinder, the6heat transferred ,to other fuel elements,

..the walls, and ceiling via radiation is computed. Correlations-are used

to determine! the :convective heat transfer 'in the buoyant plume of hot

gases above the flameks. Provisions are also made to simply model the

effect of the walls and ceiling, thermal barriers within :the& room,- and

the layer of hot gases, accumulating near the ceiling as thermal sources.

The times required for--:ignition of the .fuel elements are then computed

and an element is considered ignited if its ignition time is less than

the time that the element has been exposed to the heat source. Time is

incremented, and the process starts over, with newly ignited fuel elements adding their contributions` to the total rate of heat release. A

flow chart for this process is given in Figure 1.

-'3 -
Compute Burning Rate m(t)'

., --- ----- .

:Compute Heat Release Rate.Ot)

Characterize Gas

Layer Near. Ceiling

-,(Temperature' Depth, etc.)

L4

Compute Rate of Heat

Transfer to Fuel

Elements Within Room

4 .1

Determine Gas. Flow

Rate Into Compartment I I

Determine If Non-Burning.

Fuel Elements Ignite.

FIncrement Time

Figure 1: 'Flow Chart for Computational Model

" 4" -*
The models used to describe these individual phenomena which comprise

the overall model are, with the exception of the ceiling. gas layer model, algebraically quite simple. The main role of the computer in

COMPBRN is to keep track of the heat input and output of each of the,

many fuel elements used tomodel a complex fuel bed.

2.2 ADDITIONAL MODELING DURING APPLICATION

Because there are a number of assumptions implicit in COMPBRN's various

models, further modeling of. the compartment is required before an actual

fire can be analyzed.

In general, a compartment fire will involve a fuel bed of arbitrary

geometry, location, and orientation. Furthermore,. the fuel bed may be

composed of a number of different fuel types., the relative amounts of

which may.vary with location.

In preparing the input data for COMPBRN, the analyst must proceed as

follows:

1. A Cartesiancoordinate system is defined with respect to

some arbitrary location:in the compartment. The z-axis

always. represents the Vertical direction. The compartment walls and ceiling are assumed to be rectangular in

shape and parallel.to the .planes defined by the coordinate axes.

2. The fuel bed is modeled with a series of homogeneous,

two-dimensional rectangular slabs, called- super modules.". The. long axis of .each slab defines its "direction," and must be parallel to one of the: coordinate

- 5-
axes. The direction normal to the slab's• face defines

its "orientation," and also must be. parallel to one of

the coordinate axes. The imposition of a regular geometry on the fuel bed is not always a severe-modeling restriction, since a physical fuel'element can be modeled as

a number of contiguous'super modules. However, COMPBRN is

presently incapable of modeling' fires progressing along

an inclined surface.

3.' The super modulest.are then further discretized (along the

direction of their long axes) into contiguous fuel packets, called "fuel cells." These fuel ce'lls should be almost square in dimension.: when modeling horizontal fuel

beds,-since COMPBRN's flame models for this.configuration

are based on results .from flames over square or round

pools. --fulel. 

4. The fuel bed is further modeled using the "porosity .factor" and the "communication matrix;" -The porosity factor

is required since a "rough" fuel element.- (e,.g. cables in

a-cable tray section)-will have more surface area available for burning than provided by 'the representing flat

super module of the:same nominal.dimensions. This factor

is defined as the actual area available''for burning per

unit base (supermodule) area.

The communication matrix is constructed, by -the analyst. It allows 'him to specify which super' modules or

fuel cells are contiguous, and which ones cannot transmit

-6-
heat to others (due to intervention by other fuel elements, barriers, etc.). If instructions to the contrary

are not provided, all fuel elements are assumed to be

able to "communicate with," i.e. transmit heat to, all

other'fuel elements.

The porosity factor and the communication matrix allow

some fine-tuning of COMPBRN's simple modeling. It must

be realized however that the memory space requirement for

the communication matrix expands as the square of the

number of fuel cells, and this is a practical limiting

factor -in the ability of COMPBRN to handle complex prob-:

lems in detail.

- 7 -
Chapter 3

EXPLANATION - PROGRAM FLOW

COMPBRN is a FORTRAN computer code which is divided into a main program,

a group of primary subroutines, and a number of auxiliary subroutines

and function subprograms. The modu-I` nature of the fcode is not necessary for execution efficiency, but is. intended to allow the analyst to

incorporate improved models as they become available. All of the subroutines incorporate variably-dimensioned arrays; if the fuel bed configuration requires a large ,number of fuel cells, only the main program

need be modified and recompiled.

As can be seen in Figure 2, COMPBRN execution begins with the MAIN

program. Fuel bed geometry, initial fire characteristics, fuel, type

character istics, communication matrix data, compartment parameters, model variation factors, and output format data are read using NAMELIST

directed input. MAIN calls subroutine INCHK, in which various -portions

of the input data can be printed :for verification,.-and then the two. initialization subroutines, INIT and TINIT.

After-,problem definition and initialization is complete, subroutine

SOURCE is used to compute the fuel burning rate, the resulting rate of

heat output, and the, flame height, for each burning fuel cell. The models vary according to the degree of compartment ventilation and whether

the fuel element is horizontal or vertical. The strength of each element as a source of heat flux is then calculated. The auxiliary. su-

-8-
Figure 2: COMPBRN Subroutine Flow

-9-
broutines BARR, DIFFUS, or CEILNG are called by SOURCE to determine the

heat fluxes from a thermal barrier, a wall, or the ceiling, respectively. Subroutine CEILNG, if called, requires the use of -function subroutine RADEQ, which solves a special version of the quartic equation to obtain an improved value for the hot gas layer temperature, and subroutine

ZDNEW, which updates the height of the neutral density plane within the

compartment given an updated hot gas layer temperature.

Subroutine TRANSF follows SOURCE, and transfers the heat fluxes from

the various sources to the target elements via convective and/or radiative pathways. The shape factors used in the flame radiation transfer

calculations are computed using the function subroutines CYLPAR andCYLPER.. The shape factors for transfer from rectangular'sources (e.g.

the ceiling) were computed earlier using SHAPE and RECT, and are stored

in the-communication matrix, Special models are used to account for the

feedback radiation of a flame over a vertical slab back to the:slab, and

for the conduction of heat to contiguous (and non-burning) fuel elements.

Subroutine IGNIT is used to determine if non-burning fuel elements

have absorbed enough heat to ignite, and subroutine OUTPUT prints data

desired by the analyst for each-time step.

1-10 -
Chapter 4

OUTPUTS

COMPBRN prints two sections of output data. The first is used to verify

the problem input data, while the second lists user-specified variables

during computations.

The input.data check is regulated using the indicator variable

INCHCK. If INCHCK equals 0, only heading material and the model variation parameters stored in the array FCTR are output. These variation

parameters can be used to multiply the computed results of a particular

correlation, if individual model uncertainties are to be propagated

through the code. A list of the models which can be modified with FCTR

is given in Section 5.12. If INCHCK equals 1, the physical parameters

characterizing the fuel are also output. If INCHCK is not equal to 0 or

1, all input data, with one exception, are printed. In the case of the

super module coordinates, the super module midpoints input (SMX, SMY,

and SMZ in Section 5.4) are not printed when INCHCK is not equal to O0or

1. Instead, the code outputs the coordinates of the first fuel cell in

each super module, where the cells are numbered as described in Section

5.4.

The output data listed each time step is specified with the indicator

variables IOUTPT and NSMOUT, and the arrays MOUTPT and MSMOUT. IOUTPT

determines the number of output variables to be printed, and the elements of MOUTPT specify which variables are to be printed (see Section

- 11 -
5.13). Similarly, NSMOUT determines the number of super modules for

which output data is to be printed, and the elements of MSMOUT specify

which particular super modules are included in this number. .If IOUTPT

and NSMOUT equal 0, only the elapsed time since initiation and the total

mass burning rate in the compartment are printed. If IOUTPT equals 0

and NSMOUT is not equal to 0, the elapsed time, the burning rate, and an

indicator variable array,,. $BURN, which specifies which.fuel ýcells (if

any) in the specified super module are burning, are printed.

--:. ii12 -
Chapter 5

INPUTS

COMPBRN requires a large amount of input data to specify each particular

problem. This data is read using NAMELIST format, and the NAMELISTs

must be input in the order presented.

5.1 OVERALL JOB PARAMETERS (&STRT)

This data block, defines parameters which are in effect for all of the

jobs to be performed. All subsequent data blocks apply only for the

particular job being executed.

Data Set Members

NJOB = The number of jobs to be run. In general, the entire

input data block must be repeated for each job,.

NTIME = The maximum number of time steps for each job. Jobs may

end earlier if certain iterative, calculations do not

converge sufficiently quickly or if all fuel which can be

ignited has been.exhausted.

NREAD The logical unit number off which remaining data is to

be read (default =5).

NWRITE = The logical unit number on which output is to be written

(default,= 6).

DELT = The time step increment, in seconds.

- 13 -
5.2 TITLE CARD

The first 80 characters on this card will be'reproduced as the job title

in the code output.

5.3, INDIVIDUAL JOB PARAMETERS (&SIZE)

The limits given in this data block are for the COMPBRN version listed

in Appendix B.

Data Set Members

NSM = The number of super modules. (NSM : 12)

NFUEL = The maximum number of'fuel types. (NFUEL • 5)

NCOM = The number of entries for the construction of the adjacency portion of the communication matrix ICOMM. See

Section 5.8. (NCOM .: 3.0)

NNCOM The number of non-communicating entries for the construction of ICOMM. See Section 5.9. (NNCOM : 200)

NPILOT The number of fuel cells initially on fire.

(NPILOT : 10)

IROOM = The indicator variable which determines if compartment

data is to-be. included or not. This latter option is

useful when the fire is sufficiently small with respect

to the room such that enclosure effects can be neglected.

IROOM=l =>compartment data is 'input

IROOM01:=> no compartment data input

INITG = The indicator variable which shows if the hot gas layer

has some initial characteristics and/or if the initial

external heat fluxes are-non-zero. Input ,only if

-' 14 -
IROOM =1.

INITG=I => Initial values for gas.layer and/or

Sheat fluxes are 'input

'INITG#1 => TG (Gas Layer Temperature) = 298 0K

DG (Gas Layer Thickness) = 0.0 m

QEXT (External Heat Fluxes

For All Fuel Elements) = 0.0 W/m2

5.4 FUEL BED DEFINITION (&FUELB)

This data block must- be entered once for each super module, i.e. NSM

times. The first entry corresponds to super module number 1, the second

to module 2, etc.

Data Set Members

SMX = x-coordinate midpoint of super module (m).

SMY = y-coord.inate.midpoint of super module (m).

SMZ = z-coordinate-midpoint of super module (m).

SLNG = Super module length (m).

SWID Super module width (m).

SDEP = Super module depth ((m).

SMASS = Super module mass (kg).

SPOR = Super module porosity factor f (Dimensionless). P

SLOSS = Undefined variable.

NFCL = Number of fuel cells in super module.

IORNT = Orientation-of super module.

IORNT=I => Normal axis is in the x-direction.

- 15 -
IORNT=2 => Normal-axis is: in the y-direction..

IORNT-3-> Normal axis is-in the z-directi6n. -

IDIREC Direction of.:super module.

IDIREC=l => Long-axis is in the x-direction. "

IDIREC=2 > Long axis is in the y-direction.

IDIREC=3 =>Long~axis is in the z-direction.

IFTYP = Super module fue-l name/identification number.

IFTYP = Super module fuel type location! identifier. This parameter identifies the particular location in the array

IFUEL of. NAMELIST.block FUELT (Section 5-6) which correspconds to the fuel for the super module described.

Example: "IFTYP=31' means that the super,,module fuel is type,

3,:where fuel type 3 is defined by position 3 of array IFUEL. Note that IFTYP < NFUEL.

COMPBRN will automatically compute the 'locations of the fuel cells in

each super module. It.will also number the fuel: cells, the- numbers in-.

creasing withi'distance in nthe direction. of-:. the long axis,.. Thus, i-fthe

super module is directed along the x-axis, a fuel cell at x=1 will have

a lower number.than a fuel cell atý x=2.

5.5 PILOT FUEL PARAMETERS (,PILOT)

COMPBRN assumes that a compartment fire starts :on top of a :small number

of selected fuel cells 'in the fuel bed. These. initial fires, called pilot fires, may involve materials not. comprisingthe main fuel bed. For

example, a segment of cable tray may -be overlaid by an amount of exter.-.

- 16 -..
nal fuel whidh .s'used' to initiate the compartment fire. The location

of this. pilot, fire -is, however, treated as being identical with ýthe; location of. the'. underlying. fuel'..cell.

Each item in this data block.is an array, with.the first. entry corresponding to the first -pilot fuel cell, the second with the second

cell, and so forth., NPILOT entries should be made for each. array.

Data Set Members'.

IPIL = Array listing the super modules of the pilot fuels.

JPIL = Array-listing, the particular fuel-cells (within the

given super module) which are assigned'pilot fires.

Example: If NPILOT=2, IPIL=(3,2), and JPIL=(4.5), this means'

that:there are two pilot fires, the first being located onthe top of-the 4th fuel cell of the 3rd super module, and the second.being ontop of the 5th

fuel cell of the 2nd super module.

IPFUEL = Pilot:fuel type array. This entry corresponds to IFTYP

'defined in Section 5.4. Fuel types are numbered consecutively. Thus, if the 3rd entry of IPFUEL is '5', this

means that the 3rd pilot fire involves fuel type 5.

PMASS Pilot fuel mass array (kg).

- 17 -
5.6" PHYSICAL'FUEL"PARAMETERS (&FUELT).

Each of the items in this data "block represeents an array -of dimension

NFUEL. The :first entry -in each array corresponds t61:"el.type 1, the

second entry corresponds to fuel type 2, and so forth.

Data Set -Members

IFUELI = Array of idendification numbers .(names) assigned to

each fuel type.

0 < IFUEL(I) < 0 =>.'Super module is combustible-fuel.

10 <IFUEL(I): < 20':=> Super module is a detector (use20

.30

40

IFUEL(I)

5 IFUEL(I)

• IFUEL(I)

ful for measuring heat fluxes).

30 => Super -module is a wall.

40 => Super module is a thermal barrier.

50 => Super module is the room ceiling.

Note that. the room ceiling can only be modeled with one super module,

and that those super modules modeling the ceiling and the walls may only

have one fuel cel.l each.

In the following list, entries must be provided for all fuel types,

but realistic entries are optional for the ranges of- IFUEL noted below

(where for notational simplicity, we let K IFUEL(I)).

=3

DENS = Density (kg/imn). (10 - K < 20 and 30 < K < 40)

SPHT = Specific heat (J/kg0 K). (1i0:< K- < 20 and 30 < K < 40,)

THK = Thermal conductivity (W/m°K). (10 < K < 20)

HEAT = Heating, value (J/kg). (10 < K < 50)

'FIGTP = Piloted ignition,.temperature ( 0 K). (10 < K < 50)

FIGTS :..-Spontaneous ignition temperature ( 0 K). (10 < K < 50)

BRATV = Ventilation controlled burning-rate

- 18
constant C (Dimensionless). (10 : K < 50) V

BRATSO = Surface controlled burning rate-constant

.." (kg/m2 s). (10< K< 50) 0 -

:BRATS1 Surface controlled burning rate constant

• C (kg/m2 s). (10 • K < 50) S

GAMMA = Fraction of. flame heat released `as

radiation s.(Dimensionless). (10 l K < 50)

FABSRP = Absorption coefficient for flame gases

(m-'). (10 5K`< 50)

REFL = Reflectivity (Dimensionless). (0 S K'< 20)

5.7 MISCELLANEOUS DATA (&MISC)

Data Set Members

RTEMP = Room-temperature (OK). (Default = 298.)

FLCF = Heat transfer coefficient for heat transfer:in a

flame (W/m2 0 K).

5.8 COMMUNICATION DATA (&SEE)

This data block is used to input information defining 'if a given fuel,

cell touches another fuel cell. "This affects the heat transfer models

used, transfer between contiguous cells being greater than transfer between disconnected cells.

- 19 -
COMPBRN .automatically establishes adjacency for consecutive fuel

cells in each super module. Furthermore, by default, all fuel cells can

transmit heat to' all' other fuel ,cells (i.e. they "communicate"). The

purpose of this data block is to.define some exceptions to these general

rules (non-communication. ishandled in Section 5.9). NCOM entries must

be.imade, one for each eexception.

Data Set Members,

IV = 4 element -array containing adjacency data. If fuel cel

jf.of super. module i is adjacent to fuel cell lI6f super

module kjIV would be (i,j,k,l). If super modules i and

k are parallel, have the same number. of fuel cells, are

of the same lengthi and are, adjacent, IV would be

(i,999,.k,0).. The '999',tells COMPBRN that the jth cell

of super; module i is adjacent to the jth cell of'super

module: k, i., j <.SNFCL. Note that the adjacency relationship is .symmetrical. Two separate entries, (i,j,k,l)

and (k,l,i,j) are not. required; only enter one *of the two..

5.9 NON-COMMUNICATION DATA. (&NSEE)

This data block defines which fuel cells cannot transmit heat to certain

other fuel cells. It is treated similarly to the preceding data block,

NNCOM- entries-being required.

Data Set Members

NV -=4 element array containing non-communication data. If

20
fuel cell j of super module-i cannot transmit heat to,

fuel cell 1 of super'mOdule'k,:NV would be (i,j,k,l)."

If all of the fuel cells in super 'module i cannot transmit heat to all of the, fuel"' cells, in super' module k, NV

would be ý(i,999,k,0). -If all of the fuel :cells in super

module i cannot transmit heat'to al:l.,of the fuel cells

in super-modules j through j', NV would be (i,888,j,j').

Note that the non-communication, relationship is asymmetrical.

5.10 ROOM DATA (&ROOM)

This data block should be enteredonly if IROOM •1. DefaUlt Values are

assigned for some of the parameters.

Data Set Members

DCF Coefficient ofdischarge for doorway', (Dimensionless).

(Default.= 1.0)':

-DHGT Height,:of doorway (m)-.

DWID Width of doorway (m).ý_

FC Fraction of forced-ventilation inflow which enters the

gas layer (Dimensionless).

FH = Fraction of forced ventilation outflow which leaves the

gas layer (Dimensionless).

GABSRP = Absorption coefficient for ceiling gas layer (m-').

(Default = 1.4)

HCEIL = Heat transfer coefficient for ceiling (W/m20 K).

- 21-
PLCF1: = Buoyant plume ,entrainment coefficient, flame tip is

above neutral plane (Dimensionless).

(Default =i •0)

PLCF2 =.Buoyant plume entrainment coefficient, flame tip is

below neutral plane .(Dimensionless).,

(Default = 05),,

-THETA = Convergence acceleration parameter used in calculation

of neutral plane height ZD.

(Default = 0.5)

VFV = Forced ventilation volumetric airflow (m 3 /s).

5.11 INITIALIZATION-DATA, (&GINIT)

This -data block should.be.entered only if IROOM 1 and INITG = 1.

Data Set Members

TG = Initial gas layer temperature (*K).

DG = Initial gas layer thickness (m).

QEXT = Two dimensional array containing initial external heat

fluxes impinging on each fuel cell, the (i,j) entry

corresponding to fuel cell j of super module i.

- 22 -
5.12 MODEL VARIATION FACTORS (&MODVAR) "

This data. block allows the user, to multiply.the results, of various models by a specified factor. Unless otherwise specified, all factors .are

set to 1.0.

Data Set Members

FCTR 15 element array containing the modification Ifactors for

15 models:

1) Modifies room-averaged value for C (ventilation con- V

trolled burning rate constant).

2) Modifies surface controlled burning rate of :fuel.

3) Modifies

4).Modifies

5) Modifies

6) Modifies

7) Modifies.

'faces in.

8) Modifies

surfaces

9) Modifies

10) Modifies

bed back

11)'Modifies

flame height over horizontal fuel slabs.

flame height over vertical fuel slabs.

radiative heat flux received by target element.

buoyant plume temperature. above flame.

heat.transfer' coefficient for .vertical sur-

.the flame.

heat transfer coefficient for horizontal

in the .flame.

convective heat flux received by'each element.

'heat transfer,.from flame over vertical fuel

to fuel bed.

conductive heat flux received.

12)

13)

14)

15)

Undefined.

Undefined.

Modifies critical time for spontaneous ignition.

Modifies critical time for piloted ignition.

- 23 -
'5.13

This

(see

OUTPUT DEFINITiON :(&OUTF) -

data block Aallows the user to, specify -the printed outputdesired

.Section 4).'.

Data Set Members

INCHCK '= 'Indicator-..variable controlling the amount of input data'

printed.'-

INCHCK=O => Only heading material and FCTR.

INCHCK=I=>. Above, plus fuel parameters

(Section 5.6).

INCHCK=2 > All input datai -

IOUTPT = Number of output variables to be printed (IOUTPT • 11).

MOUTPT •' Array whore it'h,•element is the identification number

of the output- variable desired'. .OUTPT -elements should

be input.

I."D. -Numbers: l)_,Total mass burning.rate:.TMDOT (kg/s).

2):..Total heat release rate TQDOT (W)'.

3) Gas layer temperature TG (OK).

"4) Gas layer-thickness DG (M).

5) Fuel cell indicator variable for burning

$BURN; $BURN(ij)=.TRUE. =>Jfuel cell j of

super module i is' burning.

6) Fuel cell source heat'flux QDOT2P (W/m2 ).

7) Fuel cell external heat flux"QEXT (W/m2 ).

8).Fuel cell mass: FMASS (kg)'..

9) "Flame 'height over fuel cell FJh"T (m) .

-.10): Fuel cell-.flame temperature FLTEMP (°K).

•"24 -" .
11)' Fuel cell surface temperature TEMP -(OK)

(useful only for walls, barriers, and ceiling).

NSMOUT = Number of super modules for which variables (I.D. numbers

5 through 11) are to be printed.

MSMOUT Array whose ith element. is the identification number of

-the super module.to be documented.

For example', if the input card is:

&OUTPT IOUTPT=-3, MOUTPT=-l,6,7, NSMOUT=-2, MSMOUT=-h,2, &END

COMPBRN will print the total mass burning rate for the

compartment, and the source and external heat fluxes for

each fuel cell in super modules 1 and 2.

Notes: 1) IOUTPT=-O and NSMOUT=O => Only TMDOT, the total mass

burning:rate, is printed.

2) IOUTPT=O and NSMOUTO0 => TMDOT and $BURN, an array

'which shows which fuel cells are burning, .are printed.

- 25 -
Chapter 6

A SAMPLE PROBLEM

In this section, we use COMPBRN to predict the behavior of a wood crib

fire in a compartment. The wood crib consists of three layers of three

sticks, as shown in Figure 3,. and is located on the floor in the middle

of the room. The room itself is shown in Figure 4. The room walls and

ceiling are concrete, and have no openings..other than the 2.4 m high by

1.0 m wide doorway in one wall.. There is no forced ventilation airflow

into the room. The fire is initiated on the central stick of the bottom

layer in the crib.

We model the crib with 9 super modules, as seen in Figure 5. Each

super module consists of 5 fuel cells. To improve our analysis of the

fire, we make 3 additional modeling assumptions:

1. The porosity factor of each fuel super module is set

equal to 3.0. This accounts for the'surface area of the

wood stick sides which are available for burning.

2. The feedback coefficient, BRATS1 (C in Reference'(1)),

S

is set equal to 0 for the bottom two layers. Our modeling thus prohibits enhanced *burning within the crib,

since the crib latticework restricts the amount of fresh

.air which can reach the fire.

3. The bottom super :modules (numbers 1, 2, and 3) cannot

transmit large amounts of heat directly to the uppermost

- 26 -
12

17-

-

20 

-

Figure 3: Wod Crib

I
0z

x X- 13.1 ft

/ Y- 19.7 ft

Z - 13.1 ft

IXv -3.30 ft

SZD 7.90 ft

Figure 4: Compartment Housing Crib
I0

T

I'

zz

"4" 

- '-, -"

20 

I

IFigure 5: Model of .ood Crib
super modules (numbers 7,-8, and 9) because the intermediate wood sticks intervene. Furthermore, each layer of

sticks cannot transmit significant amounts of heat to any

underlying layers.

These assumptions have proven useful in the earlier wood crib

;fire simulations described in Reference (1). More importantly,

they illustrate how the porosity-..factor.SPOR, the communication

matrix ICOMM, and the artificial definition of. various fuel

types can be used to improve our modeling.

The enclosure is modeled with two super modules: one represents the ceiling';and the, other represents one of the walls. To

measure the heat flux contribution from the wall (to see if more

walls should be added to .the model), we also place a heat flux

detector near the wall, and use the communication matrix to prevent the :.detector from seeing any heat. sources other than the

wall.

The problem input and output .'are respectively given in Ta-

,:ble 1 and Table 2. A graph of the crib's total mass burning

rate is presented in Figure 6, andjthe heat flux from the wall

to the detector is plotted in Figure 7. Finally, the average

hot gas layer temperature'is plotted in Figure 8.

-'30 -
TABLE 1

COMPBRN Sample Problem Input

&STRT NJOB=l, NTIME=30,1 NREAD=5, NWRITE=6, DELT-60., &END

&SIZE NSM=12, NFUEL=5, NPILOT=-l, NCOMO0, NNCOM=29,. IROOM"1,

INITG=O, &END

S********* TEST CASE: WOOD CRIB IN ENCLOSURE **** * **

&FUELB SMX=2.0, SMY -2.7968, SMZ=0.0508, SLNG=O.508,

SWID=. 1016, SDEP=0.1016, SMASS=4.20, SPOR=3.0O,

SLOSS=l., NFCL=5, 'IORNT=-3,; IDIREC=l, IFTYP=l, &END

&FUELB SMY=3.0, &END;-:,

&FUELB SMY=3.2032, &END .

&FUELB SMX=l.7968, SMY=3.0, SMZ"0.1524, IDIREC=2,/ &END '

&FUELB SMX=2.0, &END,,.,.

&FUELB SMX=2.2032, &END

&FUELB SMX=2.0, SMY=2.7968, SMZ=0.2540, tIFTYP=2, IDIREC=l, &END

&FUELB SMY=3.0, &END

&FUELB SMY=3.2032, &END

&FUELB SMX=2.0, SMY=0.0, SMZ=2.0, SLNG=4.0, SWID=4.0,

SDEP=0.1, SMASS=1.0, SPOR=l.0, SLOSS=I.0,

NFCL=l, IORNT=-2, IDIREC=l, IFTYP=3,, &END

&FUELB SMY=3.0, SMZ=4.0, SLNG=6.0,

IORNT=3, IDIr-EC=2, IFTYP=4, &END

&FUELB SMY=0.1, SMZ=2.0, SLNG=0.1, SWID=0.1,

SSLOSS=l.0, IORNT=2, IDIREC=l, IFTYP=5, &END

&PILOT IPIL=2, JPIL=2, IPFUEL=l,; PMASS=O.l, &END

&FUELT IFUEL=1,2,20,40,10, DENS=2*800. ,3*"000., SPHT=2*1045. ,3*2000.,

TK=2*0.092,3*I., BRATV=3*0., FABSRP=2*l.4,3*0., REFL=5*.I,

BRATS1=0.,5.53E-7,3*0., GAMMA=2*.4,3*0.,

BRATSO=2*..0062,3*0,,

HEAT=2*I.86E7,3*0.,

FIGTP=2*480. ,3*0.,

FIGTS=2*530. ,.3:*0. , &END

&MISC RTEMP=298. ,:FLCF=22., &END

&NSEE NV=l,999,3,0, &END

&NSEE NV=,,999,9,0, &END

&NSEE NV=1,888-,1,12, &END

&NSEE NV=21,888,10,12, &END

&NSEE NV=3,,999,l,0, &END

&NSEE NV=3,999,7,0, &END

&NSEE NV=3,888,1.0-,12, &END

&NSEE NV=4,888,1,3,' &END,

&NSEE NV=4,999,10,0i &END

&NSEE NV=4,999,12,0, &END

&NSEE NV=5,888,1,3, &END

&NSEE NV=5,999,10,0, &END

&NSEE NV=5,999,12,0, &END

&NSEE NV=6,888,1,3, &END

&NSEE NV=6,999,10,0, &END

&NSEE NV=6,999,12,0, &END

&NSEE NV'7,888,i,6, &END

&NSEE NV=7,'.999,12,0, &END

- 31 -
&NSEENV=8,888,1,6, &END

&NSEE NV=8,)999,10,0, &END

&NSEE NV=8,999,12,0, &END

&NSEE NV=9,888,1,6, &END

&NSEE NV=9,999,10,0, &END

&NSEE NV9,999,12,0, &END

&NSEE NV=I0,888,2,6, &END

&NSEE NV=10,888,8,9, &END

&NSEE NV=1I;888,1.,6,i &END-

&NSEE NV=11,999,12,0, &END

&NSEE NV=12,888,'1,11 &END

&ROOM.DCF=1.0, DHGT=2.4, DWID=1.,0 FC=0.0,'FH=0.0, GABSRP=1.4,

HCEIL=3.0, PLCF1=1.O, .PLCF2=1.5, THETA=0.5, VFV=0.0, &END

&MODVAR FCTR=15*1., &END,

&OUTF INCHCK=2, I0UTPT=7, MOUTPT=1,2,3,4,5,6,7,

NSMOUT=12, MSMOUT-I,2,3,4,5,6,7),8,9, 10,1 1 ,12, &END

- 32 -
TABLE 2

COMPBRN Sample Problem Output

PROGRAM COMPBRN - A DETERMINISTIC CODE TO COMPUTE THE PROGRESS OF FIRE

OVER A GIVEN FUEL ARRAY WITHIN ENCLOSING,.

BOUNDARIES. ALL UNITS ARE IN THE MKS SYSTEM.

*~~******** TEST CASE: WOOD CRIB IN ENCLOSURE *

INPUT DATA:

JOB 1 OF 1 JOBS

- 331'.
VARIABILITY FACTORS FOR FIRE MODELS:

VENTILATION CONTROLLED BURNING RATE

FUEL-SURFACE CONTROLLED BURNING RATE

FLAME HEIGHT FOR HORIZONTAL FUEL

FLAME HEIGHT FOR VERTICAL FUEL

RADIATIVE HEAT FLUX INTERCHANGE

BUOYANT PLUME TEMPERATURE

CONVECTIVE HEAT TRANSFER COEFFICIENT FOR FLAME

CONVECTIVE HEAT TRANSFER COEFFICIENT FOR PLUME

CONVECTIVE HEAT FLUX

HEAT TRANSFER TO SELF FOR VERTICAL FUEL

HEAT TRANSFER TO ADJACENT FUEL

HEAT FLUX FROM CEILING HOT GAS LAYER

HEAT FLUX FROM REFLECTIONS OFF. WALLS

CRITICAL TIME TO SPONTANEOUS IGNITION

CRITICAL TIME TO PILOTED IGNITION

1.0000

1.0000

1.0000

1.,0000

1.0000

1.0000

1.0000

1.0000

1.0000

1.0000

1.0000

1.0000

1.0000

1.0000

- l•34 .. -I• 
NUMBER OF FUEL TYPES: 5

DATA FOR FUEL TYPE 1:

FUEL'TYPE: 1

DENSITY (KG/M-*3):

SPECIFIC HEAT (J/KG-K):

THERMAL CONDUCTIVITY (W/M-K):

HEAT OF COMBUSTION (J/KG):ý PILOTED IGNITION TEMPERATURE (OEG. K):

SPONTANEOUS IGNITION TEMPERATURE (DEG. K):

VENTILATION CONTROLLED BURNING RATE FACTOR (KG/M--2.5-S):

SURFACE CONTROLLED SPECIFICOBURNING RATE (KG/M**2-S):

SPECIFIC BURNING RATE RADIATION AUGMENTATION (KG/J-M**2):

FRACTION OF HEAT RELEASED AS RADIATION:

SMOKE ATTENUATION FACTOR (M**-l):

REFLECTIVITY:

DATA FOR FUEL TYPE 2:

FUEL TYPE: 2

DENSITY (KG/M**3):

SPECIFIC HEAT (J/KG-K):

THERMAL CONDUCTIVITY (W/M-K):

HEAT OF COMBUSTION (J/KG):

PILOTED IGNITION TEMPERATURE (DEG. K):

SPONTANEOUS IGNITION TEMPERATURE (DEG. K):

VENTILATION CONTROLLED BURNING RATE FACTOR (KG/M**2.5-S): SURFACE CONTROLLED SPECIFIC BURNING RATE (KG/M**2-S):

SPECIFIC BURNING RATE RADIATION AUGMENTATION (KG/J-M**2):

FRACTION OF HEAT RELEASED AS RADIATION:

SMOKE ATTENUATION FACTOR (M**-l):

REFLECTIVITY:

DATA FOR FUEL TYPE 3:

FUEL TYPE (WALL): 20

THERMAL CONDUCTIVITY (W/M-K):

THERMAL DIFFUSIVITY (M**2/S):

REFLECTIVITY:

DATA FOR FUEL TYPE 4:

FUEL TYPE (CEILING): 40

THERMAL CONDUCTIVITY (W/M-K):

THERMAL DIFFUSIVITY (M**2/S):

REFLECTIVITY:

DATA FOR FUEL TYPE 5:

FUEL TYPE (DETECTOR): 10

800.

1045.

0.092

.1860E+08

480.

530.

.0

.6200E-02

.0

0.400

1.400

0.100

800.

1045.

0.092

.1860E+08

480.

530.

.0

.6200E-02

.5530E-06

0.400

1.400

0.100

1.000

.5000E-06

0.100

1.000

.5000E-06

0.100

- 35 -
ROOM PARAMETERS:

CEILING LENGTH.'WIDTH, HEIGHT (M):

DOOR HEIGHT, WIDTH (M), ORIFICE COEFFICIENT:

FORCED.VENTILATION.(M**3/S):

FORCED VENTILATION CONSTANTS (FH AND FC):

PLUME ENTRAINMENT CONSTANTS (PLCF1 AND PLCF2): HOT GAS ABSORPTION COEFFICIENT (M**-1):

NUMERICAL ACCELERATION PARAMETER (THETA):

CEILING HEAT TRANSFER COEFFICIENT (W/M**27DEG.K):

6.00 *-4.00 4.00

2.40 -,1.00 1.00

0.0

0.0 0.0

1.00 1.50

1.4000

0.50

.3.000

-36-
PILOT FIRE DATA:

NUMBER OF PILOT FIRES: 1

DATA FOR PILOT. FIRE I.:

LOCATION-OF PILOT FIRE (FUEL ARRAYFUEL CELL): 2, 2

FUEL TYPE: 1

MASS (KG): 0.100.

- 37- -
NUMBER OF FUEL ARRAYS: 12

DATA FOR FUEL ARRAY 1:

.INITIAL FUEL CELL COORDINATES (XY,Z) (M):

DIMENSIONS (LENGTH,WIDTH,DEPTH) (M):.

NUMBER OF FUEL CELLS:

DIRECTION OF AXIS, ORIENTATION:

MASS (KG):

POROSITY FACTOR (DIMENSIONLESS):

HEAT LOSS FACTOR (UNDEFINED):

FUEL TYPE:

DATA FOR FUEL ARRAY 2:

INITIAL FUEL CELL COORDINATES (XY,Z) (M):

DIMENSIONS (LENGTH,WIDTH,DEPTH) (M):

NUMBER OF FUEL CELLS:

DIRECTION OF AXIS, ORIENTATION:

MASS (KG):

POROSITY FACTOR (DIMENSIONLESS):

HEAT LOSS FACTOR (UNDEFINED):

FUEL TYPE:

DATA FOR FUEL ARRAY 3:

INITIAL FUEL CELL COORDINATES (X,YZ) (M):

DIMENSIONS (LENGTH,WIOTH,DEPTH) (M):

NUMBER OF FUEL CELLS:

DIRECTION OF AXIS, ORIENTATION:

MASS (KG):

POROSITY FACTOR (DIMENSIONLFSS):

HEAT LOSS FACTOR (UNDEFINED):

FUEL TYPE:

DATA FOR FUEL ARRAY 4:

INITIAL FUEL CELL COORDINATES (X,YZ) (M):

DIMENSIONS (LENGTH,WIOTHDEPTH) (M):

NUMBER OF FUEL CELLS:

DIRECTION OF AXIS, ORIENTATION:

MASS (KG):

POROSITY FACTOR (DIMENSIONLESS):

HEAT: LOSS FACTOR (UNDEFINED):

FUEL TYPE:

DATA FOR FUEL ARRAY 5:

INITIAL FUEL CELL COORDINATES (X,Y,Z) (M):

DIMENSIONS (LENGTH.WIDTH,DEPTH) (M):

NUMBER OF FUEL CELLS:

DIRECTION OF AXIS, ORIENTATION:

MASS (KG):

POROSITY FACTOR (DIMENSIONLESS):

HEAT LOSS FACTOR (UNDEFINED):

FUEL TYPE:

DATA FOR FUEL ARRAY 6:.

INITIAL FUEL CELL COORDINATES (X,ýyY,Z) (M):

DIMENSIONS (LENGTH,WIDTH,DEPTH) (M):

NUMBER OF FUEL CELLS:

DIRECTION OF AXIS, ORIENTATION:

MASS (KG):

POROSITY FACTOR (DIMENSIONLESS):

HEAT LOSS FACTOR (UNDEFINED):

1.80 2.80

0.51 0.10

•X

1.80

0.51

3.00

0.10

X

1.80 3.20

0.51 0.10

X

1.80 2.80

0.51 0.10

Y

2.00 2.80

0.51 0.10

Y

2.20 2.80

0.51 0.10

Y

0.05

0.10

5

z 0.84.O0

3.0000

1.0001

0.05

0.10

5

z 0.8400

3.0000

1.0001

0.05

0.10

5

z 0.81400

3.0000

1 .000

:: 1

0.15

0.10

5

z 0.8400

3.0000

1.000

1

0.15

0.10

5

z 0.8400

3. 0000

1.0001

0.15

0.10

5

z 0.8400

3.0000

1.000

-7,38 -
FUEL TYPE:

DATA FOR FUEL ARRAY 7:

INITIAL FUEL CELL COORDINATES (X,Y,Z) (M):.

DIMENSIONS (LENGTH,WIDTH,DEPTH) (M):

NUMBER OF FUEL CELLS:

DIRECTION OF AXIS. ORIENTATION:

MASS (KG):

POROSITY FACTOR (DIMENSIONLESS):

HEAT LOSS FACTOR (UNDEFINED):

FUEL TYPE:

DATA FOR FUEL ARRAY 8:

INITIAL FUEL CELL COORD.INATES (X.Y,Z) (M):

DIMENSIONS (LENGTH,WIDTH,DEPTH) (M):

NUMBER OF FUEL CELLS:.

DIRECTION OF AXIS, ORIENTATION:

MASS (KG):

POROSITY FACTOR (DIMENSIONLESS):

HEAT LOSS FACTOR (UNDEFINED)::

FUEL TYPE:

DATA FOR FUEL ARRAY 9:

INITIAL FUEL CELL COORDINATES (X,Y,Z) (M):

DIMENSIONS (LENGTH,WIDTH,DEPTH) (M):

NUMBER OF FUEL CELLS:

DIRECTION OF AXIS, ORIENTATION:

MASS (KG):

POROSITY FACTOR (DIMENSIONLESS):

HEAT LOSS FACTOR (UNDEFINED)-

FUEL TYPE:

DATA FOR FUEL ARRAY 10:

INITIAL FUEL CELL COORDINATES (X,YZ) (M):

DIMENSIONS (LENGTH.WIDTH,DEPTH) (M):

NUMBER OF FUEL CELLS:

DIRECTION OF AXIS, ORIENTATION:

MASS (KG)::

POROSITY FACTOR (DIMENSIONLESS):

HEAT LOSS FACTOR (UNDEFINED): -

FUEL TYPE:

DATA FOR FUEL ARRAY 11:

INITIAL FUEL CELL COORDINATES (X,Y,Z) (M):

DIMENSIONS (LENGTH,WIDTH,DEPTH) (M):

NUMBER OF FUEL CELLS:

DIRECTION OF AXIS, ORIENTATION: "MASS (KG):

POROSITY FACTOR (DIMENSIONLESS):

HEAT LOSS FACTOR (UNDEFINED):

FUEL TYPE:

DATA FOR FUEL ARRAY 12:

INITIAL FUEL CELL COORDINATES (X,Y,Z) (M):

DIMENSIONS (LENGTH,WIDTH,DEPTH) (1):

NUMBER OF FUEL CELLS:

DIRECTION OF AXIS, ORIENTATION:

MASS (KG):

POROSITY FACTOR (DIMENSIONLESS):

HEAT LOSS FACTOR (UNDEFINED):

1

1.80 2.80

0.51 0.10

X

1.80 0. 51 3.00

0.10

X

1.80 3.20

0.51 0.10

-X

0.25

0.10

5

z 0.8400

3.0000

1.000

2

0.25

0.10

5

z 0.8400

3.0000

1.000

2

0.25

0.10

5

z 6.8400

3.0000

1.000

2

2.00

0.10

1

Y

1.0000

1.0000

1.000

3

4.00

0.10

1

z 1.0000

1.0000

1.000: 4

2.00

0.10

1

Y

1.0000

1.0000

1.000

2.00 0.0

4.00 4.00

X

2.00 3.00

6.00 4.00

Y

2.00 0.10

0.10 0.10

X

- 39 -
FUEL TYPE:

- 40 -
COMMUNICATION DATA:

NON-COMMUNICATION DATA:

NUMBER OF NON-COMMUNICATING FUEL CELLS: 29

FUEL CELL PAIR

1

2,

3

4

5

6

7

8

9

1011

12

13

14

15

16

17

18

19

20

21

22

23

24

25

26

27

28

29

(FUEL ARRAYFUEL CELL) (FUEL ARRAYFUEL CELL)

1 999 3 '0

1 999 9 0

1 888 11.- 12

2 888 10 12

3 999 1 0

3 - ,.ý999 7 0

3 888 .10 12

4 888. 1 3

4 999 10 0

4 999 12 0

5 888 1 3

5 999 10 0

5 999 12 0

6 888 1. 3

6 999 10 0

6 999 12 0

7 888 1 6

7 999 12 0

8 888 1 6

8 '999 10. . 0

8 999 12 0

9 888 1 6

9 %999 10 0

9 .999 12 0

10 888. 2 6

10 888 8 9

11 888 1 6

11 999 12 l0

12 .888 1 11

-41 -
MISCELLANEOUS:

ROOM TEMPERATURE (DEC. K): 298.

CONVECTIVE HEAT TRANSFER COEFFICIENT FOR FLAME (W/M*-2-K): 22.00

TIME INCREMENT (S): 60.00

NUMBER OF TIME STEPS FOR JOB: 30

- 42 -
*********** TEST CASE: WOOD CRIB IN ENCLOSURE "

TIME (SEC).: 0.

TOTAL MASS BURNING RATE (KG/S):

TOTAL HEAT RELEASE RATE (W)

HOT GAS LAYER TEMPERATJRE (K)

HOT GAS LAYER THICKNESS (M)

. 1920E-03 3571.

298.1

1.600

MODULE:

''MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE: SMODULE:

1

2

3

4

5

6

7

8

9

10

11

12

1

2

3

4

5

6

7

8

9

10

11

12

BURNING 7

BURNING ?

BURNING ?

BURNING ?

BURNING ?

BURNING ?

BURNING ?

BURNING ?

BURNING ?

BURNING 7

BURNING ?

BURNING ?

F

F

F

F

F

F

F

F

F

F

F

F

F

T

F

F

F

F

F

F

F

F

F

F

F

F

F

F

F

F

F

F

F

F

F

F

F

F

F

F

F

F

F

F

F

F

SRCE

SRCE

SRCE

SRCE

SRCE

SRCE

SRCE

SRCE

SRCE

SRCE

SRCE

SRCE

FLUX

FLUX

FLUX

FLUX

FLUX

FLUX

FLUX

FLUX

FLUX

FLUX

FLUX

FLUX

.0

.0

.0

.0

.0

.0

.0

..0

.o0

403.

448.

.0

.0

. 306E+05

.0

.0

.0

.0

.0

.0

.0

.0

.0

.0

.0

.0

.. 0

.0

.0

.0

.0

.0

.0

.0

.0

.0

.0

.0

.0

-.0

.0,

.0

"0., .

-.0

.0

.0

.0.0"

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

1

-2

3

4

5

6

7

8

9

10

11

12

-EXT

EXT

EXT

EXT

,EXT

EXT

EXT

EXT

EXT

EXT

EXT

EXT

FLUX

FLUX

FLUX

FLUX

FLUX

FLUX

FLUX

FLUX

FLUX

FLUX

FLUX

FLUX

j,135E404, 126 E+05

130E+04

109E+04

516.

.128E+04

.. 317E+04

.124E+0O4

98.1

46.2

402.

.153E+0*4

,.0

1.148E+0*4

.201 E+014

.201 E+014

6142. ý

.1148E+014

.308E+05

.143E+0k4

.135E*O*4

.126E+05

... 130E+014,

.300 E+0*4

.300E+04

699r.I

.128E+04,

.317E40*4

.12*4E+04

.103E+04

.148E+04

978.

.201E+04

.201E+04

642.

947.

1143 E+04

.898.

757.: 907.

707.

.109E+04

.109E+04

516.

691.

828.

642.

- 43 -
TEST. CASE: WOOD CRIB IN ENCLOSURE ****.*.****,****. *

TIME (SEC)': 60.

TOTAL MASS BURNING RATE (KG/S):

TOTAL HEAT RELEASE RATE (W)

HOT GAS LAYER TEMPERATURE (K)

HOT GAS LAYER THICKNESS (M)

*.1920E-03 .3571.

298.1

1.600

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:,

MODULE:

MODULE:•

MODULE:.

MODULE:-

.MODULE:

MODULE:

MODULE::

MODULE: -

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

1

2

3

4

5

6

7

8

9

10

11.

12

1

2

3

4.5,

6

7

8

9

10

11

12

BURNING ?

BURNING ?

BURNING 7

BURNING 7

BURNING ?

BURNING ?

BURNING ?

BURNING ?

BURNING ?

BURNING ?

BURNING ?

BURNING ?

F

T

F

F

F

F

F

F

F

F

F

F

F

T

F

F

F

F

F

F

F

T

F

F

F

F

F

F

'F

F

'F" F"

F

F.

F

F

F

F

F

F

F

F

F

F

F.

F'

SRCE

SRCE

S SRCE

_SRCE

SRCE

SRCE

SRCE

SRCE

SRCE

SRCE

SRCE

•SRCE

FLUX

FLUX

FLUX

F LUX

FLUX

FLUX

FLUX

FLUX

FLUX

FLUX

FLUX

FLUX

.0

.0"

.0

.0

.0

.0

.0

.0

.0

414.

448..'

..0 ,

.306E+O5;.."- .0

.0 .0

.0 ..0

.0 .0

.0 .0

.0 .0

.0 .0

0. .0

.0

.0

.0

.0

.0

.0

.0.0

.0

.0

'.0

.0.0

.0

.0

.0

.0

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

1 EXT

2 EXT

3" EXT

4. EXT

5 EXT

6, EXT

7 , ,IEXT

8 'EXT

9 EXT

10 EXTi

11. EXT

12' EXT

FLUX

'FLUX.

'FLUX'

FLUX FLUX

FLUX

FLUX .FLUX

FLUX..

FLUX

FLUX

FLUX

.135E+04'

126E+05',ý

.130E+04

.109E+04

.109E+04

_516.

.- 129E+04 .. 317E+04

124E404

,ý98.+1 ' '47.1

413.

153E+04

.148E+04

201E+04

201E+04!

642.

'.148E+04,

.308E+05!

143E+04,

.135E+04 .103E+04

.126E+05•148E+04

.130E+04 978.

.300E+04 .201E+04

.300E04ýý .201E+04

699. 642.

.129E+041 948.

.317E+04 .143E+04

1240E+04-898.'

759.

907.

707.

109E+O4

109E+04

516.

692.`

828.

642...-

- •44,-.
******** TEST CASE: WOOD CRIB IN ENCLOSURE *************************

TIME (SEC): 120.

TOTAL MASS BURNING RATE (KG/S):

TOTAL HEAT RELEASE RATE (W).:

HOT GAS LAYER TEMPERATURE (K)

HOT GAS LAYER THICKNESS (M)

.1295E-02 -

.2409E+05

299.0

1.600

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:'

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

1

2

3

4

5

6

7

8

9

10

11

12

1

-2

3

4

5

6

7

8

9

10

11

12

BURNING

BURNING

BURNING

BURNING

BURNING

BURNING

BURNING

BURNING

BURNING

BURNING

BURNING

BURNING

.7

?

7

?

?

7

7

7

.?

?

7

?

F

T

F

F

F

F

F

T

F

F

F

F

F

T

F

F

F

F

F

T

F

WF

T

Tý

F.

T.

F.

F

T

F

F

F

F

F

F

F

F

F

F

F

SRCE

SRCE

SRCE

SRCE

SRCE

SRCE

SRCE

SRCE

SRCE

SRCE

SRCE

SRCE

FLUXFLUX

FLUX

FLUX

FLUX

FLUX

FLUXFLUX

FLUX

FLUX

FLUX

FLUX

.0

: 306E+05

.0

.0

.0

.0

.0

.0

.0

415.

453.

.0

";381E+04

148E+04

: 376E+04

.3114E+04

.314E+04

.164E+04

.568E+04

, 569 E+05,-

.563E+04

99.2

.538E+04

413.

.0

,. 306E+05

.0

.0

.0

..0

,0

,533E+05

.0

.0.306 E+05

.0

.0

.0

.0

.0

.381E+04

148 E+04-

.376E+04:-

* 349E+05

.349E+05

.241E+.04

lm,.,.568E+04

:569E+05 I563E+04

.0

.0

.0

.0

.0

.0.

.0,

o.0 0,

.0

,0

.0

.0"

.0:

.0,

.0o

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

'MODULE:

1

2

3.

5

6

7

8

9

10

11

12

EXT FLUX

EXT FLUX

EXT FLUX

"EXT FLUX

EXT FLUX

!EXT FLUX

EXT.FLUX

EXT FLUX

EXT FLUX

EXT FLUX

EXT FLUX EXT FLUX

.413E+0£4

.0

.408 E+04

.610E+0£4

'. 610E+0£4

.214E+04

.630 E+0£4

:0'184E+05-

..625 E+04

304 E*04

,.;150E+05

'299E+04

.;610E+04

:-f10E+04`

.214E+04

-438E+04

1.778E+04.

.433E+04

225E+04

299 E+04

ý.220 E+04

31£4E+04'

314E+04

.164E+04

'3210-04

.'426E+.04

3 V6 E+04-

- 45 -
*.*****~* TEST CASE: WOOD CRIB IN ENCLOSURE *********************

TIME (SEC): ISO.-

/,

TOTAL MASS BURNING RATE (KG/S):

TOTAL HEAT RELEASE RATE (W).

HOT GAS LAYER TEMPERATURE (K)

HOT GAS LAYER THICKNESS (M)

MODULE: 1 BURNING ?

MODULE: 2 BURNING 7

MODULE: 3 BURNING ?

MODULE: 4 BURNING ?

MODULE: 5 BURNING ?

MODULE: 6 BURNING ?

MODULE: 7 BURNING. 7

MODULE: 8 BURNING ?

MODULE: 9 BURNING ?

MODULE: 10 BURNING ?

MODULE: 11 BURNING ?

MODULE: 12 BURNING ?

.4317E-02

.8029.E+05

302.4

2.292

F

T

F'iF

T

F

T

T

F

F

F

F

T

F

T

T

F

T.

T

T

F

T

F

T

T

T

T

T

T

F

Tý

F

T

T

F

T

T

T

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

MODULE:

2

3

4

5

6

7

8

9 *10

11

12

1

2

3.

4

5

6"

7;,

8

9

10

11

12

SRCE

SRCE

SRCE

SRCE

SRCE

SRCE

SRCE

SRCE

SRCE

SRCE

SRCE

SRCE

FLUX

FLUX

FLUX

FLUX

FLUX

FLUX

FLUX

FLUX

F LUX

FLUX

FLUX

FLUX

.0

.,306E+05

.0

.0

.0

.0.0

. 649E+05

.0

416.

474.

.0

.452E+04

*.238E+04

446E+04

.611E+04

:.669E+04

.422E+04ý

,.144E+05

.354E+05

.143E+05

115.

.134E+05

428.

.0

.306E+05

.0

.0

.0

ý.0

.0

.462E+05

.0.

.511E+04

.148E+04

.506E+04

;206E+05

.220E+05 ;604E+04

.160E+05

.276E+05

.:159E+05

.;0

.306E+05

.0

.306E+05

.306E+05

.0

.0

.649E+05

.0

.511E+04

.148E+04

.506E+04

.187E+05

:210E+05:

.?750E+04

:157E+05 .,365E+05

.156E+05

.0

306E+05

.0

.0

.0

.0

.0

384E+05

.0

.452E+04

238E+04

:"-,..446E+04

.206E+05

.220E+05

.604E+04

136E+05 .,..255E+05

.135E+05

F

T

F

F

T

F

T

T

T

.0

.0

.0-

.0

.010

.0

.0

.355E+04

.:156E+05

.350E+04

,.611E+04

.669E+04

.422E+04

.106E+05

.293E+05

.105E+05

EXT FLUX

EXT FLUX

EXT FLUX

EXT FLUX

EXT FLUX

,EXT FLUX

.EXT FLUX

EXT FLUX

EXT FLUX

EXT FLUX

EXT. FLUX

EXT FLUX

1- 46,--.
TEST CASE: WOOD 'CRIB IN ENCLOSURE " *

TIME (SEC):.. 240.

TOTALMASS BURNING RATE (KG/S): .1061E-0i

TOTAL HEAT RELEASE :RATE (W) -: .1974E+06

HOT GAS LAYER TEMPERATURE (K) : 544.6

HOT:GAS LAYER THICKNESS (M) : 4.060

MODULE: I :.BURNING ? T T T T F

MODULE: 2 -BURNING7 T? T T T T

MODULE: 3 BURNING ? T T T T FMODULE: £4 BURNING ? T T T T T

MODULE: 5 BURNING? T T T T T

MODULE: 6 -BURNING ? ,T T T T T

MODULE: 7 BURNING-? T T T T T

MODULE: 8 'BURNING ? .T T T T T

MODULE: 9 BURNING? T T T T T

MODULE: 10 BURNINGr? F

MODULE: 11 BURNING ? F

MODULE: 12 BURNING? F.

MODULE: 1 SRCE FLUX .0 .0 .0 .0 .0

MODULE: 2 ýSRCE FLUX .306E+O5 .306E+05 .306E+05 .306E+05 .306E+05

MODULE: 3 SRCE FLUX .0 .0 .0 .0 .0

MODULE: 4 SRCE FLUX .0 .. 306E+05 .306E+05 .306E+05 .0

MODULE: 5 SRCE FLUX .306E+05 -,306E+05 .306E+05 .306E+05 .306E+05

MODULE: 6 SRCE FLUX .0 .0 .306E+05 .0ý- 0

.MODULE: 7 SRCE FLUX .435E+05 .446E+05 .44•LE+05 ..429E+05 .4O6E+05

MODULE: 8 SRCE FLUX .556E+05 .;516E+05 .562E+05 .505E+05 .525E+05

-MODULE: 9 SRCE_FLUX .434E+05 .445E+05- .443E+05 .429E+05 .O6E+05

iMODULE: 10 SRCE FLUX 419..

MODULE: 11 SRCE FLUX .499E+04

MODULE: 12 SRCE FLUX .0

.MODULE: 1 EXT FLUX .577E+04 :-.656E+04 .683E+O04 .656E+04 .577E+04

MODULE: 2 EXT FLUX .299E+014 <.238E+04 .295E+04 .238E+04 .299E+04

MODULE: 3 EXT FLUX .498E+I04 i576E+O1 .604E+04 .576E+O4 .498E+O4

MODULE: 4 EXT FLUX .245E+05 .150E+05 .2143E+05 .150E+05 .2L45E+05

MODULE: 5 EXT FLUX .120E+05 .166E+05 .293E+05 .166E+05 .120E+05

MODULE: 6- EXT FLUX .120E+05 .272E+05 .255E+05 .272E+05 .120E+05

MODULE: 7 EXT FLUX, .347E+05 .38OE+05 .489E+05 .365E+05 .317E+05

MODULE: 8 EXT FLUX .646E+05 .610E+05 .721E+05 .583E+05 .593E+05

MODULE: 9 EXT FLUX .339E+05 .372E+05 .481E+05 .357E+05 .309E÷05

MODULE:' 10 EXT FLUX .260E+04

MODULE: .11 EXT FLUX .495E+05

MODULE: 12 EXT FLUX .116E+04

- 47 °
0.04

0.03

.,

S0

00

o 0.02

m 0,

C0

m•

a 0.01 -

0 * 0 * , I I

0 10 20 30

Time - min

Figure 6: Crib Mass Burning Rate
1.0"1

0

0.e ,I.

0

0

S

*000

0

S

0.6 - a N

E

S

| •

N

0

-1.

0.4

0.2-

0

S

.

0

0

I 1*0e

0 r-4 m |

0 10

I

20

ff O

30

Time - min

Figure 7: Heat Flux from Wall to Detector
0

8o0

700-

600-

0

9

41

0 0

SO. .. ~...000 0000

0: 0

0

I..

Ut

0

0 500w

.3

0~.

E

.3

I0

500.

400-

3004 .0 0

V

0 10 V 0 20 V

30

Time min

Figure 8: Average Hot Gas Layer Temperature
Chapter 7

REFERENCES

1. N. 0. Siu, "Probabilistic Models for the Behavior of Compartment

Fires," NUREG/CR-2269,UCLA-ENG-8090, University of California at

Los Angeles, August 1981.

2. G. Apostolakis, M. Kazarians, and D. C. Bley, "Methodology for

Assessing the Risk from Cable Fires," Nuclear ::Safety, 23(4).:

391-407(July-August 1982)

3. M. Kazarians and G. Apostolakis, "Fire Risk Analysis for Nuclear

Power Plants," NUREG/CR-2258,UCLA-ENG-8102, University of California at Los Angeles, September 1981.

4. N. 0. Siu, "Physicil Models for Compartment Fires," Reliability

Engineering, 3: 229-252(1982)

5. N. 0. Siu and G. Apostolakis, "Probabilistic Models for Cable

Tray Fires," Reliability Engineering, 3: 213-227(1982),

6. B. Carnahan, H. Luther, and J. Wilkes, Applied Numerical Methods, John Wiley and Sons, Inc., New York, 1969.

- 51 -
Appendix A

_SOME DETAILS ON THE MODELS FOR WALLS AND

THERMAL BARRIERS

In our analysis of compartment fires, we' are often interested in the be-,

havior of slabs of, non-combustible materials as they absorb heat from

,the flames. .We. are particularly interested in. their strength as heat

flux sources, since the heat from these slabs may enhance the growth

rate of a fire.. If we are •interested in--the heat flux-radiating from

the back (unexposed) side of the slab, we treat: the slab as a "barrier."

If the heat flux re-radiated and reflected from the front face of the

slab is of interest, we treat the slab'as a "wall.""

A.1 THERMAL:• BARRIER MODEL

Consider. aslab of width L,::with an isotropic heat flux.of strength..q

impinging on the left face (x=O). 'The.- time-dependent diffusion equation, with appropriate boundary conditions, then governs the temperature

profile within the barrier:

-T 3 2 T (A.1)

t 0: T(x,O) = T

0

x 0: E% = -[T(Ot)-T + (0,t)-Tj AE[T k xTJ (A.2)

x= L: -k [L =. h L[T(Lt )-TA1. + o0[TV(L,t)--1] (A.3)

where a = barrier thermal diffusivity (m2 /s)

-52-
= barrier emissivity.

a =Stefan-Boltzmann' constant (5.669,7,x iJ O W/m2 0

K')

h = barrier surface heat transfer coefficient (W/m2 OK)

k = barrier thermal conductivity (W/mOK)

TA ambient temperature.,( 0 K)

If we assume that steady-state conditions• exist, the solut:ion of

Equation (A-1) is a linear temperature profile,_ i.e.

T(x) T(O) + fT(L) -T(O)](x/L) 

(A.4)

and. the boundary. conditions become (assuminig a constant h):

E h[T(O)-T] + oCE[T4 (0)-T"1 - [T(L)ý-T(0)](k/t-) (A.5) A 

"

- k[T(L)-T(O)] = h[T(L).-TA] + oaE: T(L)-TI L. A (A.6) A

Equations (A.5) and (A.6) can be solved numerically for _T(O) and T(L).

The radiative heat flux emerging form the back side of the barrier is

then

'qaE o T' ( ) ..... -- (A . 7:)

The convective component is relatively small for high values of T(L),

and will not impact on distant objects.

To compute the net heat flux impinging on an object maintainedat ambient temperature, we use

I,

q OE= oE[T(L)_.- T4 ] - A (A.8)

The values of q computed from Equations (A.7) and (A,8)ý vary only slightly for significant levels of thermal, radiation.

-. 53.-
We note that the steady-state assumption :"is -conservative, since a

thermal wave takes a finite amount' of time to penetrate the-slab. The

time scale for this penetration is.of the order

and may 'be very long for dense, poorly.-conducting materials. For example, the&-time constant-for two inches of concrete is aboutf 14.,5 hours (by

which time, the fire may be- extinguished).

A.2 WALL MODEL

The governing -equations- for: the wall temperature are also Equations

(A.l), (A.2)", and (A..,3). However, since the walls tend to be thick and

are often composed of- concrete, we cannot ignore their transient behavior under thermal loads. We must therefore solve the governing equations numerically.

In COMPBRN, we Choose to model the wall with a very coarse spatial

grid consisting only'of I10 mesh points:. These mesh points are lumped

near the wall's exposed surface;, since the penetration of the thermal

wave is very slow. In fact, the time scale of many fires is sufficiently short that the actual wall thickness has little impact on the calculations; we can easily fix T(L,t) at ambient levels, and ignore Equation (A.3) without significantly changing our results. The actual

numerical solution is obtained using a Crank-Nicolson scheme described

in Reference (6).

Once the temperature profile within the wall is found for any particular time step, the' ,radiative heat flux. source strength of the wall's

exposed face. is given by

- -(A.9

• ,54
where

p = wall reflectivity

As in the barrier model, the target object's re-radiation can be accounted for by subtracting a term proportional to T4. We also neglect A

the convective heat flux from the wall in COMPBRN, although suitable models can be incorporated easily by the user, if desired.

- 55 -
Appendix B

CODE LISTING

56 -
C *

C COMPBRN, - 'A COMPUTER PROGRAM FOR ýMODELING. THE BEHAVIOR" *

C OF COMPARTMENT FIRES *

C ... .*

C AUTHOR "- NATHAN 0. SIU *

C *

C VERSION- IBM 3033 (7/82) .. *

C, C -.- A ** AA*

C

C.

C MAIN PROGRAM (READS',-INPUT DATA', CALLS: SUBROUTINES)i. *

C " 7 .

C

C

C

.HUEL.lI, ARAUIi .CEL 

FUEL CELL ARRAYS

DIMENSI

1

2

3

4

5

ON BRAT(12,5.), FLHT(12.,5); FLTEMP(12,5)., FMASS(12,5),-

FX(12,5), FY(12,5), FZ(12,5)., .GAMN(12,5), QDOT2P(12;5), Q(12,5.), QDOTC(12,5), QEXT (12•:5)7 , QEXTO (12-,5) _i-::

.SUMQSQ(12,i5 , j TE•MP(12,5,10),, ICOUNT(-12 i5),. ..- "

$BURN (12,5),-"$FMASS (12,5), .$PIGN(121,5), $STRTO(12,5),

$TOP(12,5)

C

C

C

C

C

C

SUPER-MODULE VECTORS

DIMENSION AREA(12), DELS(12), DEP(12), .FLNG(1-2), FLOSS(12),

1 FMASSO(12), RAD(12), POR(12),. WID(12),.

2 IDIR(12),, IORT(12), ITYP(12), MSMOUT(12),- NFC(12)-,

3 $B(12), $C(12), $D(.12), $W(12)

FUEL TYPE VECTORS . ."

DIMENSION BRATSO(5), BRATS1(5), BRATV(5), DENS(5), DIFF(5),

1 FABSRP(5)., FIGTP(5); FIGTS(5), GAMMA(5);" HEAT(5),

2 SPHT(5), REFL(5), THK(5),

3 IFUEL.(5),

4 $BF(5), $CF(5), $DF(5), $WF(5)

PILOT FIRE VECTORS.

DIMENSION PBRATV(10), PBRTSO(10), PBRTSI(1O)., PHEAT(1O),

1 PMASS(10), PMDOTS(10),:

2 IPIL(10), IPFUEL(10), JPIL(10),

3 $PMASS(10)

COMMUNICATION ARRAYS AND VECTORS

DIMENSION IV(4), NV(4), ICOMM(60,60), IAD(30,4).i NAD(200;-4)

C

C

C

C

C

C

C•

- 57 -
C COMMON STATEMENTS,,:

C

COMMON/ALL/! ACEIL, AWALL, DELT, ýFCTR(15), FLCF, RTEMP, TIME,

1 TITLE(20),

2 ICEIL, IJOB, INCHCK, INITG, IOUTPT, MOUTPT(11),

3 $DOOR, SEND, $ROOM

COMMON /GAS/ A4, CDOOR, CF, CFV, DCF, DEN, DG, DHGT, DWID,

1 FC, FH, GABSRP, HCEIL, PLCF1, PLCF2,

2 TG, TGO, THETA, VFV, WFL, WFVC, WFVH, ZD, ZN

DATA IDIM1/12/, IDIM2/5/, IDIM3/60/, IDIM4/30/, IDIM5/200/,

1 NREAD/5/, NWRITE/6/

NAMELIST /STRT! NJOB, NTIME, NREAD, NWRITE, DELT

NAMELIST /SIZE/ NSM, NCOM, NFUEL, NNCOM, NPILOT, IROOM, INITG

NAMELIST /FUELB/ SMX, SMY, SMZ, SLNG, SWID, SDEP, SMASS, SPOR,

1 SLOSS. NFCL,ý IORNT,. IDIREC, IFTYP

NAMELIST /PILOT/ IPIL, JPIL, IPFUEL, PMASS

NAMELIST /FUELT/ IFUEL, DENS, -SPHT, THK, HEAT, -FIGTP,

1 FIGTS, BRATV, BRATSO, BRATS1, :GAMMA,. FABSRP, REFL

NAMELIST /MISC/ RTEMP,-FLCF

NAMELIST /SEE/ IV

NAMELIST /NSEE/ NV

NAMELIST /ROOM/ DWID, *DHGT, DCF, ,.FC, FH,: GABSRP,. HCEIL, PLCF1,

I PLCF2, ,THETA, VFV

NAMELIST /GINIT/i' TG, DG, QEXT,.

NAMELIST /MODVAR/ FCTR

NAMELIST /OUTF/ INCHCK, IOUTPT, MOUTPT, NSMOUT, MSMOUT

READ (NREAD,STRT)

IF (NJOB.EQ.0) NJOB = 1

"DO 8000 IJOB=1,NJOB'

C

C INPUT PROBLEM DATA

C

READ (NREADSIZE) "

READ- (NREAD,9100)-TITLE.

C

C READ FUEL BED CHARACTERISTICS

C

DO.1000 I=1,NSM ..

READ (NREAD,FUELB)

XNFCL = FLOAT(NFCL)..

DELS(I) = SLNG/XNFCL

AREA(I) = DELS(I)*SWID ,- .

RAD(I) = SQRT(AREA(I)/3.14159)

IF (IORNT.EQ.3) RAD(I) = SWID/2.

FMASSO(I) = SMASS/XNFCL"

WID(I) = SWID

POR(I) =-SPOR

FLOSS(I) = SLOSS

FLNG(I) = SLNG :

DEP(I) = SDEP

IDIR(I) = IDIREC.

IORT(I) = IORNT

ITYP(I) = IFTYP"

NFC(I) = NFCL

- 58
$D(I) = .FALSE.

$W(I) = .FALSE.

$B(I) = .FALSE.

$C(I) = .FALSE.

CALL SETUP(IIDIREC,IDIM1,IDIM2,NFCL,SMX,SMY,SMZ,SLNG,DELS(I),

1 FX,FY,FZ)

DO 500 J=INFCL

FMASS(I,J) = FMASSO(I)

500 CONTINUE

1000 CONTINUE

C

C READ PILOT FIRE CHARACTERISTICS

C

READ (NREAD,PILOT)

C

C READ FUEL TYPE CHARACTERISTICS

C

READ (NREAD,FUELT)

DO 2900 I=lNFUEL

DIFF(I) = THK(I)/(DENS(I)*SPHT(I))

$DF(I) = .FALSE.

$WF(I) = .FALSE.

$BF(I) = .FALSE.

$CF(!) = .FALSE.

ITYPE = IFUEL(I)

IF (ITYPE.LT.10) GO TO 2900,

IF ((ITYPE.GE.10).AND. (ITYPE.LT.20))

IF ((ITYPE.GE.20).AND.(ITYPE.LT.30))

IF ((ITYPE.GE.30).AND.(ITYPE.LT.40))

IF ((ITYPE.GE.40).AND.(ITYPE.LT.50))

DO 2500 J=1,NSM

IF (ITYP(J).NE.I) GO TO 2500

IF ($DF(I)) $D(J) = .TRUE.

IF ($WF(I)) $W(J) = .TRUE.

IF ($BF(I)) $B(J),= .TRUE.

IF (.NOT.$CF(I)) GO TO 2500

$C(J) = .TRUE.

-ICEIL = J

2500 CONTINUE

2900 CONTINUE

READ (NREAD,MISC)

DO 3100 I=1,IDIM3

DO 3000 J=1,IDIM3

ICOMM(I,J) = 1

3000-CONTINUE

3100 CONTINUE

IF (NCOM.EQ.0) GO TO 3500

DO 3300 I=1,NCOM

READ (NREAD,SEE)

DO 3200 J=1,4

IAD(I,J) = IV(J)

3200 CONTINUE

3300 CONTINUE

3500 IF (NNCOM.EQ.0) GO TO 3900

$DF(I) = .TRUE.

$WF((I) = .TRUE.

$BF(I) := "TRUE.

$CF(I)= :.TRUE.

- 59 -
C

C

DO 3800 .I=INNCOM

..READ (NREAD,NSEE),

DO 3700 J=l,4

NAD(I,J) .= NV(J)

3700 CONTINUE

3800 CONTINUE

3900 CALL COMM(NCOM,NNCOM,NSM,IDIM1, IDIM3,IDIM4,IDIM5,NFC, iAD,NAD,

1 ICOMM)

$ROOM = (IROOM.EQ. 1)

IF (.NOT.$ROOM) GO TO 5000

READ ROOM CHARACTERISTICS

READ (NREAD,ROOM)

ADOOR = DWID*DHGT

$DOOR = (ADOOR.LE...E-2) CDOOR = 3.56*DCF*ADOOR*SQRT(DHGT)

CFV = 352.6*VFV

CF = 1000.*CFV

IF (INITG.EQ.1) --READ (NREAD,GINIT)

-IF (DCF.EQ.0.).DCF = 1.

IF' (GABSRP.EQ.0.) GABSRP = 1.4

IF (PLCF1.EQ.0.) PLCF1 = 1.

IF,.(PLCF2.EQ.0.) PLCF2-= .5

IF (RTEMP.EQ.0.) RTEMP = 298.

IF (THETA.EQ.0.) THETA =.5

WFVC =-CFV/RTEMP5000 READ (NREAD,MODVAR)-

DO 6000 I=1,i.15

IF (FCTR(I). LE".0;.O) FCTR(I) =.i.0 -..

6000 CONTINUE

READ (NREADOUTF)

C

C

C

C

-C

.C

C

C

C

ECHO CHECK FOR INPUT DATA

CALL INCHK (FX,FY,FZ,'QEXT,DEP,FLNGFLOSS,FMASSO,POR,WID,BRATSO,

1 BRATS1,BRATV,DENS,DIFF.,FABSRP,FIGTP, FIGTS,GAMMA;HEAT,REFL,SPHT,

2 THK,PMASSIDIR,IORT,ITYP,NFC,IFUEL,IPIL,IPFUEL,JPIL, IAD.NAD,

3 $BF,$CF,$DF,$WF,IDIM1,IDIM2,IDIM4;IIDIM5,NCOM,NFUEL,NJOB,-

4 NNCOM,NPILOT,NSM,NTIME,NWRITE)

INITIALIZE JOB: VARIABLES

CALL INIT (QEXT,SUMQSQFX,FY,FZ,TEMP, $BURN,$FMASS,AREA,DELS,,POR,

1 WID, BRATSO,BRATS1,BRATV,HEATPBRATV,PBRTSO,PBRTS1I,PHEAT,PMASS,

2 IDIR, IORT,NFC,IPIL, IPFUEL,JPIL,$PMASS,ICOMM,$B,$C,$W,IDIM1,

3 IDIM2,IDIM39,NFUEL,NPILOT,NSM)

TIME = 0.0

WDO 7000 ITIME=I,NTIME

TIME-STEP INITIALIZATIONS

CALL TINIT(FLHTFLTEMPQ,QDOT2P,QEXT,QEXT0,$PIGN,$STRT0,NFC,

1 IDIM1IDIM2,NSM).

- 60 -
C

C COMPUTE SOURCE HEAT FLUXES

CALL SOURCE (BRAT,FLHT,FLTEMP,FMAS$,FZ,GAMN,Q,QDOTC, QDOT2P,QEXT,

1 QEXT0,TEMP,AREA",DELSDEP,FMASSO,POR,RAD,WID,,PBRATXV, PBRTSO,PBRTS1,

2 PHEAT,PMASS,PMDOTS, BRATSO,BRATS1,BRATV,DIFF, GAMMA, HEATREFL,THK,

3 ICOUNT,IORT,ITYP,NFC,IPIL,,IPFUEL,JPIL,$BURN,$ FMASS, $ STRT0,$TOP,

4 $B,$C,$D,$W,$PMASS, IDIMI,,IDDIM2ý,NFUEL,NPILOT•,NSM,NWRITE)

IF ($END) GO TO 8000

C

C COMPUTE HEAT FLUX TRANSFER TO RECEIVERS

C

CALL TRANSF(FLHT,FLTEMP,FX,FY,FZ,QDOTC,QD0T2P,QEX'T, TEMP,AREA,

1 DELS,:RADWID,FABSRP,GAMMA,ICOUNTIORTITYP,NFC,'ICOMM,$SBURN,

.2 $PIGN, $B ,.$C, $D, $W, IDIM1, IDIM2, IDIM3, NFUEL",NSM, NWRITE)

IF ($END) GO TO 8000

C

C DETERMINE IF OTHER, FUEL. CELLS IGNITE

C

CALL IGNIT(QEXTSUMQSQ,TEMPDIFF, FIGTP.,FIGTS,THK,ITYP,NFC,

1 $BURN,$FMASS,$PIGN,$B,$C,$D,$W,IDIM1,IDIM2,NFUEL,NSM)

C

C OUTPUT FOR EACH .TIME STEP

C

CALL OUTPUT(FLHT,FLTEMPFMASS,QDOT2Pe, QEXT,TEMP,$BURN$,NFCMSMOUT;

1 IDIM1,IDIM2,NSM,NSMOUT,ýNWRITE)

TIME = TIME + DELT

7000 CONTINUE

8000 CONTINUE

9100 FORMAT (20A4)

9000 STOP

END

- 61 -
C •

C SUBROUTINE INC}HK (READ CHECK FOR INPUT DATA) *

C*

SUBROUTINE I NCHK. (FX,FY,FZ,QEXT,DEP,FLNG:,FLOSSFMASSO,POR,WID,

1 BRATSO,BRATS1,BRATV,DENS,DIFF,•FABSRP,FIGTP,FIGTS,GAMMA,HEAT,REFL,

2 SPHTTHK,PMASSIDIR, IORT,ITYP,NFC, IFUEL,IPIL,IPFUEL,JPIL,IAD,

3 NAD,$BF,$CF,$DF,$WF,IDIM1, IDIM2,;IDIM4,IDIM5,NCOM,NFUEL,NJOB,

4 NNCOM,NPILOT,NSM,NTIME,NWRITE)

IMPLICIT LOGICAL ($)

DIMENSION QEXT(IDIM1,IDIM2), FX(IDIM1,IDIM2), FY(IDIM1,IDIM2),

1 FZ.(IDIM1,IDIM2),t.

DIMENSION DEP(NSM'), FLNG(NSM), FLOSS(NSM), FMASSO(NSM), POR(NSM),

1 WID(NSM), IDIR(NSM),, IORT(NSM), ,ITYP(NSM), NFC(NSM)

DIMENSION BRATSO(NFUEL), BRATS1(NFUEL), BRATV(NFUEL), DENS(NFUEL),

1 DIFF(NFUEL), FABSRP(NFUEL), FIGTP(NFUEL), FIGTS(NFUEL),

2 GAMMA(NFUEL), HEAT(NFUEL), REFL(NFUEL), SPHT(NFUEL),

3 THK(NFUEL), IFUEL(NFUEL),

4 $BF(NFUEL):, $CF(NFUEL), $DF(NFUEL)', $WF(NFUEL)

DIMENSION PMASS(NPILOT),'

1 IPIL(NPILOT), IPFUEL(NPILOT), JPIL(NPILOT)

DIMENSION IAD(IDIM4,4),-NAD(IDIM5,4)

COMMON

.2

3,

COMMON

2

/ALL/,ACEIL, AWALL, DELT, FCTR(15), FLCF, RTEMP, TIME,

• TITLE'(20).

ICEIL, IJOB,'INCHCK., I-N-ITG, IOUTPT, MOUTPT(11),

- $DOOR, $END,,$ROOM

/GAS/ A4, CDOOR, CF, CFV, DCF, DEN, DG, DHGT, DWID,

FC, FH, GABSRP,.HCEIL, PLCF1, PLCF2,

TG, TGO, THETA,"VFV, 'WFL, WFVC, WFVH, ZD, ZN

DATA X/IHX/,-Y/1HY/, Z/!HZ/

WRITE (NWRITE,9000);TITLE, I

WRITE (NWRITE,9100)

WRITE (NWRITE,9101) FCTR(1)

WRITE (NWRITE,9102) FCTR(2)

WRITE (NWRITE,9103) FCTR(3)

WRITE (NWRITE,9104) FCTR(4)

WRITE (NWRITE,9105) FCTR(5)

WRITE (NWRITE,9106) FCTR(6)

WRITE (NWRITE,9107) FCTR(7)

WRITE (NWRITE,9108) FCTR(8)

WRITE (NWRITE,9109) FCTR(9)

WRITE (NWRITE,9110) FCTR(10)

WRITE -(NWRITE,9111) FCTR(11)

WRITE (NWRITE,9112) FCTR(12)

WRITE (NWRITE,9113) FCTR(13)

WRITE (NWRITE,9114) FCTR(14)

WRITE (NWRITE,9115) FCTR(15)

IF (INCHCK.EQ.0) RETURN

WRITE (NWRITE,9200) NFUEL

DO 6000 I=1,NFUEL

WRITE (6,9210) I

IF. ($DF(I))GO TO 1000

IF ($WF(I.)) GO'TO 2000

JOB, NJOB

- 62 -
IF ($BF(I)) GO TO 3000

IF ($CF(I)) GO TO 4000

WRITE (NWRITE,9220) IFUEL(I)

WRITE (NWRITE,9230) DENS(I)

WRITE (NWRITE,,9240) SPHT(I)

WRITE (NWRITE,9250) THK(I)

WRITE (NWRITE,9255) HEAT(I)

WRITE,: (NWRIT1'E,9260) FIGTP(I)

WRITE 9!-(NWRITE, 9270) FIGTS(I)

WRITE (NWRITE,9280) BRATV(I)

WRITE (NWRITE,9290),BRATSO(I)

WRITE (NWRITE,9300) BRATS1(I)

WRITE (NWRITE'1;9310) GAMMA(I) ,

WRITE (NWRITE.9320) FABSRP(I)

WRITE (NWRITE,9330) REFL(I)

GO,. TO 6000

1000'WRITE (NWRITE,9340) IFUEL(I)

GO TO 6000

" 2000 WRITE. (NWRITE,9350) IFUEL(I)

GO TO 5000

3000 WRITE (NWRITE,9360) IFUEL(I)

GO TO 5000

4000 WRITE (NWRITE,9370) IFUEL(I)

5000 WRITE (.NWRITE,9250) THK(I)

WRITE (NWRITE,9400) DIFF(I)

WRITE (NWRITE,9330) REFL(I)

6000 CONTINUE

IF (INCHCK.EQ.1) RETURN

IF (.NOT.$ROOM) GO TO 7000

-WRITE (NWRITE,9500)

WRITE' (NWRITE,9501) FLNG-(ICEIL), WID(ICEIL), FZ(ICEIL,1)

'WRITEý (NWRITE,9510) DHGT,. DWID, DCF

WRITE O(NWRITE,9520) VFV

WRITE,(NWRITE,9525) FH, FC

WRITE (NWRITE,9530) PLCF1, PLCF2

WRITE4 (NWRITE, 9540) GABSRP

WRITE, (NWRITE,9550) THETA..

WRITE (NWRITE,9570). HCEIL

IF (INITG.EQ.0) GO TO 7000.

WRITE (NWRITE,9580) TG, DG

WRITE: (NWRITE,9585):..

DO 6900 I=I,NSM

NFCI = .NFC(I)

WRITE. (NWRITE,9590) (QEXT.(I.,J),J=1,NFCI)

6900 CONTINUE

.7000 WRITE (NWRITE,9600) NPILOT

DO 7500 I=1,NPILOT •

WRITE (NWRITE,9605) I,

WRITE (NWRITE,9610) IPIL(I), JPIL(I)

WRITE (NWRITE,9620) IPFUEL(I)

WRITE4(NWRITE,9630),PMASS(I,)

.7500 CONTINU,.

'WRITE (NWRITE,9700) NSM

DO 8000.'I=1,NSM .

- 63.-
WRITE: (NWRITE,9705) .I.

WRITE (NWRITE,9710) FX(I,1), FY(I,1), FZ(I,1)'

WRITE (NWRITE,9720) FLNG(I), WID(I), DEP(I)

WRITE (NWRITE,9730) NFC((I)

IF` (IDIR(I).EQ.1) DIREC = X

IF .(IDIR(I).EQ.2) DIREC = Y

IF (IDIR(I).EQ.3) DIREC Z

IF (IORT(I).EQ.1) ORNT= X

IF (IORT(I).EQ.2) ORNT Y,

IF (IORT(I).EQ.3) ORNT = Z

WRITE (NWRITE,9740) DIREC, ORNT

WRITE (NWRITE,9750) FMASSO(I)

WRITE (NWRITE,9760) POR(I)

WRITE (NWRITE,9770) •FLOSS(I)

WRITE (NWRITE,9780) ITYP(I)

8000 CONTINUE

IF ((NCOM.NE.0).OR.(NNCOM.NE.0)) WRITE (NWRITE,9790)

IF (NCOM.EQ.0) GO TO 8500

WRITE (NWRITE,9800) NCOM

DO 8200 I=1,NCOM

WRITE (NWRITE,9810) I, (IAD(I,J),J=1,4)

8200 CONTINUE

8500i-IF (NNCOM.EQ.0) GO TO 8900

WRITE (NWRITE,9820) NNCOM

DO 8800 I=1,NNCOM.

WRITE (NWRITE.,'9810) I,- (NAD(I,J),J=1,4).. "

8800 CONTINUE

8900 WRITE (NWRITE,9900)

.WRITE (NWRITE,99910) RTEMP :

WRITE (NWRITE,9920) FLCF

WRITE (NWRITE,9930) DELT

WRITE (NWRITE,9940) NTIME

9000 FORMAT" (IHI,'PROGRAM COMPBRN - A DETERMINISTIC CODE.TO COMPUTE',

1 ' THE PROGRESS OF FIRE' ,/,19X,'OVER A GIVEN FUEL ARRAY WITHIN',..

2 ENCLOSING'',/,19X,-'BOUNDARIES. .,ALL UNITS ARE IN THE ',

3 'MKS SYSTEM.',///,IX,20A4,

4 //,' INPUT DATA:',///,' JOB" ,13,"OF',13. JOBS')-

9100 FORMAT (1H1I, 'VARIABILITY FACTORS FOR FIRE MODELS ' )

.9101 FORMAT ..(1HO,5X,'!VENTILATION CONTROLLED BURNING -RATE'. ,T89,F7.4)

9102 FORMAT (1HO,5X, 'FUEL-SURFACE CONTROLLED BURNING RATE',T89,F7..4)

9103 FORMAT (1H0,:5X,'FLAME HEIGHT FOR-HORIZONTAL FUEL' ,T89,F7.4)

9104"FORMAT (1HO,5X,'FiAME HEIGHT- FOR VERTICAL FUEL' ,T89,F7.4). :

•9105 FORMAT (1H0,5X,'RADIATIVE HEAT FLUX INTERCHANGE'.,T89,F7.4)..

9106 FORMAT (1HO,5X,'BUOYANT PLUME.-TEMPERATURE' ,T89,F7.4)-.:

9107 FORMAT (-IH0,5X,'CONVECTIVE HEAT TRANSFER COEFFICIENT, FOR ',

1 'FLAME' ,T89,F7.4)

9108.FORMAT. (IH0,5X., .'CONVECTIVE HEAT TRANSFER: COEFFICIENT:FOR .',

I .""PLUME'.,T89,F7.4)

9109 FORMAT (1H0,5X,'CONVECTIVE HEAT FLUX' ,T89,,F7.4).

9110 FORMAT (iHO,5X,•''HEAT •TRANSFER TO:-SELF FOR VERTICAL FUEL' ,T89,F7.4)"

9111 FORMAT (1HO,5X,.'HEAT ýTRANSFER TO"ADJACENT-FUEL' ,T89 ;F7ý.4-)

9112 FORMAT (IH0,5X,'HEAT:FLUX FROM CEILING HOT GAS LAYER' ,T89,F7.4)

9113 .FORMAT i(lHOf5X, 'HEAT FLUX FROM REFLECTIONS OFF .WALLS' .,T89 ,F7.4).

9114 FORMAT (1H10,SX'CRITICAL TIME TO. SPONTANEOUS IGNITION"',T89,F7.4)

'64 -
9115 FORMAT9200"FORMAT

..9210 FORMAT

9220 FORMAT

•9230 FORMAT

9240.FORMAT

9250 FORMAT

.9255' FORMAT

9270.FORMAT

9280 FORMAT

9290 FORMAT

"1

9300 FORMAT

" .93.10 FORMAT

:i.19320 FORMAT

930FORMAT, 9340 FORMAT

9350 FORMAT

9360 FORMAT

9370 FORMAT,

ý9400 FORMAT

1

9500ý FORMAT

.. 9501 FORMAT '9510 FORMAT

S1 , 9520 FORMAT

95235 FORMAT

9530 FORMAT

1

9540 FORMAT

9550 FORMAT

1

9520 FORMAT

9530 FORMAT

1

950FRA

('IHO,5X, 'CRITICAL .TIME TO PILOTED IGNITION' ,T9,F7.4)

(1H1,5X,'NUMBER."OF.FUEL TYPES:',15) (1H0,5X,'DATA FOR FUEL TYPE',13,':')

(IHO,10X,'FUEL TYPE:',3X,13)

(IHO,10x,'DENSITY (KG/M**3)t:',T90,F6.0)

(IH.,10X,'SPECIFIC HEAT (J/KG-K):',T90,F6.0)

(1H ,10X,'THERMAL CONDUCTIVITY (W/M-K):',T90,F6.3)

(1H i10X,'HEAT OF COMBUSTION (J/KG):',T84,G12.4)

(1H ,10X,'PILOTED IGNITION TEMPERATURE (DEG. K):',T90,

F6. 0)

(1H ,10X,'SPONTANEOUS IGNITION TEMPERATURE (DEG. K):',.

T90,F6.0)

(1H ,10X,'VENTILATION .CONTROLLED BURNING RATE FACTOR',

.',(KG/M**2.5-S):',T84,G12.4)

(1H ,10X,'SURFACE CONTROLLED SPECIFIC BURNING RATE ',

'(KG/M**2-S):.',T84,GI2.4)

(1H ,10X,'SPECIFIC BURNING RATE RADIATION AUGMENTATION',

(KG/J-M**2):',T84,GI2.4)

(IH ,10X,'FRACTION OF HEAT RELEASED AS RADIATION:',

T90,F6,.3)

(1H ,1OX,'SMOKE:ATTENUATION FACTOR (M**-I):',T90,F6.3)

(1H ,10X,'REFLECTIVITY:',T90,F6.-3)

(1H",1OX,'FUEL TYPE (DETECTOR):',14)

(IH.,10X,'FUEL TYPE (WALL):',I4)-'

(1H 10X;'FUEL TYPE (BARRIER):',I4)

(1i ,10X,'FUEL TYPE (CEiLING):',I4)

(1H ,1OX,'THERMAL DIFFUSIVITY (M**2/S):',T84,

G12.4)

(.H1,,5X, 'ROOM PARAMETERS:')

:(1HO,'1OX,.'CEILING LENGTH, WIDTH, HEIGHT (M):',T78,3F6.2)

(1H ,10X,'DOOR HEIGHT, WIDTH (M), ORIFICE ',

'COEFFICIENT:',T78,3F6.2)

(1H ,10X,'FORCED VENTILATION (M**3/S):',T89,F7.2)

(1H ,10X,'FORCED VENTILATION CONSTANTS (FH AND FC):',

T84,2F6.2)

(1H ,10X,'PLUME ENTRAINMENT CONSTANTS (PLCF1 AND',

PLCF2):',T84,2F6.2)

(1H ,10X,'HOT GAS ABSORPTION COEFFICIENT (M**-I):-',

T90,F6.4)

(1H ,10X,'NUMERICAL ACCELERATION PARAMETER (THETA):',-

T89,F7.2)

(IHO,10X,'CEILING HEAT TRANSFER COEFFICIENT ',

'(W/M**2 DEG.K):',T90,F6.3)

,(1H0,//,6X,'INITIAL CEILING GAS LAYER PARAMETERS:'"//,

" .•11X,'AVERAGE GAS TEMPERATURE:',T89,F7.2,/,

liX,'GAS LAYERDEPTH:',T89,F7.3)

(1H1,'INIT.IAL HEAT FLUXES TO FUEL BED:',//)

(1H, 'MODULE:',I3,/,10D12.4)

(1H1,'PILOT FIRE DATA:',/,6X,'NUMBER OF PILOT FIRES:',I5)

(IHO,5X,'DATA FOR PILOT FIRE',I3,':'.)

(1HO,10X,'LOCATION OF PILOT FIRE (FUEL ARRAY,FUEL CELL)', ',T89,I3,',',13)

(1H ,1OX,'FUEL TYPE:',T91,I5)

(1H ,1OX,'MASS (KG):',T90,F6.3)

9585

9590

9600

9605

9610

FORMAT

FORMAT

FORMAT

FORMAT

FORMAT

1

9620 FORMAT

9630 FORMAT

- 65 -
9700 FORMAT (H.1,5X,'NUMBER. OFe FUEL.ARRAYS:',15) -

9705 FORMAT (IHO,5X,'DATA FORýFUEL ARRAY,"13,':')

9710 FORMAT (1HO,10X, 'INITIAL FUEL CELL COORDINATES (X,Y,,Z),.(M):' ,T70,

1 F8.2,lX,F8.2,1X,F8.2)

9720 FORMAT (lH , 1OX, 'DIMENSIONS -(LENGTH,WIDTH,DEPTH) (M): 'T7.0, $..."

1 F8.2,lX,F8.211X,F8.2)

9730 FORMAT (1H ,10X, 'NUMBER OF FUEL CELLS:',T91,15)

9740 FORMAT.., (IH_,10X,'DIRECTION OF AXIS, ORIENTATION:',T86,AI,8XA1).

9770 FORMAT (1H.,10X,'HEAT LOSSý FACTOR (UNDEFINED):',T90,F6'.-3)

9750 FORMAT (1H ,0OX,'MASS (KG):',T88,F8.4)

9760 FORMAT (.H ,10X,'POROSITY FACTOR (DIMENSIONLESS):' ,T89,F7.4)

9780 FORMAT (IH ,10X,'FUELTYPE:',T91,15)

9790 FORMAT (1H-1,5X, 'COMMUNICATION DATA:',/)

9800 FORMAT (IH0-,5X,'NUMBER OF ADJACENT FUEL CELLS NOT IN SAME. FUEL',

2 ' ARRAY:'-,15,//,14X,'FUEL CELL,-PAIR-',3X,

3 '(FUEL ARRAY,FUEL-CELL) , (FUEL ARRAY,FUEL CELL',/)

9810 FORMAT..(1H.,19X,13,12X,j3,8XI3,11X,13,8X,13).

9820 FORMAT (IH0,//!,6X, 'NON-COMMUNICATION DATA: ',/,11X,'NUMBER OF ',

1 'NON-COMMUNICATING FUEL CELLS: ',1i5ý,//,,

2 14x,'FUEL CELL PAIR',3X,'(FUEL ARRAY,FUEL CELL)',,

3 ' , (FUEL ARRAY,FUEL CELL)',/)

9900 FORMAT (1H1,6X,'MISCELLANEOUS:')

9910 FORMAT (1H ,10X,'ROOM TEMPERATURE (DEG. K):',T90,F6.0)

9920 FORMAT (IH ,10X,'CONVECTIVE HEAT TRANSFER COEFFICIENT FOR FLAME!',

1 ' (W/M**2-K):',T9.0,F6.2) t

9930 FORMAT (1H ,10X,'TIME INCREMENT (S):',T90,F6.2).

9940 FORMAT (1H ,10X,'NUMBER OF. TIME STEPS FOR JOB:',T91!,I5)

RETURN'

END

- 66 -
C *

C SUBROUTINE INIT (INITIALIZATION OF PROBLEM PARAMETERS)

C *

SUBROUTINE INIT(QEXT,SUMQSQ,FX,FY,FZTEMP,$BURN,$FMASS,AREA,DELS,

1 POR,WID,BRATSO,BRATS1,BRATV,HEAT,PBRATV,PBRTSO,PBRTS1,PHEAT,

2 PMASS,IDIR;IORT,NFC,IPIL, IPFUEL,JPIL,$PMASS,ICOMM,$B,$C,$W,

3 IDiM1,IDIM2,IDIM3,NFUEL,NPILOT,NSM)

IMPLICIT LOGICAL ($),

DIMENSION

1

2

DIMENSION

1

2

DIMENSION

DIMENSION

1

2

3

DIMENSION

QEXT(IDIM1,IDIM2), SUMQSQ(IDIM1,IDIM2), FX(IDIM1,IDIM2),

FY(IDIM1,IDIM2), FZ(IDIM1,IDIM2), TEMP(IDIM1,IDIM2,10),.

$BURN(IDIM1,IDIM2), $FMASS(IDIM1,IDIM2)

AREA(NSM), DELS(NSM), POR(NSM), WID(NSM),

IDIR(NSM), IORT(NSM), NFC(NSM),

$B(NSM), $C(NSM), $W(NSM)

BRATSO(NFUEL), BRATS1(NFUEL), BRATV(NFUEL), HEAT(NFUEL)

PBRATV(NPILOT), PBRTSO(NPILOT), PBRTS1(NPILOT),.

PHEAT(NPILOT), PMASS(NPILOT),

IPIL(NPILOT), IPFUEL(NPILOT), JPIL(NPILOT),

$PMASS(NPILOT)

ICOMM(IDIM3,IDIM3)

COMMON /ALL/ ACEIL, AWALL, DELT, FCTR(15), FLCF, RTEMP,-TIME,

1 " TITLE(20),

ICEIL, IJOB, INCHCK, INITG, IOUTPT, MOUTPT(11),

3 $DOOR, $END, $ROOM

COMMON /GAS/ A4, CDOOR, CF, CFV, DCF, DEN, DG, DHGT, DWID,

1 FC, FH, GABSRP, HCEIL, PLCF1,. PLCF2,

z TG, TGO, THETA, VFV, WFL, WFVC, WFVH, ZD, ZN

COMMON /VENTC/ AFUEL, AFUELO, CVA, CVA0,.FHTA, FRADA, FUELA, QMAX,

1 QTOT, QTOTO, QTOT1, QTOT2, TMDOT, .TMDOTO,

TMDOTS, TQDOT, TQDOTC, WIN, ZOA-1

3 NBURN, NBURNO,

f$DECAY, $FLUX1, $FLUX2, $VENT, $VCONT

C

C COMPUTE SHAPE FACTORS FROM WALLS, BARRIERS, AND CEILING TO OBJECTS

C

IJ = 0

DO 1900 I=1,NSM

NFCI = NFC(I)

IF ($W(I).OR.$B(I).OR.$C(I)) GO TO 1100

IJ = IJ + NFCI

GO TO 1900

1100 DO 1700 J=1,NFCI

IJ = IJ + 1

KL = 0

DO 1500 K=I,NSM

NFCK = NFC(K)

DO 1300 L=1,NFCK

KL = KL + 1

IF (I.EQ.K) ICOMM(IJ,KL) = 0

IF (ICOMM(IJ,KL).EQ.0) GO TO 1300

F12 = SAAPE(IORT(I),IORT(K),IDIR(I),FX(IJ),FY(I,J),

1 'FZ(I,J),FX(K,L),FY(KL),FZ(KL),DELS(I),WID(I))

- 67 -
ICOMM(IJ,KL) = IFIX(F12*1000.)

1300 CONTINUE

1500 CONTINUE

1700 CONTINUE

1900 CONTINUE

IF (.NOT.$ROOM) GO TO 3000

IF ((INITG.EQ.1).AND.(ZD.GT.O.)) GO TO 2500

TG = 300.

ZD = 0.2*DHGT

ZN = 2.*ZD

2500 IF (TG.EQ.0.) TG = 300.

:"CVA =O0.

TGO = TG,.

QTOT = 0.

QTOTO = 0.

QTOT1 = 0.

QTOT2 = 0.

AWALL = 2.*FZ(ICEIL,I)*(DELS(ICEIL) + WID(ICEIL))

ACEIL = AREA(ICEIL)

QMAX 4.21E6*FZ(ICEIL,1)*ACEIL

$VENT = .FALSE.

$VCONT = .FALSE.

$DECAY = .FALSE.

$FLUX1 .FALSE. '

3000 TMDOT = 0:.

TMDOTS = 0.

DO 3700 I=5,NSM

NFCI = NFC(I)

DO 3500 J=i,NFCI. El

$BURN(I,J) =._FALSE.'

$FMASS(I,J) = .FALSE. --

IF ($ROOM.AND.(INITG.EQ.I)) GO TO 3100,.'

SUMQSQ(I,J) = 0..

QEXT(I;,'J)•= 0.

GO TO 3200

3100 SUMQSQ(I,J) = QEXT(IJ)**2.

3200 DO 3300 K=1,10

TEMP(I,J,K) RTEMP

3300 CONTINUE

3500 CONTINUE

3700 CONTINUE .

$END = .FALSE._-

AFUEL = 0.

NBURN = 0

C

C INITIALIZE PILOT FIRES

C

DO 4000 I=-,NPILOT

NBURN = NBURN + 1

AFUEL = AFUEL + POR(IPIL(I))*AREA(IPIL(I)) $BURN(IPIL(I),JPIL(I)) =.TRUE,., ,•

IPTYP = IPFUEL(I)

PBRATV(I) = BRATV(IPTYP)

PBRTSO(I) =,BRATSO(IPTYP)

-68 -
PBRTSI(I) =,IBRATS1(IPTYP)

PHEAT(I) = HEAT(IPTYP)

$PmASS(I) = .TRUE.

4000 CONTINUE

RETURN

END

- 69 -
C *

C SUBROUTINE TINIT (PERFORMS TIME STEP INITIALIZATIONS) *

C *

SUBROUTINE TINIT(FLHTFLTEMP,Q,QDOT2P,QEXTQEXTO,$PIGN,$STRTO,

.INFC,IDIMI,IDIM2,NSM)

IMPLICIT LOGICAL ($)

DIMENSION FLHT(IDIM1,IDIM2), FLTEMP(IDIM1,IDIM2), Q(IDIM1,IDIM2),

I QDOT2P(IDIM1,IDIM2),-QEXT(IDIMI,IDIM2), QEXTO(IDIM1,IDIM2),

2 NFC(NSM),- $PIGN(IDIMI,IDIM2), $STRTO(IDIM1,IDIM2)

COMMONI/ALL/ ACEIL, AWALL, DELT, FCTR(15), FLCF, RTEMP, TIME,

1 TITLE(20),

2 ICEIL, !JOB, INCHCK, INITG, IOUTPT, MOUTPT(11),

3 $DOOR, $END, $ROOM

COMMON /GAS/ A4, CDOOR, CF, CFV, DCF, DEN, DG, DHGT, DWID,

1 FC, FH, GABSRP, HCEIL, PLCF1, PLCF2,

2 TG, TGO, THETA, VFV, WFL, WFVC, WFVH, ZD, ZN

COMMON /VENTC/ AFUEL, AFUELO, CVA, CVA0, FHTA, FRADA, FUELA, QMAX,

1 QTOT, QTOTO, QTOTI, QTOT2, TMDOT, TMDOTO,

2 TMDOTS, TQDOT, TQDOTC, WIN, ZOA,

3 NBURN, NBURN0,

4 $DECAY, $FLUXI, $FLUX2, $VENT, $VCONT

IF (.NOT.$ROOM) GO TO 1000

CVAO = CVA*FCTR(1)

CVA = 0.

FHTA = 0.

FUELA = 0.

TGO.= TG

ZOA = 0.

$FLUX2 = $FLUX1

$FLUX1 = .FALSE.

1000 AFUELO = AFUEL

AFUEL = 0.

TMDOT0 = TMDOT

TMDOT 0.

TMDOTS = 0.

TQDOT =0.

TQDOTC = 0.

NBURNO = NBURN

NBURN = 0

DO 2500 I=1,NSM

NFCI = NFC(I)

DO 2000 J=1,NFCI

QEXT0(I,J) = QEXT(I,J)

QEXT(I,J) = 0.

$PIGN(I,J) = .FALSE.

$STRTO(I,J) = .FALSE.

FLTEMP(I,J) = 0.

FLHT(I,J) = 0.

QDOT2P(IJ) = 0.

Q(I,J)= 0.

2000 CONTINUE

2500 CONTINUE

- 70-
RETURN

END

- 71 -
c*************:********************~*~*:*:*********************************

C *

C SUBROUTINE SOURCE .(COMPUTES HEAT FLUX SOURCE STRENGTH OF FUEL *

C ELEMENTS, WALLS, BARRIERS, AND CEILING) *

C*

SUBROUTINE SOURCE(BRAT,FLHT,FLTEMP,FMASS,FZ,GAMN,Q,QDOTC,QDOT2P,

1 QEXT,QEXTO,TEMP,AREA,DELS,DEP,FMASSO,POR,RAD,WID,PBRATV,PBRTSO,

2 PBRTS1,PHEAT,PMASSPMDOTS,BRATSO,BRATS1,BRATV,DIFFGAMMA,HEAT,.

3 REFL,THK,ICOUNT,IORT,ITYP,NFC,IPIL,IPFUEL,JPIL,

4 $BURN,$FMASS,$STRTO,$TOP,$B,$C,$D, $W,$PMASS,

5 IDIM1,IDIM2,NFUEL,NPILOT,NSM,NWRITE)

IMPLICIT LOGICAL ($)

DIMENSION BRAT(IDIM1,IDIM2), FLHT(IDIM1,IDIM2),

1 FLTEMP(IDIM1,IDIM2), FMASS(IDIM1,IDIM2), FZ(IDIM1,IDIM2),

2- GAMN(IDIM1,IDIM2), Q(IDIM1,IDIM2), QDOTC(IDIM1,IDIM2),

3 QDOT2P(IDIM1,IDIM2),. QEXT(IDIM1,IDIM2), QEXTO(IDIM1,IDIM2),

4 TEMP(IDIM1,IDIM2;10),

5 ICOUNT(IDIM1,IDIM2), $BURN(IDIM1,IDIM2),

6 $FMASS(IDIMI,IDIM2), $STRT0(IDIM1,IDIM2), $TOP(IDIMI,IDIM2)

DIMENSION AREA(NSM), DELS(NSM), DEP(NSM), FMASSO(NSM), POR(NSM),

1 RAD(NSM), WID(NSM),

2 IORT(NSM), ITYP(NSM), NFC(NSM),

3 $B(NSM), $C(NSM), $D(NSM), $W(NSM)

DIMENSION PBRATV(NPILOT), PBRTSO(NPILOT), PBRTS1(NPILOT):,

1 PHEAT(NPILOT), PMASS(NPILOT), PMDOTS(NPILOT),

2 IPIL(NPILOT), IPFUEL(NPILOT), JPIL(NPILOT),

3 $PMASS(NPILOT)

DIMENSION BRATS0(NFUEL), BRATS1(NFUEL), BRATV(NFUEL), DIFF(NFUEL),

1 GAMMA(NFUEL), HEAT(NFUEL), .REFL(NFUEL), THK(NFUEL)

COMMON /ALL/ ACEIL, AWALL, DELT, FCTR(15), FLCF, RTEMP, TIME,

1 TITLE(20),

2 ICEIL, IJOB, INCHCK, INITG, IOUTPT, MOUTPT(11),

3 $DOOR, $END, $ROOM

COMMON /VENTC/. AFUEL,.AFUELO, CVA, CVAO, FHTA, FRADA, FUELA, QMAX,

1 " QTOT, QTOTO, QTOT1, QTOT2, TMDOT,! TMDOTO,

2 TMDOTS, TQDOT, TQDOTC, WIN, ZOA,.

3 -NBURN, NBURNO,

4 $DECAY, $FLUX1, $FLUX2, $VENT,.$VCONT

C

C PILOT FIRE BURNING RATE

C

DO 1900 I=1,NPILOT

IF (.NOT.$PMASS(I)) GO. TO 1900

IPILOT = IPIL(I)

JPILOT = JPIL(I)

PAREA = AREA(IPILOT)

IF (PMASS(I).GT.0.0) GO TO 1500

PAREA = 0.0

PMASS(I) .0.0

$PMASS(I) = .FALSE.

1500 PMDOTS(I) = (PBRTS0(I) + PBRTS1(I)*QEXTO(IPILOT,J.PILOT))*

I PAREA*POR( IPILOT)*FCTR (2)

$STRTO(,IPILOT,JPILOT) =. ($PMASS (I.))

--- 72 -
1900 CONTINUE

DO 8000 I=1,NSM

I...-IF ($D(1).OR.$W(I)`. OR.-$B(I) .OR.$C(I)) GO TO 8000

FRMASS = 3*FMASS0(I)

ITYPE = ITYP(I)

NFCI = NFC(I)

-$ORNT = (IORT(I).NE.3)

DO 4000 J=1,NFCI

IF (.NOT.$BURN(I,J)) GO TO 4000

NBURN = NBURN +. 1

$STRT = $STRTO(I,J)

IF (.NOT.$STRT) GO TO 2500

DO 2000 II=1,NPILOT

IF ((I.NE.IPIL(II)).OR.(J.NE.JPIL(II))) GO TO 2000

INDEX = II

GO TO 2500

2000 CONTINUE

C

C BURNING RATES (VENTILATION AND FUEL SURFACE CONTROLLED)

C

2500 AFUEL = AFUEL + AREA(I)*POR(I)

.•FMDOTS = (BRATSO(ITYPE) + BRATS1(ITYPE)*QEXTO(I,J))*

1 'AREA(I)*POR(I)

IF ($STRT) FMDOTS = PMDOTS(INDEX)

IF (FMDOTS.GT.0.) GO TO 2700

BRAT(I,J) = 0.

$BURN(I,J) = .FALSE.

GO TO :4000

2700 FMDOT = FMDOTS

IF (.NOT.$ROOM) GO TO 2900

FBV = BRATV(ITYPE)

IF ($STRT) FBV = PBRATV(INDEX)

'CVA = CVA + FBV*AREA(I)*POR(I)

I-IF ($VENT) FMDOT = WIN*AREA(I )*POR(I )*CVAO/AFUELO** 2

2900 TMDOT = TMDOT + FMDOT

TMDOTS = TMDOTS + FMDOTS

BRAT(I,J) = FMDOT

C

C HEAT PRODUCTION RATE

C

QCOMB = HEAT(ITYPE)

IF ($STRT) QCOMB:= PHEAT(INDEX)

Q(I,J) = FMDOT*QCOMB

TQDOT = TQDOT + Q(I,J)

GAMN(I,J) = 1. -- GAMMA(ITYPE)

IF ($STRT) GAMN(I,J) = 1. GAMMA(IPFUEL(INDEX)),

C

C DECREMENT FUEL MASS

C

IF ($STRT) GO TO 3000

XMASS =. FMASS(I,J) - FMDOT*DELT

IF (XMASS.GT.FRMASS) GO TO.3500

$BURN(I,J) = .FALSE.

$FMASS(I,J) = .TRUE.

- 7.3 -
FMASS(I,J).-= FRMASS

GO TO 4000

C

C DECREMENT PILOT.MASS

C

3000 PMASS(INDEX) =PMASS(INDEX) FMDOT*DELT

GO TO 4000

3500 FMASS(I,J) = XMASS

4000 CONTINUE

C

C FLAME HEIGHTS

C

VMDOT = 0.0

KCOUNT = 0

DO 5900 J=1,NFCI

IF (.NOT.$BURN(I,J)) GO TO 5900

FMDOT = BRAT(I,J).

IF ($ORNT) GO TO 5100

C

C FLAMES OVER HORIZONTAL FUEL SLABS

C

FHT = (30.64*RAD(I)*(BRAT(I,J)/(AREA(I)*SQRT(RAD(I))))**.61)

1 *FCTR(3)

IF ($ROOM) ZOA = ZOA + FZ(I,J)*TMDOT

GO TO 5700

C

C FLAMES OVER VERTICAL FUEL SLABS

C.

5100 VMDOT VMDOT + FMDOT

KCOUNT = KCOUNT + 1

IF ((J.LT.NFCI).AND.$BURN(I,J+1)) GO TO.5500

FHT = (29.,7*(VMDOT/WID (I))9**.6667)ý*FCTR(4) + DELS(i,)e

IF ($ROOM).ZOA = ZOA + FZ(IJ-KCOUNT+1)'*VMDOT

VMDOT = 0.0"'". ."

KCOUNT = 0

$TOP(I,J) = TRUE.

GO TO 5700

5500 FHT = DELS(I)

$TOP(I,J) = .FALSE.

5700 FLHT(I,J) = FHT

IF (:NOT.$ROOM) GO TO 5900

FUELA = FUELA + AREA(I)

FHTA = FHTA + FHT*FMDOT

.5900 CONTINUE

C

C SOURCE STRENGTH OF FLAMES

C

IF ($ORNT) GO TO 7000

C

C HORIZONTAL FUEL SLABS

C

DO 6000 J=I,NFCI

IF (.NOT.$BURN(I,J)) GO TO 6000

HTFLUX = Q(I,J)/i(3.14159*RAI(I)*(RAD(I) + 2.*FLHT(I,J)))

:- J74 -
QDOT2P(I,J) = HTFLUX

FLTEMP(I,J) = SQRT(SQRT(HTFLUX/5.6697E-8))

QDOTC(I,J) = GAMN(I,J)*Q(I,J)

TQDOTC TQDOTC + QDOTC(IJ)

6000 CONTINUE

GO TO 8000

C

C VERTICAL FUEL SLABS

C

7000 IBOTP = 1

7100 JCOUNT 0

$BOT = .FALSE.

IBOT = IBOTP

DO 7300 J=IBOT,NFCI

IF (.NOT.$BURN(I,J)) GO TO 7300

IF (.NOT.$BOT) IBOTP = J

$BOT = .TRUE.

JCOUNT = JCOUNT + 1

ICOUNT(I,J) = JCOUNT

IF ($TOP(I,J)) GO TO 7500

7300 CONTINUE

GO TO 8000

7500 ITOP = IBOTP + JCOUNT - 1

TFLHT =0.0

VQDOT = 0.0

DO 7700 J=IBOTP,ITOP

TFLHT = TFLRT + FLHT(I,J)

VQDOT = VQDOT + Q(I,J)

7700 CONTINUE

QDOTC (I, ITOP) = VQDOT*GAMN(I, ITOP)

HTFLUX = VQDOT/(2.*TFLHT*WID(I))

FLAMT = SQRT(SQRT(HTFLUX/5.6697E-8))

DO 7900 J=IBOTP,ITOP

QDOT2P(I,J) = HTFLUX

FLTEMP(I,J) = FLAMT

7900 CONTINUE..

IBOTP = ITOP + 1

IF (IBOTP.LE.NFCI) GO TO 7100

8000 CONTINUE

IF (NBURNO.EQ.0).GO TO 9500

IF (.NOT.$ROOM) GO TO 8200

FRADA = SQRT(FUELA/3.14159)

FHTA = FHTA/TMDOT

ZOA = ZOA/TMDOT

8200 DO 8900 I=1,NSM

ITYPE = ITYP(I)

NFCI = NFC(I)

EPS = 1. - REFL(ITYPE)

TK = THK(ITYPE)

TDIFF = DIFF(ITYPE)

DELX = DEP(I)/10.

IF ($W(I).OR.$C(I)) GO TO 8500

IF ($B(I)) CALL BARR(IDIM1,IDIM2,I,NFCI,NWRITE,TK,EPS,

1 DEP(I),FZ,QDOT2P,QEXTO,TEMP)

- 75 -
8500

8700

8900

9000

9500

9600

IF ($END) GO TO 9000

GO TO 8900

IF ($C(I)) CALL CEILNG(EPSREFL(ITYPE),QDOT2P(ICEIL,1),

1 QEXTO(ICEIL,I),TEMP(ICEIL'1,1).,FZ(,ICEIL,1),NWRITE)

IF (.$END) GO TO 9000

DO 8700 J=1,NFCI

HL = HCEIL

CALL DIFFUS(IDIM1,IDIM2,I,J,NWRITE,DELX,EPS,FZ(I,.J),HL,

1 QEXT(I,J),QEXTO(I,J),TDIFF,TK,TEMP)

IF ($W(I))QDOT2P(I,J) = 5.6697E-8*EPS*TEMP(I,J,1)**4 +

1 REFL(ITYPE)*QEXT0(IJ)

IF ($END) GO TO 9000

CONTINUE

CONTINUE

RETURN

WRITE (NWRITE,9600)

FORMAT (lH0,'END OF JOB:

$END = .TRUE.

RETURN

END

NO FLAMES OVER FUEL BED')

- 76 -
C .

C SUBROUTINE TRANSF (COMPUTES HEAT FLUX TRANSMITTED TO ,FUEL *

C ELEMENTS, WALLS, BARRIERS, AND.CEILING)• *

C*

.SUBROUTINE TRANSF(FLHT,FLTEMP,FX,FY,FZ,QDOTC,QDOT2P,QEXT,TEMP,

1 AREA,DELS,RAD,WID,FABSRP,GAMMA,ICOUNT,IORT,ITYP,NFC,ICOMM,.

2 $BURN,$PIGN,$B,$C,$D,$W,IDIM1,IDIM2,IDIM3,NFUEL,NSMNWRITE)

IMPLICIT LOGICAL ($)

DIMENSION FLHT(IDIM1,IDIM2), FLTEMP(IDIM1,IDIM2.), FX(IDIM1,IDiM2),

1 FY(IDIM1,IDIM2), FZ(IDIM1,-IDIM2), QDOTC(IDIM1,IDIM2),

2 QDOT2P(IDIM1,IDIM2), QEXT(IDIMI,IDIM2),

3 TEMP(IDIM1,IDIM2,10), ICOUNT(IDIM1,IDIM2),

4 $BURN(IDIM1,IDIM2), $PIGN(IDIM1,IDIM2)

DIMENSION AREA(NSM), DELS(NSM), RAD(NSM), WID(NSM),

1 IORT(NSM), ITYP(NSM), NFC(NSM).,-

2 $B(NSM), $C(NSM), $D(NSM), $W(NSM)

DIMENSION FABSRP(NFUEL), GAMMA(NFUEL)

DIMENSION ICOMM(IDIM3,IDIM3)

COMMON /ALL/ ACEIL, AWALL, DELT, FCTR(15), FLCF, RTEMP, TIME,

1 TITLE(20),

2 ICEIL, IJOB, INCHCK, INITG, IOUTPT, MOUTPT(11),

3 $DOOR, $END, $ROOM

COMMON /GAS./ A4, CDOOR, CF, CFV, DCF, DEN, DG, DHGT, DWID,

1 FC, FH, GABSRP, HCEIL, PLCF1, PLCF2,

2 TG, TG0, THETA, VFV, WFL, WFVC, WFVH, ZD, ZN

COMMON /VENTC/

1

2

3

4

IJ 0

AFUEL, AFUELO, CVA, CVAO, FHTA, FRADA, FUELA, QMAX,

QTOT, QTOTO, QTOT1, QTOT2, TMDOT, TMDOTO,

TMDOTS, TQDOT,,TQDOTC, WIN,-ZOA,

NBURN, NBURNO, $DECAY, $FLUX1, $FLUX2, $VENT, $VCONT

C

C HEAT FLUX SOURCE LOOP

C

DO 5900 I=I,NSM

NFCI = NFC(I)

IF (.NOT.$D(I)) GO TO 1000

IJ = IJ + NFCI

GO TO 5900

1000 FRAD = RAD(I)

ITYPE = ITYP(I)

$WALBR = ($W(I).OR.$B(I).OR.$C(I))

$ORNT = (IORT(I).NE.3)

DO 5800 J=1,NFCI

IJ = IJ + 1

IF ($WALBR) GO TO 2000

IF (.NOT.$BURN(I,J)) GO TO 5800

C

C FLAME LOCATION

C

FLX = FX(I,J)

FLY =,FY(I,J)

- 77. -
.... FLZ1 FZ(I,J)

IF ($ORNT) FLZ1 = FLZ1 - DELS(I):/2.

FLZ2 = FLZ1 + FLHT(I,J)

HTFLUX = GAMMA(ITYPE)*QDOT2P(I,J.)

2000 KL = 0

C

C HEAT FLUX RECEIVER LOOP

C

DO 5500 K=1,NSM

NFCK = NFC(K) ý $IK = (I.EQ.K)

3000 IF ($WALBR) GO TO 3500

$ORNT2 = (IORT(K).NE.3)

$ORNT3 = ($ORNT.AND.$ORNT2)

3500 DO 5000 L=1,NFCK

KL = KL + 1

IF (ICOMM(IJ,KL).EQ.0) GO TO 5000

IF ($WALBR) GO TO4700

IF (ICOMM(IJ,KL).EQ.2) GO TO 5000

C

C

C

C

C

C

C

C

RECEIVER LOCATION (FLAME LOOP)

FZKL = FZ(K,L)

FR = SQRT((FX(K,L) - FLX)**2

DELZl = ABS(FZKL - FLZ1)

DELZ2 = ABS(FZKL - FLZ2)

IF (FR.LE.RAD(I)) GO TO 4300

+ (FY(K,L)- FLY)**2)

RADIATIVE HEAT FLUXES (WHEN FR-.GT.FRAD)

IF ($ORNT2) GO TO 4000

QI = CYLPAR(DELZ1,FR,FRAD)

Q2 = CYLPAR(DELZ2,FR,FRAD)

GOTO 4100

4000• Qi = CYLPER(DELZI,FR,FRAD)

Q2 = CYLPER(.DELZ2,FR,FRAD)

4100 IF (FZKL.GT.FLZ2) Q2 = -Q2

IF (FZKL.LT.FLZ1) QI = -QI

DELTQ = ((DELZ1*Q1- + DELZ2*Q2)/FLHT(I,J)

GO TO 4800

C

)*HTFLUX*FCTR (5)

C

C CONVECTIVE AND RADIATIVE HEAT FLUXES (WHEN FR.LE.FRAD)

C

C

4300 IF ((QDOTC(I,J).LT.100.).OR.(FZKL.LT.FLZ1).OR.

1 ((I.EQ.K).AND.$BURN(K,L))) GO TO 5000

IF (FZKL.GE.FLZ2) GO TO 4400,

$PIGN(K,L) = .TRUE.

DELTQ = GAMMA(ITYPE)*QDOT2P(I,J)

IF (.NOT.($BURN(K,L).OR.$ORNT2)) DELTQ = QDOT2P(I,J)

GO TO 4800

- 78 -
C

C PLUME TEMPERATURE USING ALPERT'S CORRELATION

ýC

4400 PLTEMP = .169*(QDOTC(I,J)**.6667)/(DELZ1**1.6667) + RTEMP

PLTEMP = AMIN1(PLTEMP,FLTEMP(I,J))*FCTR(6)

IF (,NOT.$ORNT2) GO TO 4500

C

C FOR VERTICAL OBJECTS, HEAT TRANSFER COEFFICIENT = FLAME COEFFICIENT

C

PLCOEF = FLdF*FCTR(7)

GO TO 4600

C

C HEAT TRANSFER COEFFICIENT FOR HORIZONTAL OBJECTS

C

4500 RKL = SQRT(DELS(K)*WID(K)/3.14159)/DELZ1

PLCOEF = .64

IF (RKL.LT.O'-18) PLCOEF = 2.06

PLCOEF = (PLCOEF*(QDOTC(I,J)/DELZ1)**.3333)*FCTR(8) C

C CONVECTIVE HEAT FLUX

C

4600 DELTQ = PLCOEF*(PLTEMP- TEMP(I,J,1))*FCTR(9).-

GO TO 4800

C

C WALL, BARRIER, AND CEILING SOURCES

C

4700ODELTQ = ICOMM(IJ,KL)*QDOT2P(I,J)/1000.

IF ((FZ(K,L).GT.ZD).AND.$ROOM)

1 DELTQ = DELTQ + HCEIL*(TG - TEMP(K,L,I))

4800 QEXT(K,L) = QEXT(K,L) + DELTQ

$FLUX1 = .TRUE.

5000 CONTINUE

5500 CONTINUE

_C

C HEAT TRANSFER TO SELF FOR VERTICAL FUEL CELLS

C (USING A BEST-FIT FOR THE SHAPE FACTOR FROM FLAME TO FUEL)

C

IF ((.NOT..$ORNT).OR.$VCONT.OR•$WALBR) GO TO 5800•: .1.

XCOUNT = ICOUNT(I,J) .VSHAPE = - .197 + .24*WID(I) + .02*FLHT(I,J) + .0538*DELS(I)

1 * (XCOUNT - .5) + .0256*DELS(I)*(XCOUNT- .5)/FLHT(I,J)

2 + .0724*FABSRP(I)

VSHAPE = AMIN1(1.,AMAX1(0.,VSHAPE))

DELTQ = GAMMA(ITYPE)ýQDOT2P (I,J)*VSHAPE*FCTR(10)

QEXT(I,J) = QEXT(I,J) + DELTQ

5800 CONTINUE

5900 CONTINUE

IF (.NOT.$ROOM) GO TO 6000

C

C DETERMINE IF VENTILATION CONTROLLED (USING SIMPSON'S RULE

C EXTRAPOLATION FOR HEAT-SUPPORT BY AIR INITIALLY IN ROOM.-'

C

DQTOT DELT*(-QTOTO + 8.*QTOT1 + 5.*TQDOT)/12.

QTOT = QTOT + DQTOT

79 -
QTOTO QTOT1

QTOT1 = QTOT2

QTOT2 = TQDOT

$FLUX = ($FLUX1.AND.$FLUX2)

IF ((.NOT.$VENT).OR.($VENT.AND.$DECAY.AND. $FLUX))- TMDTSO = TMDOTS

TMDOTV = CVAO*WIN/AFUELO

$VENT = ((TMDOTV.LT.TMDTSO).AND., (QTOT.GE.QMAX))

$DECAY = ($DECAY.OR.(NBURN.LT..NBURNO)).

SVCONT = ($VENT.AND..NOT.$DECAY)

C

C HEAT TRANSFER TO ADJACENT (NON-BURNING) FUEL CELLS FROM

C CONTINUOUS PROPAGATION OF FLAMES

C

6000 IJ = 0

DO 6900 I=1,NSM

NFCI = NFC(I)

IF (.NOT.($D(I).OR.SW(I).OR.$B(I).OR.$C(!))) GO TO 6050

IJ =IJ+NC

GO TO 6900

6050 ITYPE = ITYP(I)

$ORNT = (IORT(I).NE.3)

DO 6600 J=1,NFCI

IJ = IJ + 1

IF (.NOT.$BURN(I,J)) GO TO 6600

FLZ = FZ(I,J)

C

C CONDUCTIVE HEAT TERM

C

QCOND = (.01.157+5.6697D-5f*FLTEMP•(I,J))* (FLTEMP(I.,J).,TEMP.(I,J,1))*

1 WID(I)

KL= 0

DO 6400 K-I,NSM

AREAK = W!D(K)*DELS(K)

NFCK =NFC(K)

$ORNT2 = (,$0RNT.AND. (IORT(K) .1NE. 3))

DO 6300 LI1,NFCK

KL = KL + 1

IF ((ICOMM(IJ,KL). NE.2).OR.$BURN(K,L)) GO TO 6300

IF ($ORNT2)":GO TO 6100

QRAD: :=. 'GAMMA (ITYPE)*QDOT2P(.I, J)*AREA (I)

DELTQ =.(QRAD, + QCOND)/AREAK

GO TO 6200

6100 DELTQ = QCOND/AREAK

IF (-(I..EQ.K).AND.(L.GTJ)) DELTQ =..QDOT2P(I J)i-:

6200 QEXT(K,L) = QEXT(K,L) + DELTQ*FCTR(11)

$PIGN(K,L) = .TRUE.

6300 CONTINUE

6400 CONTINUE

6600 CONTINUE

6900 CONTINUE

RETURN .

END,

- 80 -
C *

C SUBROUTINE*IGNIT (DETERMINES IF EXPOSEDFUEL ELEMENTS IGNITE)

C . t . , .

SUBROUTINE IGNIT.(QEXT, SUMQSQ:,TEMP,DIFF,FIGTP,FIGTS,THK,ITYP,NFC,

1 SBURN,$FMASS,$PIGN,$B,$C,$D,$W,IDIM-1', IDIM2;NFUEL,NSM)....

IMPLICIT LOGICAL ($)

DIMENSION QEXT(IDIM1,IDIM2) , SUMQSQ(IDIM.I,IDIM2),

1 TEMP(IDIMLIDIM2,-10), $BURN(IDIM14I:DIM2), $FMASS(IDIMI,IDIM2),

2 $PIGN(IDIM1,IDIM2)

DIMENSION DIFF(NFUEL), FIGTP(NFUEL), FIGTS(NFUEL), THK(NFUEL),

S - ITYP(NSM), NFC(NSM), $B(NSM), $C(NSM)., $D(NSM), $W(NSM)

COMMON /ALL/ ACEIL, AWALL, DELT, FCTR(15),,FLCF. RTEMP, TIME,

. "TITLE(20),

2 ICEIL, IJOB, INCHCK, INITG, IOUTPTMOUTPT(11),

.3 $DOOR, SEND, $ROOM

DO 8000 I=1,NSM,

IF ($D(1).OR.$W(1) OR.$B(I )OR.$C(I)) GO TO 8000

ITYPE ITYP(I),

.,NFCI = NFC(!)

DO 7000 J=1,NFCI

IF ($BURN(I,J).OR.$FMASS(I,J)) GOTO 7000

TIG = FIGTS(ITYPE)

IF ($PIGN(I,J)) TI• = FIGTP(ITYPE)

SUMQSQ(I.,J)_ = SUMQSQ(I,J) + QEXT(I,J)**2*DELT-<

IF (SUMQSQ(I,J).NE.0.) GO TO 6000

TWRIT = 1.E30

GO TO 6100

6000 IF (TIME.NE.0.) AVGHT = SUMQSQ(I,J)/TlIME

IF (TIME.EQ.0.) AVGHT = QEXT(I:,,J)**2

TIMIG = .7854*.(THK(ITYPE)*(TIG:ý TEMP(I,J.,l)))**2/

1 (DIFF(ITYPE)*AVGHT)

TCRIT = TIMIG*FCTR(14)

IF ($PIGN(I,J)) TCRIT = TIMIG*FCTR(15).

6100 IF ((TIME.LT.TCRIT)..OR. (TIME:.EQ.0.))" GO TO 7000

$BURN(I,J) = *.TRUE"

$PIGN(I,J) = .TRUE.

7000 CONTINUE

8000 CONTINUE.

RETURN

END

- 81.-
C*

C SUBROUT-INE•0OUTPUT:• (PRINTS PROGRAM OUTPUT) *

C *

SUBROUTINE OUTPUT(FLHT;:FLTEMP,FMASS,QDOT2P,QEXT,TEMP, $BURN,.NFC,

1 MSMOUT, IDIMI ,IDIM2,NSMý,NSMOUT,NWRITE)

IMPLICIT LOGICAL ($)

:DIMENSION FLHT(IDIM1,IDIM2), FLTEMP(IDIM1,IDIM2),

1 FMASS(IDIM1;IDIM2), QDOT2P(IDIM1,IDIM2),

2 QEXT(IDIMý1,IDIM2), TEMP(IDIM1,ID.IM2,10),

3 $BURN(IDiMI,IDIM2), MSMOUT(NSMOUT), NFC(NSM);".

COMMON /ALL/' ACEIL, AWALL, DELT, FCTR(15), FLCF, RTEMP,-TIME,

1- TITLE(20),

2 ICEIL, IJOB., INCHCK, INITG, IOUTPT, MOUT-PT(11),.

3 :;$DOOR, $END, $ROOM

COMMON /GAS/ A4, CDOOR, CF, CFV, DCF, DENIDG, DHGT, DWID,

I FC, FH, GABSRP, HCEIL, PLCF1, PLCF2,

2 " TG, TGO, THETA, VFVjWFL, WFVC,,"WFVH, ZD, ZN

COMMON /VENTC/ AFUEL, AFUELO, CVA, CVAO, FHTA, ,FRADA, FUELA, QMAX,

1 QTOT, QTOTO,'QTOT1, QTOT2, TMDOT,-TMDOTO, -

2 TMDOTS, TQDOT, TQDOTC, WIN, ZOA,

3 NBURN, NBURN0,.

4 ," $DECAY, $FLUX1, $FLUX2, $VENT, $VCONT

WRITE (NWRITE,9000) TITLE, TIME

$OUT .= (NSMOUT.EQ.0)

IF (IOUTPT.NE.0) GO TO 2000

WRITE (-NWRITE,9100) TMDOT

IF ($OUT) RETURN

DO 1000 I=1,NSMOUT .

K = MSMOUT(I)

•NFCK = NFC(K)

WRITE (NWRITE,9500) K, ($BURN(K,J),J=1,NFCK)

1000 CONTINUE

RETURN

2000 DO 5000.II=i,IOUTPT

INDEX = MOUTPT(II)

GO TO (2100,2300,2500,2700,2900,2900,2900,2900,2900,2900,2900),

1 INDEX

2100 WRITE (NWRITE,9100) TMDOT

GO TO. 5000

2300 WRITE (NWRITE ,9200). TQDOT

GO TO 5000

2500 WRITE (NWRITE;9300)_TG

GO TO 5000

2700 WRITE (NWRITE,9400) DG

GO TO 5000

2900 IF ($OUT) GO TO 5000

WRITE (NWRITE,9450)

DO 4000 I=1,NSMOUT :. K = MsMOUT(i) , :...- . ..

NFCK = NFC(K)

GO TO (4000.,4000,4000,4000,3000,3100,3200,3300,34001,3500,3600),

1 INDEX

82-,
3000 WRITE (NWRITE,9500) K, ($BURN(KJ-),J=.) ,NCK)

GO TO 4000

3100 WRITE (NWRITE,9600) K, (QDOT2P:(KJ),J=1,NFCK

GO TO 4000

3200 WRITE (NWRITE, 9700) K, (QEXT(K,J),J=1, NFCK)

GO TO 4000

3300 WRITE (NWRITE,9800), K, (FMASS(K,J),J=1,NFCK)

GO TO 4000

3400 WRITE.(NWRITE,9900) K, (FLHT(K,J),J=1,NFCK)'

GO TO 4000

3500 WRITE (NWRITE,9930) K, (FLTEMP(K,J),JZ=INFCK

GO TO 4000

3600' WRITE .(NWRITE,9950) K, (TEMP (K,J,1;-),J=1-,NFCK

4000. CONTINUE

5000 CONTINUE

)i

..) i • ..

.)" " • =

9000

9100

9200

9300

9400

9450

9500

9600

9700-

-9800

9900

9930

9950

FORMAT

FORMAT

FORMAT

FORMAT

FORMAT

FORMAT

FORMAT

FORMAT

FORMAT

FORMAT

FORMAT

FORMAT

FORMAT

RETURN

END

(!Hl,20A4,/,' TIME (SEC):',F6..0,//)j:

(IH ,"TOTAL MASS BURNING RATE .(KG/S):

(IH ,'TOTAL HEAT RELEASE RATE :(W) :

(IH ,'HOT GAS LAYER TEMPERATURE (K).:

(1H ,'HOT GAS LAYER THICKNESS (M) :

(iHO)

(1H ,'MODULE:',I4,3X,'BURNING ?',1OG1

(1H ,'MODULE:',I4,3X,'SRCE 'FLUX',lOG1

(1H ,'MODULE:',14,3X,'EXT FLUX ',1OG1

(lH ,.'MODULE:',I4,3X,'FUEL MASS2,0OGI

(1H ,-'MODULE:',14,3X,'FLAME HGT',I1OG1

(1H ,'MODULE:',I4,3X,'FLAME. TMP',lOG1

(1H ,'MODULE:',I4,3X,'FUEL TEMP'.,10G1

',G12.4)

',G12.4)'

',G12.ý4)`

!,G12.4),ý

1.3)

1.3)

1.3)

1.3)

1.3)

t1.3)

1.. 3) - o .

- 83 -
C --

C

C SUBROUTINE ADJCTP FILLS IN BLOCKS OF -THE ADJACENCY MATRIX.

C IF NADJ(IND,2) EQUALS '999', MODULE NADJ(IND,1) DOES NOT

C COMMUNICATE WITH MODULE NADJ(3). IF NADJ(2) EQUALS '888',

C MODULE NADJ(IND,1) DOES NOT COMMUNICATE WITH MODULES NADJ(IND,3)

C THRU NADJ(IND,4)

C

C -----------------------------------------------------------------------

SUBROUTINE ADJCTP(IJKL,NFCL,NSTOP, IDIM3, ICOMM)

IMPLICIT LOGICAL ,($)

DIMENSION ICOMM(IDIM3, IDIM3)

IIJ = IJ

DO 2000 LI=1,NFCL

KKL =,KL

DO 1000 LL=I,NSTOP

KKL = KKL + 1

ICOMM(IIJ,KKL) = 0

1000 CONTINUE

IIJ = IIJ + 1

2000 CONTINUE

RETURN

END

- 84 -
C-,.

-- - - -/ - - -- - - - -. -- - - ... -- 

C

C SUBROUTINE BARR COMPUTES THE HEAT FLUX FROM-A THERMAL.BARRIER,

C ASSUMING THAT THE THERMAL WAVE PENETRATES THE:..BARRIER: INSTANTLY

C

C -- - - -. -. -. -. -. -. -. -. -. -. -. -. -. -.. -------.. . .- . .- ........ .. .....-- ....-- -- -

SUBROUTINE BARR(IDIM1,IDIM2,I,NFCL,NWRITE,TK,EPS,FDEP,

1 FZ,QDOT2P,QEXTO,TEMP)

IMPLICIT LOGICAL .($)

DIMENSION FZ(IDIM1,IDIM2), QDOT2P(IDIM1;, IDIM2),ý

1 QEXTO(IDIMIIDIM2),, TEMP(IDIM1,IDIM2,10)

COMMON /ALL/ ACEIL, AWALL, DELT, FCTR.(15), FLCF., RTEMP, TIME,

1 TITLE(20),

2 ICEIL, IJOB, INCHCK, INITG, IOUTPT, MOUTPT(11),

3 $DOOR, $END, $ROOM".

COMMON IGAS/ A4, CDOOR,' CF, CFV, DCF, DEN, DG, DHGT, DWID,

1 FC, FH, GABSRP, HCEIL, PLCFI, PLCF2,.

2 TG, TGO, THETA, VFV, WFL, WFVC, WFVH, ZD, ZN

DO 1000 J=1,NFCL

J1 = j

TE = RTEMP

IF ($ROOM.AND.(FZ(I,J).GT.ZD)) TE TG

BOLD = TEMP(I,J,I)

AOLD = (TEMP(I,J,10) - BOLD)/FDEP

JCOUNT = 0

500 TB1 = BOLD

TB2 = AOLD*FDEP + BOLD

F = -TK*AOLD - HCEIL*(TE - TB1) 5.6697E-8*EPS*(TE**4.. -=TB1**4)

1 - EPS*QEXTO(I,J)

G = -TK*AOLD - HCEIL*(TB2 - TE) - 5.6697E-8*EPS*(TB2**4 - TE**4)

-FA = -TK

FB= HCEIL + 2.268E-7*EPS*TB1**3

GA = -(TK + HCEIL*FD.EP + 2.268E-7*EPS*FDEP*TB2**3)

GB -(HCEIL + 2.268E-7*EPS*TB2**3)

D = FA*GB - GA*FB

ANEW = AOLD - (F*GB - G*FB)/D

BNEW = BOLD - (G*FA - F*GA)/D

TBIP = BNEW

TB2P = ANEW*FDEP + BNEW

ERR1 = ABS(TB1P - TB1)

ERR2 = ABS(TB2P - TB2)

IF ((ERR1.LE.1.).AND.(ERR2.LE.1.)) GO TO 700

IF (JCOUNT.GT.30) GO TO 8000

AOLD' = ANEW

BOLD- BNEW

JCOUNT = JCOUNT + 1

GO TO 500

700 TEMP(I,J,1) = TB1P

TEMP(I,J,10) = TB2P

QDOT2P(I,J) = 5. 6697E-8*EPS*(TB2P**4 -: TE**4)

1000 CONTINUE

GO TO 9000

8000 WRITE (NWRITE,8100) I, Ji, JCOUNT, TBIP, ERR1, TB2P, ERR2

8100 FORMAT (iHO,'ITERATION MAX EXCEEDED IN SUBROUTINE BARR:',/,

- 85 -
1.' MODULE:',!4-,3X,' FUEL CELL:'.,14,/,' NO.. OF ITERATIONS:' ,14,

2 I,' HOT WALL TEMPERATURE ::',Gl2.4,5X,'ERROR:',G12.4,/,

3 ' COLD WALL.TEMPERATURE"' ,G12-4,5X,'ERROR: ',G12.4)

$END= .TRUE.

9000 RETURN

END

- 86 .-
C

C SUBROUTINE CEILNG CALCULATES THE CEILING HOT GAS LAYER :BEHAVIOR

C USING "ANIMPLICIT DIFFERENCE SCHEME TO SOLVE THE HEAT CONDUCTION

C PROBLEM.FOR THE'CEILING, AND A STEPWISE ITERATION SCHEME FOR TG

SUBROUTINE oCEI LNG(EPS, WREFL,QCEIL,QEXTC,TCEIL,ZCEIL,NWRITE)

IMPLICIT LOGICAL :($)

DIMENSION Z(2,2)

COMMON;*,JALL/ ACEIL,- AWALL, DELT, FCTR(15), FLCF, RTEMP, TIME,

1 TITLE (20),

2 .CEIL, IJOB INCHCK,- INITG, IOUTPT,-MOUTPT0(1),

3 PDOOR, $END; $ROOM

COMMON /GAS/ A4, :CDOOR, CF, CFV, DCF, DEN, DG, DHGT, ýDWID,

I FC, FH, GABSRP, HCEIL, PLCF1, PLCF2,

2 TG, TG0, THETA-, VFV, WFL, .WFVC, WFVH, ZD., ZN:

COMMON /VENTC/ AFUEL, AFUEL0, CVA, CVAO, FHTA, FRADA, FUELA, QMAX,

1 QTOT, QTOT0, QTOT1, QTOT2, TMDOT, TMDOTO,.

2 TMDOTS, TQDOT, -TQDOTC, WIN, ZOA,

3 NBURN, NBURNO,

4 $DECAY, $FLUX1, $FLUX2, $VENT, $VCONT;

C

C CHARACTERIZE HOT GAS LAYER

C

WFL = TMDOT*(PLCF1*FHTA/FRADA + 1.)**2.5

FSIG = EPS*5.6697E-8

IF (VFV.'.LE.1.E-3) GO TO 1000

DF = FH - FCDFMX = (TQDOTC-+ ACEIL*(HCEIL*(TCEIL - ýRTEMP) + FSIG*(TCEIL**4

1 RTEMP**4)))/CF ."

A3 ="HCEIL/FSIG

A4 = -(HCEIL*TCEIL/FSIG + TCEIL**4 +(TQDOTC CF*DF)/(FSIG*

-1 ACEIL)).'

$DFLO = (DF.LE.0..)

$DFHI = (DF.GE.DFMX)

$DFLH = ($DFL0.OR. $DFHI).

IF ($DOOR.AND.$DFLH) GO TO 7000

IF ($DFLH-) GO 'TO 1000

C

C SOLUTION WHEN ZN > DHGT C "

TG = RADEQ (A3,A4)

DFR = FH/TG - FC/RTEMP

IF (DFR.LT..0.) GO TO 1000

WFVH CFV/TG

WPL =CFV*DFR

ZD (FRADA/PLCF1)*((WPL/TMDOT)**.4 - 1.) + ZOA

IF (ZD.GT.FHTA) ZD = (FRADA/PLCF2)*((WPL/WFL)**.6 1.) 1 

1 + FHTA

IF (ZD.LT.DHGT) GO TO 500

ZN DHGT.

GO TO 7999

500 WDC = WPL - DF*WFVC

- 87 -
._IF (WDC. LE. 0 GO. GO TO -:1000

ZD2 = 0.5*ZD/DHGT

Q = ZD2**2

R = ZD2**3 ;+ ,0.5*(WDC/(C*SQRT(I - RTEMPiTG)))

D R**2 +Q*3

ZN = DHGT*,((R + SQRT (D,)) *1,/3.) (R -SQRT(D-))**(1./,3..))

IF (ZN.GT.DHGT) GOT.79991

C

C SOLUTION WHEN ZN < DHGT.. "?:

C

1000 A4= (HCEIL*TCEIL/FSIG + TCEVL**4.+ (TQDOTC 4+,CF*FC)/(-FSIG*

1 ACEIL))

WPLMX = (.TQDOTC + ACEIL*(HCEIL*.(TCEIL- - .RTEMP) •+ FSIG* (TCEIL**4

1 ."- RTEMP**4)))/ (1000.*RTEMP)

IF (WFL.GE.WPLMX).%ZDMX . (FRADA/PLCF1)*((WPLMX/!TDOT-),**.4 - 1.)

+ ZOA

IF; (WFL.LT.WPLMX) ZDMX = (FRADA/PLCF2)*((WPLMX/WFL)**.6 - 1.)

+FHTA

IF (ZD:.;GE. ZDMX) ZD -0.08*ZDMX .

ZDO = 0.

ZDO1-= ZDO

ZD02- ZD0& -

ZD03 = ZDO

ZD1 0.

CALL :ZDNEW (ZDI-ZDFI,FSIGi)

TGP TG -

IF (ZDP1.LE.O.) GO TO 7800

ZD2 = AMIN1(DHGT,ý.995*ZDMX)

CALL ZDNEW (ZD2,ZDP2,FSIG,1),...

AsLOPE =,(ZDP2 ý- ZDP1)/(ZD2-"- ZD1)- .

IF- ((ASLOPE .GE. 0.,)..AND. (ZDP1. GT. ZD1-)., AND. (ZDPZ.GT0ZD2)). GO TO 7800

IF ('(ASLOPE.LE.0. ).ANDI.(ZDP2. GT. ZD2)) :GO 7900W ,

. •DO - 1500 J=1,2. . . "

CALL ZDNEW _(ZD ,ZDP,FSIG,-O)•

Z(J,1)-= ZD

Z(J,•2)= ZDP

IF (ABS(ZDP-ZD).LE.0.001), GO TO 7999 "

IF (ZDP.GE.ZD) ZD:- ZD +: THETA*(AMIN1(ZDPV DMX.),- ZD)

IF (zDP.LT.ZD) ZD = ZD + THETA (AMAX1:(ZDP,:0.) ZD) ,/:

1500 CONTINUE

DO 6000 ITR=1,20

_Z12 7"Z(',2)

Z22 =-Z(2,2),

IF (ZI2.GT.Z22) GO TO 1750

Zll-= Z(l,l)

Z21= Z(2,1)

Z(1,1) = Z21

Z(,1,2) = Z22'.

Z(2,1) = Z1l

Z(2,2) = Z12

1750.SLOPE = (Z(2,2).- Z(1,2))/(Z(2,1) - Z(II))

ZD = (Z(l;2)- SLOPE*Z(1,1))/(I. :- SLOPE)"

IF (ZD.GE.ZDMX) ZD = Z(2,1) + THETA*(ZDMX -Z(2;1)).

IF (ZD.LE.0.) ZD = (1. - THETA)*Z(1,1)

-.. ,88 -`
IF ((SLOPE.GE.0'.).AND. (ASLOPE.LE,0.).AND.(ZDP2.LE.ZD2).AND.

1 (ABS(ZD-Z(2,1)).LE.I.E-3)) ZD = ZD +. THETA*(ZDMX-ZD)

CALL ZDNEW (ZD,ZDP,FSIG,O)

DZ = ZDP -ZD

IF (ABS(DZ).LE.O.OO1) GO TO 7999

IF (Z(2,2).LT.ZD) GO TO 2000

Z21 Z(2,1)

Z22 = Z(2,2)

Z(2,1) = ZD

Z(2,2) = ZDP

..Z(I,1) = Z21

Z(1,2) = Z22

GO TO ,5000

2000 IF (Z(2,1).GT.ZD) GO TO 21W00

Zil = Z(1,1)

Z12 = Z(1,2)

Z(1,i) = ZD

-Z(1,2) = ZDP

Z(2,1) = Zll

Z(2,2) = Z12

GO TO 5000

2100 IF (DZ.LT.0.) GO TO 2200

Z(1,1) = ZD

Z(1,2) = ZDP

GO TO 5000

2200 Z(2,1) = ZD

Z(2,2) = ZDP

5000 ZD01 = ZD02

ZD02 = ZD03

ZD03 = ZD

IF (ITR.LE.3) GO TO 6000

ZDOA = (ZD01 + ZD02 + ZD03)/3.

ERRMSQ = ((ZDO1-ZDOA)**2 + (ZD02-ZDOA)**2 + (ZD03-ZDOA)**2)/3.

ERRMSQ = SQRT(ERRMSQ)

IF (ERRMSQ.LE.I.E-4) GO TO 7999

IF (ITR.EQ.20) GO TO 8000

6000 CONTINUE

GO TO 7999

7000 IF ($DFHI) GO TO.7500

TG = RADEQ (A3,A4)

ZD = 0.

ZN =0.

GO TO 7999

7500 TG = RTEMP

ZD = 0.

ZN = 0.

GO TO 7999

7800 TG = TGP

ZD = 0.

ZN = 0.

GO TO 7999

7900. ZD = DHGT

ZN = DHGT

7999 WDC = CDOOR*SQRT(RTEMP/TG)*SQRT(.. - RTEMP/TG)*(1. - ZN/DHGT)**- 5

- 89 -
WIN = WDC-+ WFVC

DG = ZCEIL -ZD. ý J , 

C DETERMINE HEAT FLUX FROM GAS LAYER

C

EX1 = EXP(-1.5*GABSRP*DG)

QCEIL = 5.6697E-8*'-(EPS*EX!*TCEIL**4

1 (1. -' EX1)*(1. + WREFL*EX1)*TG**4) + WREFL*QEXTC*EX1**2

GO TO 9000

8000 WRITE (NWRITE,8100) ZN, ZD, TG

'8100 FORMAT (1H0,"INNER ITERATION MAX EXCEEDED IN SUBROUTINE CEILNG:',

1 /, ' NEUTRAL PLANE HEIGHT:',G12.4,/,

2 ' PLANE HEIGHT AT DOOR:',G12.4,/,

2 ' GAS TEMPERATURE :.',G12.4)

$END = .TRUE.

9000.RETURN

END

"2.

,,• 90 -
C

C

C SUBROUTINE COMM CONSTRUCTS THE ADJACENCY (COMMUNICATION) MATRIX

C FOR THE INDIVIDUAL FUEL CELLS - 0 = NO COMMUNICATION,

C 1 = COMMUNICATION, 2 = ADJACENCY

C

C ----------------------------------------------------------------------

SUBROUTINE tOMM(NCOMINNCOM,NSM, IDIM1,IDIM3, IDIM4, IDIM5 ,NFC, IAD,

1 NAD,ICOMM)

IMPLICIT LOGICAL ($)

DIMENSION IAD(IDIM4,4), NAD(IDIM5,4), ICOMM(IDIM3,IDIM3),

1 NFC(IDIMI)

IF (NCOM.EQ.0) GO TO 1000

DO 900 IND=1,NCOM

IJ = 0

DO 800 I=1,NSM

NFCI = NFC(I)

IF (IAD(IND,1).EQ.I) GO TO 100

IJ = IJ + NFCI

GO TO 800

100 DO. 700 J=1,NFCI

IJ = IJ + 1

$IADJ = (IAD(IND,2).EQ.999)

$ID = ($IADJ.OR.(IAD(IND,2).EQ.J))

IF (.NOT.$ID) GO TO 700

KL = 0

DO 600 K=1,NSM

NFCK = NFC(K)

IF (IAD(IND,3).EQ.K) GO TO 200

KL = KL + NFCK

GO TO 600

200 IF ($IADJ) GO TO 400

DO 300 LI1,NFCK

KL = KL + 1

IF (IAD(IND,4).NE.L) GO TO 300

ICOMM(IJ,KL) = 2

ICOMM(KL,IJ) = 2

GO TO 900

300 CONTINUE

400 DO 500 L=1,NFCK

KL = KL + 1

ICOMM(IJKL) = 2

ICOMM(KL,IJ) = 2

IJ = IJ+ 1

500 CONTINUE

GO TO 900

600 CONTINUE

700 CONTINUE

800 CONTINUE

900 CONTINUE

1000 IF (NNCOM.EQ.0) GO TO 2000

DO 1900 'IND=1,NNCOM

IJ = 0

DO 1800 I=1,NSM

- 91 -
NFCI = NFC(I)

IF (NAD(IND,1).EQ.I) GO TO 1100

IJ = IJ + NFCI'

GO :TO 1800

1100 DO 1700 J=1,NFCI

IJ = IJ 4 1

$NADJP = (NAD(IND,2).EQ.888)

$NADJ = (NAD(IND,2).EQ.999)

$ND = ($NADJ.OR.$NADJP.OR.(NAD(IND,2).EQ.J))

IF (.NOT.$ND) GO TO 1700

KL= O

DO 1600 K=1,NSM

NFCK = NFC(K)

IF (NAD(IND,3).EQ.K) GO TO 1200

KL = KL + NFCK

GO TO 1600

1200 IF ($NADJP) GO TO 1330

IF ($NADJ) GO TO 1380

IF (NAD(IND,4).EQ.999) GO TO 1400

DO 1300 L=INFCK

KL:= KL + 1

IF (NAD(IND,4).NE.L) GO TO 1300

ICOMM(IJ,KL) = 0

GO:TO :1900

1300 CONTINUE

1330 NSTOP = 0

KL1 = NAD(IND,3)

KL2 = NAD(IND,4)

DO 1350 KLDUM=KL1,KL2

NSTOP = NSTOP + NFC(KLDUM)

1350 CONTINUE

CALL ADJCTP(IJ,KLNFCI,NSTOP,IDIM3,ICOMM)

GO TO 1900

1380 NSTOP = NFCK

CALL ADJCTP (IJKL,NFCI,NSTOP',IDIM3,ICOMM)

GO TO 1900

1400 DO 1500,L=I,NFCK

KL= KL +

ICOMM(IJ,KL) = 0

1500 CONTINUE

GO TO 1900

1600 CONTINUE

1700 CONTINUE

1800 CONTINUE

1900 CONTINUE

2000 IJ 0

DO 3000 I=1,NSM

NFC1 = NFC(I)

DO 2500 J=I,NFCI

IJ = IJ+÷1

ICQMM(IJ, IJ) = 0

KL = IJ + 1

IF (J.EQ.NFCI) GO TO 3000

ICOMf(IJ;,KL) 2=."2

- 92 -
ICOMM(KLIJ) = 2

2500 CONTINUE

3000 CONTINUE

RETURN

END

- .93 -
C

C FUNCTION CYLPAR FINDS THE HEAT FLUX FROM A CYLINDRICAL FLAME TO

C AN INFINITESIMAL SURFACE PARALLEL TO THE GROUND

C

C----------------------------------------------------------------7

FUNCTION CYLPAR(Z,FR,FRAD)

IMPLICIT LOGICAL ($)

IF (Z.GT.0.O) GO TO 100

CYLPAR = 0.0

RETURN

100 X = Z/FRAD

Y = FR/FRAD

YPLS = Y + 1.

YMIN = Y 1.

A = YPLS**2 + X**2

B = YMIN**2 + X**2

CYLPAR = ATAN(X/SQRT(Y**2 - 1.))

1 + X*((A - 2.*Y)*ATAN(SQRT(A*YMIN/(B*YPLS)))/SQRT(A*B)

2 - ATAN(SQRT(YMIN/YPLS)))

CYLPAR = CYLPAR/(3:14159*Y)

RETURN

END

- 94 -
C -. .--- ---------- - -- - - - -. -- -- - - -. ------- . . ' -- - -.. -.. --• -. ----------- . . -.. . -.. . '- -

C FUNCTION CYLPER FINDS-THE HEAT FLUX FROM .-A CYLINDRICAL"'FLAME TO

C AN INFINITESIMAL SURFACE' PERPENDICULAR:TO THE GROUND

C

C------ --------------------- -----------------------------------------------

FUNCTION CYLPER(Z,FRFRAD)

IMPLICIT LOGICAL ($)

PI = 3.14159

IF (Z.NE.0.0) GO TO 100

CYLPER = 0.0

RETURN

100 X =.Z./FRAD

Y = FR/FRAD

A = (1.+Y)**2 + X**2

B = (I.-Y)**2 + X**2

YMIN =.Y - 1.

YPLS = Y + 1.

CYLPER = (i./(PII*Y))*ATAN(X/SQRT(Y**2-1.)) + (X/PI)*

1 ((A-2.*Y)*ATAN(SQRT(A*YMIN/(B*YPLS)))/(Y*SQRT(A*B))

2 - (I./Y)*ATAN(SQRT(YMIN/YPLS)))

RETURN

END

- 95 -
C - - - -.-- - - - - - - - - - - - - ---- - - - - - - - - - - - - - -

C

C SUBROUTINE DIFFUS SOLVES THE TIME-DEPENDENT HEAT CONDUCTION

C EQUATION FOR GEOMETRY USING A CRANK-NICOLSON SCHEME

C

C-----------------------------------------------------------------------

SUBROUTINE DIFFUS(IDIMi:,IDIM2,Il,J1.,NWRITEDELXEPS,FZIJ,HL,

1 QXT,QXT0,TDIFFTK,TEMP).

IMPLICIT LOGICAL ($)

DIMENSION BETA(10), GAMMA(10):, D(10), V(10), TEMP(IDIMl,IDIM2.,10)

COMMON /ALL/ ACEIL, AWALL, DELT, FCTR(15), FLCF,',RTEMP, TIME,

1 TITLE(20),

2 ICEIL, IJOB, INCHCK, INITG, IOUTPT, MOUTPT(II),

3 $DOOR, $END, $ROOM.

COMMON /GAS/ A4, CDOOR, CF, CFV, :DCF, DEN,.DG, DHGT, DWID, .

1 FC, FH, GABSRP, HCEIL, PLCF1, PLCF2,.

2 TG, TG0, THETA, VFV,, WFL, WFVC,:WFVH,"ZD, ZN

C

C

C

INITIALIZATIONS

FSIG ="EPS*5.6697E-8

TLAM = -5*TDIFF*DELT/DELX**2

BJ = 1. + 2.*TLAN

Bl= -4.*TLAM

Cl = .5*TLAM

CB = -3.*DELX.

DO 3500 J=2,10

D(J) = TEMP(ILJl,J) + TLAM*(TEMP(Il,Jl,J-l) -2.*TEMP(Il,Jl,J)

i + TEMP(I.l,Jl,J+l))

3500 CONTINUE

VIOLD = TEMP(Il,Jl,l)

C

C

C

BEGIN ITERATIONS

DO 8000 ITR=I,20

I' T.W5: 1RlTNT fAPY

C

C

C

Al = 1. + 3.5*TLAM - TLAM*CB*(HL + FSIG*VlOLD**3)/TK

Dl = TEMP(Il,Jl,l) ±.÷.TLAM*(-3.5*TEMP(Il,Jl,l) + 4.*TEMP(I1,Ji,2)

1 -.5*TEMP(I1,J1,3))

2 + CB*HL*(-TLAM*TG + TLAM*(TEMP(Il,JI,I)- TGO))/TK

3 +.CB*(TLAM*(FSIG*(TEMP(Il,Jl,l)**4 - TGO**4) - EPS*QXTO)

4 - TLAM*(FSIG*TG**4 + EPS*QXT))/TK

C2 = -TLAM*(l. - Cl/Al)

TRANSFORMED COEFFICIENTS (BETA -AND.GAMMA)

BETA(2) = BJ + TLAM*B1/Al . GAMMA(2) = (D(2) + TLAM*Dl/Al)/BETA(2)

DO 6000 I=3-,9

BETA(I) =BJ - TLAM*TLAM/BETA(I-.)

GAMMA(I) = (D(I) + TLAM*GAMMA(.I-l))/BETA(I)

6000 CONTINUE

BETA(10) = BJ - TLAM*TLAM/BETA(9) :

GAMMA(10) = (D(10) + TLAM*(RTEMP + GAMMA(9)))/BETA(10)

.- 96
C

C SOLVE FOR TEMPERATURES

C

V(1O) = GAMMA(10)

DO 7000 1=2,8

K= 11-I

V(K) = GAMMA(K) + TLAM*V(K+1)/BETA(K)

7000 CONTINUE

V(2) = GAMMA(2) - C2*V(3)/BETA(2)

V(1) = (Dl-- BI*V(2) - CI*V(3))/AI

C

C CONVERGENCE BASED ON BOUNDARY TEMPERATURE ERROR

C

ERRL = ABS(V(1),- VlOLD)

IF (ERRL.LT.1.), GO TO 8900

VlOLD = V(!)

8000 CONTINUE

C

C TEMPERATURE DUMP IF ITERATION MAXIMUM EXCEEDED

C

WRITE (NWRITE,8100) Ii, J1, V(1), ERRL

8100 FORMAT (1Ho,'ITERATION MAX EXCEEDED IN SUBROUTINE DIFFUS:',/,

1 'MODULE:',14,3X,'FUEL CELL:',14,/,' ITERATION NO. 20',

2 /,' SURFACE TEMPERATURE:',EI5.7,3X,'ERROR:',E15.7)

$END = .TRUE.

GO TO 9000

8900 DO 895.0 I=1,10

TEMP(Il,J1,I) = V(I)

8950 CONTINUE

9000 RETURN

END

- 97 -
- - -.- , .-- 

C

C FUNCTION RADEQ SOLVES THE QUARTIC.EQUATION: X**4 + A3*X + A4 0.0 -

C----------------------- ---------------------------------------------

FUNCTION RADEQ (A3,A4)

IMPLICIT LOGICAL ($)

R =,0.5*A3**2

Q = -4.*A4/3.

D = R**2 + Q**3

Yl = (R.+ SQRT(D))*•(`I./3.) -'(SQRT(D).- R)**(-./3.)

RADEQ = .5*(SQRT(2.*SQRT(Yl**2 "-4.*A4) - Yl) -. SQRTt(Yl).:)

RETURN

END

- 98 -
C

C FUNCTITON RECT(A NEC,$P) FINDSITHE SHAPE FACTOR FROMs A DIFFERENTIAL -

C. ELEMENT TO A RECTANGLE OF DIMENSIONS A X B, ANDTSPACED C AWAY.

C FROM THE ELEMENT. z$r=. TRUE. => ELEMENT IS PARALLELITO RECTANGLE,-

C $P=.FALSE. => ELEMENT IS PERPENDICULAR -

C

C ----- ------------------------------- ------------ -------------------------

FUNCTION RECT(ZDEP,XLEN,YWID, $PARL)

!IMPLICIT LOGICAL ($)

" iF ((ZDEP'.EQ.0.DO).OR.(XLEN.EQ.O.DO).OR.(YWID.EQ.O.DO))

1 GO TO 200

'PI = 3.14159

IF (.NOT.$PARL) GO TO 100

X = XLEN/ZDEP:

Y = YWID/ZDEP

A = SQRT(1. + X**2)

"B = SQRT(l. + Y**2)

RECT = .5*(X*ATAN(Y/A)/A + Y*ATAN(X/B)/B)/PI

RETURN

100 X = XLEN/YWID

Y = ZDEP/YWID

A = SQRT(XLEN**2 + ZDEP**2)

RECT = .5*(ATAN(YWID/ZDEP) - ZDEP*ATAN(YWI.D/A)/A)/PI.

RETURN

200 IF ($PARL.AND.(ZDEP.EQ.O.DO)) RECT = .5

IF (.NOT..$PARL.AND.(ZDEP.EQ.0.DO)) RECT = .25

IF. ((XLEN.EQ.O.DO).OR.(YWID.EQ.0.DO)) RECT = 0.

RETURN

END

- 99 -
C

.C SUBROUTINE SETUP COMPUTES THE FUEL CELL LOCATIONS. FROM THEI.NPUT.

C SUPER-MODULE DATA -

C -- -- -- -- -- -- -- -- -- - - --- ------ ------ ---

SUBROUTINE SETUP(I, IDIREC, IDIM1, IDIM2,NFCL,SMX,SMY:,SMZ,SLNG,

1 "DS,FX,FYFZ).

IMPLICIT LOGICAL Cs)

DIMENSION FX(IDIM1,IDIM2), FY(IDIM1,IDIM2), FZ(IDIM1,IDIM2)`

IF (IDIREC.EQ.3).GO TO 400

IF (IDIREC.EQ.2) GO TO 200

FXO SMX -. (SLNG + DS)/2.

DO.100 J=1,NFCL

FX(I,J) = .FX0 + J*DS

ýFY(I,J) =.SMY

FZ(I,J) =. SMZ

100 CONTINUE

GO TO 600.,

200 FYO SMY - (SLNG + DS)/2..

DO 300 J=1,NFCL

FY(I,J) = FYO ,+-,J*DS

FX(I,J) = SMX.

FZ(IJ)- SMZ

300 CONTINUEGO TO 600

400 FZO " SMZ - (SLNG + DS)/2.,'-

DO 500 J=1,NFCL.-

FZ(I,J) = FZO + J*DS

FX (I;,J) = SMX

FY(I,J) = SMY

500 CONTINUE

600 CONTINUE

RETURN.

END

- 1001,
C -- . . . .. --- ----... . 11 ------ ----- --- ' -- -- - ------------

C

C FUNCTION SHAP .E.FINDS THE SHAPE.FACTOR FROM A.DIFFERE'NT.IAL FUEL

C ELEMENT TO AN ARBITRARILY LOCATED RECTANGLE, USING FUNCTION RECT

C

C - - - - - - - - - - - - - - - -- - - - - - - - - - - - ---. - - - - -

FUNCTION SHAPE(IOR,KOR,IDIRC,X1,Y1,Zl,X2,Y2,Z2,Dl,D2)

IMPLICIT LOGICAL ($)

DU = ABS(Xl - X2)

IF (IOR.EQ.2) DU =ABS(Y - Y2)

IF (IOR.EQ.3) DU = ABS(ZY - Z2)

GO TO (100,200,300), IDIRC

100 W0 = X2

WI = Xl - D!/2.

VO = Y2

VI =YI D2/2.

IF (IOR.NE.2) GO TO 400

VO = Z2

Vl = Zl - D2/2.

GO TO 400

200 WO = Y2

Wi = YI - D1/2.

VO = X2

Vl =X1 - D2/2.

IF,(IOR.NE.1) GO TC 400

VO Z2

VI = Zl - D2/2.

GO TO 400b

300 WO = Z2

Wl = Zl - DI/2.

VO = X2

Vi = Xl - D2/2.

IF (IOR.NE.1) GO TO 400

VO = Y2

Vl = Y1 - DZ/2.

400 W2 = Wl + Dl

V2 = Vi + D2

WSIGN = 1.

VSIGN = 1.

IF ((W0.GT.W2).OR.(WO.LT.Wl)) WSIGN = -1.

IF ((VO.GT.V2).OR.(VO.LT.V1)) VSIGN = -1.

DW2 = ABS(W2 -- WO)

DWI1= ABS(Wl -. WO)

DV2 = ABS(V2 - VO)

DVl = ABS(V1 - VO)

DWA = AMAXi(DW2,DWi)

DWB = AMIN1(DW2,DW1)

DVA = AMAXI(DV2,DVI)

DVB AMIN1(DV2,DV1)

$PR = (KOR.EQ.IOR)

IF (.NOT.$PR.AND.(KOR.NE.IDIRC)7) GO TO 500

SHAPE = RECT(DU,DWA,DVA,$PR) + VSIGN*RECT(DU,DWA,DVB,$PR)

1 WSIGN*RECT(DU,DWB,DVA,$PR) + VSIGN*WSIGN*RECT(DU,DWB,DVB,$PR)

RETURN

101 -
500 SHAPE = RECT(DU,DVA,!DWA,$PR) + VSIGN*RECT(DU,DVB,DWA,$PR)

1 WSIGN*RECT(DU,DVA,,DWB,$PR) + VSIGN*WSIGN*RECT(DU,DVB,DWB,$PR)

RETURN

END

- 102 -
C ------- ---------- VALUE

C

C SUBROUTINE ZDNEW FINDS A NEW VALUE FOR ZD GIVEN AN OLD VALUE

C

C------------------------------------------------------ --------------I

SUBROUTINE ZDNEW (ZIN,ZOUT,FSIG,IBRNCH),

IMPLICIT LOGICAL ($)

COMMON /ALL/ ACEIL, AWALL, DELT, FCTR(15), FLCF, RTEMP, TIME,

1 TITLE(20),

2 ICEIL, IJOB, INCHCK, INITG, IOUTPT,_MOUTPT(II),

3 $DOOR, $END, $ROOM -

COMMON /GAS/ A4, CDOOR, CF, CFV, DCF, DEN, DG, DHGT, DWID,

1 FC, FH, GABSRP, HCEIL, PLCF1, PLCF2,

2 TG, TGO, THETA, VFV, WFL, WFVC,, WFVH, ZD, ZN

COMMON /VENTC/ AFUEL, AFUELO, CVA, CVAO, FHTA, FRADA, FUELA, QMAX,

1 QTOT, QTOTO, QTOT1, QTOT2, TMDOT, TMDOTO0

2 TMDOTS, TQDOT, TQDOTC, WIN, ZOA,

3 NBURN, NBURNO,

4 $DECAY, $FLUX1, $FLUX2, $VENT, $VCONT

ZFLA = FHTA + ZOA

DO 1000 K=I1,50

IF -(ZIN.LT.ZOA) WPL = TMDOT

IF ((ZIN.GE.ZOA).AND.(ZIN.LE.ZFLA))

1 WPL = TMDOT*(PLCF1*(ZIN - ZOA)/FRADA + 1.)**2.5

IF (ZIN.GT.ZFLA).

1 WPL = WFL*(PLCF2*(ZIN -,ZFLA)/FRADA + 1.)**(5./3.i)

A3TH = (HCEIL*ACEIL + (WPL + FC*WFVC)*1000.)/(FSIG*ACEIL)

TG = RADEQ (A3TH,A4)

Cl = CDOOR*SQRT(I. - RTEMP/TG)'

C2 = CI*SQRT(RTEMP/TG)

WFVH = CFV/TG

RN = 1. - (((WPL - FH*WFVH + FC*WFVC)/C2)**2)**(1./3.)

ZN = RN*DHGT

A1ZD =3.*RN

A3ZD = 4.*((WPL + (FC FH)*WFVC)/CI - RN**3).

IF (IBRNCH.EQ.1) GO TO 2000

IF (A3ZD.GE.O.) GO TO 500

ZDO ZIN

GO0TO 2000

500 ZIN = ZDO + THETA*(ZIN - ZDO)

1000 CONTINUE

ZDO 0.

2000 R = -. 5*A3ZD- (A1ZD**3)/27.

Q =-(AIZD/3.)**2

D = R**2 + Q**3

IF (D.LT.O.) RD =RN*(2.*COS(ARCOS(R/SQRT(-Q**3))/3.) - 1.) -

IF: (D.GE.0.) RD = (ABS(R) + SQRT(D))**(l./3.)

. + ..(ABS(R) - SQRT(D))**(1./3.)

IF ((R.LT.0.).AND.(D.GE.0.)) RD =-RD

ZOUT = RD*DHGT

RETURN

END

- 103 -
p
NRC FORM.335 u." s. NUCLEAR REGULATORY I. 1. REPORT NUMEE R (AssignedubyDC) U A-81S NURE./CR-5239

BIBLIOGRAPHIC DATA SHEET, UCLA-ENG-8257

4. TITLE AND SUBTITLE'(Add Volume No.; if appropriate) 2. (Leave blank)

COMPBRN -A COMPUTER CODE-FOR MODELING :COMPARTMENT -R . T S

FIRES -3. RECIPIENT'SACCESSlON NO.

•7.)..AUTHOR(S) .5. DATE REPORT COMPLETED

Nathan,0. Siu MONTH IYEAR

-: _ __ ___ __, __ __ _ August 1982 9. PERFORMING ORGANIZATION NAME AND MAILING, ADDRESS (Include Zip Code) DATE REPORT ISSUED:

School of Engineering and Applied Science MONTH YEAR

University of California May 1983

Los Angeles, CA 90024 6. (Leave -blank)

8. (Leave blank)

12. SPONSORING ORGANIZATION NAME AND MAILING ADDRESS (Include Zip Code) ,' .. • .. :,IOPROJECT/TASK/WORKýUNIT .1. POETAS/O N NO" NO.' •

Division of Risk Analysis

Office of Nuclear Regulatory Research 11. FIN NO.

U.S,.. Nuclear • Regulatory . .. Commission . • NRC NRC FIN.... Fi N. B0493 93 ,: Washington, DC 20555.

13. TYPE OF REPORT PERIOD COVERED (Inclusive dates)

..COMPBRN Manual .. .- August 1982

15. SUPPLEMENTARY NOTES 14. (Leave blank)

16 ABTAT 20wod res

1•6. ABSTRACT (200 words or •leiss):}. %•;'•.I '. 

The computer codel COMPBRN deterministically models the behavior of fire in a

compartment. This manual presents information necessary to run COMPBRN,, including"

:descriptions of, required input and resulting output. Also included are a sample,

probl..em and a isting of the code (written in FORTRAN for an IBM 3033- compute`),.

This manual is to be used. in conjunction with NUREG/CR-2269, "Probabilistic Models

;.for, the Behavior of Compartment Fires," August 1981, which describes the fire

models employed by COMPBRN.

17. KEY WORDS AND DOCUMENT ANALYSIS. 17a. DESCRIPTORS

Fire, fire growth,

computer fires,

deterministic fire growth.

17b. IDENTIFIERS/OPEN-ENDED TERMS

18. AVAILABILITY STATEMENT 19. ,SECURITY CLASS (This report) 21. NO. OF PAGES Unclassified 103

Unlimited " 20"ECURITY CLASS (Thispage) 22. PRICE

•, .assified. S

NRC FORM 335 011-811