Text

Cobra Iiic/Mit Computer Code Manual,Vepco Version.
	ML18130A261
Person / Time
Site:	Surry, North Anna
Issue date:	08/31/1979
From:	Bowring R, Moreno P; MASSACHUSETTS INSTITUTE OF TECHNOLOGY, CAMBRIDGE
To:	;
Shared Package
ML18130A262	List: ML18130A261; ML18136A095;
References
	PROC-790831, NUDOCS 7910240550
	Download: ML18130A261 (284)
	v • d • e

I: I I .* I .,. I *1 *e I I 'I ,, I. I I I I I I I* ii .,, *.'I_ COBRA lllC/MIT *COMPUTER CODE MANUAL VEPCO VERSION AUGUST, 1979 -' . ' ' ' -I _., -

.36e/36<1 t.tr \o--1e-1-9 IOLL!()51tlt .

7910240 \..,) ., . p I I I I I '., I I I I I I I I I I I I I COMPUT::'.R COD:::

March 1975 :1ASSACHUS:::TTS DTSTITUT:S O? TECHNOLOSY 77 Massachusetts Avenue Cambridge, 02139 Principal

!nvestigato:r's:

Robert W. 3owring Pablo

?repared.

f o::s Electric Fower Research 3412 Hillview Avenue* Palo Alta, California

?raj ec-: )!anager:

3urt A. Zalotar I I ,, I I I I* I I I ' I I I I I* I I I I AUTHORS AND CONTRIBUTORS The present code COBRA IIIC/MIT was developed in two phases, the first phase being sponsored by the Electric Power Research Institute (EPRI) and the second by New England Electric System and Northeast Utilities Service Co. The first phase was carried out by Robert Bowring and in this period the bulk of the modifications introduced in COBRA IIIC were executed.

Research Assistants who made contributions to the project in this first phase include John Bartzis and Chong Chiu. Programming assistance was provided by Himagshu Bhattacharya in the Information Processing Center, M.I.T. Many useful suggestions were also made by the staff of Battelle Pacific Northwest tories (BNWL), in particular, by Donald Rowe. Their contributions are gratefully acknowledged.

In the second ohase, some new capabilities were introduced into the code and small oroblems with the code prepared in the previous phase were phase Nas carried out by Pablo Moreno under the direction of Professor Neil Todreas. ful suggestions from John Valente, Chong Chiu and Ehsan Khan are acknowledged.

I--2-I I TABLE OF CONTENTS Page I Chapter 1 -Overall Comparison of COBRA IIIC/MIT and COBRA IIIC ................................. . 10 I Chapter 2 -Outline of Improvements Incorporated in COBRA IIIC/MIT .............................. . 2.1 Reduction.in Running Time ............... . 12 I 12 2. 2 Dynamic Storage ......................... . 12 *I 2.3 Polynomial Expressions for Physical Properties

............................ . 2.4 New Correlations

........................ . 12 I 12 2.5 Analysis of BWR's in an Assembly to Assembly Basis ........................ . 13 I 2.6 Metal and Coolant Modal Average Power 13 I 2.7 Various Errors and Anomalies in COBRA IIIC were Remedied ......................... . 2.8 Input Data Presentations 13 I 13 Chapter 3 -Modifications Incorporated in COBRA IIIC/MIT ... 14 I 3.1 Reduction in Running Time ............... . 14 3 .1. 0 Introduction

..................... , . 14 I 3.1.1 Description of the Problem ....... . 3.1.2 Channel Topography 14 I 16 3.1.2.1 Channel Boundary Inter-action ................ . 16 *I 3.1.2.2 Channel Boundaries Within the Stripe ............ . 17 I 3 .1. 3 Setting the Array LOCA ........... . 3.1.4 Reducing the Size of the Striped 18 I Matrix ......................... . 20 I I I I I I I I I I I I I I I 'I I I I I TABLE OF CONTENTS -Continued 3.1.5 Setting the AAA Array in DIVERT ...... . 21 3.1.6 Processing the AAA Array in DECOMP and SOLVE ......................... . 22 3.1.7 Summary and Modifications

........... . 22 3.1.8 Tests and Results to Check the Reduc-3.2 tion in Running Time .............. . 23 3 .1. 8 .1 Test Cases . .. . . . . . .

.. . .. . . . .. 23 3.1.8.2 Timing Method .. . . . . .. . . . . ... . 24 3.1.8.3 Timing Runs ... . .. . . . . .. . . . . .

25 3.1.8.3.1 Original COBRA IIIC . . . . 25 3.1.8.3.2 New Cross-Flow Solution.

26 3 .1. 8. 3. 3 Faster "Rest of Calcu-lation ..... , .. , ... , . . . 2 7 3.1.8.4 Discussion . . . . . . . . . . . . . . . . . . .

28 3.1.8.4.1 Result of Modification..

28 3.1.8.4.2 Prediction for Reactor Case . . . . . . . . . . . . . . . . . .

28 3.1.8.4.3 Cross-Flow Solution Time. 28 3.1.8.4.4 Running Large Cases......

29 3.1.8.4.5 Scope for Further Improve-ment...................

30 3.1.8.5 Summary and Conclusions . . ... .. 30 Increase in Array Sizes..................

31 3.2.0 Introduction

.................... . 3.2.1 Method of Dynamic Storage ........ . 3.2.2 Modifications

............. , ..... . 3.2.3 Common Lists .................... . 3.2.3.1 Blank COMMON .. , ...... , .. . 32 32 33 33 34 3.3 3.4 TABLE OF CONTENTS -Continued 3.2.3.2 Named COMMON, COBRA 1 ........ 3.2.3.3 Named COMMON, COBRA 2 ........ 3.2.3.4 Named COMMON, COBRA 3 ........ 3.2.3.5 Named COMMONS, LINK 2' LINK 3 3. 2. 4 Subroutine CORE ......................

3.2.5 Input Data 34 34 34 35 35 36 37 3. 2. 6 Running Time ....................... . Polynomial Expressions for Physical 37 37 38 40 Properties

.............................. . 3. 3. 0 3.3.1 Extra 3.4.1 3.4.2 Introduction

....................... . Physical Properties

................. . Correlations

....................... . Smith Slip Ratio ................... . Baroczy Two-phase Friction Pressure Drop ............................. . 3.4.2.0 Introduction

................ . 3.4.2.1 Description of Correlation

... 3.4.2.2 Form of Correlation for Coding in COBRA IIIC ............. . 3.4.2.3 Construction of Multiplier Array ..................... . 40 40 40 42 43 43 3.4.2.3.1 Stage 1: Calculate Physical Property Index (PPI)... 43 3.4.2.3.2 Stage 2: Calculate ers (M) at G=l.O Mlb/ft2hr 4li 3.4.2.3.3 Stage 3: Calculate Multipli-46 er at all values of G ... 3. 4. 2. 3. 4 Coding Sequence .......... . 47 3.4.2.3.5 Multiplier Values ........ . 48 I I* I I I *1 I I I I I I I ,, I I I I ,, j I I I I ,, I I I ,, I I TABLE OF CONTENTS -Continued 3.4.2.4 Interpolation for G and X.. .. . . . 48 3.4.2.5 Accuracy of Coding Baroczy Functions. . . . . . . . . . . . . . . . . . . . . 4 8 3.4.2.5.1 Physical Properties.......

49 3.4.2.5.2 Multiplier at a Mass Velo-city of 1.0 MLb/ft2hr...

50 3.4.2.5.3 Mass Velocity Correction Factors.................

51 3.4.2.5.4 Interpolation.............

51 3.4.2.5.5 Overall Errors............

51 3.4.3 Thom Two-Phase Heat Transfer Coefficients

... o ...... o ***** o ** o * * * * * *

52 3.4.3.0 Introduction....................

52 3.4.3.1 The Equations of Thom et al..... 52 3.4.3.2 Use of Heat Transfer Coefficient in Code .............

,, ... o o o o o

5 3 3.4.3.2.1 Heat Flux.................

53 3 . 4. 3. 2 2 JBOIL. . . . . . . . . . . . . . . . . . . . . 5 4 3.4.3.3 Modifications...................

55 3.4.3.3.1 PROP......................

55 3.4.3.3.2 HCOOL (N, I, JJ)... ..... .. 56 3.5 Analysis of BWRs in an Assembly to Assembly Basis" ... o ** v ** Cl * , *** o ****** e **** o c ll " o o o" o o .. 56 3. 5. O Introduction............................

56 3.5.1 Iteration Theory........................

57 3.5.2 Iteration Strategy......................

58 TABLE OF CONTENTS -Continued 3.5.3 Coding................................

59 3 . 5 . 3 .

1 IND AT ........................ . 3.5.3.2 CARDS4 and CHAN .............. . 3 . 5 . 3 . 3 SC HEM:E . . . . . . . . . . . . . . . . . . . . . . . .

3. 5. 3. 4 SEPRAT ..................... .. 3.6 Fuel and Coolant Modal Average Powers ...... . 3. 6. 0 Introduction

......................... . 3.6.1 Power Input in the Original COBRA IIIC 3. 6. 2 Modifications

........................ . 3.6.2.1 IQP3 Trigger ................. . 3 . 6 . 2 . 2 QPR 3 ......................... . -3 . 6 . 2 . 3 HEAT . . . . . . . . . . . . . . . . . . . . . . . . . . 3 . 6 . 2 . 4 DIFFER ....................... . 3,7 Various Errors and Anomalies in COBRA IIIC .. 3.7.1 Negative to Heat Flux irt PROP ......... . 3.7.2 Array out of Range in DIFFER ......... . 3.7.3 Initialisation Quirk in COBRA IIIC .... 3,7,3.0 Introduction

................. . 3.7.3.1 Initialisation

............... . 3.7.3.2 Modification to COBRA IIIC .... 3.7.4 Anomalous Behavior of COBRA IIIC Using the Subcooled Void Option when Bulk Boiling at Channel Inlet .......... . 3.7.4.0 Introduction

................. . 3.7.4.1 Explanation

.................. . 59 59 60 60 61 61 62 63 64 64 65 65 65 65 66 67 67 67 69 . 69 69 70 I .1 I I I I I I I I I I I .., *I I I

I I I I I I I I I Chapter 4 -,, I 1-I I I .,. I I I TABLE OF CONTENTS -Continued 3.7.5 COBRA IIIC Modification to Prevent Occasional Overflow in Subroutine SCQUAL................................

71 3.8 ._-Input Data--for COBRA I-IIC/MIT

........ -..... -..... 71 3. 8. o Introduction. . . . . . . . . . . . . . . . . . . . . . . . . . . .

71 3.8.1 Input Data Presentation Based on the COBRA IIIC Presentation...............

73 3.8.l.l COBRA IIIC Input Data...........

73 3.8.1.2 Simplified COBRA IIIC Input Data to be Used for Assembly to Assembly Analysis of LWR......

74 3.8.2 New Input Data Presentation.............

75 Description of the Organization of COBRA IIIC/MIT Code and of its New Subroutines..................

77 4. O Introduction. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

77 4.1 Organization of the COBRA IIIC/MIT Code.......

77 4.1.1 Read in the Input Data..................

78 4.1.2 Print the Input Data....................

80 4.1.3 Thermal/Hydraulic Calculations..........

80 4.1.4 Print Out of the Results................

81 4. 2 Description of Subroutines. . . . . . . . . . . . . . . . . . . .

81 4 . 2 . 1 MAIN. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

81 4.2.2 INDAT ...................................

82 4 . 2

3 INPRIN . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

8 2 4 . 2 . 4 CALC . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

8 3 4 . 2 . 5 EXPRIN. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8 3 4 . 2 . 6 CARD 2 0 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8 3 4 . 2 . 7 ITHO ........................

. . . . . . . . . . . 8 3 4

2 . 8 CHAN . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8 4

I! 1: TABLE OF* CONTENTS -Continued I *page 4.2.9 MODEL .................................. . 84 I 4 . 2 . 10 0 PERA . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

4.2.11 FIZPRP ................................. . 84 I 84 4.2.12 TABLES ................................. . 4. 2 .13 READIN ................................. . 85 I 85 4.2.14 TIDY ................................... . 85 I 4. 2 . 15 PRECAL ................................. . 4.2.16 CARDSl ................................. . 85 I 85 4 . 2 . l 7 CARDS 4 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

402.18 CORE ..................................

I. 86 I' 86 4.2.19 SEPRAT ................................. . 86 .I 4 . 2 . 2 0 QPR 3 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 86 4.2.21 ACOL ................................... . 86 I 4.2.22 BAROC .................................. . 4. 2 . 2 3 SURTEN ......................... , ....... . 87 I 87 4. 2. 2 4 HAPROP ............................... , .. 87 I 4.2.25 Functions ROLIQ, ROVAP, HLIQ, HVAP and SATTEM" .............................. . Chapter 5 -Presentation of Results ........................... . 87 t 88 5. 0 Introduction

................................. . 88 *1 5. l COBRA IIIC ................................... . 88 5.2 COBRA IIIC/MIT with the Old Input Data Presen-I tation ...................................... . 89 5.3 COBRA IIIC/MIT with the New Input Data Presen-tation .........

o ............................. . 89 I Chapter 6 -Listing of the COBRA IIIC/MIT Code ................ . 90 1. I I I I I I I I '! I I I I I I I I ---**-*--*--

-.-***--------

---------* **----SA-Figure 1. Figure 2. Figure 3 . Figure 4 . Figure r-:J

Figure 6. Figure 7. Figure 8 . Figure 9. Figure 10. Figure 11. Figure 12. Table 1. Table 2 . Table 3. Table 4. Table 5. Table SA. Table 6. Table 7. LIST OF FIGURES Channel and Boundary Numbering Scheme Channel and Boundary Numbering Scheme Channel Layout for 10 to 128-channel cases Baroczy Two-Phase Friction Multiplier at G = 1. 0 Mlb/ft2 hr. Baroczy Mass Velocity Correction to Phase Friction Physical Property Index (PPI) vs. Pressure PPI = (µ1/µg)0.2/(p1/Pg)

Baroczy Two-Phase Friction Multiplier at G = 1.0 Mlb/ft2hr.

Baroczy Mass Velocity Correction to Phase Friction Multipliers Organization of COBRA IIIC/MIT Organization of the Reading of the Data When Card Group 20 is Selected Organization of the Reading of the Data When Card Group 20 is Not Selected Arrangement of Channels and Rods in the Example Problem LIST OF TABLES Comparison of Storage Requirements and Costs With Problem Size Using COBRA IIIC/MIT on an 93 94 . 95 96 97 98 .99 100 101 102 103 104 IBM 370/168 105 LOCA Array for Case of Figure 1 106 LOCA Array for Case of Figure 2 *107 Timing for Original COBRA IIIC (10 Axial Intervals) 108 Timing With New Cross-Flow and New nRest 11 Subroutines (10 Axial Intervals) 109 Running Times per Iteration

-in seconds 110 Coordinates of Two-Phase Pressure Drop Correlation 111 Coinparison of Actual and Calculated Values of PPI 112

..

le Table 8. Value of n in M = 1-X + Xn/PPI 113 9-. -*--

of Veloci_t_y I I Table l:J. Table 11. Factor for Interpolation Between in 9 114 Multiplier Arrays at Various Pressures 115 Comparison Between Costs of COBRA IIIC and COBRA IIIC/MIT for the Problem of Fig. 12 116

--9-APPENDICES

1. Subroutines ACOL, DIVERT, DECOMP and SOLVE 2. Input Data for 32-Channel Case 3. Timing Method 4. Faster "Rest of Hydraulic Calculations" 5. Subroutine CORE and BLOCK DATA 6. Subroutine ZIGET 7. Physical Properties Subroutines

8. Coding of Subroutine BAROC 9. SEPRAT Coding 10. Input Data Presentation based on that of COBRA IIIC 11. Simplified COBRA IIIC Input Data Presentation to be used for Assembly to Assembly Analysis of LWR 12. New Input Data Presentation

13. 14. 15. Results Results tat ion) Results tation) obtained with COBRA IIIC for obtained with COBRA IIIC/MIT for the example problem obtained with COBRA IIIC/MIT for the problem 16. Coding of Timing Subroutine PRNTIM 17. Subroutine DIFFER . 18. Listing of the COBRA IIIC/MIT Code the example problem (Old Input Data Pres en-(New Input Data Pres en-I I I Page I 117 . 124 I 127 131 I 1341 139 146 ,, 149 153 I 156 I 176 211 I 271 I 383 I 498 619 I 622 I 626 I II I. I I I I I I I ., I I ., I I I I I 1* I I .,, CHAPTER 1 OVERALL COMPARISON OF COBRA IIIC/MIT AND COBRA IIIC COBRA IIIC/MIT was developed, primarily, to be the hydraulic part of MEKIN(l), a neutronic-thermal/hydraulic code developed at M.I.T. by J. Stewart and R. Bowring. For this pur-pose, it was necessary to adapt a thermal-hydraulic code able to-analyze LWR.s using 200 channels to simulate the core. COBRA IIIC was the code selected at the starting point. The idea was to keep the same conservation equations (i.e. mass, energy and momentum) and the same basic organization of the code but make modifications in COBRA IIIC as necessary to solve such big problems.

Also be made to resolve some minor problems in COBRA IIIC and to introduce a streamlined method of Data Input. It was planned then, that COBRA IIIC/MIT should yield cal results to those of COBRA IIIC when applied to the same problems, but the spectrum of problems tfiat could be economically analyzed with COBRA IIIC/MIT would be much larger, and even for the same problems, some improvements would result. The two main characteristics of COBRA IIIC/MIT compared with COBRA IIIC, are its lower running time and the use of the Dynamic Storage Data Management scheme developed at the Savannah River Laboratory.

These two factors make possible the*simulation of cases with a rather large number of channels at a reasonable cost. As an example, in Table 1 the storage requirements and the cost of three cases analyzed using COBRA IIIC/MIT with an IBM 370/168 computer are presented.

Comparative requirements for COBRA IIIC are not included since the 30 and 101 channel cases cannot be run with this code. Some other minor characteristics are also introduced, all of which will be described in full detail in Chapter 3. They are outlined in Chapter 2 to provide the reader with the necessary general picture of the improvements.

I I I I I I I I I I I I I I I I t I I I I ' I I I I. I 'I I I I I II I I I I I CHAPTER 2 OUTLINE OF IMPROVEMENTS INCORPORATED IN COBRA IIIC/MIT 2.1 Reduction in Running Time The most important characteristic of the present code is the reduction in running time. It is possible to achieve this reduction by taking of the sparsity of the cross-flow matrix while keeping basically the same method of solution as in COBRA IIIC. This change makes it possible to analyze cases with up to 200 channels at a reasonable cost. 2.2 Dynamic Storage The dynamic storage data management scheme developed at the Savannah River Laboratory is used in this version of the code. Utilizing this scheme the code itself sets the computer space required for each particular case. 2.3 Polynomial Exnressions for Physical Prooerties While in COBRA III.C the physical properties of the wate:::' and

steam need to be in COBRA IIIC/MIT such properties are calculated within the code from polynomials at pressures specified by the Input Data. 2.4 New Correlations Three new correlations were introduced.

These are the Smith Slip Ratio, Baroczy Two-Phase Pressure Drop and the Thom Two Phase Heat Transfer Coefficient. 2.5 Analysis of BWRs. in an Assembly to Assembly Basis A different iteration strategy was programmed for BWRs since for them, there is no interconnection between fuel assemblies.

2.6 Metal and Coolant Modal Average Power COBRA IIIC does not allow different axial profile fluxes for different channels.

This limitation was removed in COBRA IIIC/ MIT where it is possible to Input different profiles for ent channels by giving the heat generated in each node of the core either for steady state or transient conditions.

2.7 Various Errors and Anomalies in COBRA IIIC were Remedied 2.8 Input Data Presentation In order to economically analyze cases with a large number of channels, some modifications are needed in the original Input Data Presentation of COBRA IIIC since it requires too many cards to describe the problem. In COBRA IIIC/MIT, two general options were made available.

One is actually that of COBRA IIIC but-with a new alternative available for assembly to assembly analysis of LWRs. The other is designed for problems where a rather large number of channels are used. It gives a. large reduction in Input cards. I I ' I I I t I 1* ,, I I I *I l I 'I I. I I I I I* :1 I ,, I I ,, I ., I I' I I I I I CHAPTER 3 MODIFICATIONS INCORPORATED IN COBRA III-CC.MIT) 3.1 Reduction in Running Time 3.1.0 Introduction It had been suggested by BNWL personnel that the running time of COBRA IIICC.2) increased as the cube of the number of channels.

This was due to the time spent in subroutines DIVERT, DECOMP and SOLVE in triple nested Do Loops to solve the set of cross-flow simultaneous equations.

Improvements had been made in COBRA IV by (a) modifying the method of constructing the coeffi-cient matrix and (b) using an iterative method of solving the equations instead of partial decomposition.

This gave a running time increasing more nearly linearly with the number of channels.

As a result of this, it became clear that the original COBRA IIIC would have a prohibitive running time for reactor cases of (say) 200 channels and hence the code would need to be modified.

.. The problem was examined at M.I.T. and the modifications are based on a verbal description by BNWL of their more efficient way of constructing the matrix and (b) an M.I.T. technique of solving the equations using basically the same method as before but taking advantage of the sparsity of the matrix. (3) 3.1.1 Descriotion of the Problem The hydraulic equations may be written down in such a way that for an axial interval, 11 guessing 11 the lateral pressure difference between channels, the cross-flows may be expressed as a set of N simultaneous equations, one for each boundary in turn, b = k N=l l=l when W. = cross-flow across jth boundary J akl = coefficients.

In the original COBRA IIIC, the coefficients akl and bk are calculated in DIVERT. The matrix akl is decomposed with partial pivoting in DECOMP to give all zeros below the diagonal.

The equations are solved in SOLVE to give the values of Wj. The original method is wasteful in .running time because no account is taken of the sparsity of the matrix ak 1. In a big problem, the vast majority of the values of akl would be zero because there is assumed to be no direct interaction between channels sufficiently separated.

However, these values of akl are calculated to be zero by the computer whereas it would be more economical to set them to zero. Furthermore in decomposing and -solving the matrix, every element has to be tested to see whether , it is finite or zero. It is this multiplicity of operations which gives the excessive running time. In the new method, firstly the boundaries influencing a cular boundary are identified and values of akl only calculated for these; the remainder are set to zero. Secondly, the boundary I I I I I I I ' I I .I I I I I I I I I I I ' I* I I I *1 I , I I I I I 'I. I I I I numbering scheme is organized so that finite coefficients of akl lie in a diagonal stripe, all *coefficients outside the stripe being zero. This striping" technique was suggested by Kent Hansen. ( 3) 3.1.2 Channel Tocograohy 3.1.2.1 Channel Boundary Interaction.

The cross-flow equation for boundary i-j involves the conditions in the channels i and j forming the boundary.

The conditions in channel i are determined by the local effects and by interaction across the four boundaries (square lattice assumed) i-j, i-i 1 , i-i 2 , i-i 3 with the four adjacent channels, j,.i 1-i 3. Similarly the conditions in channel j are determined by the local effects and interactions across the four boundaries j-i, j-j 1 , j-j 2 , j-j 3. Thus the equation for w 1 j contains terms involving the other six boundaries i-j 1 , i-j 2 , i-j 3 , j-i 1 , j-i 2 , j-i 3; no other boundaries influence the equation for wlj" This is illustrated in Figure 1. Boundary 9 is between channels 6 and 7. Channel 6 is bordered by boundaries 5, 8 and 12 and channel 7 by 6, 10, 13. Thus the equation for boundary 9 contains terms from boundaries 5, 6, 8, 10, 12, and 13. For a square lattice geometry, the maximum number of "other boundaries" is 6 and it may be smaller for edge channels (e.g. boundary 5 is influenced by only five others, namely, 1, 2, 8, 9, 12). When generating the matrix akl' only those values of l appropriate to k are used to calculate a; all other values of a are set to zero. The information relating k and l is stored in the array LOCA. 3.1.2.2 Channel Boundaries Within the Strine. The positions of the non-zero elements of the array akl depends upon the numbering system used for the boundaries and thus is under the control of the User. It is then necessary to find a system which will minimize the width of the diagonal stripe in the matrix. It is considered that the best method is that illustrated in Figure 1, i.e. numbering left to right, top to bottom. Boundary 9 is influenced by boundaries between 5 and 13 (i.e. 5, 6, 8, 10, 12, . 13) so that the stripe width would be 8 (= 13-5). However, boun-dary 12 is influenced by boundaries between 5 and 19, i.e., a stripe width of 14, and this would be the maximum value. In general, if there are N channels arranged as a square, there would be IN channels in a row and the stripe width would be 4/N -2. In decomposing the matrix and solving the equations, it is only necessary to operate within the stripe width. One would then expect the number of operations to be proportional to the product of the number of boundaries and the stripe width, i.e. approximately to N3 12 instead of between N 2 and N3 if the whole matrix must be searched.

For up to 10 channels, one would not expect any gain as the stripe width is comparable to the number of boundaries, but I I I I I I I I I I I I I I I ' I I I I I ,. I I I I* I I I I I I I I I I I I there should be a considerable gain in large problems.

The COBRA Input Data has been modified(4) so that the channel boundaries would be ordered as in Figure 1. 3.1.3 Setting the Array Loca The array LOCA contains two sets of information.

Firstly, it relates a particular boundary to the up-to 6 secondary boundaries which can influence it (as described in Section 3.1.2.1).

Secondly, it gives the sign of the cross-flow term. The convention used in COBRA is that for two channels I and J, the cross-flow direction I-J is taken to be positive when I < J and negative when I > J. Consider the coefficient akl where k is the principal boundary I-J and 1 is the secondary boundary J-M. After computing the value of akl it is multiplied by S where S = 1.0

{ S = -1.0 when I > J otherwise S = 1.0. to rulei: > M or I < J < M; The secondary boundaries in the LOCA array are given the sign of S. This is a quicker way of setting S in the construction of the akl matrix than calculating it ab initio each time. The LOCA array is set in the subroutine ACOL (coding in Appendix 1) immediately after the channel boundary Input Data are read in. It is called from INDAT, CARDS 4 or CHAN depending upon whether the original, simplified or new Input Data presentation is being used. LOCA (NK,8), where NK =number of channel boundaries is set as follows: ( LOCA(K,l) = K (i.e. the principal boundary) (LOCA (K,L), L = 2,7) = secondary boundaries; if none, set to zero. Sign as discussed above. LOCA(K,8)

=Total number of boundaries.

The boundaries of LOCA are illustrated in Table 2 by consider-ing some of the boundaries in Figure l. The method of setting LOCA is as follows: Do Loop K = 1, NK (number of boundaries)

LOCA(K,l)

= K (LOCA (K,L), L = 2,7) = 0 Identify II and JJ, the channels defining boundary K. Test all boundaries in inner DO Loop by identifying their channels III and JJJ. If either III or JJJ coincide with II, then a secondary boundary has been located and is stored in LOCA. It is given the sign of (II -LL)/(II -JJ) where the three channels are LL, II, JJ. Repeat inner Do Loop, reversing II and JJ, testing against old JJ. Set LDCA(K,8) to number of boundaries.

End Do Loop. The array name and size were chosen to make them the same as used in COBRA IV. The size could have been (NK,6) since LOCA(K,l)

= K I I I I ,, j' I I I ,,, I I I ', I' I 1 I' 'I I I .I I I I I ,, 1. I I ,. I *I I I ,. I I I is superfluous and the information stored in LOCA(K,8) could be obtained by testing for the number of non-zero elements.

However, the COBRA IV size has been used to simplify any future tion of the code into MEKIN. The subroutine ACOL was programmed without reference to the COBRA IV coding. COBRA IV probably does not use the negative sign convention in LOCA but calls the function S(K,I) to find the sign of akl' Subroutine ACOL ends by (a) calculating MS, the stripe width and (b) checking that enough storage has been allocated for the array AAA. In cases where IPILE = O, which imply that we may have more than four channels surrounding any single one, the LOCA array was changed to allow for the account of these channels.

The new LOCA array is now for these cases (LOCA (K,L), L = 1,14). This change will permit the simulation of cases like the one described in Figure 2 for which the LOCA array is that of Table 3. It is clear -that with this array a wider number of cases is covered but that number is limited and in some problems the split of channel 4 (in several channels) may be needed in order to keep the array inside the 14 spaces. 3.1.4 Reducing the Size of the Strioed Matrix With the boundaries numbered as described in Section 3.1.2.2, the non-zero elements of the matrix of akl coefficients would be within a diagonal stripe of width MS (calculated in ACOL). It is only necessary to store the numbers within the stripe and these may be contained in an array of size (NK, MS); for most cases this would be smaller than (NK, NK). To use the smaller array, an element AAA (K,L) in the square array is located at (K,MID-K+L) in the smaller array where MID is the mid point of the stripe, i.e., MID= (MS+l)/2.

Without dynamic storage, the full benefit of the reduced size of the AAA array cannot be achieved.

The array is set in MAIN as a single array and brought to DIVERT via SCHEME through argument lists. In DIVERT, it is set as a double array of size (NK,MS) and carried through to DECOMP and SOLVE via argument lists. This gives some flexibility in that the amount of the array used (and thus running time) and its shape can be problem-dependent.

Note that a check was made in ACOL that the product NK*MS did not exceed the available storage. 3.1.5 Setting the AAA Array in DIVERT The coefficients akl in the AAA array are set in subroutine DIVERT. It has been modified primarily by replacing the triple nested Do Loop (DO 80 K = 1, NKK) by a single Do Loop (DO 310 K=l, NK) -see Appendix 1. The procedure to calculate AAA and B is the same except (a) only the relevant boundaries, as stored in LOCA, are scanned, (b) for convenience the statement Function ABIT is used .and (c) the diagonal stripe is set into a smaller array as described in Section 3.1.4. Other modifications to DIVERT are I I I I I I I I I I I I I I I 1 I I I I I a I I I I I I I I I I I 1 I I I I relatively trivial. 3.1.6 Processing the AAA Array in DECOMP and SOLVE The AAA array is decomposed with partial pivoting in DECOMP. The subroutine has been modified by (a) removing the because of the strong diagonal dominance, (b) only scanning within the diagonal stripe (this saves running time) and (c) the diagonal stripe is set into the smaller array as described in Section 3.1.4. In SOLVE, the modifications are similar; only the diagonal stripe is scanned and the smaller array is used. 3.1.7 Summary of Modifications MAIN: (NOW INDAT) CARDS 4 OR CHAN: ACOL: SCHEME: DIVERT: DECOMP: Set AAA as a single dimension array of size MA (temporarily MGxMG) Set MA Call ACOL from Card Group 4 coding Call SCHEME with AAA in argument list Order boundaries Call ACOL Set LOCA from boundaries.

Calculate MS AAA in argument list Call DIVERT transferring AAA through argument Rewritten.

AAA transferred through argument, NK through COMMON and.MS through COMMON/ MEKIN. Set AAA (NK,MS) using LOCA, thus avoiding unnecessary calculation of zero elements Rewritten.

AAA transferred through argument (thus eliminating COM..MON/BUL).

Matrix composed, only scanning within diagonal stripe to save running time. I I SOLVE: Rewritten.

AAA transferred through Equation solved only scanning within diagonal shape to save running time. 3.1.8 Tests and Results to Check the Reduction in Running Time 3.1.8.l Test Cases. Test cases with 10, 16, 32, 64 and 128 channels were constructed.

The channels were typical of a PWR assembly with a flow area of 37.43 in 2 and containing 264 rods. They were identical in geometry but their radial power factors varied between 1.003 and 1.182 and were distributed about lines of symmetry (Figure 3). Thus all cases should give identical results for the 10 basic channels (lettered A-J in Figure 3) -this was checked -and should follow the same calculation path; this gives a valid basis of comparison between the running times for various cases. For economy, only 10 axial intervals were used.* A steady state calculation, requiring one hydraulic iteration, was followed by a transient calculation requiring seven iterations.

The Input Data for the 32-channel case are reproduced in Appendix 2 ; note that the simplified presentation (Appendix

11) of the original COBRA IIIC input Data Presentation is used for Card Group 1 (fluid properties) and for Card Groups 4, 7, and 8 (channel parameters).

The channel arrangement.

is in a 4 x 8 rectangular matrix with four symmetrical quadrants (Figure 3). I I I I I I I I I I I I I I I. I I I I I I I I ,. I I I I I I I I I I I 3.1.8.2 Timing Method. The timing method is detailed in Appendix 3, Briefly, a new subroutine PRNTIM (I) was called from strategic points in the program, the calling position being identified by the value of I from 0 to 8. The subroutine called TIMING (ITIM) to obtain the CPU time in centiseconds.

Using PRNTIM, the following were printed: (a) the time between entering and leaving SCHEME (i.e. the total time per hydraulic iteration) and (b) the sum of the times spent between entering and leaving DIVERT (i.e. the total time per iteration spent in solving the cross-flow equations).

The times given in Tables 4 and 5 are 11 Cross-flow" = (b) and "Rest"= (a) -(b). The following precautions were taken to obtain representative times.* Printing Time: The times at successive calls to PRNTIM were stored and processed within the subroutine but not printed until after the conclusion of the hydraulic calculations.

This avoided confusing the lengthy printing time with the time being measured.

Comoiler:

The results given in Tables 4 and 5 are for times compiled using the H Compiler.

It was found that its code optimisation capability reduced the running time of the cross-flow subroutines by nearly a factor of 2 (0.52 to 0.28 sec. for the 32-channel case) compared to the cheaper G Compiler which was_ used during modification development. Variation Between Runs. Nominally identical runs gave slightly different times when repeated.

The variation between individual sections of computation was random up to 10%, but the variation of the average time per iteration was smaller, up to (say) 5%; for example, the average times for the 32-channel case in successive runs were 0.275, 0.275, 0.278, 0.285 and 0.289 seconds. When several results were available for the same case, those presented in Tables 4 and 5 are for the median run of the set. Because of its importance in extrapolation to reactor conditions, the 128-channel case in Table 5 was repeated.

The two times agreed within 1%. Load Module Size. For economy, a "small" load module was used for the 10-, 16-and 32-channel runs and a "large" module for the 64-and 128-channel runs. As a check, the 32-channel case was also run with the large load module; the times agreed within 2% of the values with the small module. 3.1.8.3 Timing Runs. 3.1.8,3.1 Original COBRA IIIC: Timing runs were made with the original COBRA IIIC for 10, 16 and 32 channels; cases with more channels were prohibitively expensive.

The results are given in Table 4, which also (a) compares the times with values calculated from the.expression below and I I ,, .I 'I I I I I f I I I I I I I I I. I I I I I I I 1. I I I I I I I I I I I I (b) gives the estimated times for a PWR case with 193 channels and 10 axial intervals.

The IBM 370/168(MIT) computer was used. The results were fitted by the expressions below. The times (seconds) for 10 axial intervals and per hydraulic iteration are given for (a) the solution of the cross-flow equations (Subroutines DIVERT, DECOMP and SOLVE) and (b) the remainder of the hydraulic calculation (i.e. the time spent in subroutine SCHEME minus the cross-flow solution):

NC =No. of channels:

NK =No. of gap boundaries Cross-Flow:

t = 0.08 + 0.000235(NK) 2*8 "Rest": t = 0.0027(NC)l.75 3.1.8.3.2 New Cross-Flow Solution:

As said before subroutines DIVERT, DECOMP and SOLVE were replaced by versions with improved coding taking advantage of the sparsity of the matrix. The modifications consisted of: (a) DIVERT: Constructing the coefficient matrix AAA using only relevant gap boundaries, setting the remaining terms to zero; previously the zeros were calculated. (b) DECOMP: Operating a diagonal stripe of width MS in the matrix; outside the stripe, all elements are zeros hence need not be accessed.

Previously all zeros were operated upon. (c) DIVERT, DECOMP and SOLVE: The diagonal stripe was con-tained within an array NK

MS. This does not reduce the running time but saves storage. Table 5 gives the timing results obtained with the new flow solution subroutines.

By considering the three stages in the solution (see Section 3.1.8.4.3), one may obtain an expression for the time as follows: New Cross-Flow:

t = 0.00282NK

+ 0.00000837 NK

MS 2 + 0.0000127NK 2 3 .1. 8. 3. 3 Faster "Rest of Calculation":

It was noted that the running time for the Rest" in the original COBRA IIIC increased as (NC)1*75 (see Section 3.1.8.3.i) and was not pro-portional to NC as one might have expected.

This suggested scope for further improvement and subroutines DIFFER and HEAT seemed the most likely candidates.

The modifications are described in Appendix 4. Briefly, DIFFER was modified by changing the order of calculation of the summed contributions from neighboring channels in the mass balance, heat and pressure drop equations.

HEAT was modified for PWR and BWR cases only (in which there is a one-to-one correspondence between rods and channels) by avoiding the unnecessary examination of the effects of other channels on each rod heat flux and power. Table 5 gives the results of the timing runs. The times were closely proportional to the number of channels (the number of daries is roughly proportional to the number of channels) and were fitted by the expression:

New 11 Rest": t = 0.01 NC I I I I I I ,, I I I I I I I I I I I I I I I I I I I I I I I I I I I I *1 I .I 3.1.8.4 Discussion.

3.1.8.4.1 Result of Modification.

In the original COBRA IIIC, the running time increased nearly as NC3 where NC = number of channels -Table 4. In the modified version, it increases between Nc 1*5 and NC 2. This represents a reduction in running time by a factor 3 for NC = 10, 25 for NC = 32 and 220 for a reactor case of NC= 193 (Tables 4 and 5). 3.1.8.4.2 Prediction for Reactor Case. For a typical PWR with 20 intervals, the results in Tables 4 and 5 for 10 intervals should be doubled. Hence for the original COBRA IIIC, the estimated time per iteration would be about 1.8 hours9.259259e-5 days 0.00222 hours 1.322751e-5 weeks 3.044e-6 months whereas with the new version, it is reduced to about 30 seconds. Assuming (say) 6 iterations, the time might be about 3 costing ---$36.--*-----3.1.8.4.3 Cross-Flow Solution Time. One would expect the time to solve the cross-flow equation to be made up as follows: (a) Constructing the matrix: t NK. It is constructed by operating on the elements of the array LOCA (NK,8). (b) Decomposing the matrix: t NK*(MS)2. It is decomposed within a triple-nested Do Loop on NK, MS and MS, i.e. for each row of the matrix, the elements are traversed across and down the stripe width. In the original version of COBRA the Do Loop limits were NK, NK and NK giving a running time nearly NK3. (c) Solving the equations.

The number of operations is (n + 1) in solving for the (NK-n)th boundary, hence the total number of operations is proportional to NK(NK + 1) i.e. for NK large, approximately NK 2. Based on this notion, it was found that the running time could be expressed in the form a NK + b NK MS 2 + c NK 2 where a = 0.00282, b = 0.00000837, c = 0.0000127.

For the reactor 193-channel case in Table 5, the magnitude of the three therms is 1.01, 10.37 and 1.6. seconds; total = 12.99 seconds. Thus the middle term, which one may tentatively associate with the time required to decompose the matrix, predominates.

k Since approximately MS x NC 2 and MK x NC, the two dominant terms give the solution time as NC 2 , as was found in practice (see Section 3.1.8.4.1).

Earlier predictions that the time should be proportional to Nc 1*5 were based upon the misconcep-tion that the diagonal stripe was only traversed across, and not down, when decomposing the matrix. 3.1.8.4.4 Running Large Cases. It is ble that the 128-channel case is the largest ever run with COBRA IIIC; the running time using an unmodified version of the code would be high, 15 min. per iteration.

It augurs well for the use of COBRA in MEKIN that no stability problems associated with size were encountered and the number of iterations and answers were the same as for the 1/16 symmetry, 10-channel case. I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I 3.1.8.4.5 Scope for further imorovement.

The majority of the running time in large problems is taken in the cross-flow solution (see Table 5). In decomposing the matrix, operating only within the diagonal stripe reduced the number of operations but it is still relatively inefficient.

For example, in the reactor case there were 356 (i.e. NK) element.s per row in the original matrix of which up to 7 were non-zero.

Use of the striped matrix reduced the number of elements scanned from 356 to 59 (i.e. MS), This is a big saving, but of these 59, there are still at least 52 zeros which are scanned unnecessarily.

An iterative method, such as used in COBRA IV, which only operated on the finite elements in the matrix might be considerably faster even though several iterations would be required.

On the other hand, there might be problems of numerical stability.

Further investigation would be required to resolve these points. 3.1.8.5 Summary and Conclusions. (a) Timing runs have been made for the original and modified versions of COBRA IIIC. The modifications reduced the running time by a factor of 25 for a 32-channel case and by an expected factor of 220 for a PWR case of 193 channels.

For a PWR case with 20 axial intervals, the predicted running times are 1.8 hours9.259259e-5 days 0.00222 hours 1.322751e-5 weeks 3.044e-6 months per iteration with the original COBRA IIIC and 30 sec. per iteration with the modified version. (b) The subroutines modified were: Cross-Flow Solution:

DIVERT, DECOMP and SOLVE Hydraulic Calculation:

DIFFER, HEAT The modifications to the cross-flow solution comprised (a) economies in constructing the coefficient matrix (b) operating within a diagonal stripe when decomposing the matrix. The modifications to the rest of the hydraulic calculation involved summing the interactions between channels using a Do Loop over the number of boundaries instead of the number of channels; this avoided IF tests of whether boundaries were relevant.

Also, for PWR and BWR cases, the zero contribution of a rod to vant channels was not calculated. (c) It is predicted that for a PWR calculation with 20 axial intervals using the modified version of COBRA IIIC, the running time per iteration would be 26 sec. to solve the cross-flow equation and 4 sec. for the rest of the calculation.

There would be little to gain by further improving the "rest" but substantial gains might be possible if the cross-flow solution could be improved.

It is probably that the majority of time is spent in decomposing the matrix, 20.7 seconds out of 26, and since the time for this tion increases as NC 2 this accounts for the overall running time dependence.

3.2 Increase in Array Sizes I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I *1 I I I I 3.2.0 Introduction COBRA IIIC /MIT was written -using the Dynamic ( 5) Storage Data Management scheme developed at the Savannah River Laboratory; the storage allocation for each of the problem-dependent arrays is set according to the space required for the particular run being made. Then the same version of the code can be used for small and large problems.

3.2.1 Method of Dynamic Storage The technique, as applied to COBRA IIIC(2) is described briefly below. Problem dependent arrays had earlier been identified(6). These were then combined into a single array DATA( ), the position in the array corresponding to the variable considered.

Thus the channel areas, which in the original COBRA IIIC would be written as (A(I), I=l,MC), became (DATA($A+I),I=l,MC), where MC =number of channels and ZA is an integer. Thus, if ZA were present for the run to (say) 3412 and MC=30 then the channel areas would be stored in DATA(3413) to DATA(3442).

Double arrays were stored in a similar way, e.g., the channel flow rates ((F(I,J),J=l,MX),I=l,MC) would be stored as DATA(ZF+I+MC*(J-1)), etc. The convention was adopted that all names beginning with the Z sign were integers.

Also all problem-dependent variables in the original COBRA were prefixed by the Z sign to identify their position in the array DATA; for 6-character variables the last was dropped, e.g., HPERIM became ZHPERI. At the start of a run, the addresses in the DATA array were set in subroutine CORE (called from INDAT-Ref.

7), the maximum number of channels (MC), gap connections (MG), etc., being given as Input Data. 3.2.2 Modifications Every statement in COBRA containing a problem dependent array was identified (the H-compiler list of arrays was a help in this chore) and rewritten in the new form. Non-integer quantities were contained in that array DATA, integers in IDAT and logical quanti-ties in LDAT were DATA(l), IDAT(l), and LDAT(l) were equivalenced.

Examples of the modifications are given below. Old Array Size New A(I) (MC) DATA($A+I)

NTYPE(I) (MC) IDAT(ZNTYPE+I)

FDIV(K) (MG) LDAT(ZFDIV+K)

LENGTH(I) (MC) DATA($LENGT+I)

F(I,J) (MC, MX) DATA(ZF+I+MC*(J-1))

TROD(L,N.J) (MN, MR, MX)

MR*(J-1)))

In COBRA IIIC, in the above examples, FDIV was declared Logical and LENGTH Real, hence the use of LDAT and DATA respectively.

The problem dependent variables are identified in the Named Common COBRA3, described below. An example of dynamic storage programming is given in Appendix 17 for subroutine DIFFER. 3.2.3 Common Lists All COM..MON variables and Own Arrays in the original COBRA IIIC I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I were collected together into four main areas (see Appendix 17) with two other smaller Commons LINK2 and LINK3. 3.2.3.1 Blank COMMON. This contains only the array DATA, to which are equivalenced the* arrays LDAT (declared Logical) and IDAT (declared Integer).

3.2.3.2 Named COMMON, COBRAl. This contains all Scalar quantities in the original COBRA IIIC Commons. 3.2.3.3 Named COMMON, COBRA2. This contains all fixed arrays in COBRA IIIC, e.g. PP(30) (the values of pressure at which physical properties are tabulated).

3.2.3.4 Named COMMON, COBRA3. This contains the problem-dependent information.

Firstly the values of MA, MC, MG, MN, MR, MS, MX are given. MA has been previously set to NK*MS where NK =actual number of gap connections

(<MG); MC, MG, MN, MR, MX are the maximum values of the number of channels, gap connections, fuel nodes, rods and axial intervals respectively; MS is the stripe width of the AAA array in DECOMP (Section 3.1.4). Next follows 333, an integer denoting the number of positions, currently 96, in the DATA array. This is followed by the positions themselves, e.g., SA (which gives the start of the storage area for the channel areas) etc. 3.2.3.5 Named COMMONS LINK2, LINK3. These contain the variables which were required to link together subroutines INDAT, INPRIN, CALC and EXPRIN when they were separated in rewriting MAIN(?). LINK2 contains fixed variables (of COBRA2) and LINK3 scalar quantities (of COBRAl). The original LINKl which contained problem-dependent arrays is not required as its information is contained in COBRA3. They have been kept separate from the other Commons as they are not required in the thermal-hydraulic calculations but only in the Input Data setting. Their contents are given in Reference

7. 3.2.4 Subroutine CORE Subroutine CORE (Appendix

5) computes the origins within the DATA array of the old, fixed-dimensioned arrays A, NTYPE, etc. CORE is first called at the beginning of INDAT, immediately after a state-ment setting the problem size, i.e., READ MC,MG,MN,MR,MX.

In Subroutine CORE the integer array ZLX(I) is used as the length of the old arrays A,NTYPE, etc. Values of ZLX(I) are comouted for I=l,2, ... zzz, e.g. ZLX(l) = length of the A array = MC ZLX(8) -length of the CCHAN array = MR*MX In CORE the elements of the array ZORG(I) are effectively equiva-lenced to the origins in the COBRA3 common list, i.e. Z0RG(l) = ZA = origin of the A array Z0RG ( 2) = ZC_CHAN = origin of the CCHAN array I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I The elements of this array are initially computed as Z0RG(l) = 1 Z0RG(I) = Z0RG(I-l)+SLX(I-l)

I=2,ZZZ ZIGET is a subroutine written in assembly language (Appendix 6). It is called for dynamic allocation of core space for the DATA array (see Reference 5). The previously computed origins are incremented by the integer displacement returned to CORE from ZIGET, i.e. Z0RG(I) = S0RG(I)+KS-l; I=l,ZZZ In the first entry to CORE the length of the old array AAA has not been computed.

A second entry to CORE is made from ACOL(S) to establish the length and origin of the array AAA, i.e. ZLX(2)=MA and $0RG(2) = ZAAA-$0RG(l)+MA.

If at e_itper entry to CORE, sufficient storage is unavailable an error message is written to indicate how much addition storage is required for the problem. 3.2.5 Input Data When using the Dynamic Storage version of COBRA IIIC, an extra Input Data card must be given before the normal. Input Deck, namely Read MC, MG, MN, MR, MX (515). when MC > Number of channels in problem -' MG > Number of gap connections (for BWR cores give MG at least 3 for the s*pecial use of the array SP during iteration)

MN > Number of fuel nodes in problem (at least 1) MR > Number of rods MX > Number of axial intervals.

Then follow the Input Data. 3.2.6 Running Time The Dynamic Storage version of COBRA IIIC was checked against the earlier version for various cases. The agreement with the printed values of enthalpy and mass velocity was exact and the values of cross-flow agreed to five significant figures, i.e. the last digit printed might be different.

The effect in timing was then investigated.

The running times for 16-, 32-, 64-and 128-channel cases were compared with those previously obtained(9)(Section 3.1.8). The results are given in Table 5A. The dynamic storage appears to have increased the running time slightly, from 14% for the 16 channel case to 4% for the 128 channel case. This increase is considered to be acceptable.

3.3 Polynomial Expressions for Physical Properties 3.3.0 Introduction The physical properties of water and steam are calculated

-. in COBRA IIIC by interpolating between values read in as Input Data at a set of discrete values of pressure.

In COBRA IIIC/MIT the physical properties, at some fixed pressure, are calculated from the polynomial expressions used in the subchannel code HAMBO and that have been programmed now in CIERA IIIC/MIT.

These expres-sions have been used simply to generate the physical properties at fixed pressures instead of giving the values on cards, thus retaining the interpolation method within-the program. I I I I I I I I I I I I I I I I I I I I I I I I I -I I I I I I I I I I I I I 3.3.1 Physical Properties The properties are coded (see Appendix 7) in a number of functions and subroutines as follows: Liquid Density (lb/ft3) Vapour Density (lb/ft3) Liquid Enthalpy (Btu/lb) Vapour Enthalpy (Btu/lb) Saturation Temperature

(°F) Liquid Specific Heat (Btu/lb°F)

Liquid Viscosity (lb/ft hr) Liquid Conductivity (Btu/ft 2 hr°F) Surface Tension (lb/ft) Where p = pressure (psia) Function " " " " ) ROLIQ (P) ROVAP (P) HLIQ (P) HVAP (P) SATTEM (P) ) Subroutine HAPROP ) ) ( P , H , C P , XMU , XK ) ) Subroutine SURTEN(P,RL,RG,ST)

H = enthalpy (Btu/lb) -from HLIQ(P) CP = specific heat (Btu/lb op) XMU -viscosity (lb/ft.hr)

XK = thermal conductivity (Btu/ft 2 hr°F) RL = liquid density (lb/ft3) from ROLIQ(P) RG = vapour density (lb/ft3) from ROVAP(P) ST = surface tension (lb/ft) For numerical checking, values computed from the polynomials are tabulated below Pressure 15 200 1000 2000 2500 Temp. 213.15 381.99 544.68 635.87 668.14 Liq. Dens. 59,7887 54.3483 46.2835 38.8453 34.7882 Vap. Dens. 0.03803 0.43726 2.24609 5,33925 7.66218 Liq. Enth. 181. 207 355,705 542.631 672.346 730.884 Lat. Heat 969.594 842.484 648.896 461.815 359,708 Vap. Enth. 1150.801 1198.189 1191.527 1134.161 1090.592 Sp. Heat 1.00962 1. 07091 1.26524 1.64743 2.08775 Vise. 0.65919 0.33974 0.22546 0.18350 0.16945 Th. Cond. 0.39028 0.38524 0.32675 0.27480 0.25004 Surf. Ten. 0.004024 0.002666 0.001222 0.000464 0.000218 I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I* I I I I 3.4 Extra Correlations The extra correlations included are: 3.4.1 Smith Slip Ratio The Smith (10) slip ratio correlation is written as S = 0.4 + 0.6 {o.4 + X (pL/pg o.4 + o.6x 0. 4) }1/2 where S =slip ratio; X =quality; p 1 ,pg = liquid and vapor densities.

It is called by setting J3=2 in the Input Data thus adding to the existing options of J3=0 for S=l, J3=1 for the Armand/Massena slip ratio. Subroutine BVOID was modified by calculating the void fraction from *the quality and the Smith slip ratio when J3=2. 3.4.2 Baroczy Two-Phase Friction Pressure Droo 3.4.2.0 Introduction.

The Baroczy correlation (11) is complicated to program involving (a) a physical property index as a function of pressure and (b) interpolation in two graphs for the property index, quality and mass velocity.

In COBRA the system pressure remains constant during the hydrau-lie iterations but may change during a transient.

To minimize running time, an array is set at the start of the calculations representing the multiplier at a given pressure, for various values of quality and mass velocity.

This array is recalculated at the start of each time step when the pressure changes. During the hydraulic calculation, the multiplier is found by interpolation within the array. The coding for the correlation is contained within a new routin BAROC (IPART,P,Q,GWV,FMULT,PPI) where IPART=l or 2; P,Q, GWV =pressure (psia), quality (O<Q<l) and mass velocity (lb/ft 2 hr); FMULT =

multiplier, and PPI = physical property index. The coding of the correlation and its development is described below. It is written so that a call to BAROC with IPART=l and P = system pressure sets the array. In subsequent calls to BAROC with IPART=2, the FMULT, is obtained for the given values of GWV and W. The Baroczy multiplier is called by setting J4=2 as Input Data, thus adding to the earlier options of J4=0 for I I I I I I I I I I I IF(J4.EQ.2)

CALL I homogeneous theory, J4=1 for Armand,* etc. CALC in COBRA IIIC was modified by a statement:

BAROC(l,PREF,O.O,O.O,RUB,PPI) immediately before the start of the iteration loop DO 430 NN=l NTRIES VOID in COBRA IIIC was modified by a call to BAROC if J4=2, calculating the mass velocity from the flow rate DATA(ZF), thus GWV = 3600.0*DATA(ZF+I+MC*(J-l))/DATA(ZA+I)

IF(J4.EQ.2)

CALL BAROC(2,PREF,XP,GWV,DATA($PHI+I),PPI)

I I I I The call is made with the other J4 IF tests e.g. between VOID0520 and VOID0530.

It returns with the interpolated multiplier DATA(ZPHI+I).

1. I I I I I I I I I I I I I I I I I I I I 3.4.2.1 Description of Correlation.

Baroczy developed his correlation from experimental data for a number of fluids, e.g. boiling water, air-water, mercury-nitrogen, diesel oil-air, etc. He correlated the various fluids in terms of a Physical Property Index (PPI) defines as (µ1/µg)0*2;(p 1/pg) where µ = viscotity, p = density and subscripts l,g = liquid and gas. He constructed a family of curves (Figure 4) for a mass city of 1.0 Mlb/ft 2 hr giving the two phase friction multiplier (M) versus PPI for various constant qualitites.

He also tabulated values of M vs. PPI and quality (X) --Table 6. A second group of curves (Figure 5) gave a correction factor for four mass velocities (G), namely G = 0.25, 0.5, 2.0, 3.0 Mlb/ ft 2 hr, versus PPI for various constant qualities.

A two phase friction multiplier is calculated knowing the fluid, pressure (P), quality (X) and mass velocity (G) as follows: (1) From the pressure and physical properties, calculate PPI = (µl/µg)0.2/(pl/pg).

(2) From PPI and X interpolate in Figure 4 to find the plier at G = 1.0 Mlb/ft 2 hr. (3) From PPI, X and G, interpolate in Figure 5 to find the multiplier correction factor. (4) The required multiplier is the product of the value at G = 1.0 (stage 2) and the correction factor (stage 3). 3.4.2.2 Form of Correlation for Coding in COBRA IIIC. In COBRA, the two phase multiplier is calculated for each channel at each axial interval and for each hydraulic iteration.

It is therefore important to use an efficient coding technique to minimize the computer running time. In the COBRA application, only water need be considered and the pressure remains constant for each time step. With these simplifications, the calculation of the Baroczy multiplier may be made in two steps. (a) Once only, at the beginning of each time step, construct an array giving the corrected two phase multiplier at discrete values of quality discrete values of mass velocity.

The 14 values of quality are O, 0.1, 1.0, 3.5, 5.0, 7.5, 10, 15, 20, 30, 40, 60, 80, and 100%. The 7 values of mass velocity are O, 0.25, 0.5, 1.0, 2.0, 3.0, and 1000 Mlb/ft 2 hr. (b) Each time the multiplier is required during the hydraulic iteration, interpolate for G and X within the already constructed array. 3.4.2.3 Construction, of Multiolier Array. The stages in the calculation are as follows: I I I I I I I I I I I I I I I I 3.4.2.3.1 Stage 1: Calculate Physical Property Index (PPI). PPI is defined as (µ1/µa)0*2;(p 1/p ) which is a function I -0 g of pressure alone. .Examination of experimental data (Table 7) I I I I I I I I I I I I I* I I I I I 1* I I showed that PPI could be expressed as a functi9n of pressure as below. These values of PPI are compared with the values calculated from the physical properties in Section 3.4.2.5.1 and plotted irt Figure 6. PR = p/3204 (p in psia) p < 1429.5 osia -4 6 a = 2.46896 x 10 , b = 0.195508, c = -0.03141 3 d = 0.264363 p > 1429.5 osia a = 0.220112, b = -0.299745, c = o.440706 d = -0.325823 A = PPI 2 . 3 pR(a+bpR+cpR

+dpR ) = A/[ 1. 0-pR+A] (Al) 3.4.2.3.2 Stage 2: Calculate Multipliers (M) at G=l.O Mlb/ft 2 hr. Baroczy gives values of the multiplier (M) for various qualities and PPI, i.e., (µ1/µg)0*2;(p 1/pg) as in Table 6. Each point in the Table may be characterized by a value of the index n in the expression n M = 1 -X+X /(PPI) This expression tends to the correct limits of M=l at a quality X=O and M=l/PPI at The index varies between fairly narrow limits as shown in Table 8. For a given value of PPI, i.e., pressure, the value of Mat each quality is calculated from the expression above with n ob-tained from Table 8 by interpolating for ln n versus ln PPI (over most of the range of PPI) at each of the required qualities.

The detailed method of calculating the mulitplier aat a given value of PPI and at each value of Xis as follows: (a) If X=O.O, set M=l.O and if X=l.O, set M=l/PPI. In-terpolate in Table 6 for each value of X between O and 1.0 by the following method. (b) If PPI < 0.001, find M by log-log interpolation in Table 6

= 0.0001 and 0.001. (c) If PPI > 0.3, find M by log-log interpolation in Table 6 between PPI = 0.3 and 1.0 . . (d) For 0.004 < PPI < 0.03 and X < 0.10; quadratic interpola-tion for n in Table 8. At required X, calculate n at PPI = 0.004, 0.01, 0.03. Express n as a quadratic in PPI, i.e., ln(ni)=a.+b.lnPPI+c(lnPPI) 2 where l. l. xi a. bi Ci l. 0.001 1. 2621 0.6749 0.073 0.01 1.9551 1.0043 0.1097 1.4985 0.8408 . 0.035 0.0971 0.05 0.7965 0.5531 0.0673 0.075 0.771 0.5638 0.0713 0.1 0.4838 0.4793 0.0657 I I I I I I I I I I I I I I I I I I. I

I I I I I I I I I I I I I I I I I I I (e) Calculate M at X from value of n. (f) For all remaining points, calculate n at tabulated values of PPI above and below required value as* in Table 8. Find value of n by log-log interpolation on PPI. Calculate M at X from value of n. 3.4.2.3.3 Stage 3: Calculate Multiplier at all value of G. The curves of Correction Factor (F) vs. ln PPI as shown in Figure 5 were approximated by four straight lines. The intersection was rounded off as illustrated below. F. 1 J.-F. ]. P. 1 J.-c I ---:::-1' .,,,. ' ,.. I .,,,. _.,,,. c ........ I ---. ........ ..__ , ___ I ------L I P. ]. E ln PPI The two lines AC' and C'E intersect at C' where PPI=P.. Two ]. values of PPI (P. 1 , P.+1) are selected such that lnP. 1 = lnP.-0.15, J.-]. J.-]. lnPi+l = lnPi+0.15 (i.e., points Band D) and the corresponding values of the correction Fi-l'Fi+l calculated.

Between Pi-l and Pi+l' the correction is taken as the mean of the values on the lines BC'D and BD. Thus the values of correction calculated follow the line ABCDE. The correction factor (Fe) is taken to be linear with PPI between points: F = zl at PPI = 0.00026 c z2 0.0038 Z3 0.057 Z4 0.198 1.0 1.0 where the values of z 1 for each quality are given in Table 9. 3.4.2.3.4 Coding Sequence.

The construction of the multiplier array CORAB is coded in subroutine BAROC (Appendix 8). The subroutine is entered with IPART=l and the coding sequence is: (a) From the pressure (P), calculate PPI as in Section 3.4.2.3.1. (b) Calculate the multiplier at G=l.O Mlb/ft 2 hr at qualities o, 0.001, 0.01, 0.035, 0.05, 0.075, 0.1, 0.15, 0.2, 0.3, 0.4, 0.6, 0.8, 1.0, as described in Section 3.4.2.3.2.

The values of multiplier are set in the array (CORAB (I,4), I=l, 14). I I I I I I I I I I I I I I I (c) The correction factors are calculated within an outer DO I Loop on quality and an inner DO Loop on mass velocity.

The correction factor is 1.0 at both X=O and X=l.O for I I. I I I I I I I I I I I I I I I I I I I I all mass velocities

--cf. Figure 5. . 2 For G=lOOO Mlb/ft hr, the correction are all set to 1.0 as indicated by the trend of Figure 5 i.e., (CORAB (I,7), I=l, 14) = 1.0. The correction factors for G=0.25, 0.5, 2.0, and 3.0 Mlb/ft 2 hr are calculated as in Section 3.4.2.3.3.

The multipliers at these mass velocities are found by multiplying the value at G=l.O by the appropriate correction factor. The multipliers are then set in the array (CORAB (I,J), I=l,14) where J=2-6 for G=0.25, 0.5, 1.0, 2.0, 3.0. Finally, the multipliers at zero mass velocity are taken to be the same as at G=0.25 and are set in (CORAB (I,l), I=l,14). 3.4.2.3.5 Multiolier Values. Calculated Values of the Multipliers are given in Table 10 pressures of 500, 1000, and 2250 psia. 3.4.2.4 Interpolation for G and X. Cross-plots were made of the multiplier versus quality at various constant G, p and of the multiplier versus mass velocity at various constant G,. p. These showed. that it was reasonable to interpolate in the multiplier array linearly with quality, linearly with mass velocity below G = 1.0 Mlb/ft 2 hr. and linearly with reciprocal mass velocity above G = 1.0 Mlb/ft 2 hr. This is accomplished in the second part of the coding of BAROC (Appendix

8) using the statement function ZRECT. 3.4.2.5 Accuracy of Coding Baroczy Functions. 3.4.2.5.1 Physical Prooerties.

Values of PL' pg' u 1 , ug were obtained from References 12 and 13 Table 7 gives (a) the values taken for the physical properties, (b) the values of PPI calculated from theproperties and from the polynomial in Section 3.4.2.3.1 and (c) the error in the polynomial expression.

Table 7 shows that the error in the calculated value of PPI is less than 1% except below 250 psia and above 2500 psia. This is adequate for the reactor range of interest.

The polynomial values may also be compared with those used by Baroczy. In his Table 1 he quotes the values of PPI he used for steam-water data at various pressures.

These are tabulated below and compared with the values calculated from physical properties.

o psia PPI (Baroczy)

PPI (properties) 590 0.033 0.037 1000 0.057 0.066 1400 0.083 0.102 2000 0.148 0.173 It appears that Baroczy used values some 10-20% lower than the values calculated from the physical properties.

It is not clear why this should be so; an error of 10% corresponds to a 10% error in the density ratio or 60% in the viscosity ratio. The discrepancy seems too high to be explained in terms of the expected differences between experimental data. If the discrepancy is genuine, the values of multipliers obtained from the Baroczy correlation at a given pressure would be somewhat lower than the values Baroczy himself used in developing the correlation (see Figure 4). I I I I I I I I I I I I I I I I I I. I I I I I I I I I I I I I I I 3.4.2.5.2 Multiolier at a Mass Velocity of 1.0 Mlb/ft 2 hr. The calculated multipliers for the conditions in Table 6 (except X = 0.005, 0.02) are correct since they were coded as the data to interpolate between. The qualities 0.005, 0.02 were not given as data and the values calculated are tabulated below. p PPI Multiplier at G = 1. 0 Mlb/ft 2 hr. psi a From Table 6 Calculated x = 0.005 0.02 0.005 0.02 --0.001 5.80 16.o 5.31 16.1 'Vll 0.001 5.60 16.8 5.11 14.4 60 0.004 4.90 11.9 4.62 11.2 160 0.01 3.30 7.00 3.02 6.72 500 0.03 1. 55 2.57 1. 43 2.47 1400 0.1 1.12 1. 48 1.12 1.44 2600 0.3 1. 02 1.13 1. 03 1.14 The interpolation method results in small discrepancies at X = 0.005 and 0.02 between Baroczy values (Table 6) and those I calculated.

However, since the calculated values are correct at I I I I the neighboring values X = O, 0.1 and 0.035, the integrated pressure drop over a range of qualities would be only marginally in error. The values of the multiplier at intermediate values of PPI were checked by using the computer graph-plotter to reconstruct Figure 4 to the same scale using subroutine BAROC --see Figure 7. imposing the two Figures and holding them to the light, the lines coincided exactly (i.e. within, say, 3%) exceot for the small kinks at PPI = 0.03 in Figure 7. This and calculated values, confirmed the programming of the multiplier at G = 1.0 Mlb/ft 2 hr. 3.4.2.5.3 Mass Velocity Correction Factors. The coding of the mass velocity correction factors was also checked using the computer graph-plotter.

The Baroczy curves given as Figure 5 were obtained by photocopying and magnifying.

The computer was used to plot to the same scale and the two sets of errors imposed and held to the light. After some trial-and-error fitting, I I I I I I I I I I the two sets of curves agreed except over some of the regions I where two lines are paired one into the other (e.g. G = 0.25 Mlb/ft 2 hr., X = 10% to 20%, PPI = 0.057 in Figures 5, 8). Over most of the I curves, the agreement was within about 3.4.2.5.4 Interoolation.

Comparing the values of the multipliers obtained using the interpolation scheme and the values from a curve joining several points, the maximum difference was about 3%. 3.4.2.5.5 Overall Errors. For each stage in the calculation of the multipliers, the error in the computer value is comparable with that of reading from the published Baroczy curves. This is considered acceptable.

I I I I I I. I I I I I I I 3.4.3 Thom Two Phase Heat Transfer Coefficients 3.4.3.0 Introduction.

The fuel surface temperature in COBRA IIIC/MIT is calculated from Thorn's(l 4) modification of the Jens and Lottes(l 5) superheat equation.

This correlation was preferred to the Chen(l 6) boiling heat transfer correlation.

We describe here the modifications to subroutines PROP and I HCOOL in COBRA IIIC. I I I I I I I I I I I I 3.4.3.1 The Equations of Thom et al. Based on experi-mental data, Thorn et al correlated the heated surface temperature (Tw) under single phase and subcooled boiling conditions as follows: where Single Phase hsp = o.134(K/D)

Re 0*65 Pr 0*4 Tw = Tb + 0/hsp Subcooled Boiling Tw = Tsat + 0.07200.5/ep/1260 h = (T -Tb)/0 tp w D =equivalent hsp' htp = single phase and two phase coefficient (Btu/ft2hr°F) heat transfer k = thermal conductivity (Btu/fthr°F) p = pressure (psia) Pr, Re = Prandtl and Reynolds Numbers = Wall temperature

(°F) Tb, Tsat = bulk liquid and saturation temperatures

(°F) 0 = heat flux (Btu/ft 2 hr) These expressions were modified versions of the eralier cor-relations of: Dittus-Boelter(l 7): -Jens-Lottes(l 5): hsp = o.023(K/D)

Re 0*8 Pr 0*4 Tw = Tb + 1.900.25/ep/900 I I I I I I I 3.4.3.2 Use of Heat Transfer Coefficient in Code. The I heat transfer coefficient is used twice in the COBRA calculations, namely to: (1) calculate the rod temperature and heat flux, (2) set JBOIL, indicating hte onset of subcooled boiling., 3.4.3.2.1 Heat Flux. In COBRA IIIC, a forward differencing scheme is used to obtain some of the quantities required in the heat balance equation, i.e., they are based upon values obtained for the previous interval.

One of these quantities is the heat transfer coefficient which is used to calculate the new rod temperature profile and the new heat flux. These calculations are made in subroutine HEAT called from SCHEME, before the heat balance equation is solved, HEAT calls HCOOL to calculate the heat transfer coefficient using values, e.g., heat flux, flow rate, etc., for the previous interval.

I I I I I I I I I I. I I I I I I I I I I I I I I I I I I I I HEAT then calls TEMP to calculate the rod temperature profile from the heat transfer coefficient and the new heat generation within the fuel. Next, the heat flux is calculated from the heat transfer coefficient and the wall-fluid temperature difference; the flux is then used in the heat balance equation to find the interval-exit enthalpy.

The COBRA has been retained although there is an inconsistency.

The heat transfer coefficient is calculated from the heat flux from the previous interval and is then used to find the heat flux for the current interval.

Ideally, one should iterate so that the heat fluxes become consistent.

However, this would only be necessary during boiling (when the heat transfer coefficient is a function of heat flux). Under those conditions, the surface temperature is only slightly above the saturation value and the error would only marginally affect the calculated rod temperature profile. 3.4.3.2.2 JBOIL. The indicator JBOIL is tially set to zero and then reset in PROP to the number of the axial position (i.e., J) at which subcooled boiling could start. This is determined as the point at which the surface temperature calculated for single phase conditions would exceed the surface temperature if boiling were occurring Tb + 0/h > T + sp sat sup In the original PROP, hsD is calculated using the Boelter(l7) equation and the superheat ) from the sup Lottes(l 5) correlation.

These are programmed in PROP. However, it seemed inconsistent to use different equations to set JBOIL from those used to calculate the heat flux. PROP has therefore been modified so that it now calls HCOOL to obtain the same set of equations.

PROP is entered after the heat balance equation, hence in calculating JBOIL the updated values of flow rate, heat flux, etc., are used to calculate the wall temperatures.

It should be noted that JBOIL merely indicates the position at which subcooled boiling could start and does not appear to have any effect upon the hydraulic calculations (it is used in CHF2 to cal-culate critical heat flux). Steam generation is calculated to start using a different criterion, namely the subcooled boiling correlation, and this is not necessarily consistent with JBOIL. For this reason, positive quality, rather than JBOIL, is used as the trigger in HCOOL to determine whether the single or two-phase heat transfer correlations should be used. 3.4.3.3 Modifications.

3.4.3.3.1 PROP. Two cards in PROP are changed, namely the calculation of DATA(ZHFILM+I) and TLBOIL. Card PROP0900 becomes:.

DATA(ZHFILM+I)

= HCOOL(-1,I,J)

I I I I I I I I I I I I I I I I I I. I I I I ,, I I I I I I I I I I I I I II I Card PROP0940 becomes: TLBOIL = TF -DTWALL + HCOOL(-2,I,J)

For the first, HCOOL is entered with the trigger N=-1 and returns with HCOOL = single phase heat transfer coefficient.

For the second, HCOOL is entered with N=-2 and returns with HCOOL =wall superheat.

3.4.3.3.2 HCOOL(N,I,JJ).

The function HCOOL replaces the dummy HCOOL in the original COBRA IIIC which set the heat transfer coefficient to 5000/3600.

The index N refers to the rod number when called from HEAT and is a trigger when called from PROP. The indices I and JJ are channel number and axial position, respectively; JJ=J-1 from HEAT and J from PROP -see Section 3.4.3.2. For N=-1 (i.e., called from PROP), HCOOL is set to the single phase heat transfer coefficient (hsp). For N=-2 (i.e., called from PROP), the two-phase wall superheat (Tw-Tb) is returned.

For N > 0 (i.e., called from HEAT), hsp is returned for zero quality and htp for positive.

quality. 3,5 Analysis of BWRs in an Assembly to Assembly Basis 3.5.0 Introduction COBRA IIIC is properly used for interconnected channels e.g. subchannels in a cluster or fuel assemblies in a PWR. The inlet flow to each channel is specified and the equations solved to give a constant exit pressure for each channel. Because of the constraint of fixed inlet flows, the inlet pressures to each channel are not necessarily the same. It is found that for PWR's the solution is not sensitive to the inlet flow distribution assumed because the cross-flow largely obliterates its effect after a few calculation intervals.

However, for a BWR, there is no cross-flow or interaction between assemblies and therefore errors in the given inlet flow I I I I I I I could cause large errors in the solution obtained.

It is necessary I to change the calculation strategy so that the channel flow rates are iterated to give the same pressure drop across all channels.

This note describes the modifications made to COBRA IIIC for BWR calculations.

3.5.1 Iteration Theory For channel i, a flow rate of m. lb/hr gives a pressure drop l of pi psi. It is required to find a new flow rate Cm 1 + om 1) such that the pressure drop shall be p , where p is the same for all 0 0 channels.

This implies p 0 = pi. + (op.lorn. )om. l ]. ]. ( 1) I I I I I I I I I I I I I I I I I I I From the mass conservation equation, Eomi = O, whence by summing the set of equations for all channels:

= E(v.p.) ]. ]. Substituting p form (2) into (1) for each channel 0 3,5.2 Iteration Strategy ( 2) ( 3 ) Equation (3) shows that if (op/om) is known for every channel, I the correction to the flow in each channel can be calculated direct-1 ly. Since there is a non-linear relationship between p and m, a I I I I I I I I I number of iterations would be required.

The value of v. for each channel is set to 6m./6P. where 6p., ]. ]. . ]. ]. 6mi are the changes in pi and mi for each channel between successive iterations.

For the first iteration, an approximate value of vi is obtained from the following argument.

The pressure drop across a channel is made up of friction, acceleration, elevation and grid components.

Of these, all but the elevation component (6p ) can be written as e a power function of m: whence = (l/n)m Vi D -6p . e The value of n is approximately 1.8 in single phase, decreasing considerably in the two phase region. A compromise value of n is 1.4 whence l/n = 0.7. The expression for v does not have to be accurate; it only has to start the iteration off in the right direction.

Thus we have: ( 4) Iteration 2 et seq.: vi= ( 5 ) where = elevation pressure drop component = changes in p and m between iterations.

3.5.3 Coding 3.5.3.1 In INDAT (nart of the old MAIN). BWR tions are called by setting IPILE=2 as Input Data. IPILE is a new parameter and is given on the Control Card (see Appendix 11). The variable J7 is set equal to IPILE for transfer between subroutines via the Common List; J7 was already in Common and was not used for any other purpose. The variable (WP(K),K=l,MG) is initialized to zero before reading in the Input Data. This is w', the turbulent (fluctuating) cross flow, which is zero for separated channels.

It would be reset subsequently in PWR cases. 3.5.3.2 In CARDS 4 and CHAN. The Input Data are organized in this subroutine.

For BWR's (IPILE=2), the interaction between channels must be suppressed.

Firstly, each channel has no adjacent channels, i.e. ((LC(I,L),I=l,NCHANL),L=l,4)=0.

Secondly I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I NSCBC, NBBC, J5, ABETA, BBETA, GK are all set to zero; these are respectively subcooled mixing option, two phase mixing option, thermal conduction option, two coefficients in the mixing expression and the geometry factor for conduction.

3.5.3.3 In SCHEME. This subroutine contains the Do Loop over the channel intervals and organizes the hydraulic calcula-tion. The modifications are: (a) Set IPILE = J7 If IPILE = 2 (b) skip the call to MIX at card SCHM0820 (c) Call SEPRAT (l,J,JUMP) after the heat balance equation; thence skip the calculation of cross flow by jumping to the pressure drop calculation at card SCHM1300. (d) call SEPRAT (2,J,JUMP) immediately before RETURN. 3. 5. 3. 4 In SEPRAT. This subroutine

.(see Appendix 9) is new and contains the bulk of the modifications required for BWR calculations.

In the first call to SEPRAT after the heat balance equation, the flow rate P(I,J) is set taking transient effects into account; this is equivalent to card SCHM1190 but with DFDX=O. Then follows a call to DIFFER (3,J) to calculate the pressure gradient dp/dx with zero cross-flow.

The second call to SEPRAT comes after the channel exit pressure has been calculated.

The calculation has converged if the fractional difference between the minimum and maximum values of pressure drop is less than the criterion FERROR set as Input Data (nominal value = 0.001). If the criterion is satisfied, the calculation returns to SCHEME and thence to MAIN with JUMP=2 indicating convergence.

If not, JUMP is set to 1 and new inlet flows calculated before Returning.

The new inlet flows are calculated as follows. For the first iteration (a) the individual flow rates are summed to give FTOT, (b) v. for each channel is calculated from equation 4. In subsequent

]. iterations, v. is calculated from equation 5. Next, rv.p. and rv. ]. ]. ]. ]. are summed to give omi from equation 3 and thus the new mi. Finally, the values of F(I,l) i.e., m., are corrected for accumulated

]. rounding *off errors by summing and comparing with the original value FTOT. The array SP(MG,MX) is used for storing values between inter-vals; its PWR use is to store pressure differentials between channels to calculate the cross flow but for BWR's, this use of the array is not required.

SP(I,l) is used for vi and SP(I,2), SP(I,3) for storing the values of mi and pi from the previous iteration.

3.6 Fuel and Coolant Modal Average Powers 3.6.0 Introduction In the original COBRA IIIC, the total power is given via the Input Data as a function of time; the axial power profile is the same for each rod. A second option has now been added to give the I I I I I I I I I I I I I I I I -1 1. I I I I I I . I I I I I I I I I I I I I I fuel and/or coolant heat generation rates as a function of time for each nodal volume in the reactor; this implies that the axial power profile may be different for each rod. The following sections describes the modifications to COBRA IIIC and their testing . 3.6.1 Power Innut in the Original COBRA III In COBRA IIIC, the sequence of the thermal calculation is (a) Give power parameters via Input Data. (b) From SCHEME first call HEAT to calculate the heat flux (FLUX (N,J)), fuel temperature and heat addition per unit length (QPRIM(I))

where I -Channel, N =rod and J = axial station. (c) From SCHEME, next call DIFFER to calculate the enthalpy gradient DHDX (I) using QPRIM(I). (d) In SCHEME, calculate the enthalpy at J from the heat balance equation (card SCHM0900) giving void fraction, mean fluid density, etc. FLUX(N,J) is the important parameter that has to be calculated.

It is used to calculate QPRIM(I), (for DHDX(I) calculation), cooled boiling, CHF etc. In the original version of COBRA IIIC it is calculated by another route but thereafter the calculations follow the same course (except DHDX(I) includes the coolant heat generation where required).

In the old method FLUX(N,J) is calculated in subroutine

. as a function of the Input Data parameter average heat flux (AFLUX), rod power factor (RADIAL) axial orofile (AXIAL) and the time factor (POWER) (see HEAT 0520). If fuel temperatures are not required, FLUX(N,J) is used to calculate QPRIM(I) etc. However, if fuel temperatures are required, then FLUX(N,J) is first calculated as above, i.e. from the input data, then it is used to calculate heat generation rate in subroutine TEMP (TEMP 0370), in which fuel temperatures are calculated.

FLUX(N,J) is then overwritten and recalculated as a function of the temperature difference between wall tempearture and fluid temoerature (HEAT 0790). It is worth pointing out that in the input data the total power is given in the units of heat flux, but this is equivalent to either (a) Power generation rate within the fuel if temperature calculations are to be performed. (b) Average heat flux when temperature calculations are not to be performed.

3.6.2 Modifications The modifications to COBRA IIIC involve: (a) Set a trigger (IQP3) to indicate the calculation route and initialize QF and QC to zero. (b) Preparation of a new subroutine (QPR3) to read in the fuel (QF) and coolant (QC) heat generative rates. I I I I I I I I I I I I *1 I I I I I I I I I I I I I I I I I I I I I I I I I (c) Modify HEAT to set FLUX and QPRIM. (d) Modify DIFFER to include QC in the coolant heat balance equation.

3.6.2.1 IQP3 Trigger. The trigger IQP3 is initialized to 2 and the arrays of QF,QC to zero before reading in any Input Data. IQP3 may be reset in Card C8 through the value of Nl, i.e. the number of points along the channel at which the axial heat flux profile is specified.

If Nl > 1: If Nl < 1: IQP3 remains equal to 2. IQP3 is reset equal to Nl. The meaning of IQP3 in subsequent calculations is: IQP3 = O. Fuel heat generation (QF) given. Coolant deposition not given and remains set to zero. IQP3 = 1. Both fuel (QF) and coolant (QC) heat tion given. IQP3 = 2. Normal Input Data and calculations.

3.6.2.2 QPR3. This new subroutine (QPR3) is called if IQP3 = 0 or 1. This subroutine will read first ZM. If ZM is tive, the values of QF and QC (which will be given after ZM) are all multiplied by ZM to convert to Btu/sec. If ZM = -1.0, the multiplier ZM is reset to 3413.0/3.6 so that QF and QC may be read in MW and automatically converted to Btu/sec. If ZM = -2.0, it is rest to 1000.0/3.6, i.e., the conversion factor for MBtu/hr to MBtu/sec.

after ZM, the values of QF are read in and if IQP3 = 1, the corresponding values of QC are also read.in. The values are then printed. 3. 6. 2*. 3 HEAT. The modifications to HEAT are: (1) Skip the calculation of QAX from CURVE (HEAT 0480) if IQP3 = 0,1; QAX is used to calculate fLUX and this route is not (2) Skip the calculation of FLUX(N,J) from QAX (HEAT 0520) if IQP3 = 0,1; instead calculate FLUX (N,J) from QF (N,J). 3.6.2.4 DIFFER. Modify the calculation of DHDX (DIFF 0500) to include QC in the enthalpy rise equation.

It should be noted that DIFFER (l,JJ) is called from SCHEME to calculate the enthalpy rise with JJ = J-1, i.e. for values of the paz:ameters appropriate to the beginning of the axial interval being considered.

The values of QC (I,J) read in pond to the interval J-1 to J, hence in DIFFER it* is to use QC (I,JJ+l) in order to obtain the correct values since JJ+l = J. 3.7 Various Errors and Anomalies in COBRA IIIC 3.7.1 Negative Heat Flux in PROP In a particular run, a failure occurred in which a negative number was raised to a non-integer power. Thi.s was traced to subroutine PROP, card PROP0940, in which QPRIM(I) was found to be negative; QPRIM is the heat addition per unit length. I I I I I I I I I 'I I I I I I I I I. I I I I I The physical picture was that in a transient, the fluid temperature was specified, via the Input Data, to increase fairly rapidly. The fluid surface temperature increased with the heat I addition but less rapidly than the fluid temperature.

The net I I I " I I I I I I I effect was a negative heat flux, i.e. from the fluid to the fuel. The failure occurred in PROP when using the Jens and Lottes equation to determine whether subcooled boiling had started (this contains q"**0.25, where qn =heat flux). For negative heat flux, this test is clearly superfluous since subcooled boiling cannot occur if the surface is colder than the fluid. The cure was to bypass the statement when the heat flux is negative:

. --* --* Between Cards 2ROP0930 and PROP09&0, insert IF (QPRIM(I).LT.0.0)

GO TO 110. This modification leaves the physical model intact. It allows negative heat flux to occur (which it could in reality) but when it does, suppresses the test to see if subcooled boiling has started. 3.7.2 Array Out of Range in DIFFER At the start of the calculation, when finding the inlet flow split between channels, subroutine SPLIT calls DIFFER (3,J) with J = 1 (i.e. at inlet). In card DIFF0860, FLOWSQ is calculated I from F(I,JMl) where JMl = J-1, i.e. zero. Normally it does not matter what rubbish is calculated for FLOWSQ at inlet by this card I I I because by the next statement, if J = 1, it is recalculated. However, on rare occasions trying to fetch F(I,O) causes run failure. The cure is to complete what the original programmer probably intended, namely to bypass the normal statement when J = 1. Card DIFF0860 should be replaced by: 380 IF (JMl.GT.O)

FLOWSQ = ABS(F(I,JMl)*F(I,JMl))

3,7,3 Initialisation Quirk in COBRA IIIC 3,7,3.0 Introduction.

Two identical cases run secutively (in the same Run) using COBRA IIIC may give somewhat different answers. This is clearly undesirable and was found to be due to incomplete reinitialisation between cases. This quirk was discovered when checking that the ordering of the channel connections in the cr9ss-flow matrix had no effect upon the solution of the set of simultaneous equations.

The order was not significant, but, in demonstrating it, the quirk in COBRA IIIC at confused the issue. The effect of the lack of initialisation in COBRA IIIC is fairly trivial and simple to put right. It is reported here in the interests of complete documentation.

3.7.3.1 Initialisation.

Quantities are initialised in COBRA IIIC before reading in the Input Data for the first case. For subsequent cases, the Data are taken to be the same unless changed, hence a complete reinitialisation before each case. would I I I I I I I I I I I I I I I I I I. I I I I I I I I I I I I I I I I I I I I be undesirable.

Since most non-Input parameters are recalculated, the original COBRA IIIC contains no reinitialisation between cases. Two identical consecutive cases gave marginally different answers and this was traced to the lack of reinitialisation of the two arrays ((SP(K,J), W(K,J), J = 1, MX), K = 1, MG) = 0.0. Tests showed that if both arrays were reinitialised, the answers were identical but if either were not, the answers were different.

The array SP is the "guessed" differential pressure between channels *and is updated between successive hydraulic iterations.

The initialised values of zero represent the first guess hence if not initialised this represents an alternative first guess and is not necessarily 11 wrong 11* It is even arguable that if the second case is similar to the first, it might be better to use the results of the first as the initial guess for the second. The differences in answers for the two first-guesses were found to be within the convergence tolerance.

Thus it would probably not be worthwhile since W must be reinitialised, it is considered better to alise SP too in the interests of consistency.

The array Wis the cross-flow between channels.

This is a calculated quantity and it is not obvious why its initialisation should affect the answer. The only preset value is zero at inlet and this would not change during calculation.

Time did not allow a detailed follow-up of why W required i I reinitialisation.

The results were consistent when SP and W were initialised to zero between cases and this allowed the more important work to continue.

3.7.3.2 Modification to COBRA IIIC. (a) Remove the following cards from their original position DO 915 K = 1, MG

W(K,J) = 0. WOLD (K,J) = O. 915 SP(K,J) = 0. (b) Replace these cards, preceded by DO 915 J = 1, MX MAIN 5860 MAIN 5870 MAIN 5880 MAIN 5890 between cards MAIN 6070 and MAIN 6080 (now in INDAT) i.e. immediately before CALL DOY(DATE).

It may be noted that in the interests of coding simplicity, WOLD is unnecessarily reinitialised.

3.7.4 Anomalous Behavoir of COBRA IIIC Using the Subcooled Void Option when Bulk Boiling at Channel Inlet 3.7.4.0 Introduction.

A COBRA IIIC case with bulk boiling at inlet was run using a deck of cards from a previous case which called for the Levy subcooled void option. This option should clearly have no effect for bulk boiling at inlet. In fact, the void fraction calculation was anomalous.

At low qualities, the calculated values were correct, then they decreased with increasing quality, reached a minimum and then rose again. Correct values were obtained by not invoking the subcooled void option. I I I I I I I I I' I I I I I I I I I I I I I I I\ I* I I I I *1 I I I I I I The reason for this behavior is given below. No modifications were considered necessary to COBRA IIIC since it is only necessary to omit the subcooled void option when a case with bulk boiling at inlet is run. 3.7.4.1 Explanation.

In subroutine PROP, the single phase heat transfer coefficient HFILM is calculated by card PROP0900_

only when the liquid enthalpy is below the saturation value. At the onset of bulk boiling, it is not recalculated but remains at the last single phase value. However, for bulk boiling at inlet, it is never calculated hence takes the default value of 1.0E78. The heat transfer coefficient is used in subroutine SCQUAL to calculate DELTAT (Card SCQ0370) and thence the "quality" XD at which subcooled void starts. Because HFILM is taken to be large, XD becomes positive.

When the equilibrium quality Xe is less than XD, the void fraction is calculated from Xe. When it is larger than Xd, the void fraction is calculated from (Xe-Xd). This gives a drop in void fraction which then increases again with increasing Xe. Thus with bulk boiling at inlet and HFILM not set, the void fraction is calculated correctly up to some quality XD and then drops to a lower value. If HFILM had been set to some nominal value, XD would have been negative and the effect of subcooled void small but finite. Because the effect should be zero, it is better not to use the subcooled void option. 3.7.5 COBRA IIIC Modification to Prevent Occasional Underflow In Subroutine SCQUAL The Levy subcooled void correlation is programmed in COBRA IIIC in subroutine SCQUAL. Subcooled void starts when the equilibrium quality (Xe) reaches a calculated negative value Xd. As Xe increases further, the subcooled non-equilibrium quality (X) in9reases from zero at Xd and tends exponentially to Xe according to the equation For Xd negative but close to zero and Xe large, the argument of the exponential can become large and negative when the computer will give an "Underflow" error signal. Under these conditions, the value of the exponential is, for practical purposes, zero. Hence the cure is to set it to zero when the argument of the exponential function is less than some number (say) -15.0. This makes x = Xe. The modification to SCQUAL is After the calculation of ARG and between Cards SCQL0450 and SCQL0460, insert IF (ARG.LT.-15.0)

GO TO 140. 3.8 Input Data for IIIC/MIT 3.8.0 Introduction In the present version of the code two general options are I I I I I 1* I *I I I I I I I I I I 1. I I I I I I I I I I I I I I I I I I I I available to input the data. The first option is based on the original version of COBRA IIIC, and has two different lities, one for the analysis of a small number of subchannels and the other for the analysis of the whole core of PWR's and BWR 1 s, each radial node being an assembly.

The first possibility uses exactly the same Input Data Presentation as COBRA IIIC, but for the second, card groups 4, 7 and 8 were lumped together to generate a new card group 4. This modification is possible because in such analysis many channels have identical characteristics and only the description of one channel is needed in order to describe all the channels identical to this one. Also in this analysis, every channel is assumed to have a rod that generates the power introduced by the fission reaction in the channel. This implies that the description of the rods and its power may be given in a more simplified way. Other improvements are made here which consist of the tion of the cards describing the physical properties of the coolant that are calculated in COBRA IIIC/MIT by the code and then are not required to be inputed. The second option of given the Input Data is very useful in the study of PWR's and BWR's, either in an assembly to assembly basis or a more detailed description of the core that allows for the representation of subchannels together with lumped regions and will make possible the study of the whole core in only one run. The arrangement of the Input Data is totally different in this option which also has some other capabilities that the previous option does not. They are, for instance, the printing of the image of each card immediately after it is read, the use of the Baroczy theory as a new option for the new-phase friction model, the possibility of introducing coupling parameter in the mixing term of the energy equation and some other characteristics that provide a large range of flexibility in the division in channels of the whole core. With respect to the reduction in the number of cards required to describe the Input Data, this second option is even better than the second possibility of the first option. The conclusion is then that the Input Data selected will depend on the case that is being analyzed.

If this is a case with a small number of subchannels the original version of COBRA IIIC may be selected.

If the case is the analysis of the core in an assembly to assembly basis either the second possibility of the first option or the second option can be chosen but we recommend the second option because of its wider range of alternatives.

And finally for a more detailed analysis of the core, the second option should be selected.

3.8.l Input Data Presentation Based on the COBRA IIIC Presentation 3.8.l.l COBRA IIIC Inout Data. The code keeps the option of inputing the data as in COBRA IIIC. This option will I I' I I I I I I I I I II I I I *I* -1 1. I I I I I I I I I I I I I I I I I I I I be helpful in analysis of a few subchannels where the fluid is not water or where some special characteristic, which in the other options were deleted, are needed.* The cards required are described in Appendix 10. 3.8.1.2 Simplified COBRA IIIC Input Data to be Used for Assembly to Assembly Analysis of LWR. As the title suggests, this option is based on COBRA IIIC but has some simplifications when analyzing that PWR and BWR in an assembly basis are possible and very convenient because of the huge reduction in cards required to describe the problem. It is known that reactor channels (assemblies) tend to have similar geometries and differ only in their powers and heat flux profiles.

Each rod is associated with one channel only (in a non assembly to assembly case, a rod may be bordered by several subchannels);

thus rod data may be given as part of the channel data. Also, one may remove the option for reading in wire-wrap spacer data. The Input Data presentation may be simplified by combining the data given in Card Groups 4, 7 and 8 of the COBRA IIIC Input Data presentation as follows: (a) define several (up to 15) channel types in terms of their geometry, rod dimensions and grids, (b) specify which channels are of which types, (c) list the radial power factors of each channel, (d) specify the interconnections between channels. The physical property data presentation may also be simplified.

In the original COBRA IIIC, the properties at various pressures are read from cards. In the. new version, polynomial expressions are used to generate, within the program, the data which were originally read from cards. The Input Data changes then with respect to the original version of COBRA IIIC and can be found in Appendix 11. 3.8.2 New Input Data Presentation This last option of giving the Input Data was developed for MEKIN, and then was further extended, because.of its simplicity, to be used in COBRA IIIC/MIT.

It allows the description of problems such as analysis of-LWR in an assembly to assembly basis or in a much more detailed pattern of channels (i.e., combining lumped channels outside of the hot assembly with a fine mesh of subchannel inside of the hot assembly, that will allow the analysis of the whole core in only one run. It has the same capabilities as the Input presentation described in Appendix 3.8.2, except for area and gap variations that were deleted because of their small use in LWR, plus some others such as the printing of the images of each card immediately after it is read, the possibility of giving the power per node as input data, the use of coupling coefficients in the mixing term of the energy equation and a new model (i.e., the Baroczy model) to calculate the two-phase friction factor. For I I I. I l-1* I I I I. I I I I I I I I I I I I I I I I I ' I I I I I I I I I I these reasons we strongly recommend the use of this option in the cases described above. This option is described in Appendix 12. CHAPTER 4 DESCRIPTION OF THE ORGANIZATION OF THE CODE AND OF ITS NEW SUBROUTINES 4.0 Introduction The organization of COBRA IIIC consists of a very long MAIN routine which combines several functions, namely (a) Read Input Data, (b) Print Input Data, (c} Initialize Variables, (d) Control Hydraulic Calculation, and (e) Print Results. The organization of COBRA IIIC/MIT is more.complex because MAIN was split in order to distribute the functions listed above among several subroutines.

Also, due to the new Input Data tion and the new capabilities incorporated to the code, some new subroutines were required.

The objective of the present chapter is to describe the organization of the code in order to facilitate an understanding of the program in cases where it is desired to incorporate some new improvement.

Also, in the last part of the Chapter the description of the task that each new subroutine executes is given. 4.1 Organization of the COBRA IIIC/MIT Code The new MAIN is very short, only 7 statements, and its only task is to guide the calculations.

It successively calls INDAT (i.e. to read in Input Data) INPRIN (i.e. to print Input Data) I I I I I *1 I I ' I I I I I I I I I I I I I I I I I I I .I I I I I I I I I I and CALC (i.e. to start the hydraulic calculations).

From CALC itself a call is made to EXPRIN which prints out the results and returns the control of the code to MAIN to start a new cycle in cases where more than one case is analyzed.

In Figures 9, 10 and 11, the order in which the successive subroutines are being called is given. These figures will be of great help in the understanding of the following description.

4.1.1 Read in the Input Data As said in Section 3.8, there are two options to Input the Data, one when NGROUP is equal to 20 (new Input Data tion) and the other when NGROUP takes successive values of 1 through 12 (original COBRA IIIC Input Data Presentation).

The last one uses subroutine CARDSl to read in card group number 1 (i.e. physical properties of the coolant) while the remaining of the Input Data is read from subroutine INDAT except when IPILE = 1 or 2 (analysis of PWRs or BWRs in an assembly to assembly basis) at which time card group number 4 (i.e., tion of channels and rods) is read in from subroutine CARDS4. Subroutine INDAT initializes also all the variables that require such initialization.

From Figure 11 it can be seen that some other subroutines are called from INDAT in order to fix the dimensions of the arrays and to determine the array LOCA (which is done in time ACOL). When card group 20 is selected (Figure 10) the new Input Data Presentation has to be used. In this case, subroutine INDAT initializes the variables, reads in control cards Cl, C2, C3 and C4 and makes a call to subroutine CORE to establish the dimension of the arrays. Then the control is transferred to subroutine CARD 20 which reads in cards defining the mapping of channels and the average heat flux (i.e. cards C5, C6, C7 and C8). Subroutine READIN is from CARD 20 to read in the axial and radial peaking factors (i.e. cards C9 and ClO) and CARD 20 itself reads in the number of axial steps and time steps (i.e. Card Cll). Then the control is transferred to subroutine ITHO which calls successively to subroutines CHAN, MODEL, OPERA, FIZPRP and TABLES which read in respectively the description of the channels and the rods (i.e. Cards from Tl to T7), the lic model (i.e. cards from TB to Tl9), the steady state and transient conditions, inlet enthalpy, mass flow at the inlet and outlet pressure (i.e. cards from T20 to T25a). Now subroutine FIZPRP will calculate the physical properties of the coolant making successive calls to HAPROP, HLIQ, HVAP, ROLIQ, ROVAP, SATTEM, and SURTEN. Finally, subroutine TABLES identifies the channels, rods and nodes for which final results are to be printed (i.e. cards from T26 to T30). The control returns again to ITHO which calls to CORE 3 to check whether or not enough storage was given to solve the problem and prints the size of those variables whose size is problem dependent.

The control is then transferred again I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I to OPERA, CHAN, MODEL and FIZPRP successively to print out the Input Data and the physical properties of the coolant. Once the reading of the Input Data is ended, a call is made to PRECAL to establish the initial conditions.

4.1.2 Print the Input Data When the reading of the Input Data is ended, control is transferred to INPRIN via MAIN to print out, if desired, the Input Data. The printing is selected with variable Nl (card C4). If card group 20 is chosen, the printing of the Input Data described here can be skipped since a similar print out was already obtained.

4.1.3 Thermal/Hydraulic Calculations The scheme used in COBRA IIIC/MIT is the same as that used in COBRA IIIC. First, CALC (part of the old MAIN) establishes the boundary conditions for the steady state problem and starts an iteration loop which sweeps the calculation through the core. Subroutine SCHEME is called in each iteration to perform the calculation for every axial step. Mass flow convergence is checked at the end of each axial step and when all values are inside the convergence factor the calculations are ended. wise a new iteration is performed.

For the transient calculation, the steps described above are performed for every time step until the end of the transient is reached. 4.1.4 Print Out of the Results Once the calculations are ended, subroutine EXPRIN is called from CALC, which will print out the results. Several options are available in the INPUT DATA to allow for the tion of results which are desired to print. Once the printing is over, the control of the program is transferred to MAIN which calls to INDAT to start a new case if there is one. wise the program stops. 4.2 Description of Subroutines Several of the subroutines used in COBRA IIIC/MIT are new, some as consequences of the new features incorporated into the program, and others because of the different organization of the code. The rest of the subroutines remain almost unchanged with only small modifications to allow for the introduction of the new features.

Here the description of the new subroutines is given. 4.2.1 MAIN The new MAIN guides the calculation:

INIT = 1 2 CALL INDAT (INIT, NOPRIN) IF (NOPRIN. EQ.O) CALL INPRIN CALL CALC (which calls EXPRIN) INIT = 2 GO TO 2 END I I I I I I I I I I I I I I I I I I I I I' I I I I I I (' I I I I I I I I I I First a call to INDAT is made that will read in the Input Data and initiate the variables.

After this step is executed, a call to INPRIN is made to print out the Input Data only if NOPRIN in Card C4 is set equal to zero. Then CALC is called to do the hydraulic calculations and CALC itself will call EXPRIN to print out the results. The control returns then to MAIN which after making INIT = 2, will call again to INDAT to read in a new case if there is one. Otherwise INDAT will stop the program. 4.2.2 INDAT If card group 20 in card C4 is not called, either the nal or the simplified Input Data Presentation (which combines card groups 4, 7 and 8 in only one) are used. This subroutine then, as explained in 4.1.1,"will read in the Input Data. But if card group 20 is called, INDAT will only read in the four first cards and then will call subroutines which pript as well as read Input Data. In ari.'y case, INDAT initiates the variables.

INDAT has cards from MAIN 5365 to MAIN 8830 of the old MAIN. 4.2.3 INPRIN As indicated before, this subroutine is only called when NOPRIN is set to zero, and then will print out the Input Data in the same fashion that COBRA IIIC does. INPRIN has cards from MADI 8840 -to MAIN 0350 of the old lVl..AIN. 4.2.4 CALC This subroutine performs the hydraulic calculations.

CALC has the cards that in the old MAIN were named MAIN 0360 through MAIN 1820 and from MAIN 2340 through MAIN 2405. 4.2.5 EXPRIN With this subroutine, the output is printed and has cards from MAIN 1822 through MAIN 2331 of the old MAIN. 4.2.6 CARD 20 This subroutine is called from INDAT where card group 20 is selected.

a) b) c) d) Its purpose is to read and set the following data: The channel map, i.e., the array NTHBOX (I,J). Channel powers, i.e., the average heat flux AFLUX: the axial heat flux profile NAX, Y( ), AXIAL( ); the radial power DATA (ZRADIA+I).

Channel length (z) and number of axial intervals (NDX). Total time of transient (TTIME) and number of time intervals (NDT). This information is passed to ITHO via the argument list and COMMON. 4.2.7 ITHO This subroutine is called from Card 20. Its purpose is to control the Input Data reading and printing.

It calls CHAN, MODEL, OPERA, FIZPRP, TABLES to read in Input Data. It theh calls OPERA, CHAN, MODEL, FIZPRP again to print out the Input. I I I I I I I I I I' I I I I I I I I I I I I I I I I I I I I I I I I I I I I 4.2.8 CHAN This subroutine is similar to CARDS4. It reads the channels and rod parameters (e.g., area, grids, rod diameters, etc.) for various types of channels and specifies which channels are of which type. 4.2:9 MODEL This subroutine sets the hydraulic model. A preset model is obtained by submitting a blank card representing a number of indicators all set to ZERO. Alternatively, various parts of the model may be changed by setting individual indicators of ZERO. For example, the mixing might be preset to S =

if Nl = 0 but changed by setting Nl = 1 and giving an extra card with a new value of S. Correlations for two-phase friction, subcooled void, slip ratio may be preset but changed if required.

4.2.10 OPERA subroutine sets the operating conditions, for example mass velocity, transient forcing functions, etc. 4,2.11 FIZPRP This subroutine sets the physical properties.

By giving a blank card, values would be calculated automatically to span the problem range. 4.2.12 TABLES This subroutine reads in and sets the variables which control the channels, rods, and nodes that will be printed on EXPRIN. 4.2.13 READIN This is an auxiliary subroutine for reading in certain quantities within a DO Loop. 4.2.14 TIDY Its purpose is to set variables which are used when printing 0 but not during the hydraulic calculations.

It is only entered if NOPRIN = O, i.e., if the original COBRA IIIC printing (in INPRIN) is required in addition to that controlled by ITHO. 4.2.15 PRECAL Its purpose is to set certain variables which would otherwise be set in INPRIN in preparation for the hydraulic calculation (e.g., the flow split at inlet). PRECAL written as a separate subroutine in order to keep intact the other Input Data route via INDAT and INPRIN. 4.2.16 CARDS 1 Its purpose is to read from cards or calculate from polynomials the elements of the physical property arrays. It is called from INDAT when NGROUP = 1 and it calls to HLIQ, HVAP, ROLIQ, ROVAP, SATTEM, HAPROP and SURTEN successively when Nl < O to calculate the physical properties from polynomials.

I I I I I I I I I I I I I I I I I I. I I I I I I I I I I I I I I I I I I I' I 4.2.17 CARDS 4 Its purpose is to read and process channel, grid and rod data when Card Group 20 is not chosen and IPILE = 1 or 2, in other words when the simplified version of the original Input Data, which combines card groups 4, 7, and 8 in only one, is used. It is called from INDAT when NGROUP = 4 and IPILE = 1 or 2 and it calls to ACOL to set array LOCA for use in setting up cross-flow matrix. 4.2.18 CORE This subroutine is described in full detail in 3.2.4. It computes the origins within the DATA ARRAY of the old, fixed arrays A, NTYPE, etc. 4.2.19 SEPRAT This subroutine is described in 3,5,3.4. It contains the bulk of the modifications required for BWR assemlby to assemlby calculations.

4.2.20 QPR 3 Its purpose is to read in QF and QC, fuel and coolant heat generation rates for each model volume. For more details see 3. 6. 4.2.21 ACOL This subroutine sets the LOCA ARRAY. It is called either from INDAT, CARDS4, or CHAN depending on the Input Data tion being used. This subroutine is described in full detail in 3.1.3. 4.2.22 BAROC This subroutine is described in 3.4.2. It calculates the two-phase friction pressure drop following the Baroczy correla-tion. 4.2.23 SURTEN This subroutine calculates the Surface Tension (lb/ft). 4.2.24 HAPROP This subroutine calculates Liquid Specific Heat (Btu/lg °F), Liquid Viscosity (lb/ft hr) and Liquid Conductivity (Btu/ft 2 hr °F). 4.2.25 Functions ROLIQ, ROVAP, HLIQ, HVAP and SATTEM Their function is to calculate Liquid Density (lb/ft3), Vapor Density (lb/ft3), Liquid Enthalpy (Btu/lb), Vapor Enthalpy (Btu/lb) and Saturation Temperature

(°F), respectively.

The subroutines described above are those which are new in COBRA IIIC/MIT with respect to COBRA IIIC. The rest of the subroutines used in the MIT version of COBRA are basically unchanged with respect to COBRA IIIC and only small modifica-tions have been incorporated in some to reflect the introduction of new subroutines in the overall code. These modifications were described in Chapter 3. I I I I I I I I ., I I I I I I I I I I I I I I I I I I I I I I I I I I I I I 5.0 Introduction CHAPTER 5 PRESENTATION OF RESULTS It is the intention of this chapter to show two different points, first that both codes, COBRA IIIC and COBRA IIIC/MIT, yield identical results when applied to the same problem and second to present the new capabilities of the code with respect to Input Data Presentation and running time. To prove the first point we had to select a case which fit within the tions of COBRA IIIC. A case with only 10 cahnnels and 13 rods was chosen. To show as much as possible of the second point, the 10 channel by 13 rod case was arranged into a rather rare pattern of channels and rods (Figure 12). The representation of Figure 12 can be as a reactor core in which channel 7 is an actual subchannel while the rest of the channels sent lumping of subchannels.

Three different runs were made of this sample problem, one with COBRA IIIC, one with COBRA IIIC/MIT using the same Input Data Presentation as before, and the final one with COBRA IIIC/MIT using the new Input Data Presentation.

The results and conclusions are presented below. 5.1 COBRA IIIC In Appendix 13 the listing of the cards required to be inputed in order to describe the problem, and the results obtained, are given. The case, like the two following, were run on the IBM 370/168 computer at MIT and in Table 11 their cost is given. 5.2 COBRA IIIC/MIT With the Old Inout Data Presentation As can be seen by comparing Appendix 14 with 13, the results obtained are identical, except for some values of the D.N.B. Ratio, to those given by COBRA IIIC. The slight difference in the D.N.B. R value is due to the fact that in COBRA IIIC/MIT the effect of the spacers upon taken into account in the W-3 correlation by taking TDC=0.019 as recommended in reference 18 as in COBRA IIIC the parameter TDC is taken as zero. 5.3 COBRA IIIC/MIT With the New Inout Data Presentation From the listing of the cards required to describe the problem, (Appendix 15), it is clear that the reduction in Input Cards is quite substantial when this option of inputing the data is selected.

Also, from Appendix 15 it can.be ed that the results are almost identical to those of the vious cases, the slight differences being due to variations in the physical properties of the coolant between their calculation from polymonial expressions, in this case versus cases 5.1 and 5.2 where they are inputed. I I I I I I ___ I I I I I I I I I I I I I I I I CHAPTER 6 I LISTING OF THE COBRA IIIC/MIT CODE I The listing of the code is given in Appendix 18 in the following order: *I 1. BLOCK DATA 17. GAUSS 33. ROVAP 49. PHNTIM I 2 . AREA 18. HA PROP 34. s 50. CHF 3 . BAROC 19. HCOOL 35. SATTEM 51. CHF2 I 4. BVOID 20. HEAT 36. SCHEME 52. DIFFER 5. CALC 21. HLIQ 37. SCQUAL 53. DIVERT I 6. CARDSl 22. HVAP 38. SE PRAT 54. IND AT I 7 . CARD20 23. INPRIN 39. SOLVE 55. MODEL 8. CHFl 24. ITHO 40. SPLIT 56. TABLES I 9. CIJ 25. MAIN 41. SUR TEN 57. ZIGET 10. CURVE 26. MIX 42. TEMP 58. ZIFREE I 11. DE COMP 27. OPERA 43. TIDY I 12. DOY 28. PRECAL 44. VOID 13. ELAP 29. PROP 45. ACOL I 14. EXPRIN 30. QRP3 46. CARDS4 15. FIZPRP 31. READ IN 47. CHAN I 16. FORCE 32. RO LIQ 48. CORE I I I I References

1. "MEKIN: MIT-EPRI Nuclear Reactor Core Kinetic Code", CCM-1, R. Bowring; J. Stewart; R. Shober; R. Sims. 2. D. S. Rowe. "COBRA IIIC: A digital computer program for steady state and transient thermal hydraulic analysis of rod bundle nuclear fuel elements." BNWL-1695 (1973). 3. K. Hansen. Personal Communication.

4. R. W. Bowring. "Simplified Input Data for COBRA IIIC when used for PWR and BWR Cases". MEK-20 (1974). 5. J. W. Stewart. "MEKIN Data Management".

MEK-9 (August 1974). 6. R. W. Bowring. "An Introduction to COBRA IIIC Coding". (November 197 4) . 7. R. W. Bowring. "COBRA IIIC Modifications to Subdivide MAIN". MEK-36 (April 1975). 8. R. W. Bowring. "The coding for a Faster COBRA IIIC Cross-Flow Soluti.on".

MEK-28 (January 1975). 9. R. W. Bowring. "Timing Runs for the Original and Modified COBRA IIIC". MEK-31 (January 1975). 10. S. L. Smith, "Void Fractions in two-phase flow; a correlation based upon an equal velocity head model", Proc. I.M.E. Vol. 1, Pt. 1, No. 38, page 647, 1969-70. 11. 12. 13. 14. C. J. Baroczy. rr A sy_stematic

correlation for two-phase pressure drop" .. AICHE-ASME H.T. Conf., Los Angeles, Aug. 1965, Preprint 37. J. H. Keenan, F. G. Keyes, P. G. Hill and Joan G. Moore. Tables". pub. John Wiley and Sons (1969). Electrical Research Association.

St. Martins Press, N.Y. (1968). "1967 Steam Charts". "Steam Pub. I I I I I I I I I I I I I I I I J. R. S. Thom, W. M. Walker, T. A. Fallon and G. F. Reising. "Boiling in subcooled water during flow up heated tubes or I anriuli".

Proc. Inst. Mech. Meg., Vol. 180 IIIC, p. 225 (1965-6).

15. W. H. Jens and P. A. Lattes, ANL4627 (1951). I

I I I I I I I I I I I I I I I I I I I 16. J. C. Chen. "A Correlation for boiling heat transfer to saturated fluids in convective flow 11* ASME Paper 63-HT-34 (1963). 17. F. W. Dittus and L. M. K. Boelter, Univ. California Publs. Eng. 44, (1930). 18. J. Weisman and R. Bowring. "Methods for Detailed Thermal and Hydraulic Analysis of Water-Cooled Reactors".

Nuclear Science and Engineering 57, 255-276 (1975).

I I I I 1 1 2 2 3 3 4 I 4 6 7 I 5 8 6 7 10 8 I I 11 2 13 11..I I 9 1 10 1 11 1 12 I 8 19 20 I I 13 2 14 23 15 21..i 16 I I I Figure 1. Channel and Boundary Numbering Scheme I I I. I I I I I I I I I I I I I I I I I I I I ..... l 1 2 ; 3 J '+"""' 4 / ., 5 9 r J i o,_ c.9 3 6 1 }1 8 10 "'12 l3J 9 '-l I+ l 5..r* "'l 6 11 ,..., l ., 12 13 .-1 8 Figure 2. Channel and Boundary Numbering Scheme NC = 10 D c c G B F A E A E N F c G D c Figure 3. / / A B c /D , E F / C?f c / H J: F B / 1 ... 'J" H E A NC = 16 B A A B c D D c B A F E E F G c c G F E I H H I F B B F I H A B c E F / c{ / H /I F / 4 H E J H E H I F E F G A B c NG = 32 A B c E F G H I F y' c B A A B c D D c B / I I I I I I I I I I. H J J H E A A E H J J H E A I I I H J J H E A A E H J J H E A I F B H H I F B B F I E E F G c c G F A A B c D D c B ._ ----NC = 128 H H E E A A ..... NC = 64 I F F G B c B c D I I I I Channel Layout for 10 to 128-channel cases. All channels were identical except for their radial power factors which were arranged symmetrically as indicated by letters A-J. I. I I I I I I I I I I I I I I I I ,. I I _. IlE:\T THA:--:SFEH.-LOS A:--!GELES E.:\GI:\EE!n\G l'IWGIU:ss

****-/MJO +---"'<--'<

Q) --<.---------4----------+-----------

1 !CC.0 ...... ..-l 0... ...... µ ..-l ;:::l c:: 0 *r-i µ u *r-i 11. SOOIUM ("fl 1200 ,oL I 1400 POTASSIUM

("Fl 1000 I 0.1 1600 1800 200J 1200 1600 1eoo 2CCO AUfllOiUM FREON-22 (°Fl 40 MERCURY ("Fl

__

WATER ("F)212 328 .. 67 j i J so ea iOO 140 _ Figure 4. Baroczy Two-Phase

?riction l1ul tip lier at G = 1.0 Mlb/£t 2 hr.

Figure 5. 16 1$ , .. 13 . 2: LZ ... .. . cl * *o-l I.I OS 09 0.6

PROPERTY INDEX, I Baroczy Mass Velocity Correction to Friction Multiplier I I I I .I I I I I I I I I I I I I I I I I I I I I """"" H c... c... ..._, >. I +.) ""' Cl) 0.. 0 ""' I 0... I I I I I I I Figure 6. I I I (!).. >-98-s 4 Polynomial I

'3 . ! Experimental

/rJ'tJ Jao .Soa /cr;a Pressure (psia) Physical Property Index (PPI) vs. Pressure PPI = (µL/µg)0.2/(pL/pg) i... Cl) *.-! *.-! +.l ;j i:: 0 *.-! +J u *.-! i... Cl) Ul -..c:: 0.. I 0 ;3:" E-< Figure 7. ic oou.*__..-------,----

1 Quality % ______ 3 ! L:_ I ' ' ---

' I I i I I I I I -;C c'I I

- --* -****--*-***---*

-

i i i ! I I I I I I l' I ! I I I a.0001 o.ool a. a I Baroczy Two-Phase Friction Multiplier at G = 1.0 Mlb/ft 2 hr. o.; CODE VALUES (Plotted by computer from subroutine BAROC) I I I I 1-1 I I I I I I I I I I I I I I I Figure 8. -100-I '-/. 'S /. 3 i. "Z. /.I /, 0 c.q 0 ....-; (J,7 ."J ...,. (I., CJ.5 a.;, o.ccc. /. 8 /.1 'J J. i, ;5 /,J., "' /.3 ;.z " """ /.I

'{ (t.) \ \ /, 0 -ac: /J. <1 .( l o-7 "* o-5 o.coal a . .:01 o.cr 0.1 Baroczy Mass Velocity Correction to Two-Phase Friction Multipliers CODE VALUES (Plotted by computer from subroutine BAROC)

MAIN -101-INDAT -(Read in the INPUT DATA, see Figures 10 and 11) INPRIN CALC-QPR3 PROP BAROC FORCE AREA PROP-HCOOL -E VOID VOID CQUAL BAROC HEAT -TEMP -GAUSS MIX Only for the first axial station SCHEME ---...+-DIFFER (Enthalpy Increase)

FORCE AREA PROP-HCOOL t BVOID VOID-SCQUAL_ BAROC DIFFER (Pressure Drop without Crossflows)

SEPRAT 1soLVE DIFFER (Mass Flow Increase) -DIFFER (Pressure Drop with Crossflows)

SE PRAT EXPRIN -CHF -CHFl or CHF2 Figure 9. Organization of COBRA IIIC/MIT I I I I I I I I I I. I I I I I I I I. I I I I I I I I I I I I I I I I I I I I INDAT --102-CORE -ZIGET CHAN ACOL DOY MODEL HAPROP TOD READ IN OPERA HLIQ FIZPRP HVAP CARD20 ITHO TABLES ROLIQ CORE3 ROVAP OPERA SATTEM CHAN SURTEN TIDY MODEL [HAPROP FIZPRP SUR TEN --f PROP PRECAL CURVE SPLIT Figure 10. Organization of the Reading of the Data When Card Group 20 is Selected CORE2 IND AT -103-CORE -ZIGET DOY ---+-TOD CARDSl ---;:;,. -[HAPROP

CARDS4 -ACOL -CORE2 (if IPILE 1,2) CORE3 ACOL -CORE2 CORE3 PROP CURVE SPLIT Figure 11. Organi-zation of the Reading of the Data When Card Group 20 is Not Selected I I I I I I I I I I I I I I I I I I I I I -104-I I I I s I 5 J I 8 I "8" I l 0 I 10 I I Channel: 5 I Rod: 5 I Figure 12. Arrangement of Channels anq Rods in the Example Problem I I I I

-105-COST($) CASE SIZE STORAGE 101 101 30 30 10 10 READING CALCULATIONS PRINTING Channels 500K . 25 22.09 Rods Channels 275K .16 6.51 Rods Channels 210K .10 1. 08 Rods Comparison of Storage Requirements and Costs With Problem Size Using COBRA I'IIC/MIT on an IBM 370/168 TABLE 1 8.53 3,34 1.83 I I I I I I I I I I I I 11 I I I I I I I I. I I I I I I I I I I I I I I I I I I I I K 1 2 3 4 5 6 7 8 9 l 4 2 -1 3 -2 4 1 5 -1 6 -2 '7 -3 I ,. 8 -4 9 106-LOCA (K,L), L = 1, 8 5 0 0 0 4 5 . 6 0 0 5 6 -7 0 0 0 4 11 0 0 0 4 2 8 12 0 6 3 9 13 0 6 10 -14 0 0 0 4 11 5 12 0 6 -8 12 6 13 7 etc. LOCA Array for Case of Figure 1 ; TABLE 2 K 1 1 2 -3 2 2 1 4 3 3 -1 5 4 4 14 2 5 5 4 6 6 4 7 7 4 -107-LOCA(K,L), L = 1, 14 0 0 0 0 0 0 0 0 0 6 9 -13 15 0 0 0 -16 0 0 0 0 0 0 0 0 6 9 15 0 0 0 6 7 9 13 15 3 -16 0 0 5 7 9 13 15 -8 0 0 0 5 +6 9 13 15 -10 0 0 0 LOCA Array for Case of Figure 2 TABLE 3 0 3 0 9 0 4 0 9 0 10 0 9 0 9 I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I -108-NC 10 16 32 193. NK 12 24 52 356 Steady State (Iteration)

Seconds Cross-Flow 0.31 1. 73 15.02 -Rest . 0.21 0.39 1.21 -Transient (7 Iterations)

Cross-Flow 2.31 12.64 105.42 -' 2.48 7.84 Rest 1. 01 -Calculated (per Iteration)

Cross-Flow 0.33 1.80 15.07 3274 Rest 0.15 0.35 1.16 13 Notes: Averaged and Calculated Times are for 10 axial intervals and per hydraulic iteration.

All results with H Compiler.

NC=l93 represents extrapolation to PWR but values should be doubled for 20 axial intervals.

Calculated values from Cross-Flow:

t = 0.08 + 0.000235(NK) 2*8 Rest: t = 0.0027(NC)l.75 , Timing for Original COBRA IIIC (10 Axial Intervals)

TABLE 4

-109-NC 10 16 32 64 128 193 NK 12 24 52 112 232 356 MS 11 15 15 31 63 59 Stead;y State ( 1 Iteration)

Seconds Cross-Flow 0.06 0.10 0.27 1.37 9.16 -Rest 0.11 0.21 0.41 0.77 1.59 -l Transient (7 Iterations) i n. Cro.ss-Flow 0.32 0.81 1.95 9.65 63.09 -l Rest 0.71 1.12 2.51 4.19 8.35 -Average(per Iteration)

Cross-Flow 0.05 0.11 0.28 1.38 9.03 -Rest 0.10 0.17 0.37 0.62 1.24 -Calculated(per Iteration)

Cross-Flow 0.05 0.12 0.28 1.38 9.05 12.99 Rest 0.10 0.16 0.32 o.64 1.28 1.93 Notes: Average and Calculated Times are for 10 axial intervals and. per hydraulic iteration.

All results with H compiler.

NC=l93 represents extrapolation to PWR but values should be doubled for 20 axial intervals.

Calculated values from: Cross-Flow:

Rest: t = 0.00282NK

+ 0.00000837NK*MS 2 + 0.0000127NK 2 t = O.OlNC Timing With New Cross-Flow and New "Rest" Subroutines (10 Axial Intervals)

TABLE 5 I I I I I I I I I I. ,-1 I I I I I I I I I I I I I I I I I I I I I I I I I I -110-No. Channels 16 32 64 128 Previous Values Cross-Flow 0.11 0.28 1. 38 9.03 Rest 0.17 0.37 0.62 1.24 Total 0.28 0.65 2.00 10.27 New Values Cross-Flow 0.13 0.32 1. 51 9.27 Rest 0.19 0.36 0.71 1.39 Total 0.32 0.68 2.22 10.66 Increase in Running Time 14% 5% 11% 4% Notes: "Previous Values" are taken from the Average Values in Table 5. The "New Values" are directly comparable and are given for (a) the solution of the cross-flow equations and (b) the rest of the Hydraulic calculation.

Running Times per Iteration

-in seconds TABLE 5A Property index .

0.1 0.5 1 . 0.0001 2.20 5.80 9.20 0.001 2.15 5.60 8.80 0.0011 2.08 4.90 7.80 ' 0.01 1.59 3.30 11. 80 0.03 1.12 1.55 1.81 0.1 1.04 1.12 L22 0.3 1.01 1.02 1.06 Two-phase multiplier

[For G = lx10 6 lb/(hr.)(sq.ft.)]

0 Quality, (%) 2 3.5 5 7.5 10 15 16.0 26.5 47.0 99.0 163. 376. 14.8 22.8 34.2 48.2 70.0 108 11.9 16.3 22.8 29.0 36.0 49.5 7.00 9.60 12.4 16.o 20.0 27.0 2.57 3.115 4.7 6.10 7.90 11.0 1.48 1. 78 2.05 2.50 2.80 3.60 1.13 1.26 1.36 1.50 1.59 1. 77 20 30 630 1,300 1118 240 63.0 86.o 33.5 43.5 13.2 17.3 4.20 5.50 1.93 2.25 For a property index equal to 1, all two-phase multipliers have a value of 1. Coordinates of Two-Phase Pressure Drop Correlation TABLE 6 : 40 60 2,050 4,300 330 538 110 155 53.0 69.0 21.2 26.0 6.50 8.00 2.48 2.86 80 6,600 760 203 85.0 30.0 9.10 3.20 100 10,000 1,000 250 100 33.3 10.0 3.33 I I-' I-' I-' I -------------------

p T psia oF 211. 97 240 119. 2 280 117.9 340 247.i 400 466.3 460 7113. 5 510 1002 545 1225 570 1541 600 1784 620 2057 640 2529 670 2892 690 II * "Prop * "Cale": µg µL io-6 lb/ft.sec 8.23 160.0 8.73 134.o 9,55 104.24 10.42 87 .111 . 11.18 75 .116 12.01 67,73 12.61 60.38 13.61 59,58 14.28 56.80 14.87 52.94 16.05 119.58 18.06 43.95 20.58 35.04 PPI = (µL/µg)0.2/(pL/pg)

PPI from polynomial l/pg l/pL ft 3;1b Prop. 16.327 0.0169 0.00187 8.650 0.0773 0.00345 3,792 0.0179 0.00761 1.866 0.0186 0.01524 0.9961 0.0196 0.02892 0.6153 0.0207 0.04765 0. 41149 0.0216 0.06641 0.3537 0.0224 0.08509 0.2677 0.0236 0.1162 0 .. 2209 0.02117 0.1441 0.1805 0.0259 0.1798 0.1278 0.0288 0.2692 0.09113 0.0325 0.3833 Comparison of Actual Calculated Values of PPI rrABLE 7 PPI Cale. .00175 .00324 .00741 .01526 .02885 .04723 .06640 .08518 .1160 .143'7 .1808 .2683 ,3901 Error % -9.3 -6.2 -2.7 0.1 -0.2 -0.9 0.0 +O.l -0.2 *-0.3 +0.6 -0.3 +l. 7 I I-' I-' f\) I x % PPI=0.0001 0.1 1. 307 0.5 1. 442 1 1. 543 2 1. 662 3.5 1. 781 5 1. 796 7.5 1. 785 10 1.790 15 1. 731 20 1. 719 30 1. 695 40 1-730 60 1.652 80 1.862 -113-n 0.001 0.004 0.01 0.03 0.980 0.788 0.743 0.813 1. 015 0.785 0.712 0.773 . 1. 054 0.782 0.710 0.805 1.094 0.800 0.718 0.778 1.141 0.833 0.731 0.774 1.136 0.814 0.723 0.729 1.178 0.844 0.730 0.719 1.161 0.853 0.719 o.678 '1.177 0.863 0.707 0.627 1.190 0.864 0.695 0.614 1.188 0.893 0.705 0.579 1.212 0.902 0.705 0.525 1.215 0.941 0.738 0.517 1.231 0.938 0.739 0.502 Value of n in M = 1-X + Xn/PPI TABLE 8 0.1 0.796 0.827 0.819 0.766 0.748 0.737 0.714 0.721 0.680 0.610 0.610 0.576 0.537 0.522 0.3 0.827 0.923 0.839 0.793 0.723 0.700 o.678 o.684 0.679 0.672 0.636 0.625 0.595 0.472 I I I I I I I I I I I I I I I I I I I

x 0.00026 G = 0.001 1.669 0.010 1.669 0.035 1. 626 0.050 1.600 0.075 1.590 0.100 1.580 0.150 l.580 0.200 1.580 0.300 1.534 0.400 1.492 0.600 1.362 0.800 1.178 x 0.001 0.010 0.035 0.050 0.075 0.100 0.150 o. 200 0.300 0.400 0.600 0.800 CORRECTION FACTOR AT PPI = 0.00380 0.05700 0.19800 0.00026 0.00380 0.05700 0.19800 0.00026 0.25 MLB/SQFT.

HR G = 0.50 MLB/SQFT.

HR G = 1.160 1.220 1.110 1.300 1.130 1.100 1.078 1.000 1.158 1. 307 1.166 1.330 1.250 1.15() 1.086 1.000 1. 059 1.355 1.420 1. 311 1.170 1.150 1.232 1.000 1.000 1.384 1. 572 1.300 1.120 1.214 1.320 1.000 1.210 1.502 1. 695 1. 300 1.148 1.210 1.334 1.000 1.420 1.360 1.818 1.300 1. 276 1.219 1.460 1.000 1. /120 1.360 1.818 1.304 1. 256 1. 223 1.472 1.000 1.420 1.360 1.818 1.308 1. 236 1.240 1.596 1. 000 1. 324 1.330 1. 619 1. 284 1.195 1.235 1.457 1.000 1.234 1.340 1.445 1.260 1.153 1.230 1.318 1.000 1.139 1.162 1.204 1. 200 1.110 1.130 1.164 1.000 1.103 1.086 1. 070 1.100 i.070 1.084 1.061 1.000 CORRECTION FACTOR AT PPI = 0.0026 0.00380 0.05700 0.19800 0.0026 0.00380 0.05700 0.19800 G = 2.00 MLB/SQFT.

HR G = 3.00 MLB/SQFT.

HR 0. 750 o. 740 0. 749 0. 75!1 0.752 0. 750 0.736 0. 722 0. 746 o. 770 0.820 0.910 0.864 0. 905 0.970 0.630 0. 780 0.865 0.660 o. 880 0.912 0.610 0.484 0.810 1. 676 1. 829 0.817 0.625 0.501 0. 741 0.686 o. 798 0.760 0.634 0.512 0. 700 o. 704 0.805 0. 730 0.634 0.551 0. 701 0. 721 0.812 o. 700 0.634 0.590 0.702 0.746 0.788 0.665 0.606 0.605 0.673 0.750 o. 764 0.630 0.598 0.620 0.643 0.788 0.730 0.602 0.624 0.667 0.593 0.806 0.696 0.574 0.650 0. 714 0. 5112 0.860 0.705 0.574 o. 718 0.782 0.542 o. 932 0.820 0. 700 0.836 0.880 0.690 Values of the Mass Velocity Correction Factor for Between in Figure 9 TABLE 9 0.9J7 0.884 0.769 0. 700 0.671 0.642 0.587 0.540 O.l193 0.454 0.454 0.580 0.00380 0.05700 0.19800 1.00 MLB/SQFT.

HR 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1. 000 1.000 QUALITY o.o o.o 1.000 0.001 1.347 0.010 2.270 0.035 4.350 0.050 5.940 0.075 8.525 0.100 10.517 0.150 14.624 0.200 17.544 0.300 22.256 Q.400 26.999 0.600 29.117 0.800 31. 647 1. 000 32.262 .... o.o -*l. 000 0.001 1.277 0.010 1. 735 o. 035 2.942 0.050 3.690 0.075 5. 020 0.100 5.477 0.150 7 .292 0.200 8. 630 0.300 10.973 0.400 13.058 0.600 13. 905 0.800 14.804 1.000 15.096 o.o 1.000 0.001 1.123 0.010 1.266 0.035 1. 915 0.050 2.322 0.075 2.837 0.100 3.267 0.150 3.796 0.200 4.238 0.300 4.607 0.400 4.680 0.600 4. 665 0.800 4.724 1.000 4. 716 -115-MASS VELOCITY MLB/SQFT.

HR 0.25 a.so 1.00 2.00 3.00 PRESSURE = 500.0 PSIA (PPI = 0. 0310) 1.000 1. 000 1.000 1.000 1.000 1.347 1.236 1.117 1.000 0.944 2.270 2.090 1.782 1.480 1.313 4.350 3. 898 3.376 2.683 2.319 5.940 5.460 4.577 3.537 3.010 8.525 7.099 5.936 4.643 3.961 . 10. 517 . 9.432 7.657 6.061 5.182 14.624 13.101 10.648 8.290 7.003 17.544 15.827 12. 773 9. 719 8.147 22.256 20.536 16.751 12.447 10.212 26.999 24. 877 20.514 14.785 11. 912 29.117 28.329 25.170 18.622 15.001 31. 647 31.386 29. 039 24.543 21.278 32.262 32.262 32.262 32.262 32. 262 PRESSURE = 1000.0 PSIA (PPI = 0.0662) 1.000 1.000 -1.000 1.000 1.000 1.277 1.161 1.058 0.966 0.925 1. 735 -* 1.537 1.345 1.189 1.102 2.942 2.504 2.159 1. 786 1. 607 3.690 3.218 2.623 2.081 1.836 5.020 4.031 3.291 2.619 2 .295 5.477 4.830 3.870 3.090 2.689 7.292 6.457 5.153 3. 984 3.414 8.630 7.823 6.098 4.560 3.845 10.973 10.144 8. 039 5.745 4.670 13. 058 11. 976 9.653 6.576 5.130 13. 905 13. 512. 11. 914 8.211 6.311 14.804 14.765 13.656 11. 000 9.241 15.096 15.096 15.096 15.096 15.096 PRESSURE = 2250.0 PSIA (PPI = 0. 2121) 1.000 1.000 1.000 1.000 1.000 1.123 1.091 1.016 o. 986 0.954 1.266 1.182 1.092 1.001 0.971 1. 915 1. 671 1.372 1.135 1.072 2.322 1. 964 1.510 1.169 1.082 2.837 2.254 1. 716 1.280 1.183 3.267 2.648 1.851 1.329 'L226 3.796 3.100 2.150 1.474. 1.314 4.238 3.738 2.400 1.566 1.361 4.607 4.164 2.915 1.825 1.523 4.680 4.285 3.299 1. 977 1. 601 4.665 4.518 3.914 2.346 1. 900 4. 724 4.687 4.430 3.182 2.679 4. 716 4.716 4. 716 4. 716 4. 716 Multiplier Arrays at Various Pressures TABLE 10 1000.00 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 LOOO 1.000 1.000 1.000 LOOO 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000 I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I -116-READING OF PRINTING OF CALCULATIONS INPUT AND INPUT CARDS OUTPUT E COBRA IIIC 0.20 7.17 5.68 COBRA IIIC/MIT (Old Input Data 0.20 3.27 5,83 Pres.enta COBRA IIIC/MIT 6.50(l) (New Input Data 0.08 3,38 Presentation)

(1) In this case, two print outs of the Input Data were obtained, while in fact, one is enough. Comparison Between Costs of COBRA IIIC and COBRA IJIC/MIT for the Problem of Figure 12 TABLE 11 I I J -117-I I APPENDIX 1 I Subroutines ACOL, DIVERT, DE COMP and SOLVE I I I I I I I I I I . ,, I I I I I I* (" -!" ,_ ,, c . -I I I* II I I I I I I (" ., I -,_ I* I I I I -118-IJ o 0 UT I 1\! F: ll C 0 L ( T F 0 '-1

p<

J K , KM A X , L 0 C A '

MS , N K , G

I P I LE ) nr'<'IFNSinN Ir<<1l .JK(ll ACOL SET LOC!*

= 1, CALLEQ

= 2. FqOM MAIN COLD LOCAC><,Ll*L=2*7 SPECIFIES UP TO llO.JAC=:"H TO \,-IAN\IELS

J '1 ,.. ><: = l ' ,.., -< IP (rcru::.GT.f))

GO T'.1 1117 D') 103 L=2*13 tn3 LOC6.

='1 r:.O TO 110 107 DO 3L=2,7 1 LOC.1(><',U=n lJ. 0 .. j: 1 LOC A ( K, l) = r< II = IKIK> JJ = JK(Kl .:.. 00 7 t<.r<= 1 ':>..JI<:

III = IKIK:><l IF (!II.GT.fl>

GO Tt/ 7 JJJ = JK (!-<'<) ! F ( < I I *

1

I I Il

0 P

c !I

E

J .J J > l 3 0 T 0 6 GO TO 7 IF ( (I!l+JJJ -II-JJ) .::a. Ol GO TO 7 N = l\J+ l LL = IT I IF (I!.E!).III>

LL=JJJ = FLOAT<II-Lll/FLOATCII-JJ)

LOCA (K

= IF (wV.LT.().0)

U1CA(K,N)=-K:-<

7 C'>"JT rr.1 uE IF GO TO ina LOCA (K

Gt') TO 109 1 *1 ;; L *) <: A ( "<

R l = \J l!iq IF (II .GF:.JJl Gn TO B r r = JK ( r< i J:J = !K(Kl GO TO t... ,:; c 0 Il'JU*

STRIPE FOR AAA MATRIX )!VERT NIAX = 0 !)0 10 K:l,\JK N=LOCA(K'.tRl IF (!iJILE.GT.Ql GO TO 111 C ( < , l t. l 111 n0 10 L=2*'l =

J =

P<-Ll(L) IF" (J.LT.M!iXJ GO TO 10 :-1AX = J KMA)C = K 10 CONT I l\iUE MS = 2*MAX + 1 COOE21MS,NKl RE TUR'! E!-lO c c c c c c c c c c c I, c c I -119-,,

IJTVF'Qr(Jl r-irs ;:ioor.i:::D1J;:lF' It\!<:;

COWA()N ANIJ TYPE STATE'-"El\ITS SY THEnVRTl'IO,.

'-' A J () P SU 8 0 U T Pl i::: S () F" C 0 o a. -l I I C *

r') V T f"l 0 ' , T >.i o L I C I T I T Er; f. q ( 'S l COMMQr-i

DIA .DT ,ox l .Ge *GK .HSURF tHF * :::> 1 t. c:; 7 2 1 4 5

,1-1r, .r2* .!3 tJ2 ' J3 *. J4 ,J5 ,J6 ,J7 tKOE'3UGtKF' ti<LJ

,\JAXL .N88C tNCHAN oNCHF tNOX tNF , "IGAP5 *"J(";P!IJ oNGP!f'lT,\JGTYPEtNGX,L

NK tN011ESF'tNPROP , NRAMP .NPOD .NSCRi. .\JV ,PITCM ,POWER ,PREF

QA1 .RHOG oSL tTF ,THICK , 11F ,llF' .VFG ,VG ,7 C0*.4vQN ICOQPD.::>/

AA(t...}

AF(7}

AFACT(lOtlO>

AV(7)

AXIALCJO},

8XC3fll* CC(4),

CFUELC2l*

nFUF'LC2l*

GACYL(l!il*

GF'AC_T(qol(lle HGAOC2l* HH'°(30l*

HHl;(3Q), TGPto(ln).

1<CLtOC?l.

-<FtJELC2lt Kl(=-('3QJ, "ICH<lOlt Nc;Ap(oJ, PPC)f)). .ssIGMA(30)'

TCLA0(2)'

UUF'(30)'

VVF(JOl. l/VGC'30l, :t.n11AL(30l, YC30l

TTC30) LOGICAL Pf::AL GC!D KTJ, KF* KKF, KCLAO, KFUEL COMMON 1coqPAJ1 ,MG .MR ,Ms ,. 1 tiA ,$AAA o$AC t'SALPHA,'§AN , l q;r,'.) *

$CO\J * $COND , SCP * $0 , $1')C , 'liOF"DX , ;:i 1iDHOX

,$OPDX ,;;QPK: ,'§OUR ,$nR *SF

1 4 .iH c;

SLC *SLEN'3TtSLOCA tSLP , , 7 t'SPWRF ,$QC ,$QF' <<;QUlll

,5RHOOL,$SP

ST .':llTOUMY,$TTNLE,'STROD , q .iU5AVi;::.'5USfAP,$V .svrsc ,$VTSC:\IJ9$VP tSVPA

A

,$WP t'SXCROS COMMON 0AT4(}) Ll1G!CAL*LOAT!ll PHEGEq IrH TC 1 l E t:l U T V AL E l'J C E < I) A T A ( 1 l , I n A T ( l l , L D A T ( 1 > > A8ITCIZ.z1.z2,73,z4,z5,z6l

=

-Z2 + Dx/0Tl/Z3

+ l Z4*ABS<ZS+7-6>*0Xl

!PILE =J7 "-!Kl< = "'IK JMl = J-1 SLOX =

DTGC =

DXGC = DX*GC I I I I I I. I I I I C USTAR I I I I c c I ,, I* I I* I I I I I* I' c c c ,.. ,, '-I I I no c:; K= 1 .1'!'<!< t !=!OAT JJ=IDAT('liJ'<+Kl 0 A T A ( '5 I J S A 1/ c:-+ '< J = f) A T A ( <li I J c; T '1 + "\ l -120-I) AT A ( ';IJ 5 TA ::i + i<) = 11

S* ( IJ J\ TA (SU+ I T l +DAT A ( liU + JJ J l c::; CONT I NUF.: SF:T AAA AOOAY IJSTNG LnCA <SET IN ACrr_. BASED ON DATA) LMA x = '*l!IJ = ("4S+ll/2 D 0 l O = 1

1'1 :< [)') L=l *L*'-1A:t

?911 l=0. I T=IDAT ($!!<+;<)

l

?

SAVi:=ABIT(l .OATj,($1..J*ITl .JATA($USTAR+"\J ,QATA($A+!Il .DATA($('1Pl<+!Tl, 1 0 A T '1 ( F

T T + 1 .. 1 C * ( JM l -l l J , 0 A T A ( 'Ii F + I ! + MC * ( J -l ) ) ) :2 + AB I T ( 1 0 AT A < 'f lJ + J. J J * ) AT A ( <;; IJ ST A i:l * -<: I

11 AT A ( 3 A* J J l , 0 AT A ( $I") PK + .J J ) , 1 IF <IOILE.GT.l)J GO TO GO TO 7214 7213 7214 DO l (Ll-ll l IF ILL.F.:0.1>

GO TO 29C:: IZ = 1 II=" <L.LT.n> !7=-1 l = IABS<Ll TJ = JJ IF< <I.!.EO.T!):H(<=;IK'+ll l .'.)R. 1 IJ=II SAVE = ABIT<I?.DATA C<iill+TJJ ,QATA('1iA+IJJ

l DATA (<i;l')OK+IJl ,QATA (<t;r+IJ+MC*(JMl-lJ l ,DATA ($F+TJ+MC*<J-ll J)

L = MI n -K'

L 300 I=

310 Cf')NTINIJE IF < J 6. LT. l l GO T n l IJ S MOl')!FY SIMULTANEO!JS TO ACCOUNT FOR SPECIFIFD VALUE<:; OF C Q 0 S SF l 0 1:J r; I V 1\1 I i::; UP i:i n I J T T E 0 PC E 00 go 1<=l*'lt<

GO rn '10 D O P 5 l = l * !"'. LL=MIO-K*L GO TO IF< LL. GT. L 1.411 X. nq. LL.LT. 1 l GO TO 85 IFCLDAT(5Fi)IV*Lll l

l* l 0ATA(5W+l+MG*(J-l))

85 CONTINUE 911 CONT PJU!: 1)0 100 K=l.NK 1")VRT040C i)VRT04lr:

r")VQT()42C r')VqT 043 i l)VRT044i i)V'H045C nVRT064C nVRT064l r)VRT065C i)VRT0651

")VRT0652 r)VRT06S:3 r)VRT065<-

nVRl'.0670

!1VRT06<:;lQ DVRT0700 l)VRT0710 DVRT0720 rJVRT0730 OVRT0740 OVRT0790 DVRT0800 f"JVRT0810

-121-tF(.'JOT.LilH(:lii='l')T\/*l<'.l l 30 TO 100 nn gs L=l *LoAA'/.

= oJ.O LL= tw1AX0(1.(l+K-1*:1(1)ll i_L=M!llJfl MOICU=M!0*"'-Ll.

g f) A TA ( ;: A *LI_ +Ni<* ( MD I CI_ 1-1 l l = 0 * (1 ion 1r)5 2 l l () 150 1011() 1 0 AT A ( A A+ I'( *

'-> ( M r ')-l l ) = l

0 CONTI IF T.] l :;o rn l 1 !') PP. I T 2 * ( < n A T A C :;; A A A + k' + l'\J <: ( L -1 l l , l = 1 , L X l , DA T A ( $ + K l

K = l , *'I K K ) F"l')RMAT(lHI).

lP7F="l5.L..) . CALL l . 0ATA($3*1>

i\JK) IF<!ERPnR.r;r.11 GO TO 1()00 CALL ,JATA('bANS'.NE+l} .DATA($i:l+l) ,NI<} rJO 150

()AT A ( 51.i

K + 1.1 G * ( J-1 l l :DAT A (

\IS \v F. +..:; l p i:rr_JR "J PP.! "lT 1 i='QonR IN i'Ji='Cnv.P, IF.:C(PQP

= 3 P.;:TURN E: f\JI') I .I i I I I I I I I I I I I ,g c IC *1 'I* I I* c I I I I 1* I I I I I -122-5 tJ RP. 0 l ) T I i*J f n E c n M 0 ( N "'

p::* :.i R 0 R ' L ..., A x * '-1 I) ' lJ L

x

R ' N K ) n I "-1f:1'J<::;

r ON ! IL ( j\Jf(. l ) * 'i ( 1 l * ( l l SIMPLIFIED vEoc:;tON OF" WITH ?!VOTING S T () ;::;> E f) ! AG '11\1 ti. L 8 M*I (1 0 'r 6 A A '-' A T q I x * :i 0 S I T I 0 N P< , l l I N S QI I A P E 4 R 4 Y 8 EC 0 *AF" S ( i< * ( M I i) -:.< + Ll ) T t\J "J * 'n' <1 q A Y

l\J =

-

PFTUR".1 NM! = l\J-1 iJO l 7 = 1

PIVnT =

KPl = K+l LIMIT = "I>.Jn (N. C'<"+MI*1-l) no 16 T =

l<K' = .*.qn+k'-!

f:M = -1.JL ( ! *"'"<) /P!VOT !IL ( T ,t<i<) = -::M IF 2n 00 21 J:KPl*lIYTT JI=t..A!f)-!

+J JK = MTG-K*J 21 IJL C I oJ ! l = I Jl CT, JI) + f-*..ti>IJL ( K .JK l 1 f:> CIJ"lTI NUE 1 7 cor-n I NtJE IF (UL ("-I *MT[)) ) l Cl , l 8

l 9 18 PC INT 112 100 113* CCllL(-<'.,L),L=l*NN>,K=l*NNl 111 112 Mt.TOIX IN nECOMPJSE.

ZERO DIVIDE IN SOLVE. !E 0 qop = l/ l Cl E"ID DC0Mf'l060 DCOMOJOO OC0M0460 OCOMOSlO OCOM0570 OCOMOS80 DCOM0610 DCOM0620 DC0!'-40630 l!')COM0640 DCOM0630 DC0M0660 r)COl-10670 c c c -123-s lJ 8 c 0 u T i: c; 0 L I I F:: P*l "' ' L M <l 'I.

M El ' u '-* x

R * \j "( ) Ii r s r 0 l\j ! IL ( "1 k"

1 ) * ;( ( l ) * ( l ) STOCE IJIAG(11\IAL

i1.1'JI')

0'" AAA "iATPIX.

o<.u SQlJAQF.:

AP.PAY i=!ECO*"'ES (K'. ("!!i)-o<+L) l Ii\J l\JE:"" = !\JN* I F' ( *,1

F Q

l l Ci 0 TI') S NPl =

X(ll = 8(1> rJ o 2 r. = 2 * = I-1 = O.O JM!N: MA.Xl'l(l.(!-M!l)+ll l . C IJOIJBLE PRECI5!1')1'1 MAY 8E Pf'.:'QJiREO FOR !N"ER LOOP. 00 1 J = JMTl\Jo!Ml JJ =

l SUM = SI.JM + UL< T. JJ)

X LJ l ? x < I l = 8 < ! l -c;1 I'" XC"ll = X(l\J) /Ul.(NoMT'"ll 00 4 !BACK = 2*"1 I = MPl-I8AC;( C I GIJE'S l * ***

l !Pl = I+l = O.O C DOllALE Pl=ilECISJ<11\J MAY 8E: l:IS::IJJIRED F"OR IN\JER LOOP.

=

l DO 3 J = rct.JMAX JJ = MIO-I*J 1 SUM= SUM ... UL(!.JJl*X(Jl 4 X<tl = (X(Tl-5UMl/UL(!,,,,.IDl s x ( l ) = 8 ( l ) /UL ( l

M l:r)) P.F'.TUPN E"JO I I I c;OLVn0'30 SOLVoo**

SOLVOlf 50LV0110 SOLVOll I SOLV0160 i;OLVOll SOLV0200 50LV02j SOLV02 c;OLV02" SOLV02SO I c;OLV021*

SOLV031 I 'I I. I I I I I I I I I I I ., I I I I I I I I -124-APPENDIX 2 Input Data for 32-Channel Case COf.'.IT*1ENT CARDS HAVE BEm.J ADDED FOR GREATER CLARITY, THEY ARE NOT IJ: .. TCLUDED nr THE I:'JPUT DATA. 2000 1 1 0 PWR NC=32 STRIPE HUMBERIHG CARD GROUP 1 (FLUID PROPERTY POLY. OPTIOH) 1 -1 1 1000.0 2600.0 0 CARD GROUP 2 2 1 0 0.184 -0.2 0 CARD GROUP 3 3 21 O.D 0.144 0.050.309

0. 251. 45 o.*3oi.54 o.351.59 0,551,40 0.601.28 0.651.14 0.850.450 0.900.312 0.950.200 0 30 0 0.100.783 o.4oi.6o 0.700.95 1.000.094 0.151.050 0.451.56 0.201.280 0.501.50 0,750,775 0.800.588 CARD GROUP L!. 7. 8 (SIMPLIFIED DATA PFE.SE:TTAT!ON) . -I I I I I *I I I I I. I , .. 0 1 1 32 7 1. 0 3 4 1 37.43 0 1 348.o 310.2 0 .122 o. o 0.374 265.0 I 4. 0110

9 7 80

5 5 5 l.1281.1521.1521.0581.0901.1651.1821.1521.0031.1641

1651.1521.0631.0031.0901.123

1. 06 31. 00 31. 0901. 12 81. 00 31. 16 41.1651, 1521. 0901. 1651. 1821.1521.1281.1521.1521.168 0.05 1.16 2.32 3.49 3.66 3.87 2.99 1 4 8 2.0 .08 640 ** 3228 8.80 .076405.0.02241000.

c CARD GROUP 9 9 1 1 0 . 5 0. 144. 0. 10 1 1. 20 c CARD GROUP 10 10 1 0 0 1 0.02 I I I I I I I I I *1 I I I I 1* t 1* I I c c CARD 11 1 2250. o.o 1. 0 o.o 1. 0 0.1.000 0.1.000 CARD 12 3 132 132 1 4 5 -126-GROUP 11 1 2 2 2 2 557.3 2.66 0.189 1. 1. 056 1.1.15 1. 0.77 1.1.000 GROUP 12 2 2 3

-127-Appendix 3 Timing !Iethod I I I I. I I I: 1* ., I I I. I I I I *1* *t *1 I I I *a I I I I 'I I ,. ** I I I *1 I ,. I I I IN --0 -1 2 3 4 5 6 7 8 -128-Successive calls to the nei:*r. subroutine PP1J':1Ir.T (r:t) -coding in Appendix 16 -were made from and SCHEI7E in COBTI.A IIIC, as tabulated below. .. Called After Description from Card Ml'\IN (old 1304 Def ore outer DO Loop on time steps MAIN (old 1306 Within and at start of DO Loop on .. _ .. tiMr. st-M.'\IN (old) 1610 and at start of DO Loop en i terc:i. tic SCHEME 0482 Within and at start of DO Loop on axial steps . SCHEME 1080 Irn.r.lediately before CALL DIVERT " SCHEME 1090 Immediately after return f=oM DIVERT MAIN (old) 1650 Immediately after return fr on SCHEME MAIN (old) 1730 Just after end of DO Loop on iteration MAIN (old) 2360 Just after end of DO Loop on steps The calculation followed the sequence.

0: Start timing run 1: Start steady state calculation 2: Start first iteration over whole channel 3. Start interval calculation to set inlet conditions 3: Start interval calculation for first interval 4-5:* Solve cross-flow equations (3,4,5): interval calculations for successive a.:<ial steps 6: return to MAIN having completed first channel iteration PRINT RESULTS FROM PRNTIM 2: Start next iteration 3-6: second iteration, as for first

-129-PRIUT RESULTS FROM PrtlJTI::-1 (2-6): Successive iterations, PRIIJT Ii'Ror1 PR11TIM at each 6. 7: iteration closed Print COBRA output for steady state 1: Start transient calculation (2-6): Successive iterations as above. P?..PJT PROM PR:TTII1 at each 6. 7: iteration closed print COBRA outout for transient 8: Calculation finished PRIHT RE.SULTS FROM PRNTIM During each hydraulic iteration, i.e. 3-6, the total tine from the first 3 to 6 was stored and also the sumnation of the times fror. each 4 to 5. Also the number of passes 3-3 and 4-5 were counted. At the close of each iteration triggered by IN=6, the following were printed: (a) the total time from 3 to 6 (b) the summed time 4 to 5 (c) the lc..st times at TtThich stations () to 6 *;rere passed. This informc..tion enabled the fo llmring to be deduced printed: (1) time to solve the cross-flow equation, i.e. (b) (2) time in rest of hydraulic calculations, i.e. (a) -(b) I 1* I ,. I I ,. I I I I I I I ., I I I I I I. I ' I I 1 I I I I I I I t I ,. I -130-(3) time for COBRA printing + time measurement printing = 1 (transient)

-last 6 steady state) and 8 (final) -last 6 (transient)

(4) total time for run = 8 -O The results given in Tables 1 and 2 are (1) and ( 2) above. The other results were used to check the internal consistency of runs. Load Modules: A symbolic source module was written with the dimensions given non-numerically, e.g. *&MC to represent the number of channels; this can be seen in the COMMON list of DIFFER in Appendix l:r.. Two com.piled load modules, a "small" and a "large", were created from the source setting the array sizes through the Job cards. The small module was used for the 10, 16 and 32 channel cases and the large for the 64 and 128. A check run with 32 channels on the large load module showed that the running time was not module dependent.

The load module dimensions are tabulated below. Symbol i.e Max. No. of Small Large &MC Channels 32 130 &MA Array AAA 3136 15K &MG Gap Connections 56 250 &MN Fuel Nodes 10 10 &MR Rods 32 1_30 &MX Axial Intervals 31 31

-131-APPENDIX 4 Faster "Rest of Hydraulic Calculations" I I' :1 I II I ,. I I I I t I I ' I' ., I I I .,, I I ,, I ,, I I I *1 I I I* 'I 1 I' I I -132-Subroutines DIFFER and HEAT were modified to speed-up the hydraulic calculation

-coding for DIFFER in Appendix 17, that for HEAT is relatively trivial. DIFFER calculates for each channel, the local axial enthalpy gradient DHDX, the flow gradient DFDX and the pressure gradient DPDX, without and with cross-flow.

Each of these terms involves the summation of the effects of neighboring channels.

In the original COBRA IIIC, the coding method was to use an outer DO Loop on the channels and then, for each channel, scan every boundary in turn in an inner DO Loop and to sum the effects. Determination of which boundaries were irrelevant required lengthy IF tests each time. The modifications to DIFFER consisted of using an outer DO Loop on the boundaries, assessing the interaction across the boundary and then apportioning this appropriately between the two* channels concerned.

With this technique there is no inner DO Loop and no IF test to determine boundary relevancy.

This gave a considerable reduction in running time. Care was taken when coding, to avoid adding successive small numbers to a large one with loss of accuracy due to round-off.

For example, in the enthalpy calculation, the small cross-flow, mixing and conduction components are summed before being added to the channel power. The COBRA answers obtained were identical in every respect to those obtained from an unmodified DIFFER -even to the last digit of the transient cross-flow.

-133-HEAT was also modified.

This calculates the heat flux (FLUX) and channel power from the rods (QPRIM). In the general case, several rods may heat a channel and the effect requires a DO Loop summation.

H.owever, for PWRs and BWRs there is a one-to-one correspondence between rods and channels and the summation is superfluous.

HEAT was modified by skipping the summation when IPILE = 1 or 2 (i.e. PWR, BWR) but not when IPILE = O (normal COBRA). The effect of the modification is to make a minute change in the COBRA answers because of round-off.

For example, by card HEAT0690, TFLUID is calculated as (0.0 + T*PHI)/PHI and this .can be different from T (the short-cut modification) due to round-off, particularly in the reactor case when PHI (the number of rods) is large, e.g. 264. The effect accumulates during the calculation to give a detectable but barely significant difference in the COBRA answers. I I I I I' I I I I I I I I *I I' I I I. I I I 1. I I I I ,. I 1. I I I I I I I I I -134-APPENDIX 5 Subroutine Core and Block Data c -135-<::: 1 J .:; :;; 0 U T T *

t: C ;) . .; =

I':r=_:,-,,:::..;

(-5)

) "'r.::. ( l l r: 0 , *" "'} 1 / c:-, ...J *::. 3 / *.., * *..j c , ,*;1 ,-; * *1 *\J , "" * "" s

M x * 'b s '£ , s n R G < q 6 > T*1Ti:1,;::.; 'f ( .:;.; } C 'l *) IC r; ...Ju Q. .::; I *-::I\: AM f S

5 L X. , S T Y P ,, ! '. f i\j FH T I :.; E ( :..; i;., }

-.:-"' S':t,*iE°S

'A;\ : l '*" = 1 T F ( Wj

L.:: * 'l ) ,\.11; = 1 Ti="

'-1\l=l

IJ l c n r = 1 * <;;; ., s 1 rir1 ..i

.1_ x < I J :r.1c <t LY. < 2 l :i) ( -'.-) =**:G :;; L .'f ( 7 } = "'"; .;; I_ X ( ) : , X <:'. ! _ '( ( C'.) = '*! 1--: -=* c:; -;: :_ i I i :1 ) :

'.1 ,(

-1; i_ >'. ( 1 ... ) = '".; -; 1_ x r 2 n J = :"i c ;:* ,... :::1_ 'i ( 23} ::"'\{] ( 24) ::'1.:> x r 25 l ::**1c ;:--.1 x 'hLX ( 2i:. l T f_ .< ( ;:'. 7 } : 1-.j (:, i:: L :t < 2 q l = '"1;::; "" X :;;1_ X ( 3 l) :MC;.'."IX

<:: !_ .X ( J,... ) =Mr;'*' S x ( 1:3 l

.;.1_ 't ( ],;., ) =*-ir; :::Lx ( 37)

I __ ( ( ] ,;: )

... '.( ':. I_ '( I ::. l ) : ,\.t !, M X -l:;t_ '.( ( ... .:.. ) =*-1;.; * *-;;1_ 't ( :r11r, :::-1_ x < 411 l ::"1r: -;,1_:, ( 4>-:i l :MG -.; I_ '( ( :+ *:) ) : M <:: -::* .:..

c;:l_X(5ll=MG'-l-l.:.. <52J =*"IR'.l-6 ,_ x ( 5 3 ) = :wi ;(

<5.:..l

.::1_;r ( 55) :"-1.'.( 'ii_ X. ( 5q l

.:: L x < 6 2 l = .*P -;:. ,::.. -.::Lx ( 6t:.l :*.i..;i :l3 l x ( 6 5 ) = \l <<:!_X. (67) :M(i:*"'1P X (

'AC'.l-"'1X X ! (:;q) :i'"1CnA X L X ( 7 2 l

= q <<ii_:< ( 73) : 0'1(iHAX '5l X ( 74. l ::*"'lC*MX o<;l_.X ( 75J =*'-'1G'-'"'1X p 0 0 'I rn t: s;:; Ac E F"' 0:;; s p r 3 'fl R IT E c:;i AT r 0 \J

I I I I .I I I I) I I I I I 1* I I I I I I I I I ,, I I I I I I I *1 I I I I I 1 1 :1 -136-T :;' ( ; i_ ,I ( i :.: )

LT * 'i:* ,.,, c ) x ( 7:; ) = 3 *Mr, : !_ f ( ., 7 : = .; . : *-. ! -/ ( 7 .* j } = **1 :'.; :1 :.. -..:-. .tj " *: ! __ *' { j ) = ' ':"! * .( ( .:. 1 ) :: .::: ::; '.(I = *!,.( :::;_'/ (".:,""' =

':: I_ f { -:_ 5 ) -: :'. i:*

M* I :,."';:;*-: ( l) =l .:::_, i 'i = !j :-.r Jlr i=*1

r.1_ f 'i : 7,1 .X i,, .. ;:.1. ( T ) '. ( T

f

l l :>i

1:.;. () ( I l :: ; Ii,;) r; < I -1 ) + i L X < I -1 l "'**" t:. y = 1 V ;( s + "*' '} ! I ( ,."" + l

2 ) ,:-( "'*'.:.

L i . L /,, x l

-<=i<::,,rn;;

lnn j\!\T,)("')

= n.c, :-i r) l ... :; t' ,. ; = 1 * .., ;;; 5 :.. 'In <;: !'. c;. :-:, { '** ) :: ( *'i ) + K' 5 -1 r-r (" ::l () ] :::.n? ". " : T ::; f c n ,.., ;: ._, ' * *1 c ;:;

t\ o J t-..:/.

r-l_ j ( ? ) = 'f . ..:. c:;n:.;.,-;

< ! J = ( l *.1>L x ( > ..., I_ 'i x = 'I I_ x ,'I, .. ' I_ x ( ? ) r ( '*i.:.

LT * *.; L;; X ) ':. :*; j T ;: ( * -' i ; ij I ) c: T *'\ ::> *' .:, r T *::: I -:.

3 Ii n :;l ) c: T . *' f p y ( .'! .,. ,=-< GO TO 902 i:-::;n'" I,.,..,, .:-(,'-i Q.;; !1'1T ,, _; T T :_ ( "-,

l r, r, \ 1 ) :.1 ' 1'<1 C * '-1 S * "1 N * '-1

S ' "1 )( ., ,; T T::: ( ":, * '..-il r. !" ) ( :.J * *.:;1'1.t..

Mi:-.::; ( \j) ' -bl x ( f\j) '1i 0 G ( l\j) ' s Ty p E ( N)

N= 1 * $$Cl; ) " ,;i r ,.. r:.* ( -* ) -} n r. ) v* '.A ti x I

= ....

.. : T 7 c:* ( . .,

l *"'. r* j S L X 1., , 1 _ 1*, .. ; ::'. r< , ') _,., " >="...,::;,

r , / / 1 * .1 '{ D. *i r r ., ;:.; )"( .:. Y s i z Es ' , / , .. 1 A = ,

r s * / , = ;, = '* TS. /. = 1, I. -=> -= *

l -* ' * \.., """ -'

r s , /

A = ' 9 r :::; 1'1".::. :-".l 7 1/. J 1*;.i*-1 IC

= '* Il4* l**OROS*, 11, ::* ": i, 1 .*, , c. ! / ; , ) "i . .1 r: L r: t -< iJ co 1 JLG >i ll v :: 1 EE N RE: n.1 ! c En 8 Y * , r 4 * -,nr:r; =--.,;*.,:--:'.,._

"*"'"'"*'f(

()F COP;::: 3JT *,no.* 1 i111iOOS*/)

I( I )

-137-<:11'\J OF COO* F"ATLED*///)

1 i) (\ ;;) ;:-*1 Q.; ;:. T I I ) ,") v.; A *.1 [ c J. LI -I)\. :\ T I 0 *\j 0 F c 0 RE 3 0 T 0 NL y '

I 1 () ' ' '*"Ci R !} 5 ' I I **ii

'"F

.;lfr;11 i Ht:;:) F!lR T:-! rs PR0t3U:.'.M I'S ,

I 10 I/ I) 1 I// .. ,

TJS.*-----

LE"'GTrl ------TYPEt/ ----------*/ I I I *1 I I I I :1 I I I I I *1 I I I I I I a 1 I *t I I I I I I I I -138-T"AP!_Ir.IT Of:1'Ll:*,Q n I "1 r:-\I<:: T m) *::: L '( ( *::i .:. )

f. IJ I JT V E "* C;::: ( <t *.1 F. <:: ( 1. l * <:: \J lv-1:: 1 ( l ) l * ( \J A "'1F 5 ( 4 7) d 1\1A*.A:'.2 ( l ) )

(4i=,)

o8rlAAA ,8HAC o8HALPHA .AHCO .RHC0N0 o8Hn .RHOC ,8HnFOX .RHnPOX l._qµn,11;) .CIHf)q .qH=-*RHFllCTOR .RHFDIV ,8HFil\!LET r::::;qµr:-r_ux

.;:i..;F*.1ULT

.,,H;-OL'.) .RHFSr=' .RHFc:;Pt_IT

.8HFXF"L('llri

': q H !'i !I P , .>:: -i *3 o. P '-J

3 *A P c:;

R H H

8 H HF I 1. "1 , 8 H i..i I N l F T 7RYHl"'LO ,K:*iI)F'tJEL ,AHIOGAO .8HTK Pr:"AL*P ,Ri..iLC QqHMrHFP

!1,Ci-iO ci .c 141_ r)RHTf"1UIYIY r,9HU<:;T FP,'-iw I';

I '°'1

P. HP'AI ,,.;.,qo.or AL * .;:i..;TI\JLET

.::, *-w *

1 * .c '"IR

0;;,,1..4TP()J oRHllP I *8HL::'JGTH , Ri"fl-.Jf YDE ,EJHO-!l .8HQC ,l:\HP-i00Lf) .SHU

FP-iV I sew oRHWSAVE

,8HNWRflO

R.HPq I 1-.ITC ,8HQF ,RHSP .AHUH ,8HVP .8HLR *

.8HPR!NTq , 8HQPP p.4 . .8HT .8HtJSAVF e8HVPA ,8HXCPOSc:;

31 /

l ' * ' ' I * * * * ' *

-139-APPENDIX 6 Subroutine Ziget I I I I I I I I I I I 'I I I I I 'I I. I I I I I. I I I I I \I, r:I I I I I I I I I I. -140-Subroutine ZIGET was developed at Savannah River atory to provide the capability for dynamic core allocation.

The subroutine is written in assembly language and is restricted.

to use on IBM Systems 360 and 370 computers.

A listing of the assembly subroutine is included in this appendix.

The use of ZIGET is described below: It is frequently desirable that the size of one or more FORTRAN arrays be established at program execution time rather than at program compile time. The FORTRAN language does not directly support this facility.

However, an array can be dimensioned at length one and dynamically "extended" using ZIGET. The call to ZIGET is programmed as CALL ZIGET(A(l), KMAX, LWORD,INCR,NPESRV,*)

where A(l) KMAX LWORD INCR NRESRV * = First element of the array to be dynamically allocated

= Number of words to be dynamically allocated

= Word-length of the array A (i.e. LWORD=4 means four-bytes per word). = Incremental address from A(l) to first word of the dynamically allocated area of core. = Number of four-byte words to be reserved for other uses. = Error return Use of ZIGET is illustrated in the coding below:

C MAIN" PROGRAM DIMENSION A( 1) READ(6,1000)

KMAX LWORD=4 NRESRV=lOOO

-141-CALL ZIGET(A(l), KMAX,LWORD,INCR,NRESRV,&901)

IF(INCR.LE.O)

GO TO 902 DO 100 I=l,KMAX 100 A(INCR+I)=O CALL ZIFREE(DATA(INCR),&903)

STOP 901 Error stop. 902 Error stop. 903 Error stop. END Failure in ZIGET. Not enough core. Failure in ZIFREE. In the example above the call to ZIGET requests that KMAX four-byte words be allocated to the array A. A minimum of 1000 words will be reserved in core for other uses. If a fatal error occurs in ZIGET return is made to statement number 901. If there is not enough core to satisfy the request, INCR is returned as a negative integer value and KMAX contains the number of words which are actually available.

If the request for core can be satisfied, INCR is returned as a positive integer address increment such that the Ith element of A can be addressed as A(INCR+I). (NB. FORTRAN statement 100.) Also illustrated is a call to entry point ZIFREE which deallocates the core previously allocated by ZIGET. I I I I I I I I I I I I I I I I I I I I I *1 I I I I I I I I I I I I 7-TGF"T !:? 11 Rl CJ? ql ::iu Re: qk C7 Qq oo Po\ qo P.r . Cf"\ ci:-Qc:-Q 1 0. cq q12 R13 RlS -142-SE*-:itJE"C:E 7IGF:T CA.1'.lf)P.1.Lt"JGTH.I5U3.*l A= ntlM"1Y nT.-1F""lc::;Tn1'1::D CALLING PQQGCA"1 i,ir, t M 5 = r 0 T .-. 1 _ D r "' F" "1 c: r o :\j c; r 7 F. o i=-a. *

= FiqST r C:IJR = P,1COF..:-AEr>.iT Al_ c::;J8SCP ror VAL JE (

l\JOR"'4Al RETURN l l="QP lq"! ISlJR *= ** -AU. 'llJ<,4q:'.;.( OF AVAILA8LE I<:;Ug = _ogg, OY\JA"ATr STJCfAGE NOT ACCES<:;IoLE IF THF:

l\JO. oi:-AP.GU"1ft\JTS AP:: 0 ASSED FR0M THE CALI !l'\IG PR!JGRAIA.

A FC)RTqAf'.1 I<=i 3IVF:N (lHC230Tl

* * * * * * * * * * * * * * * * * * * * * * * * *

r:Sf::CT l="Gi I I)

J 1 :::-:111 2 ::::-Q!J 3 :::-'.)u ,_.. :::-*1ti c; i:-11; " i:")lJ 7 ;:-Giu c:-QU C::QIJ FQU C::QU i=: ') t ! ::::-'1! I C:-)f,I :::"'JlJ C:-QU C:QIJ i:" r) I.I i=lC nc rJC SPA U=! ;; 9 l l'l 11 12 13 1 c:; 111 11 1? 11 l 'i ls.12c1c:;i

'I. I 7 f CL7 1 ZTGET 1£.*12*12(13>

12. 15 15*15 CiT U-1 T'-1 t:'.NZ '-L SLL rlPEN 7.7 7 .1_i:_:NriT"1 2.c;.r)Cll P<ll *X*8!i* PPESET [F NO. t:;.;d6CPll Rk .,z--c;ET (OIJTPl.tTl l MIN!f'.AUM LENGiH LOAD OF ARGUMENTS TEST NIJMR>:R nF APGS. OF' APGc::;. 'gc:;i Al'ICH VALUE IM ""::ms UNITS FOP DEBUGGING

---------Gn GET ALL r:? l 0

Al_LG;:'.T

--------MAKE SUP.F: CiYSTEM OVHi'l PEGlUF:ST IS J'.j A 2"\ 80lJ"l0Aqy LA c::;pl P"ioll SLL

-143-*---------COMPUTE L

SP oi::,

A c:;a.Anl')P ST qo,:>Fc::.\nf)

ST 8AL .? l 11, F;::;iE>:OVl.-ff)

SNAP

,?QATA=QEGS nF TO TO USER RAL c:;T 1:n.O(Q*=>3l

CLOSE FOR OERUGG!NG PPESET c;vE P C\ITPL P,11=" vvc !S1'J+2(2)

2fl4} 1 .* A l*PLT5T l ll.*15 I_ 7.11(3) M 6.0(4) I PP 7. 7 5T 2.oEG:::>

LA 7.15(1'),7)

SRL 7.1 SLL 7.1 ST 7.LENGT'-H:4 C::PACE 1 PICK UP I\JTERl\JAL STAT*MENT NllMRER TO POUTINE LOAD Nn. OF WOPOS PEnUESTED.

MULTIPLY 8Y LENGTH R7=

OF RYTES TO RF:QUEST "!AKE SURE IS ON 8 8YTF: 80lJNOAPY AND ADO 8 BYTES FOR CONTROL WORDS.

LENGTH (MAY) *---------FIQST.

GET ALL AVA!LIA8LE aPEA '3 A I_ R 1 0 , A l_f_ GET l.IJNTl FOU

C R7.LENA SPACE 1 COMPARE LGTH WITH REQ, LGTH CALLER TOO MUCH. GO AnJUST ------GET

(t;l)UEST 6NTO TK LA R7.2047(0.o7>

<;RL R7.11 SLL 07,11 R7.0PGLGTH A R7,AOOC( ST R7.RESA!)Q L

c:; R?.ORGLGTH C R7o7-EQO CNH C'\JT=>LINF" ST cn.Q*SLGTl-l SP.ACE l ACEA COMPUTE A9EA SAVE VALJE LOAO OF 08TA!NEn SU8 LGTH AREA CHECK SYSTEM AVAILIA8L5 IF' NOT, 5YSTEM AREA 8AL

<:;PACE l . -ao-------:-SET UP CONTROL !NF()::lMATION INTO 03TA!NEO APEA L Ri::..AQOP MVC L A "

A ( 0 * .&, l GET AOOP OF CORE ALLOCATED STORE CJNTP.OL WORDS AT ADO 8 BYTES FOR WOROc; """""'*'"""'

1*-l."""\t"'"',...rt.

I I I I I I I I \I I I I I I I I. ., I 11 ',

I ..... -.. ...... . ...... . I I I I I I \I .I I I I I I I I I *I EPR()C i::'.PRliP 1 FINIS!-' c::T C::I:! CIM xR RAL *AVC I_ A c:;T LA i:;. ;:i.

q l (I

I) TI/ rI'JI51-!

,, '"'*S) .;::cnnE lSe4(0*n>

l ":;.

7e0(0ei:;)

Q 1 l

r F"=n:r:

<;P.ACF l -144-5UBTQGCT ADDRESS OF T F " E G A TI V E

R E T 1J R l'\J

'JOE-STORAGE

,._,r}T Ar,CESS IBLE ERROR (RF.TURN ll FUU * *---------CALLER Tnn MUCH CORE, OBTAINED AREA ANn o PAc:;c:: RAC!<' VALIJI=" JF "1AX Tr-IAT CAN 8E ORTAINED FREE"44PI V*A=AD.QP.<;c:O LA SET CODE L o 1 f1

C" T 11 1 I 51-i l_ c"' .1_ F'. Ill .\ 5 D"-* 5Y5tWHn 5 AO JUST CONTRL INF"OQMATION 8NP CC:: TALL :i OT\/ SPACE 2 0:=== OQUTINF.:5 FOP ZIGET =============================8=============8=6

<;Pl\CE l TO GWT ALL AVATLABLE GETMAIN l TP c:q5,;:ns c:iz A 8 EN I) 111 ? 0 * !J '1 'Ac C::P.:'.Cf l *--------POUTPIF.:

TO F"qF'.F.:

5YSTF'."'

ECU 0 SDt.<:E 1 o--------ADJU<;T LF"N4 CONTCOL L

S P 7

P.ESUiTH ST 8R PlO *---------AfJJUST RE*1 CONTRL INFO RF:T ALI_ L p,c.

l r:r,1 A .s o---------CQMCUTE Qi:'

DTV 5ROA D n.n.(4) RR RlO c:;PACE '2 *====

POINT

=============================

CALL ZIFCEE CA,o)

A = AOOPE5<; 0F" MAIN TO FREEO * = ERQQR RETURN ZIFqF.:t:

qc OC X*7*

.;> -!>

PC"TIJRM E!1J *"I A Pl)Cf3 Pt_ r c:;T I c::; l\J LF:NGTH OOGLGTH AOOR .LENA qc:-G2 DEG3 ECOOE Tr::'.ST IPEC! ZERO 0'-IE OF:SAOI) CF.:Sl_GTH GF.:TALL SYSOIJHO DC STM nRnP l.P USTl\JG SR TM L SLL QCT CL7'7fF:JEF:

I 1? 11 ol c:;

l ic;.is !)(ll*"'q0' f:OCTC 2.n(Q.l> 1.1 6*L00P -145-SAVE RE3ISTERS.

TEST RIGHT NO. OF LOAD OF AREA TO 8E FREEn TPY R7 BY LA 15.4<'1,0l ERQOR (PF.:TURl\J l> 8 TPANSF'Fc TO F'OPTOAt-r TRACEBACK R::>UTINE MVC IC:::N+2t2l

2(l4l PTCK UP STATEMENT NUM8ER l *PLIC:::T L

8Al.q 14'15 L 3.4(0.2l R3= NO ::>F BYTES TO RELEASfD STM ?.J

FREE MA I 1\1 V

A= A 0 0 P

o: 0 oPnJ

<OUTPIJTl)

FOR SNAP DC8=9JO.OOC8.J0=2*S.JATA=(C8l

,?OATA=REGS OEBUGG CLOSE FOR DEBUGGING PUoPOSF: F.QU

f:QU

L ll..d2(13l l Rl.ADOR LM

'-1VI l?(lJJ.x*n:* 15.}4 PETURN l)C8 CIC oc !JC oc DC f"JC 115 ns ns DC DC DC DC QC DC DC DC DC oc r.ic END OONAME=SNAPCACn FOO DEBUGGING PUPPOSF A ( TSN l F 11'1' P:' n ' F*nr P:' n ' F' (l I l F lF lF FHP Ft-qqqt F' l'1 ' V (

F*O' F ' l ' F I I) ' FI 0' F

n

F' l 63'l4() ()I F'l(')240' INTEqNAi_:

STATE"!EMT JF STORAGE REQUESTED. (MIN>

OF STORAGE OF 5TOPAGE ALLOCATED JF STORAGE ALLOCATED AOOR. OF DUMMY VARIA8LECAPGll ADDRESS OF AREA TO BF RELEASFD AODR. OF DYNAMIC ALLOC ERROR CODE TEST wo=<D AOORESS OF FQRTRAN TPACEBACK ROUTINE I I I I I I II *,1 I I I I I' I I I I I I I I I I I I I I I I I I I I I I I I I -146-APPENDIX 7 Physical Properties Subroutines c 2 2 c 1 2 c 2 l="UNC".TUIN

!-il_TQ 101 NF.:w. AUGUST 1074 U=ALOG(Pl GO TO 2 LJ:IJ-7.0

... 147-1-11_!'1=( ( ( ( (

lG00*U*O.

l l490AllDOl>

l.l*0.741S344A001)

U

-tJ*O.

lJP.g1384 1)02l *U*0.37492429002)

II 2*0. l607'll5Anl')J)

PET URN ( ( ( ( ( ( (-0.477ln-o4*rJ*0.846l8f')-03l*U-O.S339?60-02l*U P.ETlJRN ein

FUNCTION ...ivA>:>(O)

M*KIN NEW.

1974 U=ALOG<Pl

!F(o.u=:.4511.0l GO T'J ? U=lJ-7.1)

HVAP:( ( ( ( (0.371704l61J0l*J-0.9lll8126)0ll*U-0.244478l002l*IJ 2*0.

l*N AP= ( ( ( ( ( ( -0. 3i:i 740-04*1!-0.

3) 0 U+ 0. 4350 7598n-02)

U i -o

l 4 s 3 so 4 on -o 1 > 'Jo'* J. o

2 2 7 7 s 91 9 o-o 1 > *t 1 RF.: TURN END SU8POUTINE MF.:i<IN NE;.i. AIJGUST 1974 X=0.001*""1 X3=X*X*X CO=l.O/C?

!F' CH-qr).0) l *2*:::?

X=X-0.25 RETURN ENO FUNCTION ROL!Q(Pl 1.1;:1<IN NE1..i. AlJGUST 1974 U=ALOGCPl GO T0 2 U=IJ-7.0 VLIO:((( ((-0.263Rl0-03*11*0.l426780-02>*U*0.212520-02)*U l

0. 119?2 70-0 2) *' J+ 0. 1974:?10-0

2) '->U + 0. 40 46Q60-0 2) *I.I ROLIO=l.O/ilLIQ VLI'1=! ( < ( ( ( co*.£i.6RD-08*11-0.7470-07>*tJ*O.J96960-06)oU QOLIQ=l.OIVLIQ RETURN ENIJ FUNCTION QOVAPtP) I I I I ,. I I I I I 1* I I I I I I I. I F'. !< I N NE 1..;

A 1 Jr, i J c; T J en 4 ll=ALOG(PJ GI") TO ? U=l.1-7. :J -148-DVG= C ( ( ( ( *1. 4 752110 l ';' 1-0. 639 l 3524n0 l l *U-0. ??430A05002 l *11 1-() * ?. 796 7 0 541"l 0 2) *' J-1). s 111('I72..j2[l0

2) OIJ-0. 6151469 l D 0 2) ou Rf)VAP=P ;Pvr-. C<F'.T!JRN CVG:( ( ( ( ( (

1-0

3 f) 71 3 qi)-n? ) J-0

f, 3 1 l ? 6 o-a 2 )

u + 0

6 0 0 0 1 6 2 9D-0 1 ) * ! I 2 + o. 11o3qJ1sno1

> *'* 1+11. 1c?.S7 40 loo 2 1 *IJ +a. J 336 o o 56G a 3 R()l/AP:P/Pifr; RETURN F.:NO 5'1TTF:*A(OJ l\Jf:\"*

1074 U=Al_OG (Pl IF(D.Lr:.i:+.sn.ni GO TO 2 U=IJ-7. I) 1 +O. l46S778Jll02l

'J*O. l::?4115875D03l 0 11+Q.5SS99496fJ03 . 5ATTEM:((I( (

l*2.3907D-3l*U

... l'l.4346l>:in-02l*U+0.1736300400l*IJ+0.2?A0814900ll*!J em c:;uPTF'.\I (P,R!_ .oG,STJ MEK IN NEi.v. AU.GUST 1q7!.. x =RL .X=0 .. 00000l'!l'l(**4 5 T=S Tl:-f,.

Rf TURN ENO

-149-APPENDIX 8 Coding of Subroutine BAROC I I I I I I I I I I I I I I I I I I I I I I I I I ., I I I I I I I c I \, I c I g. c I c c I -150-n t M N 5 ! 0 A l ( 4 l ' A 2 ( 4 l

c (1 :a a ( l 4

7 ) ' c C) :: ::" ( l 2

8 } ' n A T ( 1 2 ' 5

s ) ' 't ( c; .> 1.

GG ( 7 l , 'IQ ( l 4 l

cc < 8 l , ZNN ( 3. 6 l n T z l

2 6 2 l , o

6 7 4. q

n

o 7 3

1

gs 5 1 , 1

o o 4 3

o

1 o q 7

l 4 g A 5

n

A 4 o * .

?Q.(li:-57/

OATA n4TA GG/O.n.0.25*0.S*l.0*2.0,J.O,lOOO.O/

DATA Q010.o.o.on1.o.01.n.03s,o.os.o.01s,o.1.o.1s.n.2.

OATA l 4300.C * .'J600.0, ?

1 4

s 7

QATA l l.492.1.362.1.178.

2 4

s 6 l.J,l,J3,l.3lltla3tle3tl*3,l.304,l.308,l.284,l.2Atl*2*l*l*

9 4

r G

- ** 1s ** ** Sn1 ** s12 ** ** 59,.60S **

- 8a,

- 7 ** 7Ql ** .J

- 642,.587,.540,.493

- 454,.454 **

OATA naTA

= l { (Y2-YY)*(Zll*<X2-XXl

-C'fl-YYl*<Zl2*CX2-X'()

Z22-i>(XX-Xlll l 1 /((Yl-Y2l""CX1-X2l) V!LUE OF YB AT XS, LINEARLY 8ETWFEN CXA.YA) A "-10 ( X C

Y C l IS OF Z AT rXX,YY>*

INTERPOLATED 8ETWEEN Zll AT (Xl.Yll* ZJ2 AT cx1.y2), Z21 AT cx2.Yll ANO Z2? AT cx2.Y2) !CAPT = 1, ANO ARRAY CORAg !CART = 2.

WITW MASS VELOCITY AND INTEqPOLATE IN C0RA8 TO OBTAIN IF' (!CART GO TO 41 SET PROPERTY FROM PRESSURE.

!F'((P.LT.ll.42ql.OR.(P.r.T.3204.0ll lOOl*P !F'(P.GT.1429.Sl GO c YY=Al(£+J no 2 I=l*J L=4-I ?

Ill ::1'.( = yy TJ l?. q C01<.JT I rin 10 I=l*l L=t:..-I 10 YY=YY*C/3?n4*A?(L)

OX : oor = AL0G(PX) 1, CO"'!Til\JUE i'.MA:X=l4 IF(OX.LT.PO(l))

>X = PPCl> J=l iu. ri=-(ox.u

.Po1Jll c;o ro ,J=J

l GO TO lt.. -151-c c:;ET MULTIOL[=:;i AT G = 1.0 c c t6 22 I=l.TvAx IFCI.EQ.ll I F ( ( I * ;: Q

l l

0 R * ( I

E Q

I 1-l l G 0 T 0 2 2 "4=!-1 IF"(J.GT.2J GO T'1 15 1..J V = 7 L I 1'1 i: ( AL() G ( Po Cl l l , Al(') G ( C :>EI=' ( , ll l ' A. l 0 G ( PP ( 2 J l , Al 0 G ( C 0 E F ( M , 2 J l , l PP I)

GO TO Jc:;

GO Tll 17 !FCCJ.LT.41 (J.GT.511 GO TO 17 GO TO 19 17 IF (J.LF.7> GO T:J 18 t./V = ZLINEIAL0G(OP(7)

J eAL0GCCOEF'P."'*7l

>, O.Q,Q.Q,PoIJ COPA8(T.4l:fXP(WVJ SO TO ?2 !'3 J=2 lO 22 30 32 34 Zt\11 = .ALOGC

-l.O

QQ(!) l*PP(J-1>

)/AL_OG<QQ(!J)

ZN2 = AL0G((COFF(M,J l -l.O

QQ(Ill*PP(j ZN=

l *ALOG(ZNl> .ALvG(PP(J)) ,PPJ) Z"'-1 = F.:XP ( z,,) r. f) p A 8 ( T

4 l = 1

0 -Q I) ( T ) + ( Q Q ( I ) ...

z \j ) Ip '-<. CONT I"!tJE I SET A3 '"'1.0. nn x tJS PJ<; l..1 VELOCITY I ON F 1'. l\ff) 1 = 1

I) Fl!T=O.lc::;

G() TO 32 I NO l =I \lO 1

1 t;Q TO 30 P-102=0. 0 DO 34 K':2.4 T I=' < C PP T

GT * ( X ( I<'. l -'3 I T l l

A \I 0 * ( o D I

L T * ( X ( K > + 8 I Tl l l I "' 0 2 = K f)O 38 I=l-!MAX "l=I-1 I I I I I I I I I I I I I I I I I I I I I I. I I I I I "I ., I I I <: c -c -152-I i:-( ( I * '1

l l

A '*Jn

C J

LT

7 l l .; 0 T 0 3 5 IF ( ( !. ;:_: 11 .r ... ox )

A \In. ( J. '** T. 7) ) r;n TI') J 3 '-1= J-1 TF(..J.i:-,1.ll

.. =J Tl:"(..J.F*l.7l GO T1 .17 y y = 7 L Ti'J i:: ( x ( r !\j I! 1_ )

DA T ( "'

T \j D l ' ""'

x ( r '.\J f) l + l )

I) A T ( N

I \In l + 1

M l. l

oi:i I l ! F C P.I

E *1

I'\

0 l G 0 T n 3 X l = 'i ( P.1 0 2 ) -):I I T 't 2=X ( PJf)2 l +;;IT Y l = 7 L I l\J ( X C T 'Ir/ 2-*1 l

0 T ( l\I

I *'-J) 2-l

l

X ( I \J 0 2 l , [')AT ( N , T \J 0 2

l , X l l V?.=7LPJE(:.<(T'rn2l oDAT('J.T\J02oM) oX{I!\Jfl2+1> .OAT(NdNl'.'2+1,M)

,'t2) y Y = n , c; * ( Z L T 1\1 E C x l , Y l , x? , Y 2

P;:;) I l + Y Y l 'iO T() JS YY=l.O 1A t;IJ TO 37 COPA8rI.Jl=l.O CrlN!

RF::TUP"t IN aooAY ro F"INn MJLTIPLIEP.

41

!F<G.GF.lonn.01G

= inno.o T"-101=1 42 l GO T0 I"IDl=INOl*l GO TO 4? 44 <:ONTil\JUE T"IC2= I 48 !'ID2=I"lf!2+

l r:n TO 4.6 G?=GG (

l G l =GG C T'-102-1 l Gl=G T F ( G

I_ E

l , r, ) G Ii TI') c:; n r..1 = l

n /i. l . (';2=1. r.:i /(i2 G"3=1.Q/G3 1".:0NTINUE 211 =

112 =

> 7.21 = CORAl=lCINnl ,p.Jl');:>-ll Z22 = ) Xl = QQCI\IGl-ll X2 = (jl)(I"'1'1l xx = FMULT =

PE'T\IRN 1!')I')1 F"QPMAT ('

= I.

1 oJTsroE vALIO oANGE OF 11.43 ro l 1;:!')4 OSIA*) E"t0

-153-I I APPENDIX 9 I SEPRAT Coding I I I I I I I I I I I I I I I I I IC r, r. I i: c I I le I I I c I I I I c I -154-CHANNELS (EG 8WPl CALLED SP usi:o F00 (l> :JM;l)P c2i J..., (3> o? IMPt_ICIT

!t-.JTEGl=.:::J

('5) 1cn=loA1.1 ARETA .AFLux ,ATOTA*-:*8BETA ,ou ,or .ox , I t:LEI/ .Ge *GK .GRID ,Hc;uRF tHF' ' HFG ,12 .IJ

ITERAT,Jl ,J2 , 1 J3 ,J4 .JS ,J6 ,J7 tKOEBUGtKF' ,KtJ , 4 ,NR8C *NCHANLtNCHF tNDX ,NF , 5 NGAPS .MGPin tNODESFtNPPOP , ,l\IPQI) .NSCC1C ,NVISCNtPI ,POWER ,PREF , 7 ,SL tTF tTFLU!DtTHETA tTHICK , 11F .vF .vFr, .vG .z CQMMQN /C()QoA?/

AA(4). AF(7), AFACT(l0,10), Al/(7), AX!AltJQ), 1 tiXL<l'1l.

iqc4l

RX()r.l, CCC4l,

OFUFLf2>, '.2 GAPXLClti), c;FACrcg.iri>.

GRIDXL(l()), HGAPC2). H!-IFCJO)t HHG(JQ), 1 IGPIDCll'll

K'CLAf)(2l, -<FUEL<2>, ><"¥.'"C30l

NCHClOl, 4 r:>P(30>,

TCLADC2l, UUF'C30lt S VVF(JOl, VVGC30l * ;u)1JAL(30), Y(JOl

TTC30) LOGICAL Kf.J, KF. KKF, KCLAD,

/CQgPA]/ MA ,MG ,MN ,MR ,MX , l

,$A .SAAA .$AC

,$AMSWF,$8 , l .ICON .§CONO .scP .so ,$OC

' <:;0!40X ,q;f")J..iYI")

,$DPOX ,$QPK t'bDUR ,$(')R ,SF

1

.5FtNLE.iFLUX t5FSP 4 .BH

,$HPEQT,$J(')ARE, t; *;;roF1_1E,c;toG-oC,'£IK

,$Jr30ILe1iJK tSLC ,'l;LENr,T,$L(')CA

,$LR , 7 $PHI t$PWRF ,$QC ,$QF R $0UAL 9 g <fiU ,1.tlJ!"'f .svrsc t$V!SCWt$VP

,$VPA ' A ,,WnLn ,$WP ,$WSAVE.SX t$XCROS DATACll LOGICAL L')ATCll PHEGER Il)t.T Cl l EQUTVO.LENCF.:

COATACll ,[nATCl> *LDAT<ll l !F' SO TO 111 on 2 I= 1 'NC:i-IAl\jl_

OATAC50FDX*tl=O.O .

1 CALL DIFFE:::>cJ,Jl RF.:TIJRN ll'l PMIN=looooo.o PMAX=-1000.0 ) .

r) 0 l 2 I = l * "IC' i-1 A l<J '11 V =OAT A (<l5P*T) IF

JMIN=WV IF (*.1'V.GT.OMAX)

JMA'(:WV

-155-12 CONTINUE" RETURN Jll'-1P= l IFCifERAT.GT.ll 30 Tn FT0T=O.O DO 14 I=l*"ICHMJL f)ATACf:SO+IJ

= '1.l*DATACc;;SP+Il I F (Q T A ( c; 5 P + I I * 'IE

0

Ii l 3 0 i 0 1 4 l

14 GO TO 20 QO 18 lR 1 2*M() ) 20 5!.J"'llJ=O.O

!J 0 2 2 I = l

11..1 <: !-! A l L 22 SUMF=n.o 00 24

?4 00 26 I= 1 * !'!Ci-IANL P*TURN ENO I I I I I SEPRATI sE=RATI I I SEPRATI I SEPRAT I I I I I I I I _J I I I I I I I I I I I I I I I I I I I -156-APPENDIX 10 Input Data Preservation Based on that of COBRA IIIC

-157-card(s) Type Cl Required *to be present.*

Problem Array Size Always FORTRAN READ list:

FORTRAN FORMAT: Read. from Subroutine:

Variable Columns MC 1-S MG 6-10 MN 11-lS MR 16-20 MX 21-2S Notes MC MG MN MR MX 10!5 INDAT Format Description CG IS IS -IS IS I.S > No.

blem. cards COBRA of channels (NCHANL) in NCHANL is set from NTHBOX on CS-C7, or in the original format, in Card Group 4. No. of gap interconnections (NK) between channels in problem. If this is not known, MG=Z*MC is ally adequate but should be checked later. For a BWR, MG may be given as zero, when it is reset to 1 in CORE. No. of fuel nodal points in problem. This should be on Card Tl. If MN is given as zero, it is reset to 1 in CORE. No. of rods (NRCD) in problem. For PWR and BWR, NROD=NCHANL, hence MR may be given=MC. No. of axial stations in problem. It may be given as NDX (Card Cll) as it is increased by 1 immediately after reading in. (1) MC to MX are used to set the array sizes in the dynamic r,ge, hence they should be set too big r*ther than too small. (2) Note that MC to MX are given in alphabetical order. I I I I I I I I I J I I I I I I I I I I I .. ' I I I I I I I I I I 1-1 I I -_) I I Card(s) Type C2 Required *to b*e present FORTRk.'1 READ list: FORTRAN FORMAT: Read. from Subroutine:

-158-Maximum Running Time Always MAXT IS, 6E12.6 INDAT Variable Columns Format Description MAXT 1-5 IS Maximum Running Time, Nominal value is 2000. CG c o* N T R 0 L

I O'\ l.f\ r-1 THE INPUT FOR A CASE REQUIRES A CASE CONTROL CARD FOLLOWED wrrrH UP TO 12 GROUPS OF INPUT INFORMATION.

EACH OF THE 12 CARD GROUPS HAS A GROUP CONTROL CARD THAT IDEWrIFIES THE GROUP NUMBER AND THE OPTIONS J\VJ\ILJ\BLE FOR THAT GROUP. 1 GO TO THE CARD GROUP SPECIFIED BY NGROUP, IF THE DATA OF A CARD GROUP THE SAME AS THE PREVIOUS CASE, THEN THAT CARD GROUP AND ITS CONTROL MAY BE OMITTED.

I I I I I I I I I I I I ,. I I I I I I Card C3 Cards (s) Type C3 Required to be present FORTRAH READ list FORTRAN FORMAT Read from subroutines Variables I PILE KASE Jl TEXT Format Il I4 I5 I7A4 -160-Case Control Card Always IPILE, KASE, Jl, TEXT Il, I4, IS, 17A4 I HD AT Columns Description 1 2 -5 6 -10 il -78 = 0 Run identification number. If > o, calculation continues O, calculation stops. Printing option for standard COBRA output -as in COBRA III-C. = 1 print entire output = 2 print only operation conditions Alphanumeric information to identify case

-161-Card Group 1 Required to be present FORTRAll READ list FORTRAN FORMAT Read from subroutine Variable Columns Format NG ROUP 1-5 I-5 6-10 I-5 ' l Always HGROUP lH I5, I5 IHDAT < 0 > 1 Description

= 1 (to select Card Group 1) calculate physical properties from polynomials I I I I I I I I I I I the physical properties' are given in the next J Cards as in the original COBRA. I I I I I I. I I I I I i I I I I I I I I I I I I I I -162-Physical Properties Required to be present when Nl(in Card Group 1) < 0 FORTRAN READ list H PH P2 iH FORTRAN FORMAT I5 Fl0.3 FlO. 3 I5 READ from subroutine CARDS l Variable Format Columns Descrintion H I5 1-5 = 1: PH defined as lm*iest pressure encountered in problem. = 2 : PH defined as lowest enthalpy encountered in prob len PH Fl0.3 s-15 Lowest pressure (psia) if Nl = 1 or lowest enthalpy (Btu/lb) if rn = 2. P2 Fl0.3 16-25 Highest pressure (psia) encou.'1tered in problem. Hl I5 26-30 Number of pressure steps generated by polynomial (maxi.mum

30) lfotes:

-163-The lowest pressure encountered in the problem is defined as that at which the lowest enthalpy would be the saturation value. For example, at 1000 psia the saturation enthalpy is 543 Btu/lb *. At an inlet subcooling of 100 Btu/lb, the enthalpy would be 443 Btu/lb and this would be the saturation value at a pressure of about 470 psia. Thus, one would require physical property data over the range 470 (or less) psia to 1000 psia in order to include data which covered the enthalpy ra...'1.ge.

To avoid translating the lowest enthalpy to pressure, the I I I I I I I I option of giving the enthalpy is included.

The program trans-I lates this value to a pressure which is safely below that required using the expression p = 6h 3 (h-1.35) I (h -0.35) when p =calculated pressure (psia), h = O.OlH, H =enthalpy (Btu/lb).

The values of p, so calculated, are given below and it may be seen that they are all less than P sat, the tabled value of pressure corresponding to H. H(Btu/lb) 181.2 300 400 500 600 700 p (psi a) 11 101 279 589 1067 1749 Psat (psi a) 15 103 319 745 1409 2236 I I I I I I I I I I I I I I I I I I I I I I I I I I I I I -164-In the original COBRA, the physical properties are read from cards into the arrays (PP(L), TT(L), etc., L = 1, In the new version, the values of (PP(L), TT(L), etc., L = 1, Nl) are generated within a Do Loop from 1 to Nl from the physical property polynomials.

With the arrays set, the subsequent use of the values is the same in both versions of the code. Note: NPROP is set to Nl for storage of the size of the arrays.

-165-Physical Properties Required to be present: :*Then :n (in Card Group 1) * > O Read from Subroutine CARDS 1 ----.---* READ IH Nl CARDS OF FLUID PROPERTY DATA. EACH CARD CONTAD-13

--SATURATIO:l PRESSURE (PSIA), TEHPERATURE(DEG-F)

LIQUID SPECIFIC VOLUME (CU-FT/LB), VAPOR SPECIFIC VOLUT.fil (CU-FT/LB)

LIQUID EHTHALPY(BTU/LB), VAPOR E1J':L'HALPY(BTU/LB), LI:':lUID VISCOSITY (LB/FT-HR), LIQUID THERI:IAL CONDUCTIVITY(BTU/HR-FT-F)

AHD SURFACE TEHSION(LB/FT), FORMAT(E5.2,F5.l,7Fl0.0).

iH MUST BE GREATER ':1HA.N ONE BUT NOT GREATER THA1l THE PARArETER MP. THIS PROPERTY TABLE I1UST HAVE PRESSURE HIGHER THAN OPERATING PRESS. AND LIQUID ENTHALPY LOWER THAN THE BUIJDLE IULET ENTllALPY.

I I I I I I I I I I I I I I I I I I I

CARD GROUP 2, FLOW CORRELA'1 1 IONS READ IN UP TO FOUR SETS OF FRICTION FACTOR COHRELJ\.TION CONSTANTS THAT CORRESPOND TO THE SUBCHANNEL TYPES, FORMAT(l2F5.3).

Nl IS THE SUBCOOLED VOID CORRELATION OPTION. Nl=O, NO SUBCOOLED VOIDS. Nl=l, LEVY SUBCOOLED VOID CORRELATION.

N2 IS THE BULK VOID CORRELATION OPTION. N2=0, HOMOGENEOUS MODEL. N2 = 1, MODIFIED ARMAND MODEL. N2 =5, READ IN SLIP RATIO, FORMAT (5X,El0.5).

N2=6, READ IN THE NUMBER OF TERMS AND COEFFICIENTS FOR UP TO A SIXTH ORDER POLYNOMIAL FUNCTION OF STEAM QUALITY, FORMAT (I5,7El0.5).

N3 IS THE TWO-PHASE FRICTION GRADIENT MULTIPLIIm OPTION. N3=0, HOMOGENEOUS.

N3=1, ARMAND. N3=5, READ IN NUMBER OF TERMS AND COEFii'ICIENTS J1 1 0R UP TO A SIXTH ORDER POLYNOMIAL fo'UNCTION Qli' QUALITY FORMAT(I5,7El0.5).

Nit IS AN OPTION TO INCLUDE A WALL VISCOSITY CORRECTION TO THE FRICTION FACTOR. IF Nll=l, IT IS INCLUDED, OTHERWISE I'l 1 IS NOT. CARD GROUP 3, AXIAL HEAT FLUX TABLE READ IN Nl PAIR OF DATA FOR THE TABLE. EACH PAIR CONSISTS OF TIIE RELATIVE POSITION (X/L) AND THE CORRESPONDING RELATIVE HEAT FLUX (LOCAL FLUX/AVERAGE FLUX). EACII CARD ACCEPTS UP TO SIX PAIR OF DATA! FORMAT(l2F5.3).

Nl MUST BE GREATER THAN ONE BUT NOT GREATER THAN THE PARAMETER MP. I I-' CJ\ CJ\ I CARD GROUP 11, SUBCHANNEL LAYOUT AND DIMENSIONS READ IN Nl CARDS OF SUBCHANNEL DATA CORRESPONDING TO THOSE SUBCHANNEL FOR WHICH DATA ARE BEING SUPPLIED.

N2 IS THE TOTAL NUMBER OF CHANNELS.

FOR EACH OF THE Nl CARDS, READ IN THE SUBCHANNEL TYPE NUMBER (IF BLANK, IT IS ASSUMMED TYPE 1), SUBCHANNEL IDENTIFICATION NUMBER, NOMINAL FLOW AREA(SQ-IN.), WETTED PERIMETER (IN.), HEATED PERIMETER(IN.)

AND UP TO FOUR SETS OF ADJACENT SUBCHANNEL CONNECTING Ilfli'ORMJ\TION, EACH SET OF CONNECTING INFORMATION INCLUDES THE ADJACENT SUBCHANNEL NUMBER (NEGATIVE IF A LINE OF SYMMETHY SPLITS A GAP AT A BOUNDARY), NOMINAL GAP SPACING AND CENTROID-TO-CENTROID DISTANCE(IN.).

IF SUBCHANNELS ARE INPUT IN ASCENDING ORDER, THEN ONLY HIGHER NUMBER SUBCHANNELS NEED TO BE IDENTili'IED AS CONNECTIONS.

CENTROID DISTANCES ARE NOT REQUIRED IF THEY ARE NOT USED IN THE MIXING COHRELATIONS.

N2 MUST BE GREATER THAN mm BUT NOT GREATER THAN THE PARAMETER MC. I I-' 0\ --1 I ---------*----------

CARD GROUP 5, SUBCHANNEL AREA VARIATION TABLE IF THERE ARE NO AREA VARIATIONS, OMIT THIS CARD GROUP. *READ N2 VALUES OF RELATIVE LOCATION(X/L)

WHERE AREA FACTORS ARE GIVEN N2 MUST BE GREATER THAN ONE BUT NOT GREATER THAN THE PARAMETER ML. READ Nl SETS OF AREA VARIATION FACTORS (LOCAL AREA/NOMINAL AREA). EACH SET CONSISTS OF SUBCHANNEL NUMBER AND N2 AREA VARIATION FACTORS, FORMAT(I5/(12F5.3)).

Nl IS LIMITED BY THE PARAMETER MA. IF Nl IS ZERO, AREA VARIATIONS ARE DELETED FOR SUCCEEDING CASES. N3 IS THE NUMBER OF ITERATIONS FOR INSERTING AREA VARATIONS.

IF N3 IS ZERO OR BLANK, N3 IS SET EQUAL TO 1. CARD GROUP 6, GAP SPACING VARIATION TABLE IF THERE ARE NO GAP VARIATIONS, OMIT IJ.'HIS CARE GROUP. READ N2 VALUES OF THE RELATIVE LOCATION(X/L)

WHERE GAP FACTORS ARE GIVEN, FORMAT(l2F5.3).

N2 MUST BE GREATER THAN ONE DUT NOT GREATER THAN THE PARAMETER ML. READ Nl SETS OF GAP SPACING FACTORS(LOCAL GAP/NOMINAL GAP). EACH SET CONSISTS OF THE ADJACENT SUBCIIANNEL NUMBERS Ji'OR THE GAP N2 GAP VAHIATION FACTORS, FORMAT(2I5/(12F5.3)).

IU IS LIMITED BY PARAMETER MS. Il" Nl IS ZERO, GAP VARIATIONS ARE DELETED FOR SUCCEEDING CASES. I I-' ()'\ co I CARD GROUP 7, SPACER DATA Ili' Nl=l, WIRE WRAP FORCED DIVERSION CROSSFLOW MIXING IS INCLUDED, OTHERWISE, IT IS OMITTED. READ ONE CARD CONTAINING THE WIRE WRAP PITCH (IN.), PIN DIAMETER AND tHRE DIAMETER (IN.), FORMAT ( 8El0. 5). IF Nl=l, N5 IS AN OPTION TO SAVE OR USE A PREVIOUSLY COMPUTED CROSSFLOW SOLUTION.

THE FLOW CONDITION MUST NOT CHANGE FOR THESE CASES NOR THE BASIC PROBLEM SETUP. THIS OPTION WOULD NORMALLY BE USED FOR CASES INVOLVING CHANGES IN POHER OR MIXING FOR NONBOILING PROBLEMS.

N5=0, CROSSFLOW SOLUTION IS COMPUTED FOR EACH CASE. N5=1, USE FIRST CASE SOLUTION FOR ALL SUCCEEDING CASES. N5=2, WRITE SOLUTION TO TAPE AND USE FOR SUCCEEDING CASES. N5=3, READ SOLUTION FROM TAPE AND USE FOR SUCCEEDING CASES. FOR EACH GAP, READ A CARD CONTAINING THE GAP NUMBER, THE EFFECTIVE FRACTION OF A PITCH FOR FORCING CROSSFLOW AND UP TO SIX RELATIVE PI'rcn LENGTHS IDENTili'YING THE LOCATION OF WRAPS CROSSING THROUGH 1\ GAP USING A POSITIVE VALUE FOR WRAPS CROSSIHG FHOM I TO J AND A NEGATIVE VALUE FOR CROSSINGS FROM J TO I WHERE I IS LESS THAN J. THE GAP NUMUERS ARE ASSIGNED IN THE ORDER THAT SUBCHANNEL PAIRS ARE IDENTIFIED IN CARD GROUP 11. -------------------

READ IN THE NUMBER OF WRAPS CONTAINED IN EACH SUBCIIANNEL AT THE START OF rr1m BUNDLE IN ASCENDING SUBCIIANNEL ORDERJ FORMNr( 10I5). USE ENOUGH CARDS TO SPECIFY ENTIRE WRAP INVENTORY.

IF Nl=2, SPACER PRESSURE LOSSES AND FORCED FLOH DIVERSION ARE INCLUDED OTHERWISE)

THEY ARE OMITTED. N2 IS THE TOTAL NUMBER OF SPACER LOCATIONS.

N3 IS TIIE NUMBER OF SPACER TYPES. Nlt IS THE NUMBER OF ITERATIONS TO INSER'l' LOSS COEFFICIEHTS OR THE WIREt*lRAP MIXING. IF N4 IS BLANK OR ZERO, ONE IS USED. READ N2 RELATIVE LOCATIONS(X/L)

WHERE SPACERS ARE LOCATED AND THE TYPE OF SPACER AT THAT LOCATION)

FORMAT(6(F5.2,I5)).

READ N3 SETS OF DATA CORRESPONDING TO EACH SPACER TYPE. EACH SET CONSISTS OF A CARD FOR EVERY SUBCIIJ\NNEL.

ON EACH CARD IS THE SUBCH NUMBER, SPACER LOSS COEFFICIENT, CONNECTION NUMBER OF GAP THROUGH WHICH FLOW IS FORCED, AND FRACTION OF FLOW DIVERTED, FORMAT(2(I5,E5.0))

IF 'l'lIE CONNECTION NUMBER IS ZERO AND, TIIE FLOW rmACTION IS ZERO, THEN rl'lIERE IS NO FORCED FLO\'! DIVERSION.

THE FORCED CROSSFLOH HAS THE SAME SIGN AS FORCED FLOW FRACTION.

I I-' -....:i 0 I CARD GROUP 8, ROD LAYOUT, DIMENSIONS AND POWER FACTORS READ IN Nl CARDS OF ROD LAYOUT DATA CORRESPONDING TO THOSE RODS FOR WHICH DATA ARE BEING SMPPLIED.

N2 IS THE TOTAL NUf1DER OF RODS. FOR EACH OF THE Nl CARDS, READ THE ROD TYPE, NUMBER, DIA.(IN.), RELATIVE ROD POWER(ROD POWER/AVERAGE HOD POllER) AND UP TO SIX SETS DATA FOR ROD-TO-SUBCHANNEL cmINECTIONS, Ji'ORMAT(Il,Ilt,IE5.2,6(I5,E5.0))

NUMBER AND FHACTION OF THE ROD POWER TO THAT SUBCHANNEL.

THE NUMBER OF FUEL ROD ARE PRESENTLY LIMlTED TO 2. H=l INDICATES ROD FUEL. N=2 INDICATES PLATE FUEL. IN EACH CASE FOR PLATE FUEL THE non DIAr.filTER( ABOVE) IS THE PLATE THICKNESS AND THE FRACTION OF POWER TO A CHANNEL IS THE FRACTION OF THE CIRCUMFERENCE REQUIRED TO SPECIFY THE PLATE WIDTH FACDIG THE SUI1CllANNEL.

N2 IS LIMITED BY THE PARAMETER MR. N3 IS rr1m NUMBER OF RADIAL FUEL NODES INCLUDING THE CLADDING.

N4 IS THE TOTAL NUMBER OF FUEL TYPES. FOR EACH FUEL TYPE, READ IN ON ONE CARD 1 THE THF.RMAL CONDUCTIVITY (B/HR-FT-F)

'* SPECIFIC HEAT (B/LB-F), DENSITY(LB(FT3), AND PELLET DIAMETER(IN.)

FOR THE FUEL 1 AND THE SAME FOR THE CLADDING EXCEPrr FOR THICKNESS (I AND r_r1rn GAP COEFFICIENT(B/HR-FT]-F).

THESE ARE ASSUMED CONSTANT.

N5 IS AN OPTION TO SELECT A CRITICAL HE./\T FLUX CORRELATION.

IF N5=0, NO CIIF CALCULATIONS ARE PERFORMED.

IF N5=1, THE BJ\W-2 CORRELATION IS USED. IF N5=2 1 THE W-3 CORRELATION IS USED AND THE USEH SHOULD VALIDATE THE TDC VALUE IN SUBROU'rINE CHF. I I--' -.:i I--' I --

W.1sr-MI ------

CARD GRQUP 9, CALCULATION VARIABLES nEAD IN DIVERSION CROSSFLOW RESISTANCE FACTOn, TURBULENT MOMENTUM FACTOR, BUNDLE LENGTH(IN.), POSITION FROM VERTICAL(DEGREES), NUMBER OF AXIAL NODES, NUMBEH OF TIME STEPS, rrOTAL TRANSIENT TIME(SECONDS)

MAXIMUM NUMBER OF ITERATIONS, ALLOWABLE FRACTION ERROR IN FLOW FORMAT CONVERGENCE AND TRANSVERSE MOMENTUM PARAMETER(S/L), FOHMAT(4E5.2,2I5,E5.2,I5,lJE5.2).

IF THE NUMBER OF ITERA'.rIONS, ALLOWABLE ERROR AND MOMENTUM PARAMETER ARE BLANK on ZERO, THE PROGRAM USES 20., l.E-3, AND .5, RESPECTIVELY.

Nl IS AN OPTION GIVING THE SPATIAL P'RINTINn INCREf.1ENT.

IF Nl=l, srrEP IS PRINTED. IF N2=-2, EVEnY OTHER STEP IS PHINTED, ETC. IF ZERO OR BLANK, rrHE PH.OGRAM SETS tll= 1. N2 IS AN OPTION GIVING 'J.1IIE TIME PHINTING INCREMENT AND IS SET UP SAME AS Ul ABOVE. N3 IS A DEBUG PRINT OPTION. IF N3=0, NO DEBUG IlHi'ORMATION IS PRINT IF N3=1 A DEBUG PRINT IS MADE FOR EArnI STEP OF Tiffi CALCULATION.

IT CAN GENEHA'1 1E A LOT OF PAPER SO IT IS NOT NORMALLY USED. I 1--' -..:i f\) I CARD GROUP 10, TURDULENT MIXING CORRELJ\:TIONS Nl IS THE OPTION FOR SUBCOOLED MIXING .CORRELATIONS.

FOR ANY SELECTED VALUE OF Nl READ IN THE CONSTANTS A AND B, FORMAT(2F5.3).

rr1m OPTIONS ARE --Nl=O, W/GS=A Nl=l, W/GS=A*RE**B Nl=2, W/GD=A*RE**B Nl=3, W/GS=D/ZIJ*A*RE**I3 UOTE THAT BETA = W/GS WHERE W IS THE TURBULEUT CROSSFLOW.

N2 IS THE OPTION FOR TUO-PIIASE MIXING. IF N2=1, Tt*IO-PIIASE MIXING IS THE SAME AS FOR SUBCOOLED CONDITIONS.

IF N2 IS GREATER THAN ONE READ IN N2 PAIR OF DATA FOR A TABLE OF TWO-PHASE MIXING DATA. EACH PAIH CONSISTS 01" THE STEAM QUALITY AND THE CORRESPONDING VALUE OF BETA. N2 IS LIMITED I3Y THE PARAMETER MP. N3 IS THE OPTION FOR THERMAL CONDUCTION MIXING. IF U3=0, NO THEHMA CONDUCTION.

IF N3=1, READ IN THE THERMAL CONDUCTION GEOMETRY FACTOR FORMAT (F5. 3). I I-' w I -------------------

CARD GROUP 11, OPERATING CONDITIONS READ IN THE OPERATING PRESSURE(PSIA), INLET ENTHALPY(BTU/LB)

' on IHLET TEMPERATURE(DEG-F)

J MASS VELOCITY(M-LD/HR-S.Q-FT)

AUD AVERAGE HEAT FLUX(M-BTU/HR-SQ-FT).

(6Fl0.0) Nl IS THE INLET ENTHALPY OPTION. IF Nl=O, INLET ENTHALPY IS GIVEN. IF Nl=l, INLET TEMPERATURE IS GIVEN. IF Nl=2, READ IN THE INDIVIDUAL SUilCHJ\.NNEL INLET ENTHALPIES, FORMArr( 12E5. 0). Ili' Nl= 3, HEAD IN THE INDIVIDUAL SUBCHANNEL INLET TEMPER/\.TURES, Ji'OHMAT(l2E5.

N2 IS THE INLET FLOW DISTHIBUTIOlJ OPTION. IF N2=0, THE SUBCHANNELS ARE GIVEN rrHE SAME MASS VELOCITY.

IF U2=1, THE INLET FLOW IS DIVID TO GIVE EQUAL PRESSURE GRADIEWr IN THE SUBCIIANNELS.

IF N2=2, READ MASS VELOCITY FACTORS FOR EACH SUTICIIANNEL, FORMAT(l2E5.0).

N3' N lj, N5 and N6 ARE OPTIONS J.i'QR TRANSIEWr FORCING FUNCTIONS.

IF ANY OI*' THESE OP'rION NUMBERS ARE ZERO OR BLANK, THE CORRESPOli!DING FORCING DATA IS NOT READ AND IS EXCLUDED FROM THE CALCULATIONS.

EACH OF THESE NUMBERS GIVE THE NUMBER OF PAIRS OF rrAJ3ULAR DNfA TO BE RE FOR EACH FUNCTION.

ALL DATA ARE READ AS OF TIME(SECONDS)

AND RELATIVE VALUE, FORMAT ( 12E5. 0)

N3 IS THE OPTION FOR REFF.HENCE PRESSURE VERSUS TIME. N lj IS THE OPTIOU FOR INLET ENTHALPY OR TEMPEHN.rURE AS A Jl'lJNCTION OP TIME DEPENDING ON THE OPTION FOR INLET ENTHALPY OR TEMPERATURE.

N5 IS TIIE OPTION FOR INLET FLOW VERSUS TIME. H6 IS rrim OPTION FOR HEAT FLUX VERSUS TIME. I I-' --1 .i:=-1 CARD GROUP 12, OUTPUT DISPLAY OPTIONS Nl IS AN OPTION FOR PRINTING ANSWERS. Nl=O, PRINT SUBCHANNEL DATA ONLY. Nl=l, PRIN'l' SUBCHANNEL DATA AND CROSSFLOWS.

Nl=2, PRINT SUBCHJ\NNEL DATA AND FUEL TEMPERATURES.

Nl=3, PRINT SUBCHANNEL DATA, CROSSFLOWS AND FUEL TEMPERATURES.

N2 IS AN OPTION FOR SUBCHANNEL DATA PRINTOUT.

IF N2=0, ALL SUBCHANNEL DATA ARE PRINTED. IF IT IS CALLED FOR BY Nl. FOR N2 GREATER THAN Z HEAD IN THE SUBCHANNEL NUMBERS FOR HHICH RESULTS ARE TO BE PRINTED FORMAT( 36I2). N3 IS AN OPTION FOR FUEL TEMPERATURE PRINTOUT.

IF N3=0, DATA FOR ALL RODS ARE PRIN'rED IF CALLED FOR BY Nl. FOR N 3 GREATER THAN ZERO, READ IN N3 ROD NUMBERS FOR 1HICII TEMPERATURES ARE TO BE PRINTED, FORMA'r(36I2).

IF CHF DATA IS CALLED FOR BY INPUT OPTION IT IS PRINTED FOR EACH SELECTED ROD PLUS A

SUMMARY

TO IDENTIFY THE ROD AND CHANNEL WITH THE MINIMUM CIIF RATIO. Nll IS AN OPTION FOR FUEL NODE PRINTOUT.

IF N'l=O, TEMPERATURES ARE PRINTED li'OR EVERY NODE. FOR N lJ GREATER THAN ZERO, READ IN N 11 NODE NUMBERS '1 1 0 BE PRINTED, FORMAT ( 36I2). TO ST AR'r A CALCULATION, READ A BLANK GROUP CONTROL CARD. TO STOP 'rHE CALCULATIONS, AFTER FINISHING A CASE, READ A I3LANK CASE * * *

END OF INPUT INSTRUCTIONS

* *

UNI'l 1 S -J\LL COMPUTATIONS ARE DONE USING FT, LB, SEC, BTU AND DEG-F. I I-' --.1 \Jl I -.. -E -, ------

I I I I I I I I I I . I I I I I I I I I Simplified used for -176-APPENDIX 11 COBRA IIIC Input Data Presentation to be Assembly to Assembly Analysis of LWR

,. -177-Card(s) Type Cl Problem Array Size Required to be present* Always FORTR!> ..... '1 READ list: FORTRAN FORMAT: Read. from Subroutine:

Variable Columns MC 1-5 MG 6-10 MN 11-15 MR 16-20 21-25 Notes MC MG MN MR MX 10!5 INDAT Format' IS IS -IS IS .. IS N . ,. *. 0.

blem. cards COBRA Description of channels (NCHANL) in NCHANL is set from NTHBOX on CS-C7, or in the original format, in Card Group 4. No. of gap interconnections (NK) between channels problem. If this is not known> MG=2*MC is ally adequate but should be checked later. For a MG may be given as zero, when it is reset to 1 in CORE. CG No. of fuel nodal points in problem. This should on Card Tl. If MN i:; gi Yen as zero, it is reset to 1 in CORE. No. of rods (NROD) in problem. For PWR and BWR, NROD=NCHANL, hence MR may be given=MC. No. of axial stations in problem. It may be given as NDX (Card Cll) as it is increased by l immediately after reading in. (1) MC to MX are used to set the array sizes in the dynamic r*ge, hence they should be set too big rather than too small. (2) Noie that MC to MX are given in alphabetical order. I I I I I I I ., I I I I I I *1 I I 1. I I I. -.: I , I I I I I I I I I 1* I I I I I I -178-Card ( s) Ty;h3 C2 Maxim.urr.

Running .,.. . ime Reqtlired

to be present Always READ list: MAXT

FORTRAN IS, 6E12.6 Read f=om Subroutine:

INDAT VariablG Columns Format Description CG MAXT 1-5 IS

Maximum Running Time, Nominal value C is 2000. O

N T R 0 L*

THE INPUT li'OH A CASE REQUIRES A CASE CONTROL CARD FOLLOWED WITH UP TO 12 GROUPS OF INPUT INFOHMNrION.

EACH OF THE 12 CAHD GROUPS HAS GROUP CONTROL CARD THAT IDENTIFIES THE GROUP NUMilER AND THE OPTIONS li'O THAT GROUP. GO 'rO THE CARD GROUP SPECili'IED BY NGROUP. IF THE DATJ\ OF J\ CARD GHOU THE SAME AS TIIE PREVIOUS CASE 1 THEN THAT CARD GROUP AND ITS CONTROL MAY BE OMI'l 1 TED. I I-' \.0 I -------------------

I I I I I I I I I I I I I I I I I I I Card(s) Type C3 Required to be present FORTRAN READ list FORTRAN FORMAT: Read from subroutine Variable Column lPILE 1 KASE Jl TEXT -180-Format Il Case Control Card Always lPILE, KASE, Jl TEXT Il, I4, I5, 17A4 IND AT Descriotion:

= 1: for PWR, with interconnected channels.

= 2: for BWR, with separated channels.

as in Appendix 10

-181-Card Group 1 Required to be present FORTRAN READ list FORTRAN FORMAT Read from Subroutine Variable Columns Format NGROUP 1-S IS Nl 6-10 IS Always NG ROUP Nl IS, I5 IND AT Description

= 1 (to select Card Group 1) < O: Calculate physical properties from polynomials.

> 1: the properties are given in the next Cards as in the original COBRA. I *I I I I I I I I *I I I I I I I I I I I 11 I l II 1\ I I I I I I I I I I I I I I I -182-Physical Properties Required to be present FORTRAH READ List FORTRAN FORHAT READ from subroutine Variable Columns 1-5 PH 6-15 P2 16-25 26-30 Format IS Fl0.3 FlO. 3 IS When Nl(in Card Group 0 N PH P2 i*H I5 FlO. 3 FlO. 3 I5 Cards l :>escriotion

= 1: PH defined as lowest pressure encountered in problem, = 2: PH defined as lowest enthalpy encountered in problen Lowest pressure (psia if Nl = 1 or lowest enthalpy (Btu/lb) if Nl = 2 Highest pressure (psia) encountered in problem Number of pressure steps generated by polynomial (maximum 30).

-183-The lowest pressure encountered in the problem is defined as that at which the lowest enthalpy would be the saturation value. For example, at 1000 psia the saturation enthalpy is 543 Btu/lb. At an inlet subcooling of 100 Btu/lb, the enthalpy would be 443 Btu/lb and this would be the saturation value at a pressure of about 470 psi a. Thus, one would require physical property data over the range 470 (or less) psia to 1000 psia in order to include data which covered the enthalpy range. To avoid translating the lowest enthalpy to pressure, the option of giving the enthalpy is included.

The program lates this value to a pressure which is safely below that required using the expression p = 6h3(h -1.35) I (h -0.35) when p =calculated pressure (psia), h = O.OlH, H = enthalphy (Btu/lb).

The values of p, so calculated, are given below and it may be seen that they are all less than P sat, the tab led value of pressure corresponding to H. H(Btu/lb) 181.2 300 400 500 600 700 p (psi a) 11 101 279 589 1067 1749 p sat (psia) 15 103 319 745 1409 2236 I I I I I I I I I I I I I I I I I I. I I I I I I I I* I I I I I I I I I I I *I -184-In the original COBRA, the physical properties are read from cards into the arrays (PP(L), TT(L), etc., L = 1, In the new version, the values of (PP(L), TT(L), etc., L = 1, Nl) generated within a Do Loop from 1 to Nl from the physical property polynomials.

Hi th the arrays set, the subsequent use of the values is the same in both versions of the code. Note: NPROP is set to Nl for storage of the size of the arrays.

*--------*-----------**-*-*--*-

Physical Properties Required to be present When Nl (in the card group 1) > 0 READ IN Nl CARDS OF FLUID PROPERTY DATA. EACH CARD COWfAINS --SATURATION PRESSURE (PSIA), rrEMPERATURE(DEG-F LIQUID SPECIFIC VOLUME(CU-FT/LB), VAPOR SPECIFIC VOLUME(CU-FT/LB), LIQUID ENTHALPY(BTU/LB), VAPOR ENTHALPY(BTU/LB), LIQUID VISCOSITY (LB/FT-HR), LIQUID THEHMAL CONDUCTIVITY(BTU/HR-Ii'T-Ji')

MID SURFACE TENSION(LB/li'T), FORMJ\T(E5.2,F5.l,7Fl0.0).

Nl MUST BE GREATER THAN mm BUT NOT GREATER THAN THE p ARAMErrER rw *. THIS PROPERTY TABLE MUST HAVE PRESSURE HIGHER THAN OPERATING PRESS. AND LIQUID ENTHALPY LOWER THAN THE BUNDLE INLET ENTHALPY * .. -> ----.. .. -.. ---* I I-' (XI Vl I -____ _J

CARD GROUP 2, FLOW CORRELATIONS READ IN UP TO FOUR SETS OF FRICTION FACTOR CORRELATION CONSTANTS THAT CORRESPOND TO THE SUBCHANNEL TYPES, FORMAT(l2F5.3).

Nl IS THE SUBCOOLED VOID CORRELATION OPTION. Nl=O, NO SUBCOOLED VOIDS. Ul=l, LEVY SUBCOOLED VOID CORRELATION.

N2 IS THE BULK VOID CORRELATION OPTION. N2=0, HOMOGENEOUS MODEL. N2 = 1, MODIFIED ARMAND MODEL. U2=5, READ IN SLIP RATIO, FORMAT (5X,El0.5).

N2=6, READ IH THE NUMBER OF TERMS AND COEFFICIENTS FOR UP TO A SIXTH ORDER POLYNOMIAL FUNCTION OF STEAM QUALITY, FORMAT (I5,7El0.5).

N3 IS THE TWO-PHASE FRICTION GRADIENT MULTIPLIER OPTION. N3=0, HOMOGENEOUS.

N3=1, ARMAND. N3=5, READ IN NUMBER OF TERMS J\ND COEFFICIENTS FOH UP TO A SIXTH ORDER POLYNOMIAL FUNCTION OF QUALPf FORMAT(I5,7El0.5).

NII IS AN OPTION TO INCLUDE A WALL VISCOSITY CORRECTION TO THE FRICTION FACTOR. IF N 4= 1 J IT IS INCLUDED, OTimmnsE rrr IS NOT. CAHD GROUP 3, AXIAL HEAT FLUX 'f ABLE READ IN Nl PAIR OF DATA FOR THE TABLE. EACH PAIR CONSIS'rs OF THE RELATIVE POSITION(X/L)

AND TIIE CORRESPONDING ffi.i:LA'fIVE IIEJ\T FLUX (Local flux/AVERAGE li'LUX), EACH CARD ACCEPTS UP TO SIX PAIR OF DATA, }i'ORMA'f(l2F5.

3). Nl r.msr Im GREA'fER THAN ONE BUT NOT GTIEA'fER THAN THE P AHAME'l'EH MP. I I-' co CJ\ I Card Group 4 Required to be present F9RTRAN READ list FORTRAH FORMAT READ from subroutine

-187-(Channel Data) Card (1) when IPILE=l or 2 NG ROUP I5 INDAT I I I I I I Variable Columns Description NGROUP 1-5 = 4 (To select Card Group 4) I NOTE: Once this card is read in the new subroutine CARDS 4 is entered for the remaining Read statements and Data processing of this Card Group 4. I I I t I I ,, I I I I I I I I I I I I 'I I I I *I ., I I I I I Required to be present FORTRAN READ list FORTRAN FORMAT RP-ad from subroutine Variable N2 NG RID NGRIDT NODESF NFUELT NCHF IMAP ITEXT o umns 1-4 5-8 9-12 13-16 17-20 21-24 25-28 29-32 33-36 -188-Card (2) when NGROUP = 4 Nl, N2, NGRID, HGRIDT, NODESF, NFUELT, NCHF, IMAP, ITEXT 9I4 CARDS4 escr1ot1on or channel types (max 15) see below) Total number of channels in problera Number of grid positions Number of types of grid Number of radial nodes on the fuel for center temp. calculo Number of fuel types = O for no CHF calculations

= 1 for B&W2 CHF correlation

= 2 for W-3 correlation

= 1 to 4 to indicate method of presenting gap interconnection data (see Cards (9) below) number of cards to be read in next which will be printed out as a message. If ITEXT=O, no message cards are read in lJOTE Channels are defined as being all of the same type if they have the same geometry, rod dimensions and grids and only differ in their power. More precisely, Cards (4) and (5) given later which define the geometry and grids must apply to all channels of the same type. In, for example, 1/4-core symmetry data, 1/4, 1/2 whole channels would be of different types.

Required to be present FORTRAN READ list FORTRAN FORI*1A'I1 Read from Subroutine Variable TEXT Columns 1-80 -189-Card (3) when ITEXT > 0 TEXT 20A4 CARDS# I I I I I I The array TEXT (20) is read I immediately printed in a DO loop from 1 to I TEXT. It is envisaged that a map of the channel nul"'.'.bering system could I be printed as an aide-memory in a large problem. ,. I I I I I I ,. I I. I I I 1* I I 1* I I *1 I I I I -190-Required to be present FORTRAH READ list FORMAT READ from subroutine Variable Columns N I FRAC AC PW PH GAPS DIST DR PHI M 1 2 -5 6-14 15-23 24-32 33-41 42-50 51-59 60-68 69-77 78-79 Card (4) Always (being !JG ROUP= 4) N,I,FRAC, AC(I), PW(I), PH(I) GAPS(I,l), DIST(I,l), DR(I), PHI(I,l.), M Il, I4, 8E9.3, I2 CARDS4 _ Description Selector for friction factor expression.

If reset to 1. Any channel number, preferably the first of the channe 1 type being described.

Factor by which AC, PW, PH should be multiplied.

Thus for 1/4 channel, one may give FRAC = 0.25 and AC, PW, PH the same as for a whole channel. Channel flow area (in 2) Channe 1 wetted pe rime te r (in) Channel heated perimeter (in) Boundary gap dimensions (in) Centroid-to-centroid channel distance (in). This is only required for a particular mixing correlation and may normally be given as zero. Rod diameter (in) Number of rods in channel Fuel type: = 1 for rod, = 2 for plate, Reset to 1 if M = O

-191-Required to be present FORTRAN READ list FORTRAH FORI*IAT Head from subroutine Variable columns CD FXF Card (5) If IWRID > 0 (CD (I.L), L=l, IWRIDT), L=l, ;-.JGRIDT 16 CARDS4 ::Jescriptions Spacer loss coefficients I Ii II 'h'x.1"'<'(')

i \-LJ 'I/ I 'II -I 11 Fraction of axial flow forced across each boundary.

It is 1* not expected that this would be used in reactor problems hence nominal value = O. O I IC If Nl (Card (2) ) is greater than one, cards describing channel and grid for channel type 2 will be given now, after these cards, the ones describing channel type 3 will be inputed and so on until the completion of the 1n channel types. I I 'I' I I I I I I, I I I I I I I, I I " I I I ., I I I -192-Required to be present FORTRAN READ list FORTRAN FORMAT Read from subroutine Variable Columns RADIAL 1-70 Notes: Card (6) Always (Radial (I), I=l, NROD) 16 E5.3 CARDS4 Description Radial power factor for rod I which is located in channel I. This is defined as the ratio of the rod power to that of the reactor average power. a) NROD is the tota1 number of rods, having set to lWHAi.'lL (total number of channels) which was its elf set to :*.r2 (Card ( 2) )

b) If all rods have the same power, RADIAL (1) alone may be given a.11d is set negative.

THis triggers setting (RADIAL (I); I=l, NROD) = 1.0

-193-Required to be present FORTRAH READ list FORTRAiI FORMAT Read from subroutine Variable Columns GRID XL I GRID Card (7) If

> 0 (GRIDXL(L), IGRID(I), (I=l, iJGRID) 8(E5.3, I5) CArt.DS4 Descritpion Relative location (z/L) wl1ere grids are located I I I I I I I I I ,I I I I I I I I I I I -19 4:-

1' ' ! I I I *I I I Card (3) Required to be present If Nl (Card(2) > 1 FORTnA:r READ list JB(I) FORTRAN FORMAT Read from subroutine.

Variable Columns Description JB 1-80 List of channels of Type 2 :'late s : The first set given is the list of channel numbers in Type 2. The list is terminated by reading in a zero (or a blank space). Hence, if the last channel number comes at the end of a card, a blanl{ card must follow in order to give the terminating zero * ..:t is safer to make a habit of punching a final zero. FollowinG Type 2, card(s) are read in for those cha."1.nels in Type 3, then Type 4 etc. UP to 111 Types. lfote that since the channel numbers for Type 1 are not read in, it is more economical to *select '::'ype 1 as that with the major-ity of channels.

An internal consistency check is made reading in JB (I)

If a set includes the channel number (I in Card ( 4) ) for Type 1 or does not include that given for its own type in Card ( 4) ' an appropriate message is printed and the run terninated.

Ii' iH = l, the JB cards above are not given.

-195-I Card (9a) Required to be present only if IPILE = 1 (If IP::CLE = 2) . BWR case, no cards are given sinl : the channels are not connected) .C. IHAP = 1 ( Card ( 2 ) ) FORTRAN READ list FORNA'I' Read from subroutine Variable ICROS ID OWN iWTES Columns 1-4 5-8 ICROSS, IDOWll 2I4 CARDS4 Description

} see notes below This option is only possible to use when the pattern of channel is recta."1.gular.

If this is the case, I CROSS is the number of columns and IDOWH the number of rows. For example, in the case represented in figure 1, ICROSS should be 4 and IDOWN 3. The maximum value for IDOWN and ICROSS is 20. The channels are sequentially numbered by the computer and the channel boundaries set in the IK, JK arrays; the

oreder is that used to illustrate the case of IMAP = 4 ( Card ( 9 d) )

I I f I 1* 11* ,. I ' ,., ti I I 1. I I -196-I. I I I I I 1 2 7 4 :; I *1. -6 7 3 J I. I 9 10 11 12 I I Figure 1. I

Rectangular Matrix of Channels I I I I I I -197-I I I I 'I -------------------..... ----. *-., I .* .. .,, .. I I -.-I --* I ,, 2 3 4 I \I * ' 5 6 7 8 9 I I I a 10 I I ----------------I I FIGURE 2 Irregular Pattern of Channels I I ,, I ---

I I I I I "' I I I I I I I ., I I I I ' I Required to be present FORTRAN READ List FORTRAN FORMAT Read from subroutine Variable IS TART IEND lJotes: Columns 1-4 5-3 -198-Card ( 9'IJ) Where IPILE = 1 and .IMA.P = 2 ISTART IEND 2I4 CARDS4 Description First channel in each row Last channel in each rm-1 One of these cards should be given for each row. !-l'ote *that this method could not be used if there were an insert blank channel in any row; for this case use IMAP = 3. The maximum value of IEND is 20 and the maximura number of rows is also 20. If less than 20 rows are to be given, a blank card (or one with two zero) should be given after the last row. The computer numbers the channels and the boundaries sequentially as illustrated in Figures 1 and 2. Examples follow:

-199-For Figure 1 the following cards should be inputed: I START 1 1 1 0 IEND 4 4 4 0 For Figure 2 the following cards should be inputed: I START IEND 3 3 3 5 :i 5 4 4 0 0 ' I I I I I 'I I I I I I I I I I I a. I

' I I I I I ' I I. I ' I 1: I I I I ' I -200-Required to be present FORTRAN READ list FORTRAN FORHAT Read from subroutine Variable Columns MAAP 1-80 Notes: Card ( 9 c) When IPILE = 1 and IHAP = 3 (NAAP(L), L= 1,20) 20 I4 CARDS4 Description The number of the cha'l.nels malcing up a row One of these cards should be inputed for each row (maximum 20 rows). The value o.f T:IAAP represents the channel number with a zero indicating no cha.'1.nel.

If less than 20 cards are to be used, the last should be all zeros (i.e., a blank card). The set of cards represents a map of the channel numbering

system, is thus under the control of the user. The boundary ordering is done by the computer.

Examples:

For pattern described in figure 1 1 5 9 0 2 6 10 0 3 7 11 0 4 8 12 0 I I. I I I I I I I I ,. I I "I I ., I I I For pattern described in 2 0 0 5 0 0 0 0 6 0 1 2 3 7 8 0 lJ 4 9

-202-Required to be present FORTRAlJ HEAD list FORTRAN F.ORMAT Card (9d) When IPILE = 1 a.'1.d H1AP = 4 (IK(L), JK(L), L = l,NK) 20 I4 I I II 'I I L.----------R_e_a

__ d __ f_r_o_m __ s_u_b __

________________

c_A_RD __ s_*_ii ____________________________ Variable ColUTILYlS Description IK } See notes be low JK ifote s: IK, JK are the channel pairs defining each boundary in turn; lJK = nur.rber of boundaries specified.

The set of numbers are read in, 20 to a card, continuing on as many cards as necessary.

They are terminated by a zero; if the final channel number is at the end of a card, the zero must be given on the next card. (lfote, the value of ?IK is not known at the time of reading in IK, JK;. it is set to the number of pairs read in). Thus, with H1AP = 4, both channel and bou.'1.dary numbering are under tne control of the user. When listing the sub channel pairs, it is preferable to give the loVTer number first; this saves the computer reversing the order. I ., I I I I, I I I I I I I I I I I 11 18 I I 11 ,, I I I 'I I I I I -1 5 2 10 2 9 2 0 2 6 For case in figure 2: 3 3 4 2 7 3 0 For Case in figure 1: 3 10 3 4 7 11 1 5 8 12 2 9 -203-Card (9d) Examples 1 8 4 9 5 6 6 7 7 8 8 Q 6 3 7 4 8 5 6 6 7 7 10 10 11 11 12 0 0

-204-Required to be present FORTRAH Read list FORTRAH FORMAT Card ( 10) When IPILE = 1 JB (L), L = 1 J 20 I4 I I I I I CARDS4 I Read from Subroutine Variable JB Notes: Columns 1-80 Description List the identification number of the channels making up each boundary", i.e. the aries that are split by a line of symmetry.

Always terminate with a zero. If there are no half boundaries, give a single card with a zero. The parameter FACTOR(K) is set to 1.0 for full boundaries and 0. 5 for "half-boundaries".

I I I I I 1-I 1* I I I I I I I. I I I I I I I I ., 1 *1 1* I ., I I I Required to be present FORTRA11 READ list FORTRAN FOill 1 TAT Read from subroutine Variable KFUEL CFUEL RFUEL DFUEL KC LAD CLAD RCLAu TC LAD HGAP Columns 1-5 6-10 11-15 16-20 21-25 26-30 31-35 36-40 41-45 -205-Card (11) When 2JODESF > 0 (K FUEL(I), CFUEL(I), fiFUEL(I), DFUEL(I), KCLAD(I), CCLAD(I), RCLAD(I), TCLAD(I), HGAP(I), I=l, NFUELT 16E5.3 CARDS4 Descriotion Fuel thermal conduction BTU ( ) Fuel specific heat BTU lb °F Fuel Density (lb/ft 3) Pellet Diameter (in) Cladding thermal conduction BTU hrft °F Cladding specific heat BTU lb °F Cladding density (lb/ft3) Cladding thickness (in) Fuel-cladding heat 2ransfer coefficient (BTU/ft hr°F)

-CAHD GROUPS 5, 6, 9, 10, 11 AND 12 ARE READ IN BY SUBROUTINE INDAT WITH THE FOLLOWING FORMAT: CARD GROUP 5, SUBCHANNEL AREA VARIATION TABLE Di' THERE ARE NO AREA VARIATIONS, OMIT THIS CARD GROUP. READ N2 VALUES OF RELATIVE LOCATION(X/L)

WHERE AREA ARE GIVEN FORMArr ( 12F5. 3)

N2 MUST BE GREATER THAN mm BUT NOT GREATER THAN THE PARAMETER ML. READ Nl SETS OF AREA VARATION FACTORS(LOCAL AREA/NOMINAL AREA). EACH SET CONSIS'l 1 S OF SUBCHANNEL NUMBER AND N2 AREA VARIATION FACTORS, FORMAT(I5/(12Ji'5.3)).

Nl IS LIMITED BY THE PARAMETER MA. IF Nl IS ZERO, AREA VARIATIONS ARE DELETED FOR SUCCEEDING CASES. N3 IS THE NUMBER OF ITERATIONS FOR INSERTING AREA VARATIONS.

IF N3 IS ZERO OR BLANK, N3 IS SET EQUAL TO 1. CARD GHOUP 6, CAP SPACING VARIATION TABLE Ill' 'rIIERE ARE NO GAP VARIATIONS, OMIT THIS CARD GROUP. READ N2 VALUES OF THE RELATIVE LOCATION(X/L)

WHERE GAP li'ACTORS ARE GIVEN, Ji'ORMAT(l2F5.3).

N2 MUST BE GREATER THAN ONE BUT NOT GREA'l'ER THAN THE PARAMETER ML. READ Nl SETS 01" GAP SPACING l"ACWORS (LOCAL GAP /NOMINAL GAP). EACH SET CONSIS'l'S OF THE ADJ AGENT SUBCHANNEL NUMBERS FOR THE GAP N2 GAP VARIATION l"ACrroRs, FORMAT( 2I5/( 12F5. 3)). Nl IS LIMITED BY THE PAHAMETEH MS. IF Nl IS ZEHO, GAP VARIATIONS ARE DELETED FOR SUCEED CASES. .. -; -.. .. .. " ........ _ ..... .. -.... .. I l\J 0 CJ\ I -

CARD GROUP 9, CALCULATION VARIABLES RirnD IN DIVERSION CHOSSFLOW RESISTANCE FACTOR, TURBULENT MOMENTUM FACTOR, BUNDLE LENGTH(IN.),, POSITION FROM VERTICAL(DEGREES), NUMBER OF AXIAL NODES, NUMBER OF TIME STEPS, TOTAL TRANSIENT Tif.lE(SECONDS)

MAXIMUM NUMBER OF ITERATIONS, ALLOWABLE FRACTION ERROH IN FLOW Ji'ORM CONVERGENCE AND THANSVERSE MOMENTUM PARAMETER(S/L), FORMAT(llE5.2,2I5,E5.2,,I5,4E5.2).

IF THE NUMBER OF ITERATIONS, ALLOWABLE ERROR AND MOMENTUM PARAMETER ARE BLANK OR ZERO, THE PROGRAM USES 20,, 1.E-3, AND .5, RESPECTIVELY.

Nl IS AN OPTION GIVING THE SPATIAL PRINTING INCREMENT.

IF Nl=l, EVERY STEP IS PRINTED. IF N2=2, EVERY OTHER STEP IS PRINTED, ETC. IF N '.':EHO OR BLANK, THE PROGRAM SI:':'S Nl=l. N2 IS AN OPTION GIVING THE TIME PRINTING INCREMENT AND IS SET UP THE SAME AS Nl ABOVE. N3 IS A DEBUG PRINT OPTION. IF N3=0, NO_DEBUG INFORMATION IS PRINTED IF N3=1 A DEBUG PRINT IS MADE FOR EACH STEP OF THE CALCULATION.

IT CAN GENERATE A LOT OF PAPER SO IT IS NOT NORMALLY USED. CARD GROUP 10, TURBULENT MIXING CORRELATIONS Nl IS THE OPTION FOR SUBCOOLED MIXING CORRELATIONS.

FOR ANY SELECTED VALUE OF Nl READ IN THE CONSTANTS A AND B, FORMAT(2F5.3).

THE OPTIONS ARE --I (\) 0 --.1 I Nl=O, W/GS=A Nl=2, W/GS=A*RE**B Nl=2, W/GD=A*RE**B Nl=3, W/GS=D/AIJ*A*RE**B NOTE THAT BETA = W/GS WHERE W IS THE TURBULENT CROSSFLOW.

N2 IS 'rIIE OPTION FOR TWO-PHASE MIXING. IF N2=1, TWO-PHASE MIXING IS THE SAME AS FOR SUBCOOLED CONDITIONS.

IF N2 IS GREATER THAN ONE READ IN N2 PAIR OF DATA FOR A TABLE OF TWO-PHASE MIXING DATA. EACH PAIR CONSISTS OF THE STEAM QUALITY AND THE CORRESPONDING VALU OF BETA. N2 IS LIMITED BY THE PARAMETER MP. N3 IS THE OPTION FOR THERMAL CONDUCTION MIXING. IF N3=0, NO THERMA CONDUCTION.

IF N3=1, READ IN THE THERMAL CONDUCTION GEOMETRY FACTO FORMAT(F5.3).

CARD GROUP 11, OPERATING CONDITIONS READ IN THE OPERATING PRESSURE(PSIA), INLET ENTHALPY(BTU/LB)

OR INLET MASS AND AVERAGE HEAT FLUX(M-BTU/HR-SQ-FT).(6Fl0.0)

Nl IS THE INLET ENTHALPY OPTION. IF Nl=O, INLET ENTHALPY IS GIVEN. IF Nl=l, INLET TEMPERATURE IS GIVEN. IF Nl=2, READ IN THE INDIVIDUAL SUBCHANNEL INLET ENTHALPIES, FORMAT(l2F5.0).

IF Nl=3, READ IN THE INDIVIDUAL SUBCHANNEL INLET TEMPERATURES, FORMAT(l2E5.0)

I ru 0 (X) I

.. --.. -*--* ... .. ** -----.. ---.J N2 IS THE INLET FLOW DISTRIBUTION OPTION. IF N2=0, THE SUBCHANNELS ARE GIVEN THE SAME MASS VELOCITY.

IF N2=1, THE INLET FLOW IS DIVIDED TO GIVE EQUAL PRESSURE GRADIENT IN THE SUBCHANNELS.

IF N2=2, READ MASS VELOCITY FACTORS FOR EACH SUBCHANNEL, FORMAT(l2E.50).

N3, N4, N5 and N6 ARE OPTIONS FOR TRANSIENT FORCING FUNCTIONS.

IF' ANY OF THESE OPTION NUMBERS ARE ZERO OR BLANK, THE CORRESPONDING FORCING DATA IS NOT READ AND IS EXCLUDE FROM THE CALCULATIONS.

EACH OF THESE NUMBERS GIVE THE NUMBER OF PAIRS OF TABULAR DATA TO BE READ FOR EACH FUNCTION.

ALL DATA ARE READ AS PAIRS OF TIME (SECONDS)

AND RELATIVE VALUE, FORMAT (12E5.0).

N3 IS THE OPTION FOR REFERENCE PRESSURE VERSUS TIME. N4 IS THE OPTION FOR INLET ENTHALPY OR TEMPERATURE AS A FUNCTION OF TIME DEPENDING ON THE OPTION FOR INLET ENTHALPY OR rrEMPERATURE.

N5 IS THE OPTION FOR INLET FLOW VERSUS TIME. N6 IS THE OPTION FOR HEAT FLUX VERSUS TIME. CARD GROUP 12, OUTPUT DISPLAY OPTIONS Nl IS AN OPTION FOR PRINTING ANSWERS. Nl=O, PRINT SUBCHANNEL DATA ONLY. Nl=l, PRINT SUBCHANNEL DATA ANC CROSSFLOWS.

Nl=2, PRINT SUBCHANNEL DATA AND FUEL TEMPERATURES . Nl=3, PRINT SUBCHANNEL DATA, CROSSFLOWS AND FUEL TEMPERATURES.

I I\) 0 \D I N2 IS AN OPTION FOR SUBCHANNEL DATA PRINTOUT.

IF N2=0, ALL SUBCHAN DATA ARE PRINTED. IF IT IS CALLED FOR BY Nl. FOR N2 GREATER THAN Z READ IN THE SUBCHANNEL NUMBERS FOR WHICH RESULTS ARE TO BE PRINTED FORMA'l 1 ( 3612). N3 IS AN OPTION FOR FUEL TEMPERATURE PRINTOUT.

IF N3=0, DATA FOR ALL RODS ARE PRINTED IF CALLED FOR BY Nl. FOR N3 GREATER THAN ZERO, READ IN N3 ROD NUMBERS FOR WHICH TEMPERATURES ARE TO BE PRINTED, ( 3612). IF CHF DA'rA IS CALLED FOR BY INPUT OPTION IT IS PRINTED FOR EACH SELECTED ROD PLUS A

SUMMARY

TO IDENTIFY THE ROD AND CHANNEL WITH THE MINIMUM CHF RATIO. N4 rs AN OPTION FOR FUEL NODE PRINTOUT.

IF N4=0, TEMPERATURES ARE PRINTED FOR EVERY NODE. FOR N4 GREATER THAN ZERO, READ IN N4 NODE NUMBERS TO BE PRINTED, FORMAT(3612).

TO START A CALCULATION, READ A BLANK GROUP CONTROL CARD. TO STOP THE CALCULATIONS, AFTER FINISHING A CASE, READ A BLANK CASE ***** END OF INPUT INSTRUCTIONS

- - - UNITS -ALL COMPUTATIONS ARE DONE USING FT, LB, SEC, BUT AND DEG-F. UNIT CHANGES FOR INPU'l 1 AND OU'l 1 PUT ARE DONE IN THE PROGRAM. I I\) I-' 0 I I I I I I I I I :I I : I ,I ,, I I I I I I I I -211-APPE?1DIX 12 New INPUT DATA Presentation I I I I I ., I I* I I I I I I I ., I I I Card(s) Type: co Required to be pr.esent:

FORTRAN READ list: FORTRAN FOIU*fAT:

Read from Subrout:i.ne:

Variable Columns NAP 1-5 NFS 6-10 -2lla-Problem Array Size Always NAP NFS 215 IND AT Format IS 15 Description Number of axial positions forinputting the heat flux shape (sets size of matrix, NAP > Nl given on card CS) Number of flux shapes input (sets size of matrix, NFS > JFS given on card .cs)

,. -212-Card(s) Type Cl Problem Array Size ., Required*

to be present Always READ list:

I FORMAT: MC MG MN MR MX lOIS INDAT -I Read from -I I. I I* I I I I* I I I I I I Variable Columns MC. 1-5 .. *. MG 6-10 MN 11-15 MR 16-20 21-25 Notes Format* IS IS *IS IS IS

blem. cards COBRA Description

.. of channels (NCHANL) in NCHANL is set from NTHEOX on CS-C7, or in the original format, in Card Group 4. No. of gap inte=connections (NK) between channels If this is net MG=Z*MC is ally adequate but should be checked later.* For a BWR .. MG mav be given as zero, it is reset l in CORE. CG No. of fuel nodal points in problem. This*should on Card Tl. If MN is given as zero, it is reset. to 1 in CORE. No. of rods (NRCTI) in problem .. For PWR and BWR, MR may be gi.ven=MC. No. of c.xial stati_cns in proble:.t.

-It may be given as NDX (Card Cll) as it is increased by l after reading in. (1) MC to MX are used to set the array sizes in the dynamic rage, hence they should be set too big rat:her than too smalls (2) Note that MC are given in order.

I I I I I I I I I I I I I I I I I I I Card(s) Type: Cla Required to be present: FORTRAN READ list: FORTRAN FORMAT: Read from Subroutine:

Variable Columns NCHF 1-5 JBXC 6-10 Jvnun 11-15 JGRID 16-20 -212a-CHF Correlation Always NCHF JBXC 4IS IND AT Format IS IS IS 15 JWRBl JGRID Description CHF correlation indicator

= 0 for no CHF calculations

= 1 for B&W correlations

= i for Westinghouse correlations B&W indicator.

= 0 for B&W-2 correlation

= 1* for BXC correlation Westinghouse indicator

= 0 for W-3 correlation

= 1 for WRB-1 correlation W-3 spacer factor indicator

= o*for no grid factor = 1 for L-grid. factor = 2 for R-grid factor I I I I 1 I I I I I I I I I I I .1 I I Card(s) Type: Clb -Required to be FORTRAN READ 1 if; t: FORTRAN FORMAT: Read from Subroutine:

Variable Columns XKS 1-10 TDC 11-20 -212b-W-3 Spacer Factor Parameters Only if NCHF = 2. JWRBl = 0 JGRID :f: 0 XKS TDC 2Fl0.S IND AT Format Ft0.5 Fl0.5 Description Axial grid spacing coefficient (k. ) s Thermal diffusion coefficient (TDC) I I I I I I I I I I I I I I I I I I I I Card(s) Type: Cle Required to be present: FORTRAH READ list:

FORMAT: Read from Subroutine:

Variable GSP PF Columns 1-10 11-20 -212c-WRB-1 Correlation Parameters Only if NCHF ,;,, 2 JWRBl = 1 GSP PF 2Fl0.5 INDAT Format no;s Fl0.5 Description Grid (in) Performance factor I I I I -*1 I I I I I I I I I I I I I I Card(s) Type: Cld Required to be FORTP.AN READ list: FORTP-AN FORHAT: Read from Subroutine:

Variable Columns IHC 1-S FLO ARE 6-15 -212d-Identification of the Hot Channel Always IHC FLOARE IS, FlO.S IND AT Fonnat IS Fl0.5 Description Channel number which.is the hot channel Nominal area of the hot channel (in )

I I I I I I I I I I I I I I I I I I 1. Card(s) Type: Cle.* Required to be pr_esent:

FORTRAN READ 1 is t : FORTRAN FORMAT: Read from Subroutine:

Variable Columns FNT 1-10

XOINCH 11-20 DAMPNG 21-30 DNBRX 31-40 HLEN 41-50 -212e-Miscellaneous Data Always FNT XOINCH DAMPNG DNBRX HLEN 5Fl0.5 IND AT Format Fl0.5 Fl0.5 Fl0.5 ,Fl0.5 Fl0.5 Description Fraction of heat generated in the fuel and cladding Axial elevation at which active fuel starts (in) Damping factor to aid convergence of flow solution (0.0 < DAMPNG < 1.0) Multiplier on all DNBRs Heated Length (ft)

I I:. --* I I I I I I I I I I I I I I I \ / I I Card(s) Type C2 Required*

to b*e present FORTRAN READ list:

FORTR.2\N FORMAT: Read. from Subroutine:

-213-Maximuia Running Time Always MAXT IS, 6El2.6 IN DAT Variable Columns Format Description

' MAXT 1-S !S Maximum Running Time, Nominal value is 2000. .. CG c 0 N T R 0 L*

-214-Card(s) Type C3 Case Control Card Required.to be present Always FORTRAN READ list: I PILE KASE Jl TEXT FORTRAN F0&."1AT:

Il, I4, IS, 17A4 Read. from Subroutine:

INDAT $" Variable Columns Format Description I PILE 1 KASE 2-5 Jl 6-10 TEXT 11-78 II = O for simplified method = 1 for PWR = 2 for BWR The value is if Card Group 20 is selected since it is overwritten on card Tl. I4 Run Identification iJumber --as in COBRA IIIC. If > O, calculation continues; if , calculation stops. I5 Printing option for standard COBRA output--as in COBRA IIIC. = O print only new input = 1 print entire input = 2 print only operating conditions

This option is only effective if NOPRIN = O, i.e., Nl =Don card C4 Alphanumeric information to identify Case. I I I I I I CG I c 0 I N T R I 0 L I I I I I *I I I 1. I I ,, I I I I I I I I *I I I I I I I I I -215-Card (s) Type C4 Select Card Group 20 Required.to be present Always FORTRAN READ list: NG ROUP Nl N2 N3 N4 NS N6 FORTRAN FORMAT: tIS Read. from Subroutine:

INDAT Variable Colur.ms Format Description NG ROUP 1-5 Nl 6-10 11-35 Notes: I5 = 20 I5 Printing trigger, NOPRIN, set to NL Nl=O, standard COBRA IIIC printing obtained as well as as "new" printout.

Nl=l, standard COBRA printing suppressed.

I5 Leave blank ( 1) If NGROUP = O, this acts as a trigger to stop reading Input Data and to start the hydraulic calculation

(.eog., after card T30). CG c 0 N T R 0 L

--... -* card(s) Type cs Required.to be present* FORTRl-.l'i READ list: . FORTRAN FORMAT: Read. from Subroutine:

-216-Channel Map parameter Always IMAP . NDlX NDZX 14IS CARDZO *1 I I I I Variable Columns Format Description CG I IMAP NDlX NDZX 1-5 6-10 11-15 If IMAP = Go to Card 1 cs IS IS IS 2 C6 J Selects method of reading channel into array NTHBOX (NDlX, NDZX). !MAP=l, Z or 3 . Size of array NTIIBOX, maximum values of each are ZS. 3 C7 I I I I I* I I -1 I I I I I I I. __,* I *1 I I I. I I .1 I I I I I I I ' I --217-The channel numbering system is in the array NTHBOX (NDlX, NDZX) with a zero for each non-channel.

The array is later used to define the interaction between adjacent channels.

Thus a channel map: l 2 3 4 5 6 7 8 9 would be represented in NTHBOX (S, 4) as 0 0 4 0 0 1 2 3' s 6 0 9 0 0 7 0 0 0 8 -0 1£ IMAP*l, there are assumed to be NDlX :* ND2X channels numbered sequentially along each row, and column by columnj to give a rectangular matrix. Thus ND1X=4, NDZX=3 gives a channel map: 1 2 s . 6 3 7 4 8 9 10 11 12 For IMAP=Z, 3 more complicated channel maps may be specified.

card(s) Type C6 Required.to be present* FORTRAN READ list: * -218-Channel Map Only if IMAP=2 !START IFIN FORTRAN FORMi'\T:

1415 NDZX Cards of this type read. Read. from Subroutine:

CARDZO Variable Columns Format* Description

!START 1-5 IS see below !FIN 6-10* IS II " A total of NDZX cards of this type are read sequentially, one for each row of the channel map. Each card gives the start and finish of the row. For example, ISTART=3, IFIN=6 would imply a row . 0 0 (N+l) (N+Z) (N+3) (N+4) o a etc. where channel N was the last channel in the previous row. For IMAP=Z, ND1X=7, NDZX=4, cards 3 3 3 6 1 7 2 2 would represent a channel map 1 2 3 4 5 6 7 8 9 10 11 12 13 I I I I I CG I I I I I I I I *I I I I I I I 1. 1-'I I I I I -219-card{s) Type C7 Channel Map Required.to be present Only if IMAP=3 FORTRAN READ lis.t: ( (NTHEOX (NDl, ND2), NDl=l, NDlX), ND2=1, ND2X) FORTRAN FORMAT: (1415) Read. from Subroutine:

CARD20 Variable Columns Format* Description NTHBOX 1-70 14!5 Channel identification number If ND1X>l4, the remaining numbers (i.e., 15-NDlX) are read on a corntinuation card. Note NDlX must not exceed 25. Each row of NTHBOX must start on a new card. CG I For IMAP=3, NDlX= 7, NDZX=4, cards. I I I I I I I I I I :Q 0 0 0 6 7 0 13 l 2 8 3 4 s 9 10 11 12 would give the same channel map as that illustrating IMAP=Z (see card C6). IMAP=3 could be used, either to specify a particular numbering system or wheri there are two channels in the same row separated by a "zero." In the simplified method, (i.e. IPILE=O) *cases as the one represented below may be required to be used. To input this kind of array only IMAP=3 is adequate*.

- lfihe cards needed are illustrated in the figur*e below ..

-220-1 2 3 4 6 7 8 9 5 10 11 12 13 14 16 17 18 19 15 20 21 22 23 24 25 26 27 28 IMAP=3, ND1X=6, ND2X=6 and 1 2 2 3 3 4 5 6 7 8 9 10 5 11 12 13 14 10 15 16 17 18 19 20 15 21 22 23 24 20 25 26 26 27 27 28 I I I I I I I I I I I I I I I I I I Card(s) Type: C8 Require& to be present: FORTRAN READ list: FORTRAN FOP.i.'1AT:

Read from Subroutine:

Variable Columns Nl 1-S JNR 6-10 JFS 11-lS AFLUX 16-25 -221-Format IS IS IS El0.5 Indicators for Axial Flux Shapes Always Nl JNR* JFS AFLUX 315, ElO. S CARD20 Description Nl=O; to read average nodal fuel powers after rest of data (Cards Cl2-14). NAX set to 0, IQP3 set to 0. Nl=l; trigger to read average nodal fuel and coolant powers after rest of Cl2-14). NAX set to 0, IQP3 set to 1. Nl_.:2; number of axial positions used for inputting the axial flux.shapes given on ca'I'd C9. Maximum value of N1=30. NAX set to Nl, IQP3 set to 2. Total number *of rods input Number of axial flux shapes input

<6) Core average heat flux in MBtu/hr-ft 2* If Nl=O or 1, the value of AFLUX is irrelevant and may be given as zero.

I I I I -I I I I I* I I I I I I I I I Card(s) Type: C8a Required to be present: FORTRAN READ 1 is t : FORTRAN FORMAT: Read from Subroutine:

Columns NRFLX 1-70 NOTES: -22la-Format 14I5 Numbering Scheme for Axial Flux Shapes Only if Nl>l and JFS>l . (NRFLX(I), I= 1, JNR) (1415) CARD20 Description Flux shape number which is aP.plicable to each rod. The flux.shape numbers are input in the same order as the rods are input.

(1) Card C9 allows the user to input from 1 to 6 *axial flux shapes. These flux shapes are numbered from l to 6, and each rod is assigned a particular flux shape via the numbering scheme provided on card C8a. The first yalue of NRFLX would be the flux shape number applicable to the first rod, the* second value of NRFLX would be flux shape number applicable to the second rod, etc. (2) If JFS = 1, then all the rods are .internally assigned the first flux shape.

I I I I I I I I I I I I I I I I I -222-Card(s) Type: C9 Axial Flux Shapes Required to be Only if 1h>l FORTRAN READ 1 is t: (YFLUX(I), (AXFLUX(J), J=l,JFS), I=l;Nl) FORTRAN FORMAT: 7Fl0.5 Read from Subrou.ti.ne:

CARD20 Variable Columns Format Description YFLUX 1-10 Fl0.5 Relative axial position, X/L (O. 0<,YFLUX<l.

0) AXFLUX 11-20 Fl0.5 Corresponding relative flux value for the first axial flux shape AXFLUX 21-30 FlO.S Corresponding relative flux value for the second axial flux shape AXFLUX 31-70 4Fl0.5 Corresponding relative flux values for the third, fifth, and sixth axial flux shapes NOTE: First relative position must be 0.0, and the last relative position must be 1.0. ---------------------

I --I 1* I I I I I I I I I I I p. I Card(s.)'

Type ClO Required.to be present* FORTRAN READ li s*t:

FORTRAN FORMAT: Read. from Subroutine:

. -223-Rod Power Factors Only required if Nl on Card (RADIAL (I), I=l, NCHANL) (14ES. O) READIN/CARDZO Variable Columns Format Description CG 1-70 14E"S ** 0 Relative Rod Power (local/average) 8 RADIAL N.CHANL No. of* channels in problem in Card Cl)

It is = set to the highest value of the channel map array NTH BOX see cards CS-C7. Note In the simplified method (IPILE7Q) some subchannels are lumped together to create .one channel, are treated as individual sub channels, (see figure below). For channel can be visualized as having only rod that generates the whole power of the channel. In order to reduce the Input Data the power given.to such a channel for its rod is here, while rods that share their power with several . channels, will be described in Card T5a. This system of *entering the Data, reduces the cards required in the old presentation (do not forget that more than 150 channels can be used and only a few of them will be real sub6hannels) and only introduce the restriction that the lumped channels need to have the same identifi-cation number as its rod .. The following example clarifies all these points:

-224-l=channel 2 3 0-l=rod 02 03 4 13 Dr 5"-9 ' 5 6 0-04 I "'\ I 9 'b 7 B" \ 7 8 .... 17 "\ T=i / t"-..°fh 10 11 12 0 10 0 11 0 12 For this case, card ClO should have the actual relative rod power for channels 1, 2, 3, 4, zero for 5, 6, 7, 8-and the actual values for 9, 10, 11, and 12. The power given to channels 5, 6, 7 and 8 from rods 13, 14, 5, 6, 7, 8, 15, 16, and 17 will be specified later in card T5a. I I I I I I I I I I I I I I I I I I "'* I ,. I I I I I I I I I I I I, i I -225-Card(s) Type Cll Required.to be present Mrscellaneous data 1 Always FORTRAN READ list:

FORT.R.:"\N FORMAT: Read from Subroutine:

Variable Colu:m.ns z 1-5 . Nnx* 6-10. NDT 11-15 TT I ME 16-20 .. Z NDX NDT TTUIE (ES.O, 2IS, lOES.O) CARD20 Format* Description ES.O IS IS ES.O r ! ' ' Channel length (in.) Number of axial intervals Number of time steps NDT=O; steady state only NDT)O; steady

+ transient Total duration of transient (sec) The length of each time step is set to TTIME/NDT

CG 9 9 9 9

** *----*--*

-**----------*---*------*.

-*-*----226-I

Card ( s) Type Tl Channel Indicators . -_,,,* Required*

to be present Always I )

READ list:

FORTRAN FOI' .... -.!AT: Read. from Subroutine:

Variable I PILE NCTYP NGRID NGRIDT NODE SF

NFXF Columns 11-lS 16-20 21-ZS 26-30 IPILE

NCTYP NGRID NGRIDT NODESF NFXF I (14I S) CHAN Format* IS IS IS IS rs IS I Description CG I Iteration trigger=O for.simplified method-* =1 for PWR, =2 for BWR ' No. of channel types to be read *in; controls reading of cards T2-T4. No. of grid positions (maximum=lO)

No. of grid types (maximum=S)

No. of fuel nodes No. of "forced flow" types. Not in use; leave blank I I I I' I I I I I

---.. I I ----t ,. I I I I I I I ,. I I I I t* f , I -227-card(s) Type T2 Required.to be present* Channel Data for Type Y Always FORTRAN READ lis-t:

FORTRAN FORMAT: Read. from Subroutine:

Variable N J FRAC GAP HNR DR If J=l: A B-C

D If J=2: A B c D Columns 6-10 11-15 16-20 21-25 26-30 31-35 36-40 41-.45 46-50 31-3S 36-40 41-45 46-50 N J . FRAC GAP HNR Dm. A B C D (215, SES.O) CHAN Format IS IS ES.O ES.O ES.O ES.O ES.O ES.O ES.*O ES. 0-ES.O ES.O ES. 0 ES.O *.J Description Friction Indicator to select tion factor for c-h..a.nnel (seeTlO).

Nominal value=l, zmaxirnum=4.

Indicator to A, B, C, D low (=l or 2) Amount by which cfuannel area, wetted and heated perimetiers and number of

heated rods are be multiplied (see below). CG 4 Effective rod gav interconnec-4 tion between chan.m.els (in.)

If IPILE=O this may be given as zero. No. of heated rod's in fuel asseJI!,bly.

8 Diameter of rods (in.) Channel area cinZ] Channel Wetted pe'r*imeter (in.) Channel heated perimeter (in.) Not used--leave hTuank No. of unheated (e.g., control) rods Diameter of rods (in.) Width of square a:ssembly (in.) Radius of corners (in.) 8 4 4 4

-,) _, -228-.. * ....... . Notes (1) In COBRA IIIC., individual cards are read for each nel and rod. For PWR and BWR smeared assemblies, considerable simplification is possible because (a) there is a one-to-one correspondence between channels and rods, hence the data may be given together, and (b) many channels have identical geometries, hence one may give a typical geometry and specify to which chan-nels it applies. (2) Channels are of the same tyP.e if they a-ire described by the same data on cards T2, T3. (3) Cards T2, T3, T4 are read sequentially in a DO Loop I=l, NCTYP. Channels making up Types 2, NCTYP are* specified on card T4. The unspecified channels are taken to be of Type 1, hence for economy, Type 1 should be defined as that which contains the majority of the channels.

(4) The channel area and perimeters may either be given directly (J=l) or calculated from the dimensicns of the assem-bly (J=2). (5) These parameters are multiplied by FRAC. Thus, if a line of symmetry divides a channel -so that it is a half-channel, the data for a whole channel may be given and FRAC set to 0.5. Alternatively, data for a single channel may be given and FRAC I I I I I I I I I I I I I I I I I I I I* , _ _ .. I I I I I I I I I I I I I* ,, I _} I I -229-set to (say) 4.0 to obtain the parameters for a smeared group of 4 channels.

If FRAC is given as zero, it is reset to l. O. (6) GAP is the "effective 1' gap between assemblies.

For no internal resistance to mixing within an assembly, GAP could be considered to be the gap between individua.l rods

the ber of gaps. This would be reduced acccTding to the internal resistance used. (7) Next card read is: NCTYP=l NGRID > 0 Card T3 NGRID = 0 Card TS NCTYP>l; I=l (i.e., first type) NGRID > 0 Card T3 NGRID = 0 Card T2 for I=2 NCTYP>l: I>l (i.e.,* subsequent types] NGRID > 0 Card T3 NGRID = 0 Card T4

-230-Card(s) Type T3 Required.to be present Grid Data for Channel Type I Only if NGRID>O :PORTRAH READ list: FORTRAN FORMAT: Read. from Subroutine:

Variable Columns CDG 1-70 (CDG(L),L=l, NGRIDT) (14ES.O) CHAN Format* 14ES. 0 Description Single phase grid coefficient for each grid type. *I I I I *I CG I 7 *1* I I I ,, I I .I I ,, .. I 1. I

, .. I --. I ,. I* I* 'I I I I* I ,. 'I I 'I I I I I card(s) Type T4 Required.to be present FORTRAN READ list: FORTRAN FORMAT: Read. from Subroutine:

-231-Channels making up Type I Only if V>l (JB (L), L=l, N) (1415) CHAN Variable Columns Format* Description JB 1-70 Notes: IS Channel Identification Number for I .. (1) The channels of Type I are listed on one or more cards. A complete card is read and the numbers up to the first zero are taken as the relevant channels.

The zero (or blank) must be given since it acts a*s a trigger, hence if the last channel number is at the end of a card, a blank ca.rd must follow to supply. the ting zero. (2) Next card read is: I

NCTYP I < NCTYP Card TS Card T2 CG

-*' Card(s) Type TS Required.to be present FORTRAN RE.AD lis*t: FORTRAN FQffi.1AT:

Read. from Subroutine:

-232-Grid Positions Only if NGRID)O (GRIDXL(I), !GRID(!), I=l, NGRID) ( 7 (ES

0 , IS) ) CHAN I I I I I Variable Columns Fonnat Description CG *I GRID XL 1-70 ES.O :-tGRID 1-70 IS Fractional distance up channel (x/L) 7 at which each grid is situated, i.e.,

' Grid Type; the coefficients for each 7 type of grid were read by T3. 1* I I 1. I I *I .'I I" I 1-1. I I. I I I I I I I I I I I I I I I I I I Card(s) Type T5a Required to be present: FORTRAN READ list: FORTRAN FORMAT Read from Subroutine:

Variables Columns Format NNll 1-5 I5 NN22 5-10 I5 NN33 10-15 I5 NN44 15-20 I5 -233-Rod Indicators only if IPILE=O NNll, NN22, NN33, NN44 (415) CHAN Description Cards of rod layout data to be read Total number of rods Number of radial fuel nodes including the cladding Total number of fuel types CG 8 8 8 8 Card(s) Type T5b Required to be present: FORTRAN READ list: FORTRAN FORMAT Read from Subroutine:

Variables Columns Format N 1 Il I 2-5 I4 DR(I) 6-10 E5.0 RADIA(I) 11-15 E5.0 ( 2 ) [LR ( I, L ) I3 PHI(I,L) E7.0 -234-I Rod layout information only if IPILE=O and NNll > 0 I N, I, DR(I), RADIA(I), (LR(I,L), PHI(I,L), ,. L=l,6) (Il, I4, 2E5.0,6(I3,E7.0))

CHAN Description Fuel rod type (1) Identification number of the rod Rod diameter (in) I I I Relative rod power (rod power/average rod pow,f{ Adjacent channel number Fraction of the rod power to that channel Then one card for every rod considered is required.

I 1* (1) N=l indicates rod fuel N=2 indicates plate fuel (2) This block is repeated 6 times (L=l,6) I I I 1* I I I I 1* I -, I II I I ,. 1* I I ,, I I I I I I I ,, "' -235-Card(s) Type T6 Required.to be present Fuel temperature data Only if NODESF)Q FORTRAN READ list: FORTR.t'\N FOR!-Lr..T:

Read. from Subroutine:

Variu.ble KF CF . RF DF KC cc RC TC HG Columns 1-5 6-10 . 11-15 *16-20 21-25 26-30 31-35 36-40 41-45 KF(I), CF(I), RF(I), DF(I), KC(I), CC(I), C(I), TC(I), HG(I), . (14ES. 0) .. CHAN Format* ES. 0

ES.O ES.O ES.O ES.O *Es. 0 ES.O ES.O ES.O Description Fuel thermal (Btu/hr ft °F) Fuel specific heat (Btu/lb Fuel density (lb/£t 3) Pellet diameter (inch) Clad thermal (Btu/hr ft "F) Clad specific heat (Btu/lb °F) Clad density (lb/ft 3) Clad thickness (inch) Fuel-to-clad heit transfer ficient (Btu/ft CG 8 s 8 a 8 8 8 8 8

.I -236-Card(s) Type T6a Required to be present: FORMAT READ list: FORTRAN FORMAT: Effective rod gap for interconnection between channels (in) Only if IPILE=O (GAPREC(I),I=l,NK) where NK is the total number of gap interconnections 14E5.0 ,, I I L---R_e_a_d

__ f_r __ om ___ s_u_b_r_o_u_t_i_n_e

__ : ___________

c_H_A_N ________________________________________

Variable GAPREC Notes Columns 1-70 Descriotion Effective rod gap for interconnection between channels (in) In order to give to each boundary its gap these gaps should be inputed in the same order as the boundaries are established.

Then a few words are required to know how the boundaries are established.

For the following case the boundaries are established for the code as follows: I I I I I ., 1* I I I I I *t I -I ,I I _, I I .1* I I I I I I I I ,. I I -237..:. 1 2 .. 4 5 3 6 7 8 9 -Boundary number . 1 2 34 5 6 7 8 9 10 11 12 13 Pair of channels making up each boundary 1-2 1-3 2-4 2-5 3-4 4-5 5-6 4-7 5-8 3-7 7-8 8-6 and in general the boundaries are established by going from lef.t to rig..11t in each row and frora top to bottom between two consecutive rows. 14 7-9 8-9 I I I: I -, I *a v* I* I I I I I I I ' I I I -237a-Card(s) Type: T6b Indicators for Gap Coefficients Require& to be present: Only if !PILE = 0 FORTRAN READ list: KXK IOIB FORTRAN Fon.MAT: 2!5 Read from Subroutine:

CH.AN Variable Co lum.T'ls Format Description K.XK 1-S IS Number .of cross flow resistance coefficients input KXB 6-10 IS Number of mixing coefficients input NOTE: KXK and KXB always equal 1 or the number gaps input on card T6a. If KXK equals 1, then the crossflow resistance coefficient input on card Tl7 is applied to all the gaps. If KXK equals the number of gaps, then individual crossflow resistance coefficients are input using card T6c. Simi+arly, if KY.B equals 1, then the mixing coefficient, ABETA, input on card T9 is applied to all the gaps. When KXB equals .the number of gaps, individual mixing coefficients can be input using card T6d.

I I I I I I I I I I I I I I I ' I I 1. Card(s) Type: T6c Required to be present: FORTRAN READ list: FORTHAN FORMAT: Read from Subroutine:

Variable Columns XFLOW 1-70 NOTE: -237b-Format 14ES.O Crossflow Resistance Coefficients Only if KXK>l (XFLOW(I) , I = 1, KXK.) (14E5.0) CHAN Description Crossflow resistance coeffic.i.ents The number of crossflow resistance coefficients input must equal the number of gaps input on card T6a. The order in which the resistance coefficients are input corresponds to the order in which the gaps are input.

I I I I I I l 1 I I ., I I I Card(s) Type: T6d Required to be present: FORTR.!\N READ list: FORTRAN FORMAT: Read from Subrouti.ne:

Variable Columns XBETA 1-70 NOTE: -237c-Format 14E5.0 Mixing Coefficients Only if KXl3 (XEETA(i), I= 1, KXB) (14ES. O) CHAN Description Mixing coefficients The number of.mixing coefficients input must equal the number of gaps input on card T6a. The order in which the resistance coefficients are input sponds to the order in which the gaps are input.

I I I I 1 1 I I ' I I I. a*:.* I I I I I I --------------

Card(s) Type T7 Required.to be present READ list:

FORTRl\N FORMAT: -238-PWR "Half-Boundaries" Only if IPILE=l (II(L), JJ(L), L=l, N) where II(N)=O (1415) Read. from Subroutine:

CHAN variable II JJ . Notes: "* Columns 1-70 1-70 Format* IS } . IS Descrintion II(L), JJ(L) are channel tification numbers define the °Lth "half-boundary." (1) A "half* bounda1*y 11 is one cut by a line of symmetry*.

In the example below the channel pairs defining the half-boundaries are

1 and 4 , 4 and 6

I *2 3 1 '4 5 I I . I CG. (2) The list of. "half-boundaries" is terminated by a zero. If the list finishes at the end of a card, a bi'ank card should follow to provide the (3) If there are no half-boundaries, give a blank card.

1* I I I I I I I I I 1. I I I I I I i I Card(s} Type TS Rcquired:to be present FORTRAN READ list:

FORT Rk'l FO Read. from Subroutine:

-239-Model Indicators Always Nl N2 N3 N4 NS N6 N7 N8 (14IS) MODEL Variable Columns Format Description Nl 1-S IS Mixing Indicator N2 6-10 IS Sin-gle Phase Friction Indicator N3 11-15 IS Two Phase Friction Indicator -N4 16-20 IS Void Indicator NS 21-25 IS Inlet Flow Indicator N6 26-30 IS Parameter Indicator N7 31-35 IS Iteration Indicator NB *36-40 . IS Physical Property Indicator "' *4 .. -. N9 '41-45 -I5 **--Coupling parameter

-in**-the_

.. _____ -------mixing term of the energy* Notes: equation (1) If all Nl-N8 given as zero blank card) a preset hydraulic model is obtained and the next card read is T20. If any are g1ven positive, the appropriate part of the model may be changed by giving extra card.(s). (Z) The preset model.is defined in the card-descriptions lowing for the Indicator=O.

CG

__ , -240-Card (s) Type T9 Mixing Model Required.to be present Only if Nl (on TS)> 0 FORTRAN READ list: . FORTRAN Read from Subroutine:

Variable Columns ABETA 1-5 .BBETA 6-10 Notes: ABETA-BBETA (14ES.O) MODEL Format* ES.O 1 ES.O Description The.mixing parameter Bis defined as SaABETA*(RE**BBETA}

where REsReynolds Number (1) If NlaO, ABETA, BBETA set to a.oz and 0.0 respectively.

-(2) Thermal Conduction between channels is suppTessed for all Nl. 1.....------------

I I I 1* I CG I 10 -I 10 I I I I I I *1 I 1 I I. I I I I ,. I I I I I I I I I I I I I t I Card(s) Type TlO Required.to be present FORTRAN READ list:

FORTR.t'\N FOR..'1l'i.T:

Read. from Subroutine:

-241-Single Phase Friction Only if NZ (on T8) "> 0 NVISCW, (A(J), B(J), C(J), J=l, 4) (15, 13ES.O) MODEL Variable Columns Format* Description CG NVISCW 1-5 IS *A 6-65 ES.O } B 6-65 ES.O c 6-65 ES.O Notes: .=l, if the wall viscosity correction to the single phase friction factor is required.

' =O, if not required.

The single phase friction factor is calculated as A*(RE**B)+C, where RE=Reynolds Number. (1) The.jriction factor defined by A(J), B(J),C(J) is applied to those channels with that value of J on card T2. If all nels have the same friction factor, J is given a5 1 on carJ T2 fot all channel types and only A(l), B(l), C(l}

on card TIO. (2) If NZ=O, NVISC is set to 0 and the smooth tube friction factor is used, i.e., A:0.184, B= -0.2 and C=O.O for a11 J=l,4. ' 2 2

'i I __, Card(s) Type Tll Required.to be .present READ list: FORTR.t"\N FORJ.'1.AT:

Read. !rom Subroutine:

-242-Two Phase Friction Model Only if N3 (on TS)>O J4 (1415) MODEL Variable Columns Forruat* Description I I I I I CG I J4 1-5 IS Two phase friction correlation trigger J4=0 =1 =2 =3,4 =s Note: Homogeneous Theory Armand Baroczy ,Not in use Polynomial in quality If N3=0, J4 is *set to O. ' I I I I I I I I I I I I I ,, I I I I 'I I I I I I I I I I j I -243-Card(s) Type Tl2 Two phase friction polynomial Required.to

be present Only if J4 (on Tll) = S READ list:

NF (AF(L), L=l, NF) FORTRAN FOPMAT: (!S, 13ES.O) Read. from Subroutine:

MODEL Variable Columns Format* Description CG NF 1-5 IS No. of terms in polynomial (max=7) 2 AF 6-40 ES.O Polynomial coefficients 2 Notes: (1) The two phase friction multipli.er*is calculated as ;'=NF* (AF(-t-)f*l) 1* . whe.re x = quality

'

.. ,.,' ) -244-Card(s) Type Tl3 Void Fraction Model Required.to be present Only if N4 (on TS)> 0 READ list: JZ J3 FORTRAN (14IS) MODEL Read. from Subroutine:

Variable J2 J3 JZ=O *=1 J3=0 =l =2 =3,4 =S =6 Note: Columns Format Description 1-5 .IS Subcooled Void Indicator 6-10 IS Slip Ratio Indicator no subcooled void Levy subcooled void correlation Slip Ratio=l Armand Slip Ratio Correlation Smith Slip Ratio Correlation Not in use. Slip ratio given (T14) Void fraction as a polynomial in quality (Tl4) If N4=0, JZ and J3 are both set to 0. \ .. -------*-*r-I I I I CG I 2 I z ,, I -, I I I I I I -* I I I I I I I I I I I I I I I I I I I 1* i I Card(s) Type Tl4 Required.to be present FORTRAJ.'l READ list: FORTRAN FORMAT: Read. from Subroutine:

-245-Slip Ratio Only if J3(on Tl3)=5 or 6 NV (AV(L), L=l, NV) (IS, 13ES.0) MODEL Variable Columns Format' Description NV AV 1-5 IS No.* of terms in polynomial 6-40

ES.O Polynomial coefficients A polynomial ty

,.I: tN r.. * (AY (dX) is calculated X=quali-For J3= 5, NV should be se.t to 1 and only one val:';.le of AV read in. The slip ratio is taken as AV(l). For up to 7 values of AV may be read in a.n.d the void tion is calculated as a polynomial in X, namely: NV (A Y-( y-) x ,.__,) "r=l l i CG 2 2 I

) Card{s) Type TIS Required.to be present*

READ lis.t: FORTRAN FOR:*lAT:

Read. from Subroutine:

-246-Inlet Flow Model Only if NS (on TS)> 0 IG (14IS) MODEL Variable Columns Format Descriot:ion IG IG = 0 IG. = 1 IG = 2 Note 1-S IS Inlet Flow Indicator

Inlet mass velocity*

same for all channels Inlet mass velocities for channels calculated give same inlet pressure gradient Inlet mass velocities (on Tl6) (1) If NS

O, IG set to O. I I I I I CG I 11 I I to I I I I I I I I .. I I I 1* 1 I I I I I I I I I I I I I I :I I Card(s) Type Tl 6

Required'to be present FORTRAN READ list:

FORTRA.i.'l FOR....,1AT:

Read. from Subroutine:

247-Inlet Flow Distribution Only if IG (on TlS) = Z (GR(I), I=l, NCHANL) (i4ES.O) READ IN/MODEL

Variable Colurins Format Description GR 1-70 ES.O . . Inlet Mass Velocity Ratio (local/ average) for all channels ' *. CG 11 I I I-I I I I I I I 1. I I I I I I I I. -248-Card(s) Type: Tl7 Parameters Required to be present: Only if N6>0 FORTRAN READ list: KIJ FTM SL THETA FORTRAN FORMAT: 14E5.0 Read from Subroutine:

MODEL Variable Columns Format Description

---KIJ 1-5 E5.0 Crossf low resistance coefficient, k FTM 6-10 E5.0 Turbulent momentum fac1:or', Ft SL 11-15 E5.0 Transverse momentum factor, s/'J.. THETA 16-20 ES.O Inclination of channel to vertical NOTE: If N6 = O; KIJ set equal to 0.5, FTM to 0.0, SL to 0.5, and THETA to 0.0 (i.e., vertical).

I I I 1* I I I I I I I I I I I I I i I Card(s) Type Tl8 Required.to be present* FORTRAN READ list: FORTRAN FORMAT: Read. from Subroutine:

-249-Convergence Criteria Only if Ni (on TS)> 0 NTRIES FERROR (IS, 13 ES.O) MODEL Variable.

Col um.Tis Format Description NT RIES 1-5 IS Maximum permis*s z.b le number draulic iterations FERROR 6-10* ES.O Flow convergence criterion Note (1) If N7=0, NTRIES set to 20 and FERROR to 0.001. CG of hy-9 9

\ J -250-card(s) Type Tl9 Physical Properties Required.to be present Only if NS (on TS)) 0 NPROP N PH P2 READ list: FORTRAN FORMAT: Read. from Subroutine:

Variable NPROP N PH P2 Notes Columns 1-5 6-10 11-15 16-20 (215, 2ES.O) MODEL Format* Description ES.O No. of pressure points in physical property table for interpolating between (Minimum=2, Maximurn=30).

IS a*l or 2 (see PH below) ES.O N=l, PH=lowest pressure (psia) in problem. N=2, PH=lowest enthalpy (Btu/ib) in problem, from which the lowest sure is calculated (see below). ES.O Highest pressure in problem. CG 1 (1) From this card, a table containing NPROP equi-spaced values of pressure from Pl (see below) to P2 is constructed g*iving relevant physical properties--calculated from polynomial expressions --at. each pressure.

Physical properties at intermediate pressures are :found by linear I I I I I I I I I I I I I (2) It is important that the table spans the physical property range of the problem. For example, with inlet subcoolir.g, the inlet en-* I enthalpy would correspond to a pressure lower than the reference value; the pressure would be that at which the enthalpy was the saturation value. Hence the first pressure in the table should be lower than I the value to the lowest steady state or transient enthalpy .encountered, so that the other physical properties at that enthalpy may be properly interpolated.

If N*l, PH is given as Pl1 I the lowest pressure in the problem and if N=Z, as the lowest py--the lowest pressure Pl is then calculated from PH. -(3) If N8=0, NPROP"is set to 30 and Pl, PZ calculated by the puter. I I I I I I I I I I I I I I I I I I I I I I Card(s) Type Tl9a Required to be present FORTRAN READ list: FORTRAN FORMAT: Read from Subroutine:

Variable Columns Format ENEH 1-70 E5.0 -251-Coupling parameters Only if N9 (on TS) > O (ENEH(K), K=l,NK) where NK=total number of boundaries 14E5.0 MODEL Description Coupling parameter introduce in the mixing term of the energy conservation equation.

Note: The order in which these coupling parameters should be entered is the same as the one described in card T6a for interconnection between channels.

-252-Card(s) Type T20 Steady State Operating C-0nditions Required.to be present Alliays FORTRAN READ list:

FORTRAN FORMAT: Read. from Subroutine:

Variable IH HIN GIN PEXIT Columns 1-5 6-10 11-15 16-20 IH HIN GIN PEXIT (IS, 13ES.O OPERA Format Description IS *Es.a ES. 0 ES.O Inlet Enthalpy Indicator IH=O: Inlet Enthalpy (Btu/lb) IH=l: Inlet Temperature

(°F) IH=2,3: HIN not used, set to zero (see T21}. Average 2 Inlet Mass Velocity .(Mlb/ f):. hY) System pressure (psia) ,;..;..., . I I I I I CG I 11 I .1..1. I I 11 11 I I 1. I .I I I I I I I I , *-" I I I I I I I I I I I I I I I I I I Card(s) Type T21 Required.to be present FORTRAN READ list: FORTRAN FORMAT: Read from Subroutine:

-253-Inlet Enthalpy Distribution Only if III = 2 or 3 (A(I), I=l, NCIIANL) (i4ES.O) READ IN/OPERA Variable Columns Format Description A 1-70 ES.O IH=2: IH=3: Inlet enthalnies for each channel (Btuilb) Inlet temneratures for each channel (° F) ....... -; -CG 11

-254-Card(s) Type T2Z Transient Indicators Required.to be present Always READ list: FORTRAN FOR.i.'1A':':

Read. from Subrnutine:

Variable Columns NP 1-5 NH 6-10 NG 11-15 NQ 16-20* Notes NP NH NG NQ (14! 5) OPERA Format* Descri?ticn IS No. of points at which pressure sient forcing function will be given (T23). Maximum=30 IS As NP but inlet enthalpy (TZ 4) . Ma-ximum=30 IS As NP but inlet flow (TZ 5) . Maxi-mum=30 IS As NP but channel power (T2 Sa) . Maximum=30 (1) NQ is only given in COBRA but not in MEKIN (leave NQ blank) as In MEKIN, the transient channel power is obtained from the .tronics.

(2) If only steady state calculations are required, T22 may be a blank card. I I I 1.-I CG I 11 I 11 I 11 I 11 I I I I .I I I I I I.

I I I I I I I I I I I I I* I I I I

I I -255-Card(s) Type T23 Pressure Transient Forcing Function Required.to be present Only if NP)l (T22) FORT.RAN READ list: FORTRAN FORi.'-lAT:

Read. from Subroutine:

Variable yp FP Notes Columns 1-70 1-70 (YP (I) , FP (I) , I= 1 , NP) (14ES. 0) READ IN/OPERA Format* ES. 0

ES.O Description Time (seconds)

Ratio of transient to steady pressure at time YP (1) YP(l), FP(l) should be given as 0.0 and 1.0 respectively.

(2) The value of FP at a time intermediate between two values of YP is found by linear interpolation.

CG 11 11

-256-card(s) Type T24 Inlet Enthalpy Transient Forcing Function Required.to be present Only if NH>l (T22) (YH(I), FH(I), I=l, NH) (14ES.O) FORTRAN READ list:

FORTR.1\N FOR..'1AT:

Read. from Subroutine:

Variable Columns YH 1-70 FH 1-70 Notes READ IN/OPERA Format ES.O ES.O Description Time (seconds)

Ratio of transient .to steady state enthalpy or temperature (depending on IH--card T20) at time Y.H. (1) As for but YH, FH instead of YP, FP. I I I I I CG I 11 I 11 I I I I I I I I I I I .,

I I I I , I I I I I I I I I I I I I I -257-card(s) Type T25 Inlet Flow Transient Forcing Function Required.to be present* Only if NG> 1 (T22) (YG(I), FG(I), I=l, NG) (14ES. O) FORTRAN READ list: FORTRAN FOTh"1AT:

Read. from Subroutine:

READ IN/OPERA Variable Columns Format Description YG 1-70 ES. 0

Time (seconds)

FG 1-70. ES.O Ratio of transient to steady average mass velocity at time Notes (1) As for card T23, but YG, FG instead of YP, FP. state YG CG 11 11

) -2 5 8---:_---------::-.:-Card(s) Type T25a Inlet Power Transient Forcing Function I I I Required-to be present* Only if NQ > 1 (T22) and IQP3=2 (CS) FORTRAN READ list: FORTRAN FORMAT: Read. from Subroutine:

Variable Columns YQ 1-70 FQ 1-10* Notes (YQ (I) , FQ (I) , (14ES. 0) READ IN/OPERA I=l, NQ) Fonnat Description

ES.O Time (seconds)

ES.O Ratio of transient

-::o steady channel power at time YQ (1) As for card T23, but YP, FQ instead of YP, FP. state I I CG I 11 I 11 I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I -259-Card ( s) Type T26 "Debug" Option Required.to be present Always FORTRAN READ list: KDE BUG FORTRA.i.'1 FORM.Z\T:

(1415) Read from Subroutine:

TABLES Variable Columns Format* Description CG KDE BUG 1-5 IS '.'Debug" option 9 =O: normal--no printing =l: "debug"--with 1test printing

) card(s) Type T27 Required.to be present FORTRAJ.'1 READ list: FORTRAN FORMAT: Read from Subroutine:

-260-Output Printing Always NSKIPX NSKIPT NOUT NPCHAN NPROD NPNODE (14IS) TABLES Variable Columns Format Description CG I I I I I I I I NSK I PX 1-5 I*S NSKIPT 6-10 IS NOUT 11-15 IS NP CHAN 16-20 *rs NP ROD 21-25 IS NPNODE 26-30 IS Axial print option =O or 1: every axial step printed )1 each (NSKIPX)th step Time step option As for NSKIPX but time (not axial) steps 9 9 =O: print channel results only 12 I =l: channel + cross flow tables =2: channel + fuel temperature tables channel + cross flow + fuel temperature tables I =O: all channels printe<l 'l: read in NPCHAN channels to be printed As for NPCHAN but rods instead of channels As for NPCHAN but radial fuel nodes instead of channels 12 12 12 I I I .I I I I I I I I I I I I I I I I I I I I I I I I I I -261-card(s) Type T28 Required.to be present Channels to be printed Only if NPCHAN (T27) .S 1 (PRINTC(I), I=l, NPCHAN) (1415) FORTfu"\N READ list:

FORTRAN Read. from Subroutine:

Variable Columns PI.UNTC 1-70 TABLES Format* IS Description CG Identification Number of channels 12 be printed. .. * . . : .. ...

Card(s) Type T29 Required.to.be present FORTRlili READ list: FORTRAN FORMAT: Read from Subroutine:

-262-Rods to be printed Only if NPROD (T27) 1 (PRINTR(I), I=l, NPROD) (14IS) TABLES Variable Column3 Format Descrioticn PRINTR 1-70 IS Identification Number of rods to be .printed. " I I I I I CG I 12 I I I I ,. I I I I I I I I

,. ' I I I I I I I I I I I I I I I i I -263-card(s) Type T30 Required.to be present Fuel nodes to be printed Only if NPNODE (T27) 1 (PRINTN(I), I=l, NPNODE) (1415) FORTR.Ai.'1 READ list: FORTRAN FORi."1AT:

Read .. from Subroutine:

Variable Columns PR I NTN 1-70 TABLES Format IS Description Radial fuel nodes be printed l=rod center, (NOTIESF + l)=outer clad surf ace ' . CG 12

) -264-Card(s) Type C4 Required.to be present End Input Data, start calculation Always READ list: BLANK CARD FORT&\N FORMAT: Read from Subroutine:

INDAT Variable Columns Format Description Note: At this point in the calculation, control returns to reading Card C4. If NGROUP 2 1-12, more Input Data are read in the original COBRA format, these** later data overwriti*ng what has already 'been read in. If NGROUP = O, calculation starts. I I I I I CG I c I 0 N T R I 0 L I I I I I I I I .. I I I 1* I 1* I I I I I I 1* I I I I I I I j I -265-Card(s) Type C12 Nodal Power tiplier Required.to be present Only if IQP3 (CS) = 1 ar 2. ZM FORTRAN READ list: FORTRAN FORMAT: (8E10.0) Read. from Subroutine:

QPR3 Variable Columns Foor.at Description ZM 1-10 ElO .*O Nodal Power Multiplier ZM= -2.0: Reset to 1000.0/3.6 ZM=

Reset to 3413.0/3.6 ZM >' 0. 0: ZM unchanged to Btu/s) (MW to Btu/ s) The nodal powers given on cards Cl3, Cl4 are all m.ml tip lied by ZM. This allows, for example, units to be converted.

CG

_,--266-Card(s) Type Cl3 Fuel Nodal Powers Required.to be present Only if IQP3 (C8) = 1 or 2 ((QF(I,J), J=l, NDX), I=l, NCHANL) (8El0.0)

READ list!

FORTRAN FORMAT: Read_ from Subroutine:

Variable Colum..""ls QF 1-80 QPR3 Format

_Descriotion Average Fuel Nodal Power for nel I, axial interval J to (J+l) The power for each channel I (I=l, NCHANL) is read in turn. Each channel-set, i.e., J=l,*NDX, starts on a new card, continuing onto the next card if NDX > 8. The uni ts of QF in the calculation.

are Btu/sec. They may be read in those units (when ZM=l.O on Cl2) or converted using ZM. NDX is read on card cii.-* -1 I* I I I CG I I I I I I I I I I I I I I I I I I I I I I I I,., I I I I I 1. I I I -267-Card(s) Type Cl4 Coolant Nodal Powers Required.to be present Only if IQP3 (CS) = 2 FORTRAN READ list:

FORTRl\N FORMAT: Read. from Subroutine:

Variable Columns QC(I,J) 1-80 ((QC(I,J), J=l, NDX), I=l, NCHANL) (8El0.0) QPR3 Format Description Average Nodal Power deposited in Coolant for channel I, axial inter,val J to J+l. As for card Cl3, but QC instead of QF. CG

) -... * .) -268-Card(s) Type Cl3 Transient Fuel Nodal Po11J.er Required.to be present Only if IQP3 = 1 or 2 and NDT> 1 ((QF(I,J), J=l, NDX), I=l, NCHANL) (8El0.0) FORTRAN READ lis.t: FORTRAN FORMAT: Read. from Subroutine:

QPR3 Variable Columns Format* Description Cards Cl3 and (if IQP3=2) Cl4 are read for the first transient time step, then both sets of cards for the next time step, etc. until data for* all time steps have been given. ., . .. I l I I I I CG I I I .I I I I I *I I I I I I 1-... ---*-----*-****-*-**.

-269-,. .........

,. 1* I I 1. I. I I I I I I I I -1 I) I. Card(s) Type Cl4 Required.to be present READ list: FORTRAN FORMAT: Read. from Subroutine:

Transient Coolant Nodal Power Only if IQP3= 2 and NDT) 1 ((QC(I,J), J*l, I=l, NCHANL) (8E10.0) QPR3 Variable Columns Format* Description See last card, "transient" Cl3. .. ' CG

) , -27f}-Card(s) Tvoe C3 Required.to be Always FORTRAN READ list: FORTRA.."I FORM}i..T:

Read from Subroutine:

Variable Columns IP ILE. KASE Jl TEXT (Il, I4, IS, 17A4) INDAT Format Description See earlier C3

Note *At the end of the calculation, control returns again to the read . statement for card C3. ' If KASE> 0; the next case is read. If KASE =O (e.g., a blank card), calculation stops. I I I I I CG I I I I I I I I *I -I I I I I NOTICE -THE ATTACHED FILES ARE OFFICIAL RECORDS OF THE DIVISION OF DOCUMENT CONTROL. THEY HAVE BEEN CHARGED TO YOU FOR A* LIMITED TIME PERIOD AND MUST BE. RETURN.ED TO THE RECORDS FACILITY BRANCH 016. PLEASE DO NOT SEND DOCUMENTS CHARGED OUT THROUGH THE MAIL. REMOVAL OF ANY PAGE(S). FROM DOCUMENT FOR REPRODUCTION MUST BE REFERRED TO FILE z:t) '.2.2i \ . DEADLINE RETURN DATE tiY l lJ r\8-t-C\ RETURN m* REACTOR FltES !fq l02l\-o5q,LL * , .** ' REcdRDS FACILITY.BRANCH

ML18130A261

Contents

Text

=============================

SUMMARY

SUMMARY

Navigation menu

ML18130A261

Text

=============================

SUMMARY

SUMMARY

Navigation menu

Search