A1110S 4bl53t

3291

M.w wimam B/ Building emulation compui ScVoO UM NO 85-32?1 V1985 C.1 NBS-PUB- Building Emulation Computer

Program for Testing of Energy

Management and Control System

Algorithms

Reference

NBS

PUBLICATIONS

William B. May, Jr. Cheol Park

sa- . NS’ Of STANDARDS & TECH FU.C

All 102461236

U.S. DEPARTMENT OF COMMERCE

National Bureau of Standards

National Engineering Laboratory

Center for Building Technology

Building Equipment Division

Gaithersburg, MD 20899

December 1985

Sponsored by:

Naval Civil Engineering Laboratory

U.S. Department of Defense

QC ce of Building and Community Systems

100 i Department of Energy

U56

85-3291

1985

NBSIR 85-3291

BUILDING EMULATION COMPUTER

PROGRAM FOR TESTING OF ENERGY

MANAGEMENT AND CONTROL SYSTEM

ALGORITHMS

NBS HZSEA2CH INFORMATION

Q C

( o

uC<<

William B. May, Jr. Cheol Park

U.S. DEPARTMENT OF COMMERCE

National Bureau of Standards

National Engineering Laboratory

Center for Building Technology

Building Equipment Division

Gaithersburg, MD 20899

December 1985

Sponsored by:

Naval Civil Engineering Laboratory

U.S. Department of Defense

Office of Building and Community Systems

U.S. Department of Energy

U.S. DEPARTMENT OF COMMERCE, Malcolm Baldrige, Secretary

NATIONAL BUREAU OF STANDARDS. Ernest Ambler. Director
'
ABSTRACT

A building emulator can be used to test energy management and control systems

(EMCS). The emulator uses a computer program to simulate the responses of a building including the equipment, building space, and building envelope to

EMCS commands. Building model software for the emulator has been developed at

the National Bureau of Standards (NBS) in an effort to assist the United

States Naval Civil Engineering Laboratory (NCEL), which is developing a sophisticated building emulator.

The concept of the building emulator and the building emulator computer

program are described in this report. The program includes the weather, the

air handling unit, the zone, and the comfort model. In addition, the energy

compilation routine is also included. The models presented here are

simplified models. With these abridged models, a single zone building with

exterior walls and a single deck air handling unit are simulated. A complete

FORTRAN source code of the building emulator computer program is appended.

KEY WORDS: air handling unit simulation model; building emulator; building

space zone model; computer simulation; EMCS algorithms; energy management and

control systems (EMCS); HVAC emulator; local equipment simulation model;

testing EMCS algorithms.

iii
ACKNOWLEDGEMENT

The authors are indebted to the Naval Civil Engineering Laboratory, Department

of Defense, and the Building Systems Division, Department of Energy, for

funding this project.

IV
TABLE OF CONTENTS

Page

ABSTRACT iii

ACKNOWLEDGEMENT iv

LIST OF FIGURES vi

LIST OF TABLES vii

CONVERSION FACTORS FROM METRIC (SI) TO ENGLISH UNITS viii

1. INTRODUCTION 1 1.1 The Concept of a Building Emulator 1 1.2 Existing and Planned Building Emulators 2

1.3 Testing of EMCS Algorithms with the NCEL/NBS Building Emulator.. 2

1.4 Software Components of the NCEL/NBS Building Emulator 3

2. BUILDING SIMULATION COMPUTER PROGRAM 5 2.1 Data Structures and Interface to Building Emulator 5 2.2 General Program Organization and Control

2.3 Initialization Routines 10

2.4 Simulation Control Routines 15

2.5 Building Component Model Routines 15

3. WEATHER MODEL 18

3.1 Assumptions and Purpose of Weather Model 18

3.2 Outside Air Temperature 18

3.3 Solar Energy 18

4. AIR HANDLING UNIT MODEL 21

4.1 Assumptions 21

4.2 Equations 21

5. BUILDING SPACE ZONE MODEL 34

5.1 Assumptions.... 34

5.2 Fundamental Equations 35

5.3 Solution of Equations.......... 38

5.4 Equations for Heat and Mass Transfer Coefficients and Terms 43

5.5 Property Equations for Moist Air 53

5.6 Local Equipment Model 53

6. OCCUPANT COMFORT MODEL 57

6.1 Assumptions and Purpose of Comfort Model 57

6.2 Comfort Envelope Equations 57

7 . COMPILATION OF ENERGY USE 58

7.1 Categories of Energy Use 58

7.2 Equations for Energy Use 59

8. SUMMARY 62

9. REFERENCES 63

APPENDIX. COMPUTER PROGRAM LISTINGS 65

v
LIST OF FIGURES

Figure 2-1. Logic diagram for emulator model program 9

Figure 4-1. Block diagram of components of generalized HVAC system

model 22

Figure 4-2. Block diagram of components of the abridged, single deck

air handling unit model 23

Figure 4-3. Air handling unit model calculation sequence 26

vi
LIST OF TABLES

Table 2-1. Emulator model state vector - Analog sensor variables 6

Table 2-2. Emulator model state variables - Energy variables 7

Table 2-3. Emulator model state variables - Command variables 7

Table 2-4. Emulator model state variables - Comfort variables 8

Table 2-5. Sample simulator configuration file 11

Table 2-6. Sample state vector initialization file 12

Table 2-7. Sample building use parameter file 13

Table 2-8. Sample building shell parameter file 14

Table 2-9. Sample HVAC system parameter file 14

Table 2-10. Sample weather data parameter file 16

Table 4-1. Nomenclature for the air handling unit model 28

Table 5-1. Nomenclature for equations in the building space zone model 36

Table 5-2. General form of the FORTRAN function for a coefficient

function 43

Table 5-3. Occupant sensible and moisture gains 46

Table 5-4. Moist air property routines 54

Table 7-1. Heating and cooling energy from emulator model 61

vix
CONVERSION FACTORS FROM METRIC (SI) TO ENGLISH UNITS

Physical

Characteristic

To Convert from To Multiply by

Length m ft 3.28

Area m

2 ft 2 10.76

Velocity m/ s mph 2.24

Temperature °C °F t F

= l*8t

c + 32

Temperature

difference °C °F 1.8

Energy J Btu 0.948 X 10‘3

Power W Btu/hr 3.41

U-value W/m2 °C Btu/hr*f

2 *°F 0.176

Thermal

resistance m

2 °C/W hr *f

2 * °F/Btu 5.678

Pressure kPa in. Hg 0.296

viii
1. INTRODUCTION

Advances in computer and communications technology have made feasible the use

of computers for controlling the heating, ventilation, and air conditioning

(HVAC) systems of buildings. Such control systems are usually called energy

management and control systems (EMCS). As for any building equipment, it is

desirable to have standard testing and rating methods for EMCS. One benefit of

such standards is that a more intelligent choice of an EMCS for a building can

be made if systems from several manufacturers are tested using standard

methods. Once an EMCS is installed, acceptance testing of the system is

possible by comparing field test results to design specifications.

1.1 The Concept of a Building Emulator

One approach to testing of EMCS systems is the use of a building emulator.

A building emulator is a special computer connected to the sensor inputs and

command outputs of an EMCS. No actual HVAC equipment or building is connected

to the EMCS. The emulator uses a computer program to simulate the responses of

a building (including the equipment, building space, and building envelope) to

EMCS commands. The EMCS then controls this simulated building as if it were an actual one. The advantages of using an emulator are that an EMCS may be

tested with any type of building if a model is available, and tests can be

repeated on different EMCS under the same conditions. Although experience with

emulators is limited, it appears that the emulator testing approach is the

most desirable in terms of the cost of testing, reproducibility of results,

and control of test conditions.

A building emulator can be used to test EMCS algorithms, which may be defined

as the software programs in the EMCS which use the EMCS sensors to collect

information, make decisions based on current and past information and

parameters, and make changes to the building through the EMCS actuators. EMCS

algorithms can either be direct digital control (DDC) to control valves and

dampers directly, for example, or supervisory to determine start/stop times

and setpoints, for example. If a building emulator is used to test EMCS

algorithms, it is not necessary to know the exact structure of the algorithm

software. Since EMCS software is often proprietary, this is a definite

advantage.

Before a building emulator can be used, it must be possible to connect an emulator to an EMCS to be tested. The recommended approach is the connection

of the emulator to the EMCS in place of the regular EMCS sensors and

actuators. Since it is likely that the EMCS sensors are electrical in nature,

they can be replaced with voltage or current sources under the control of the

emulator. The EMCS will be unable to distinguish between an actual sensor

producing an electrical signal and the emulator producing the same electrical

signal. This method allows the EMCS algorithm to be tested directly without

the effects of sensors. This implies that a building emulator cannot be used

1
to test the instrumentation on an EMCS since no sensors or actuators are

involved.

1.2 Existing and Planned Building Emulators

Two building emulator prototypes have been built in order to demonstrate the

concept of EMCS testing using an emulator. At the United States National

Bureau of Standards (NBS), a simple building emulator was constructed in 1983.

This emulator demonstrated the basic concept of a building emulator but was

not transportable, was limited in capacity, and had a model which was too

simple for quantitative algorithm testing.

A second example of a prototype emulator was developed by a contractor for the

United States Naval Civil Engineering Laboratory (NCEL) in 1983 [Wise (1984)].

This device, which was referred to as an EMCS exerciser, was designed to be

transportable and was used in field tests of actual EMCS installations.

Complex software in the exerciser allowed an operator to easily prepare for

and run tests. The tests were intended to determine if EMCS supervisory

algorithms were operating properly. The exerciser produced reports of the time

of occurrence and nature of all EMCS actions during the test. These reports

were compared to reports from the EMCS to determine if any problems existed.

The unit was not designed to allow a quantitative testing of algorithms, only

to determine if the algorithms were operating and controlling the correct

command outputs. The exerciser did work as designed and was capable of

revealing areas requiring attention in real EMCS installations.

NCEL, with the assistance of NBS, is currently developing a building emulator

of greater sophistication than the examples given above. This emulator will

have the field transport and general connection features of the NCEL exerciser

plus a detailed building simulation model and improved control and evaluation

software. NBS is responsible for developing the building model to be used in

the NCEL emulator. This report is intended to describe the building model

software which has been developed by NBS. The NBS building model does not

simulate a complete building, but simulates a single building zone with

exterior walls and an air handling unit.

The NCEL/NBS emulator will be used for either factory or field acceptance

tests of EMCS algorithms. A factory test will be performed at an EMCS

manufacturer facility to determine the capabilities of algorithms on a specific EMCS. A field acceptance test will be used to determine whether an EMCS installed at a specific site is operating in accordance with design

specifications. The initial version of the emulator will be designed to test

supervisory EMCS algorithms and will not test DDC algorithms.

• 1.3 Testing of EMCS Algorithms with the NCEL/NBS Building Emulator

Eight types of supervisory control algorithms have been selected for testing

2
with the NCEL/NBS building emulator. The concepts of the algorithms will not

be discussed in this report. Discussion of these EMCS algorithms may be found

in other reports [May (1983)* Park (1983)* Park, Kelly* Kao (1984), Park

(1984), and May (1984)]. The types of algorithms to be tested are listed

below:

1. Scheduled Start/Stop

2. Duty Cycling

3 . Demand Limiting

4. Optimum Start/Stop

5 . Day/night Setback

6. Economizer

7. Ventilation/Recirculation

8. Air Setpoint Reset

The emulator model described in this report will not currently support the

testing of the demand limiting algorithm. A model of electrical demand in a building would have to be added to the emulator simulation software.

The test procedure to be used for testing EMCS algorithms will determine the

exact characteristics that the NCEL building emulator will have. The final

form of the test procedure is not known as of the writing of this report. The

test procedure cannot be determined entirely by objective means and must

result from negotiations between building owners* equipment manufacturers,

testing personnel* and standards organizations. For the NCEL emulator it will

be assumed that the test procedure will be based on the concept of testing of

computers using a 'benchmark* computer program which exercises various

operations of any computer on which it is run. Like a program which runs at

different speeds depending on the computer on which it is executing, an EMCS

algorithm will behave differently depending on the building and climate where

it is tested. If 'benchmark* or 'reference' buildings are used for a building

emulator, then results of testing different algorithms of the same type on

these buildings may be compared.

The building simulation program in the emulator will simulate the response of

a building to EMCS control actions. The emulator will be used to provide

reproducible predictions of energy use and temperature levels for a set of

hypothetical 'reference' buildings which are being controlled by the EMCS

system under test.

1.4 Software Components of the NCEL/NBS Building Emulator

The NCEL/NBS emulator computer will employ a multi-tasking disc operating

system which will allow several software tasks or processes to effectively

appear to be running concurrently. When the emulator is being used to test an EMCS algorithm, there will be two important processes running. One of these is

the building simulation program which is -described in this report. The other

process is an input/output (I/O) program which will read from and write to the

3
hardware connected to the EMCS sensors and actuators. These two processes will

have a common data area in computer memory that both may read from and write

to.

If the EMCS sends a command by changing the state or signal level of an EMCS

command output, the I/O program detects this and changes the data area that

the simulation program uses for command input variables such as start/stop

fan. or supply air setpoint. Based on the input variables and model

parameters, the simulation program predicts values for variables such as temperatures in the building and HVAC equipment and status points. The

simulation program places these values in the common memory area. The I/O

program reads these values from the common area and sets the signal levels on

the EMCS sensor inputs.

In addition to the model and I/O software, the NCEL emulator will have control

software. An operator control program will be used to start and stop the

simulation and I/O programs, select parameter files to be used, and determine

duration of testing. At the end of tests the simulation program will deposit

accumulated energy use by the emulated building into data files. Evaluation

software will be used to examine test results and produce evaluations of

specific EMCS algorithms such as scheduled start/stop. demand limiting, or setpoint reset.

The next section of this report describes the model software for the emulator

in detail. No additional information will be given on the other software

components -of the emulator since these were not developed at the time this

report was written.

4
2. BUILDING SIMULATION COMPUTER PROGRAM

The building emulator model computer program that was completed by NBS was

designed to simulate any type of HVAC system and all possible modes of heat

transfer in a zone. However, not all subroutines were completed in order to

create a more simplified version. This simplified version will be termed the

abridged version. The extent of simplifications will be noted in the sections

following.

The building emulator building simulation computer program was written

completely in the high level computer language FORTRAN 77 and is appended in

the APPENDIX. This section will describe the actual computer program that was

developed but will not present the specific equations used to model the

building. The emphasis in this chapter will be on describing the data

structures, data interface, programming logic and component parts of the

program. Following chapters will present detailed equations for the component

models of the simulation program.

2.1 Data Structures and Interface to Building Emulator

The building emulator model program is connected to the other programs in the

emulator by three data paths. These are the common data area, parameter input

files, and report output files.* In terms of the building model, the parameter

input files control duration of the simulation, timestep, initial conditions,

and configure the building systems to be simulated. The report files are the

program output, and the common data area represents the state of the

simulation and allows input of time dependent boundary conditions.

The common data area is implemented in the emulator model program as a FORTRAN

common block. Within the model program, this common block is given the name

STATE1 and contains 100 real variables and 100 logical variables. The common

block STATEl contains the 'state vector' for the model. This means that at any

time, the values of the variables in STATEl describe the state of the

simulation. Tables 2-1 through 2-4 list all of the variables in the 6tate

vector and the engineering units used. The number on the left is the position

of the variable in the state vector array. Not all of the state vector

positions are used to allow for any future expansion of the model. The

variables in table 2-1 are all set by the emulator model program. Selected

values in the common area will be converted to electrical signals on the EMCS

inputs

The variables listed in table 2-2 are energy variables which contain

cumulative energy used by the building in various categories. This information

is written to report files at the end of an emulator test run. More

information on energy variables is given in section 7.

5
Table 2-3 lists variables vhich the emulator model program uses as input

(independent variables). The values of these variables in the common area are

set depending on vhat commands have been received from the EMCS. Some of the

variables are logical and some are real variables. The logical variables have

two states for use in controlling a two-state device such as a fan. The real

variables are used to determine values such as setpoints.

The variables in table 2-4 are set by the emulator program but are not

converted to signals for the EMCS. These are comfort variables vhich are

written to a report file at the end of a test period.

Table 2-1. Emulator model state vector - Analog sensor variables

01: TEMPERATURE of outside air (F)

02: HUMIDITY RATIO of outside air (lbv/lbda)

03: ENTHALPY of outside air (Btu/lb)

04: Windspeed (MPH)

06: TEMPERATURE of zone air (F)

07: HUMIDITY RATIO of zone air (lbv/lbda)

08: TEMPERATURE of zone wall interior (F)

09: TEMPERATURE of Zone wall surface (F)

10: TEMPERATURE of zone glass (F)

20: TEMPERATURE of return air (F)_

21: ENTHALPY of return air (BTU/lb)

22: HUMIDITY RATIO of return air (lbv/lbda)

23: TEMPERATURE of mixed air (F)

24: ENTHALPY of mixed air (BTU/lb)

25: HUMIDITY RATIO of mixed air (lbw/lbda)

26: TEMPERATURE of cooling coil discharge air (F)

27: ENTHALPY of cooling coil discharge air (BTU/lb)

28: HUMIDITY RATIO of cooling coil discharge air (lbw/lbda)

29: TEMPERATURE of heating coil discharge air (F)

30: ENTHALPY of heating coil discharge air (BTU/lb)

31: HUMIDITY RATIO of heating coil discharge air (lbw/lbda)

32: TEMPERATURE of supply air (F)

33: ENTHALPY of supply air (BTU/lb)

34: HUMIDITY RATIO of supply air (lbw/lbda)

35: ENTHALPY of return air plus ventilation air(BTU/lb)

51: MASS FLOW RATE of return air (lb/s)

52: MASS FLOW RATE of supply air (lb/s)

53: MASS FLOW RATE of ventilation air (lb/s)

60: HEATING POWER to heating coil (Btu)

61: COOLING POWER from cooling coil (Btu)

62: HEATING POWER to local zone equipment (Btu/s)

63: COOLING POWER from local zone equipment (Btu/s)

6
The input parameter files are read by the emulator model program during the

inital izat ion of the program. These files are described in section 2.3. The

input parameters files are either model parameter files or test control

parameter files.

2.2 General Program Organization and Control

The emulator building model computer program is written in the form of a FORTRAN subroutine which can be called from the main program section of the

emulator model process running in a multi-tasking system (see section 1.4).

The main logic diagram for the emulator model is shown in figure 2-1. It is

assumed that the emulator subroutine is called repeatedly without any time

delays between calls.

Table 2-2. Emulator model state variables - Energy variables

64: HEATING ENERGY (requirements) (Btu)

65: COOLING ENERGY (requirements ) (Btu)

66: HEATING ENERGY (economizer) (Btu)

67: COOLING ENERGY (economizer) (Btu)

68: HEATING ENERGY (load) (Btu)

69: COOLING ENERGY (load) (Btu)

70: HEATING ENERGY (reheat) (Btu)

71: COOLING ENERGY (recool) (Btu)

72: ELECTRICAL ENERGY to fans (kWh)

73: HEATING ENERGY to local zone equipment (Btu)

74: COOLING ENERGY from local zone equipment (Btu)

Table 2-3. Emulator model state variables - Command variables

80: CONTROL

81: CONTROL

82: CONTROL

83: CONTROL

Logical 01:

Logical 02:

Logical 03:

Logical 04:

Logical 10:

Logical 11:

POINT ADJUSTMENT supply air temperature (P)

POINT ADJUSTMENT zone temperature (F)

POINT ADJUSTMENT outside air damper position

POINT ADJUSTMENT supply air RH (2)

COMMAND air handling unit on/off

COMMAND economizer on/off

COMMAND ventilation on/off

COMMAND zone setback on/off

STATUS air handling unit on/off

STATUS zone occupied/unoccupied

(F)

7
Table 2-4. Emulator model state variables - Comfort variables

90: COMFORT time outside of comfort range (S)

91: COMFORT maximum dry bulb temperature (F)

92: COMFORT minimum dry bulb temperature (F)

93: COMFORT relative humidity at maximum temperature (%)

94: COMFORT relative humidity at minimum temperature (Z)

The emulator program consists of three parts. These are the model equations

section* the control section* and the initialization section. During the

initialization of the program* a simulation time step* for example thirty

seconds* is selected. In a conventional simulation* the control section of the

program would advance the variables which contain the time of day and cause

the model equations to be solved. The conventional simulation will execute at

the full speed of the computer. For an emulator* however* the simulation is

artificially slowed so that the time intervals as measured by a real time

clock between solutions of the model correspond to the value of the simulation

time step. The control software for an emulator must obtain the current time

of day from a real-time clock in the computer and only allow the model

equations to be solved when a real amount of time has passed.

8
Figure 2-1. Logic diagram of the emulator model program

9
In addition to monitoring a real-time clock for the passage of the time-step,

the emulator building model program control section must also determine if a command event has taken place before the time- step interval has elapsed. If a command event occurs, the program must temporarily adjust the time step to be

equal to the time interval that has passed since the most recent solution of

the model equations and cause the model equations to be solved immediately.

If this is not done, there may be unacceptable delays in the response of the

emulator to commands from the EMCS.

Each time that the emulator model subroutine is called, a determination is

made as to whether initialization has taken place. If it has not, the program

initializes. This process is described in section 2.3. If initialization has

already occurred, the program control section determines if either the

timestep interval has passed, or a command event has taken place. If neither

of these conditions has occurred, then the subroutine is exited. If either of

the conditions has occurred, then the model equations section causes the model

equations to be solved for the current time of day.

After either initialization or model equations solution, the program

calls the energy compilation subroutine. If called after equations solution,

the energy compilation subroutine determines the energy used by the emulated

building since the last model equations solution and stores a cumulative total

of energy used since the beginning of the emulator model program execution. If

called after initialization, the cumulative energy totals are reset. After the

energy compilation if finished, the emulator model subroutine is exited.

2.3 Initialization Routines

The initialization section of the emulator model program loads simulation

control parameters, model parameters, and initializes time variables.

2.3.1 Control Parameter Initialization

The control parameter input subroutine, named CONPAR, is called during the

initialization phase to load control parameters. The control parameters are previously stored in a file called the configuration file. CONPAR opens the

configuration file and reads the control parameters from it. Table 2-5 shows a sample configuration file. All reading from the file is performed using

FORTRAN list-directed read operations and therefore quotes must be placed

around the character information, but otherwise no formatting or special

delimiters are required. The entries on the right are strictly comments and

are not read by the program.

The first entry in the configuration file is the start time of the test using

the building emulator. In the example configuration file, this is 12:30. The

next entry is the duration of the test in days, hours, minutes, and seconds.

In table 2-5, the example test duration is one day and three hours, or a total

10
of 27 hours. The next four entries in the configuration file are the names of

files containing the model parameters. The contents of these files are

discussed in the following sections. The names of these files would normally

be changed to select a different set of model parameters to allow a different

building to be simulated* or the same building to be simulated in a different

climate. The seventh entry in the configuration file is the name of the file

used to initialize the state vector. The eighth configuration file entry is

the weather 'season' to use for the simulation. This will be discussed more

fully in section 2.3.2.5. The last entry in the configuration file is the

simulation time-step. The time-step is discussed in section 2.2 and section

2.4. The smaller the timestep* the greater the accuracy of the results, but

the program execution time rises quickly with decreasing time step. A larger

timestep decreases execution time up to a point. Beyond this point the

execution time does not change.

Table 2-5. Sample simulation configuration file

12 30 0

13 0 0

'HVAC.PAR'

'USE. PAR'

'CLIMATE. PAR'

'SHELL. PAR'

'WINTER. INI'

'WINTER

20 ‘

'start time: hours* minutes* seconds'

'execution time: days* hours* minutes* seconds'

'BUILDING HVAC SYSTEM PARAMETER FILE’

'BUILDING USE PARAMETER FILE'

'WEATHER DATA PARAMETER FILE*

'BUILDING SHELL PARAMETER FILE'

'STATE VECTOR INITIALIZATION FILE'

'SEASON TO USE FOR EMULATION'

'TIMESTEP IN SECONDS'

After the control parameters have been read* the CONPAR subroutine calculates

the absolute ending time for the emulator test run by adding the test duration

to the start time. The ending time is used to determine when to halt the

execution of the emulator model program and write the report files.

2.3.2 Model Parameter Initialization

After the control parameters are initialized* the initialization section of

the emulator model program calls the subroutine named PARAMA. This subroutine

is used to initialize all model parameters. The model parameters are stored in

several FORTRAN common blocks which are installed in the PARAMA subroutine and

in the various model equation solution subroutines.

2.3. 2.1 State Vector Initialization

The first function performed by the PARAMA subroutine is the initialization of

the state vector. The variables in the building emulator simulation are

located in a large array in a FORTRAN common block. This state vector contains

100 real variables and 100 logical variables. During the solution of the model

equations* some of the state variables are calculated completely from the

11
current values of other variables and do not need to be initialized. Other

state variables are calculated from the previous value of the state variable

one time step back as well as from the current value of other variables. This

latter type of variable must be given an initial value. The state vector

variables are initialized using values read from the state vector

initialization file. Table 2-6 shows a sample initialization file.

Table 2-6. Sample state vector initialization file

06 52.0 ' TEMPERATURE of zone air (F)'

07 0.008 'HUMIDITY RATIO of zone air (lbw/lbda)'

08 21.0 'TEMPERATURE of zone wall interior (F)'

09 48.0 ‘TEMPERATURE of Zone wall surface (F)’

10 73.0 ’TEMPERATURE of zone glass (F)‘

80 81.0 'CONTROL POINT ADJUSTMENT supply air temperature ( F

)

81 55.0 'CONTROL POINT ADJUSTMENT zone temperature (F)

83 50.0 'CONTROL POINT ADJUSTMENT supply air RH (%)

92 200.0 'COMFORT minimum dry bulb temperature (F)’

01 .FALSE. ’COMMAND air handling unit on/off’

02 .FALSE. 'COMMAND economizer on/off'

03 .FALSE. 'COMMAND ventilation on/off'

04 .FALSE. 'COMMAND zone setback on/off’

10 .FALSE. 'STATUS air handling unit on/off

There are three entries on each line of the initialization file. The first

entry on the line is the position in the state vector to be initialized. This

is a number from 1 to 100. It is assumed that the real numbers are initialized

first. The lines in the file must have state vector position numbers in

ascending numerical order. When a position number for a line is lower than the

number in the preceding line* it is assumed that the following lines will now

be used to initialize logical variables.

The second entry on the lines in the initialization file is the initial value

to be assigned to the state vector variable. This is a real number constant

for the real state variables and a logical constant for the logical state

variables. The third entry on the initialization file lines is a comment

identifying the name of the state variable. This is not read by the program.

2. 3. 2. 2 Building Use Parameter Initialization

The second function performed by the PARAMA subroutine is the initialization

of building use parameters. The building use parameter file contains

information on the occupancy and equipment operating times for the building.

An example of a building use parameter file is listed in table 2-7.

In the building use parameter file each line is used to determine a change in
the usage status of the building. Between the times indicated by the entries,

the usage status is assumed to remain the same. The first column indicates the

time in decimal hours that the change in building usage is to occur. The

second column indicates the level of lighting in the space in watts. The third

column determines the type of lighting used in the building. *F’ indicates

fluorescent, and 'I* indicates incandescent. The fourth column is the power

used by equipment in the zone in watts. The fifth column is the ratio of

radiated energy to total energy for the equipment. The sixth column is the

zone moisture gain rate due to the equipment in pounds of water per hour. The

last five columns are counts of the number of occupants in one of five

activity levels with Activity 5 being the most active. The activity levels are described in section 5.4. 5.2.

Table 2-7. Sample building use parameter file

HRS |LIGHT(W) | TYPE | EQUIP, SENS (W) |EQUIP,LAT(LB/H) |ACT 1 ACT2 ACT3 ACT4 ACT

0.0 100. »

f 0 . 0.5 0 . 0 0 0 0 0

8.0 500. »F» 200. 0.5 0.25 16 8 4 0 0

17.0 100. *F* 0 . 0.5 0 . 0 0 0 0 0

After the parameters have been read from the file. The PARAMA subroutine

performs some simplifying calculations. The purpose of the calculations is to

determine the sensible heat gain rate from equipment, the sensible heat gain

rate from people, the moisture gain rate from equipment, and the moisture gain

rate from people for each usage state. The equations used are described in

section 5.4.5.2. The calculated building use parameters are stored in the

common block USE in five arrays. The five arrays contain the above calculated

parameters and the time of occurrence for each usage state. A parameter in the

program determines the maximum number of usage states.

2. 3. 2. 3 Building Shell Parameter Initialization

After initializing the building use parameters, the PARAMA subroutine

initializes the building shell parameters. The building shell parameter file

contains parameters which describe the characteristics of the zone, zone

walls, and zone windows, for the building to be emulated. Parameters

describing the thermal capacity and insulation value of the zone enclosure are in this file. Also in the file are parameters for air infiltration and

moisture absorption. An example of the contents of the file is given in table

2-8. The meaning of the entries is assumed to be explained by the comments on. each line and will not be discussed here. Only the first value on each line is

read by the program. No additional calculations are performed on the building

shell parameters in subroutine PARAMA. The parameters are stored in a FORTRAN

common block with the name SHELL.

13
2. 3. 2. 4 HVAG System Parameter Initialization

After initializing the building shell parameters* the PARAMA subroutine

initializes the building HVAC system parameters. The HVAC System parameter

file contains a description of the heating and cooling equipment available in

the building to be emulated. An example of this parameter file is given in

table 2-9. The meaning of the entries is assumed to be explained by the

comments on each line. Only the first value on each line is read by the

program.

Table 2-8 . Sample bui:

1000. ' FT~ 2

4.05E-4 *BTU/S*F*FT~:

1.1IE-3 'BTU/S*F*FT~:

0.282 'FT'

0.695 'FT'

3.00E-2 'BTU/S*F*FT'

7.12E-5 ' BTU/ S*F*FT'

100.0 'BTU/F'

25000. ' FT~3

0.19 'BTU/LB*F'

130.0 *LBM/FT~3

0.22 'BTU/LB*F'

52.8 ' LBM/FT 3'

0.5 'ACH'

3.0

Table 2-9 . Sample HVA<

68.75 'ft 3/sec.'

15.0 'kW'

10.0 *kW'

'REHEAT'

05.0 'Btu/ sec .

0.40000 'Btu/sec. F

0.00010 'Btu/s2 F' 60.0 ' sec.

.FALSE.

.TRUE.

0.05 *Btu/sec.

0.10

.FALSE.

0.06 * Btu/ sec.

.FALSE.

10.0 * Btu/ sec.

15.0 'Btu/ sec .

20. *F*

'Area of the building wall exposed to the outside*

'inner wall film coeff. (STILL AIR, VERTICAL) R0.68'

'outer wall film coeff. (7. 5 MPH WIND) R0.25'

'thickness of wall outside the mass center'

'thickness of wall inside the mass center'

'wall conductivity > the mass center (BRICK)

'wall conductivity < the mass center (MASONRY)'

'thermal capacity of zone furnishings (300LB WOOD)'

'volume of air in zone (50x50x10 ft) ' 'specfic heat of wall outside the mass center'

'density of wall outside the mass center'

'specfic heat of wall inside the mass center'

'density of wall inside the mass center'

'base air infiltration rate without wind and thermal'

'furnishings moisture absorption mass relative to air'

'Air handling unit fans rated supply air volume'

'power requirement of supply fan at rated volume'

'power requirement of return fan at rated volume'

'Local equipment type'

'Local equipment capacity* heating'

'Local equipment controller proportional gain'

'Local equipment controller integral gain'

'Local equipment controller sensor time constant'

'Variable Air Volume System?'

'Return Air Fan in system?'

'Return Air Fan air heating rate'

'Minimum ventilation air, fraction of supply air mass'

'Dual Deck System?*

'Supply Air Fan air heating rate'

'Humidity Control?’

'Cold Deck Heating Coil Capacity*

'Cold Deck Cooling Coil Capacity'

'Local equipment controller setpoint setback'

14
2. 3. 2. 5 Climate Parameter Initialization

After initializing the building HVAC system parameters* the PARAMA subroutine

initializes the building climate parameters. The building climate parameter

file contains the parameters used to simulate the weather to which the

emulated building is exposed. The information in the file is for a representative day at five different times of the year for a particular

climatic area. These times of the year are termed 'seasons' and correspond to

winter* late winter* spring* early summer* and summer. During any particular

testi parameters for one season* which is selected by an entry in the

emulator model configuration file* will be used to generate a one day

variation in weather conditions. This day of weather will be repeated for any

other days that the test is run. An example climate file is given in table 2-

10. The first line in the file identifies the climate description. This is

followed by the five groupings of parameters for five seasons. The season

title is the first season information entry. The other entries are explained

by the adjacent comments. Additional information about these parameters is

given in chapter 3.

2.4 Simulation Control Routines

The simulation control section of the program is found in a subroutine given

the name ECTROL. ECTROL has two output arguments* TIMEST and GO. TIMEST is an

integer variable which returns the current time step for the simulation to

use. GO is a logical variable which is returned as true if the model equations

should be solved and false if the equations should not yet be solved.

ECTROL uses two subroutines* TIMINC and NDSTEP, to determine values for

variables TIMEST and GO. TIMINC has two logical output arguments, STEPDN and

TIMEUP. If a time interval equal to the simulation time step has elapsed since

the last solution of the model equations* STEPDN is set true, otherwise it is

set false. If the duration of the emulator test specified in the emulator

configuration file has elapsed* TIMEUP is set true* otherwise it is set false.

If TIMEUP is true* ECTROL calls REPORT* which generates the emulator model

report files. EVENT is a global logical variable which is true if the command

inputs in the state vector have changed since the last model equations

solution, or false otherwise. If either EVENT or STEPDN is true, GO is set to

TRUE and ECTROL calls NUSTEP. NUSTEP returns in the variable TIMEST the time

elapsed since the last solution of the model equations.

2.5 Building Component Model Routines

The solution of the model equations is effected by successive calls to WEATHR,

the weather simulation routine* AHU, the air handling unit simulation routine,

ZONE, the zone simulation routine* COMFRT, the comfort model routine, and

CMPILE, the energy use compilation routine. All of the models, with the

15
exception of ZONE* are steady state models not involving differential

equations.

Table 2-10. Sample weather data parameter file

f Cold humid winter* hot humid summer* moderate sun'

WINTER’

35.0 'Maximum temperature* F'

15.0 'Minimum temperature* F'

14.0 'Time of maximum temperature* hr s' 0.003 'Humidity ratio

500. 'Peak solar heat gain on horizontal

10.0 'Average windspeed* MPH'

f LATE WINTER*

40.0 'Maximum temperature* F'

20.0 'Minimum temperature* F'

14.0 'Time of maximum temperature* hr s' 0.003 'Humidity ratio

1000. 'Peak solar heat gain on horizontal

10.0 'Average windspeed, MPH'

* SPRING*

60.0 'Maximum temperature* F' 35.0 * Minimum temperature* F'

13.0 'Time of maximum temperature* hrs'

0.004 'Humidity ratio

1700. 'Peak solar heat gain on horizontal

11.0 'Average windspeed, MPH'

* EARLY SUMMER*

80.0 'Maximum temperature, F*

60.0 'Minimum temperature, F'

14.0 'Time of maximum temperature, hrs'

0.010 'Humidity ratio

2000. 'Peak solar heat gain on horizontal

9.0 'Average windspeed* MPH'

’SUMMER*

85.0 'Maximum temperature, F' 60.0 'Minimum temperature, F' 14.0 'Time of maximum temperature, hrs'

0.013 'Humidity ratio

1800. 'Peak solar heat gain on horizontal

8.0 'Average windspeed, MPH'

WEATHR determines the outside air temperature using a sinusoidal

approximation. All other values such as humidity ratio and windspeed are

16
assumed constant. Solar gain is calculated geometrically. Section 3. describes

the equations solved in the veather model.

AHU is used to determine the energy required to heat or cool a mix of return

and outdoor air to the supply air setpoint. The model has heating coils,

cooling coils, and an outside air damper. These devices are assumed to be

under perfect sequenced control. There is no local economizer control of the

dampers. Section 4 describes the equations for the air handling unit model.

The ZONE subroutine performs a solution of 5 simultaneous equations of which

three are differential equations and two are steady-state. These equations

describe the dynamics of heat and moisture transfer between perfectly mixed

zone air and a two-node wall and the outside environment. A model of local

zone heating or cooling equipment under PI control is included in subroutine

ZONE. Details of the solution technique and the equations used in subroutine

ZONE are discussed in section 5.

The routine COMFRT compares the current zone space conditions to the ASHRAE

comfort zone and determines if the zone is comfortable. The comfort model

equations are described in section 6.

The subroutine CMPILE, using values from the state vector, accumulates the

energy used by the building in various categories. The equations used to

calculate energy are described in section 7.

17
3. WEATHER MODEL 3.1

Assumptions and Purpose of Weather Model

The purpose of the weather model is to provide values for the dry bulb

temperature* humidity ratio* vindspeed* and solar radiation on a vertical

surface for each timestep in the simulation. It is assumed that tests of an

EMCS system will require the generation of one to five days of continuous

weather data. For such a time period* all of the days in the period will be

assumed to have identical weather. The humidity ratio and windspeed will be

assumed to remain constant throughout the test period. The dry bulb

temperature is assumed to vary smoothly in a sinusoidal manner with a period

of one day. Solar radiation is calculated from geometric relationships and

statistical values for direct solar radiation.

Model parameters for the weather model were listed in table 2-10. These

parameters allow the weather model to be modified to represent different

climates and seasons within climates to which a building might be exposed

within the constraints listed above. The parameters for humidity ratio and

windspeed specify the constant value that these variables will have throughout

the day. Other specific parameters will be discussed in the appropriate

sections below.

3.2

Outside Air Temperature

The model parameters which determine the variation of outside dry bulb

temperature with time for the weather model are the maximum daily temperature,

the minimum daily temperature, and the time of occurrence of the maximum

temperature. For the generation of a sinusoidal temperature profile the

amplitude, period, offset, and phase of the sinusoid must be specified. For

temperature variation the period is assumed to be 24 hours. The amplitude is

assumed to be half of the difference between the maximum and minimum daily

temperatures and the offset is assumed to be the average of maximum and

minimum temperatures. The phase can be determined from the time of occurrence

of the maximum temperature.

Rather than a mathematical sine function, the emulator model uses a table of

24 discrete points on a curve. The values in the table are normalized to lie

between -1 and +1 with a zero offset. Values between the specified points are obtained by interpolation. The interpolated value is then scaled with the

amplitude and offset to obtain the temperature.

3.3

Solar Energy

The effects of solar energy on the building are not included in the abridged

model that was produced at NBS. The prediction of solar irriadiation on the

18
building could be added to the model using the equations in the following

sections

3.3.1 Solar Radiation Intensity on the Building Wall

The total solar irriadiation on the outer surface of the building vail

(including windows) is given by:

Direct and diffusive radiation are evaluated from the equations given by

ASHRAE (1981), p. 27.2, and Walton (1983):

where CL

is a clearness number betveen 0.75 and 1.0, 8 is the solar altitude

angle, and A, B, and C are constants. Monthly values are shown in Table 1 of

ASHRAE (1981), p. 27.2. Walton (1983) provides a detailed procedure for

obtaining A, B, and C for a given Julian date if desired.

3.3.2 Sine of Solar Altitude

Sin 6 = sin(d) sin(L) + cos(d) cos(L) cos(h) (3-4)

where d : declination angle (degrees) [ASHRAE (1981), p. 27.2]

L : latitude (degrees)

h : hour angle (degrees)

3

: solar altitude angle (degrees)

The hour angle is given by:

I = Ib cos0 + I

8 + I g

(3-1)

where I : total solar radiation

Ib

: direct radiation

I

g

: diffusive radiation

Ig : ground reflective radiation (assume Ig = 0)

0 : angle of incidence of the sun's rays

Xb

= &L A eXp[-B/sin(S)] (3-2)

(3-3)

h = 15 (t

g - 12 + E/60 + TZN) - X (3-5)

where t

g

: standard time of day (h)

E : equation of time (min) [ASHRAE (1981), p. 27.2]

X : longitude (degree)

TZN : time zone number [Kusuda (1976), p. 9a]

= 4 for Atlantic

5 for Eastern

6 for Central

7 for Mountain

8 for Pacific

19
3.3.3 Cosine of Incident Angle

[ASHRAE (1981). p. 27.3]

Cos 8 = cos(B) cos(f) + Y) sin(v) + sin(s) cos(v) (3-6)

where 6 : 0 s Y :

v :

co 8 ( 3)

angle of incidence of the sun's rays

solar azimuth angle from south

angle between the normal to the vertical surface and the south

(facing angle)

tilt angle (the angle between the normal to the surface and the

normal to the ground surface)

cosine of solar altitude angle

= 1 - sin2

3 (3-7)

3.3.4 Sol-air Temperature

The abridged version of the emulator model program uses dry bulb outside air

temperature for heat transfer between the wall and the outdoors. A more

complete version would use sol-air temperature as given by the following

equation for an opaque wall [Threlkeld (1970). p.342]:

T 1 = T + sol o

a I

ow

(3-8)

where a solar absorptance of the exposed wall

solar radiation intensity (direct + diffuse)

outside air dry-bulb temperature

thermal conductance of the wall outside of the mass center

20
4. AIR HANDLING UNIT MODEL

4.1 Assumptions

The air handling unit model uses steady state air enthalpy calculations to

determine the amount of energy required for cooling and heating by the unit.

Steady state equations are used since the time step used for the air handling

unit model is similar in size to the time constants of HVAC systems. Possible

EMCS control actions that provide input to the air handling unit model are

starting and stopping the air handling unit* closing or alloving modulation of

the outside air dampers* opening or closing minimum outside air dampers* and

supply air temperature setpoint. The general model is designed to be

configurable to represent either a constant or variable air volume (VAV) unit

with either single or dual decks and heat exchangers for heating or cooling.

The abridged model (section 2.0) that was implemented in the building emulator

model computer program does not allow simulation of dual deck systems or

systems with humidity control. See subroutine AHU in APPENDIX.

4.2 Equations

A sample of the model parameters for the air handling unit simulation

subroutine was listed in table 2-9. The parameters in the table fall into two

categories* which are configuration parameters and analog parameters. The

configuration parameters are logical constants which are used to select a particular configuration for the air handling unit. There are four of these

parameters* which are DUal-DECK (DUDECK), HUMidity-CONtrol (HUMCON), ReturnAir-FAN (RAFAN), and Variable-Air-Volume (VAV), and which may be either true or false.

Figure 2-1 is a diagram of a generalized air handling unit showing the major

components of the model. The two basic configurations for the air handling

unit are dual-deck or single deck. To select dual deck* the logical parameter

DUDECK is set true. If dual deck* there will be a set of heat exchangers for

two paths and the supply fan will be located before the split. (Dual deck is

not implemented in the abridged model) If the system is single deck. DUDECK

will be false, the second set of heat exchangers will not be used, and the

supply fan will be located after the last exchanger.

The decks in the model may contain a humidifier* a cooling coil, a heating

coil or any combination of the three. The return air fan is optional and the

logical parameter RAFAN determines whether it will be included in the model.

If the parameter VAV is true* the air handling unit will be modeled as being

Variable Air Volume (VAV), and the mass flow rate of air through the unit will

be allowed to change to match the zone loading. If humidity control is to be

used* parameter HUMCON must be set to true* and the supply air will be

maintained at a setpoint for humidity as well as a dry bulb temperature

setpoint (humidity control is not implemented in the abridged model).

21
Figure 4-1. Block diagram of components of generalized HVAC system

model

22
outside air (oa)

Figure 4-2. Block diagram of components of the abridged, single deck

air handling unit model

23
Table 2-9 also shows a number of analog parameters which have to be specified

for the model. These describe the capacities of the coils and fans in the

unit. The fans in the unit will have a power requirement* a rated air volume

flow rate* and a parameter to specify how much energy was added to the air

passing through the fan by the fan waste heat. The fans for VAV will use the

rated flow rate as the maximum but will have a flow rate and power consumption

which can be lowered to match zone loads. Cooling and heating coils will have

a certain capacity in energy per unit time that can be removed or added to the

air in the unit. If a coil is given a capacity of zero* this indicates that

the coil is not present.

In addition to the parameters there are inputs in the FORTRAN common area

which contains the state vector. These are variables which represent commands

from the EMCS which is being tested by the emulator. The variable ECONOmizer

(ECONO) will be a signal to use or not use the outside air dampers for

cooling. UNIT-ON (UNITON) will be true if the unit is to be running and false

if the unit is turned off. VENtilation-OFF (VENOFF) will be false during

occupied operation* but if the EMCS closes minimum outside air dampers VENOFF

will be true.

The locations shown in figure 4-1 with two or three letter abbreviations

beside them and location descriptions in parenthesis represent states of the

air passing through the air handling unit. Each state can be characterized by

the dry bulb temperature* T* the humidity ratio* W, the flow rate m and the

enthalpy* h. The abbreviations on the diagram represent subscripts to variable

names in the equations which are used to describe the air handling unit model.

For example* the temperature of the air returning from the building zones is

symbolized by T_ a .

It will be assumed that if the capacity of the air handling unit will allow,

the condition of supply air from the the air handler will be the same as specified by the air handler supply air setpoints. For a constant volume air

handling unit* the supply air mass flow rate will be constant. If the system

is VAV, the supply air flow rate will be determined by the local equipment

model contained in the zone model (section 5.6). The return air flow rate will

be assumed to equal the VAV supply air mass flow rate. If the system is dualdeck, the individual air mass flow rates for the hot and cold deck will be

determined by the local equiment model. The return air mass flow rate will

then equal the sum of the two deck mass flow rates.

The following FORTRAN functions are assumed for psychrometric calculations

(the implementations of these functions are described in section 5.5). For the

arguments listed, T is dry bulb temperature* W is humidity ratio, P is

atmospheric pressure, RH is relative humidity, h is enthalpy, and DP is

dewpoint temperature (also see Table 5-3).
1. h = ENTHAL(T,W,P) - enthalpy of moist air

2. T = TEMP(h,W,P) - dry bulb temperature of moist air

3 . DP = DEWPT(W.P) - dew point temperature of moist air

4. W = HUMRAT(T,RH,P) - humidity ratio of moist air

5. RH = RELHUM(T,W,P) - relative humidity of moist air

6. T = TDBSAT(h,P) - dry bulb temperature of saturated air

Figure 4-2 is a flow diagram shoving the logic which determines the equations

to use for calculating the states of the air in the air handling unit. The

calculations for each time step are performed based on the return air state

which is known from the zone air model.

Since the equations used for the air handling unit model depend on the state

of several logical variables* a special format will be used to present the

equations. The construction:

IF (conditionl) is true (or false) THEN

equation 1 ELSE IF (condition2) is true (or false) THEN

equation 2

ELSE

equation 3

ENDIF

will be used to denote logical decisions about which equation to use. Indentations will be used for clarity. See Figure 4-2 for logic flow and

Table 4-1 for nomenclature of the equations which follow.

The following equations will apply if logical variable UNIT-ON is true. If it

is not. then the air handling unit is off, no air is flowing, all mass flow

rates will be zero, and all energy use by the air handler will be zero.

IF Variable-Air-Volume is true THEN

The supply air volume is assumed to be determined by the local zone VAV

equipment.

ELSE

Constant volume is assumed. The mass flow rate m

ga equals the rated air

flow rate.

ENDIF (if Variable-Air-Volume is true).

IF Return-Air-FAN is true THEN

Assume that the return air fan causes a constant amount of thermal energy

to be added to the return air per second. The return air fan heating rate

must be provided as a parameter, qrf. Then the temperature rise due to

the return air fan is:

25
Figure 4-3. Air handling unit model calculation sequence

26
( 4

- dT 1 ) rf

and T

*rf

“ra c

ra

= T

za

pa

+ dTrf ( 4

- 2 )

ELSE

T = T

ra za

( 4

- 3 )

ENDIF (if Re turn-Air-FAN is true)

W

ra

= W

za h

rfl = ENTHAL(Tra

,Wra,P)

( 4

- 4 )

( 4

- 5 )

The next calculations determine the state of the mixed air (subscript me) if

the outside air dampers are closed. It will be assumed that fya is a model

parameter which the fraction of air mass flow through the air handling unit

which is the minimum outside air flow rate for ventilation purposes. If the

EMCS cuts off ventilation air as an algorithm action* then this value will be

zero:

IF VENtilation-OFF is true THEN

^mc - h

ra W „ = w me ra

T = T me ra

ELSE

*Hnc “ ^va ^oa + ^ ^va^ ^ra

^mc “ ^va ^oa + (1 “ ^Va^ ^ra

*mc = TEMPOi^.W^.P)

ENDIF (if VENtilation-OFF is true)

IF DUa

1 -DECK is false THEN

( 4

- 6 )

( 4

- 7 )

( 4

- 8 )

( 4

- 9 )

( 4

- 10 )

( 4

- 11 )

For a single deck air handling unit* the supply fan temperature rise is

calculated as the return fan rise was calculated:

27
Table 4-1. Nomenclature for the air handling unit model

cpa specific heat of air

DP dew point

dT temperature rise

f fraction of air mass flow for the minimum outside

air for ventilation

h enthalpy

m air mass flow rate

P barometric pressure

q heating/ cooling rate

^c.max cooling capacity

^h!max heating capacity

RH relative humidity

T air temperature

Tset setpoint temperature

W humidity ratio

Subscripts:

c cooling

ca cooled air

h heating

ha heated air

ma mixed air me state of mixed air if the outside air dampers are closed

mo state of mixed air if the outside air dampers are completely open

oa outside air

ra return air

rf return air fan

sa supply air

sf supply air fan va humidified air

za zone air

dTsf

^sf

m8al c pa

(4-12)

IF HUMidity-CONtrol is false THEN

If there is no humidity control, there will only be a dry bulb

temperature setpoint. The humidity will be allowed to float to any

level.

IF ECONOmizer is true THEN

28
The EMCS allows the air handling unit to use the outside air

for conditioning supply air. No provision is made in this

version of the model for direct EMCS control of damper

position. It is assumed that local controls determine proper

outside air damper position. However* once the outside air

damper is fully open when cooling is required* it remains open

even if outside air enthalpy is high. The EMCS must take

responsibility for closing the damper in this case. To

determine the outside air damper position:

^mo ” ^oa • (4-13)

T__ = T„. and (4-14) mo oa

= W . mo oa

(4-15)

where the mo subscript indicates the state of the mixed air if

the outside air dampers are completely open.

If T

tnc < T

t - dTs £

* this indicates that the outside air

damper should be closed and that heating of the supply air will

be required. In this case:

T = T ma me

U a = Wmo and ma me

^ma “

kjac

(4-16)

(4-17)

(4-18)

However, if T__ > T..,. - dT_ f and Tmn > T„ . - dT 0<: , this

. me set . sr . mo set gr indicates that that the outside air damper should be wide open,

but additional cooling of supply air will be needed using a cooling coil. In this case:

^ma “ ^mo • ( 4- 19

• W„o aad (4-20)

^ma “ ^mo * (4-21)

But ’ if Tmc

” T set : dT sf and Tmo * Tset' dT sfr then tha

outside air dampers will be able to provide all conditioning of

the air passing through the air handler. In this case, the

outside air dampers will be at some intermediate position. The

supply air temperature will equal the setpoint but since the

fraction of outside air is not known in this case, the exact

psychrometric state of the supply air (and thus the enthalpy)

is not known. Theoretically, since the dry-bulb temperature

lines on the psychrometric chart are not parallel, an iterative

29
technique would be required to determine W

ma and hma. However,

the lines are approximately parallel and an approximate

relationship can be used since the state of the mixed air lies

on a straight line between the state of the air with the

dampers fully open and the state with the dampers closed. This

is described by

^ma o

1 Wma

- Wmnmo

or ^mc 0

1 Wmc

- Wmo

Then, since

T = ma ^set dTsf •

U —

- me mo

(Tv Set dTsf ma

- T

- W me

Tmo

W.mo

( 4

- 22 )

( 4

- 23 )

( 4

- 24 )

ELSE

= ESTHALtT^.W^.P)

If the EMCS keeps the outside air dampers closed:

( 4

- 25 )

T = T

ma me

Wma = ma meand

*\na ~ ^mc *

ENDIF (ECONOmizer is true)

( 4

- 26 )

( 4

- 27 )

( 4

- 28 )

If outside air is providing all of the conditioning, as determined

above (Tma = T

get - dT

g j

), then no heating or cooling of the mixed

air is required, and:

T — t

cal ma Wcal = Wma •

^cal “ ^ma *

^hal “ ^ma *

W

hal = Wma • and

hhal = hma *

( 4

- 29 )

( 4

- 30 )

( 4

- 31 )

( 4

-32 )

( 4

- 33 )

( 4

- 34 )

30
However, if heating of the mixed air is required (Tma < Tget - dTgf ) then the air must be heated from Tma to T^^, where:

L hal

= Tset

- dTsf (4-35)

The simplified psychrometric process for air heating is to increase

the temperature of the air at constant humidity ratio (no moisture

is added). Therefore:

Hhal

' W

ma

and

hhai = EHTHAl(Thal .W

ma

.P).

Then the energy added to the air is:

^h “ msal ^ ^hal

- ^ma ^ *

(4-36)

(4-37)

(4-38)

If q h

is greater than the heating capacity of the air handling unit,

then the supply air temperature will be below the supply air

setpoint. The supply air temperature is calculated as follows:

%•max hal

msal

hal = «h

hal •

+ Ka * and

Then the energy added to the air is:

(4-39)

(4-40)

^h “ %,max (4-41)

Alternatively, cooling of the air may be required (Tma > T

ge(. - dT

g f). Then the air must be cooled from Tma to T

ca ^, where:

cal = T

set

- dTsf (4-42)

The assumed psychrometric process for air cooling is to remove

energy from the air at constant humidity ratio until the dew point

temperature of the mixed air is reached, at which point energy and

water are removed from the air, which remains saturated, until T

ca ^ is reached.

The dewpoint temperature must be known for the calculations:

DP

ma

* DEWPT(Wma . P) (4-43)

If T ca ^ > DPma » then the supply air will not be at a saturated

condition after it is cooled and:

31
h (4-44) cal

s h ^ Tcal* W

ma* * and

Wcal, - Wmfl ma

(4-45)

However* if T

ca^ < DP then condensation of moisture from the air

will occur and the supply air will be assumed to be saturated:

RHcal = 100 Z . Wcal = HDMRAT(Tcal . SH^j. P) . and

h

cal = EHTHAl(Tcal . Wcal . P).

( 4-46

(4-47)

(4-48)

Then the cooling energy required to maintain the setpoint is:

’c * msal ( hma

‘ hcal > (4' 49)

If q c is larger than the rated capacity of the cooling coil, then

the supply air state will be above the setpoint of the air handling

unit. The supply air state and cooling energy will then have to be

calculated as follows:

. , H c,max

h

cal * >Va • and

msal

= qc,max

(4-50)

(4-51)

If > ENTHAL(DPma *

ma »P) , then the supply air at maximum

cooling capacity is not saturated:

^cal = TEMP(hcal . Wma, P>, and

W„„, = W__. cal ma

(4-52)

(4-53)

However, if h

ca j

<_ h(DPma ,W ma ,P), then the supply air at

maximum cooling capacity is saturated:

T

cal

= TDBSAT(hcal , P). and

Wcal = HUMRAT(Tcal , RH=100%, P)

The final supply air state is then:

Wsal = Wcal

T

sal ” T

cal + dT sf

(4-54)

(4-55)

(4-56)

(4-57)

32
h

8ai = ENTHAL(T8al , Wsal . P) ( 4

- 58 )

ELSE

If it is in the heating season and the humidity is too lov, a humidifier can be used to increase the humidity ratio of the supply

air. If it is in the cooling season and humidity is excessive, the

air can be dehumidified by overcooling and reheating of air.

Humidity control was not implemented in the abridged version of the

emulator model and equations will not be presented.

ENDIF (if HUMidity-CONtrol is false)

ELSE

Equations for dual deck were not implemented in the abridged version

of the emulator model and vill not be presented.

ENDIF (if DUal-DECK is false)

The states of the air passing through the air handling unit are stored in

the state vector and are made available to other model subroutines in the

building emulator model program.

33
5. BUILDING SPACE ZONE MODEL

The building zone model subroutine* unlike the other model equation

subroutines* implements a model formulated in terms of differential equations.

A dynamic rather than a steady-state model is necessary since the time step

magnitude used in the model (smaller than one minute) is much smaller than the

characteristic time constants of the zone.

5.1 Assumptions

A number of simplifying assumptions are made in order to allov a relatively

simple set of model equations. The building zone to be modeled is assumed to

be a typical perimeter zone of a large multi-story building as Kao (1983) used

in his energy conservation study for a large office building. The zone is

enclosed by an exposed vail with a glass window* three interior walls*

ceiling* and floor. Since there is one exposed wall and one glass window the

governing equations are simplified. The following assumptions are made for

the building zone model:

(1) Uniform temperature distributions on surfaces.

(2) No stratification of zone air (fully mixed air).

(3) No heat or mass flow between zones.

(4) Simple, uniform optical and thermal properties for any glass

windows

(5) The term ’wall 1 includes walls, ceiling, and floor.

(6) All wall structures (including ceiling and floor) which are not

exposed directly to the outside are treated as thermally adiabatic.

(No conductive heat flow across the wall structure is assumed).

(7) Simple two-node representation of the exterior wall with one node at

the interior wall surface and one node at the mass center point of

the wall.

(8) Convective and radiative heat transfer coefficients for the inner

surfaces of walls are considered constant.

(9) The surface heat transfer coefficient of outside wall surfaces is

considered constant and is not a function of outdoor conditions

(including wind).

(10) The wind speed is considered constant and has an influence only on infiltration.

34
(11) All static pressures in the HVAC system, zone and outdoor air are

assumed to be constant with a value of 29.92 inches of mercury.

(12) The latent heat gain from furnishings in the zone is

neglected.

(13) The heat capacity of zone furnishings is assumed to remain the same

at all times.

(14)

Radiant heat from lights (short wave radiation) is assumed to be

uniformly distributed on all inner surfaces of the zone.

(15) For all glass windows, the same constant values of transmittance

and shading coefficient are used.

(16) For all exposed opaque walls, solar absorptance is a constant.

(17) Solar radiation intensity consists of direct and diffuse radiation.

Ground reflective radiation is ignored. The solar constants A, B,

and C appearing in equations (3-2) and (3-3) are assumed constant

for a given specific month.

(18) It is assumed that solar radiation entering the zone through windows

is not absorbed by glass but is distributed uniformly on the inner

surfaces of the zone.

Based on these assumptions, it is assumed that five state variables are

sufficient to model the zone. The state variables are room air dry bulb

temperature, room humidity ratio, exterior wall interior surface temperature,

glass interior surface temperature, and exterior wall internal temperature.

The abriged model (see section 2) implemented in the building emulation

simulation program included some additional assumptions. The equations

presented in the following sections will usually be for a complete model. The

changes to equations for the abridged model will be described whenever

necessary in the following sections. The assumptions for the abridged model

are:

1) No effects of solar energy.

2) No radiative heat transfer between zone surfaces.

5.2 Fundamental Equations

Since there are five state variables in the building emulator, five equations

relating to these variables are required to solve for the values of the

variables. Table 5-1 contains nomenclature for the fundamental equations of

the zone model.
5.2.1 Zone Air Heat Balance

The first fundamental equation for the emulator zone model expresses the

balance in heat transfer between the air and furnishings in the zone and the

vail and window surfaces* entering supply air from the air handling unit*

infiltration heat transfer, and internal sources of heat [Borresen (1981),

Modest and Hacken (1977)* Bullock (1984), Park (1981)]. The equation is:

zone

dTj

dt

^supply + ^infilt + ^convec + ^internal (5-1)

The left hand side of the equation represents the storage of energy in the

zone air and furnishings. The state variable is the interior air temperature.

Of the four terms on the right side of the equation* the energy from the

supply air and the internal gains may. be considered boundary conditions. The

supply air condition is determined by the air handling unit model and the

internal gain* consisting of heat gains from occupants* lights, and equipment,

is determined by the building schedule. The other two terms for convection and

infiltration heat transfer are functions of the other state variables.

Table 5-1. Nomenclature for equations in the building space zone model

State Variables:

T^ - space dry bulb temperature

Tw

- wall mass center temperature T^

T- - wall interior surface temperature

space humidity ratio

window glass temperature

Other Parameters:

A

C

M

air

^mr

G

S gv

conduction coefficient B

- thermal mass T

- air mass for water absorption Q

- space mean radiant temperature R

- zone moisture increase rate

wall heat gain by solar radiation through

convection coefficient

temperature

heat transfer rate

radiative coefficient

windows

Subscripts:

convec - convection o - outside air g - glass

infilt - infiltration w - at the wall mass center

5.2.2 Exterior Wall Mass Heat Balance

The exterior wall is assumed to be the only surface enclosing the zone which

exchanges energy with the outdoor environment and the indoor environment. The

state variable for wall temperature is assumed to be the temperature at the

36
mass center of the wall. The mass center is defined as that point in the wall

where for a unit area of wall the total mass of the wall outside of the mass

center is equal to the total mass of the wall inside the mass center. For a masonry construction* this may result in the mass center being different from

the geometric center of the wall. The equation for the exterior wall is:

dT

’"wall “ ^inside + ^outside + ^solar ^ ^

dt

The left hand side of the equation represents the storage of energy in the

building wall. The first term on the right is the heat transfer between the

wall interior and the inner wall surface. This term is dependent on state

variables. The second two terms on the right represent the heat transfer

between the wall interior and the outdoor environment, including the

impingement of solar irradiation on the outside of the wall. The right two

terms are time-dependent boundary conditions.

5.2.3 Zone Air Moisture Balance

The air in the zone is assumed to contain moisture. The state variable used to

describe the water content of the air is the humidity ratio. The equation used

to describe the moisture balance is [Chi (1982)* Bullock (1984)]:

M

air

dWi

dt

^infilt + ^supply + ^internal (5-3)

The left hand side of this equation represents the storage of water vapor in

the air and furnishings of the zone. The right hand side of the equation has

three terms. The first term is the moisture transfer with the outdoor

environment due to infiltration, which contains a state variable. The second

term is the net moisture added to or removed from the room by the air handling

unit. This is the difference in absolute moisture between the supply and

return air from the air handling unit, and is a boundary condition. The third

term is the moisture generated internally by occupants and equipment. The

deficiency in this equation is that it assumes the air and furnishings can

absorb an infinite amount of moisture when in reality, the air can become

saturated. The emulator model zone subroutine performs a check of the relative

humidity of the air to prevent oversaturation.

37
5.2.4 Inner Wall Surface Heat Balance

For the inner wall surface of the exterior wall, the heat balance is expressed

by an algebraic equation rather than a differential equation since there is no storage of energy at the wall surface. The equation is:

V T - T* ) + B (T- - T- ) + R. (T - T- ) + S_ = 0 (5-4)w iw' w

v l iw' mr iw y gw

This equation is based on the assumption that all heat transfer in and out of

the node at the wall surface sums to zero. The first term is the conduction

heat transfer between the wall interior and the inner wall surface. ^ is a conduction coefficient. The second term is the convection heat transfer

between the wall surface and the zone mean radiant temperature. The final

term is due to any solar irradiated energy which passes through the zone

window and is absorbed by the zone inside surfaces.

5.2.5 Inner Glass Surface Heat Balance

For the inner glass surface of the zone window, if present, an algebraic

equation describes the heat transfer balance since there is no storage of

energy:

A ( T

g o

- T- ) + B (T- - T- ) + R (T lg' g

v l lg' gv mr

- T

i g

> = 0 (5-5)

This equation expresses that the sum of all heat transfer in and out of the

node at the glass surface is zero. The first term is the conduction/convection

heat transfer between the inner glass surface and the outdoor air. The

coefficient Ag is a combination convection and conduction coefficient. The

second term in the equation is the convection heat transfer between the glass

surface and the interior air. The coefficent Bg is for convection. Radiant

heat transfer between the glass surface and the mean radiant temperature of

the room is represented by the third term of the equation.

5.3 Solution of Equations

The equations described in section 5.2 are solved by the zone model subroutine

using a simple explicit forward difference numerical technique (Euler method)

to predict the values of the state variables one step forward in time. The

time step used within the zone model is different from the emulator time step

used for the other models in the emulator simulation program. Since simple

differential equation solution technique is used within the zone model, the

time steps must be small or the solution will not converge. The emulator time

step is divided into a number of smaller zone model time steps and the

differential equations are solved for each zone model time step until the

38
emulator time step has passed. However* the zone model is not constrained to

execute zone model time steps in real time.

The subroutine ZONE performs the steps required to determine the values of the

state variables at the next emulator time step. This routine first determines

the best value for the zone model time step. This process is described in

section 5.3.2. Then ZONE calls the two subroutines EULER and SURFAC in

sequence* and advances the time by one zone model time step until the emulator

time step has passed. The subroutine EULER is used to solve the three

differential equations listed above (5—1 through 5-3) for one zone model time

step. SURFAC determines the wall and glass surface temperatures for the

current time step using the algebraic equations described above (5-4 through

5-5).

5.3.1 Euler Forward Difference Technique

The subroutine EULER is a generalized routine to advance a system of ordinary

first order differential equations by one time step using the explicit

difference method. The arguments to the routine include the number of

differential equations* an array containing the dependent (state) variables*

the time step* and the name of a FORTRAN function which will evaluate the

right hand sides of the differential equations when they are expressed in the

form:

dY;

— f (t» Yi i Y« *. .. .Y_ ) * Is 1*2. . . *n

dt

(5-6)

where Y- is the ith dependent variable and n is the number of equations. The

value of Y at time t+ t can approximately be predicted by:

Y^( t+ At ) = Y^t) + At * FUNCTION (i,t,n,Y1 *Y2

* . . .YQ

) (5-7)

As long as the product of the time step At and the function value is small, a solution to the differential equations can be obtained.

5.3.2 Adaptive Time Step

The time step used by the EULER subroutine is determined adaptively to become

smaller when the zone is undergoing rapid changes in the state variables and

larger as the zone approaches steady state. Subroutine ZONE estimates a new

time step each time it is called. The estimation is based on the requirement

that the change in any state variable for one Euler time step must not exceed

a selected maximum value. For example, in the building emulator model program,

this maximum has been chosen to be 0.1 degrees F for temperature state

variables. Since the value of the function which describes the right hand side

of the differential equation (described in section 5.3.1) is the rate of

39
change of the state variable, the smallest desirable zone model time step may

be determined by dividing the maximum allowable change in the state variable

per time step by the current value of the right hand side function. Of

course, calculations must be performed to ensure that the zone model time step

is not larger than the emulator time step or smaller than a certain minimum.

In the building emulator program zone model subroutine the time step can be no smaller than 0.001 seconds.

After the zone model time step is determined, it is adjusted so that it evenly

divides the emulator time step. The quotient of this division is the number of

times that the ZONE subroutine must call the EULER and SURFAC subroutines for

each emulator model time step.

5.3.3 Differential Equation Function

The FORTRAN function ZONEFN is used by the zone model subroutine to calculate

the current value of the right hand sides of the differential equations

describing the state variables. Since there are three differential equations

for the zone model, there are three possible right hand sides of the

differential equations. The right hand side to be calculated is determined by

an index which is 1 for the zone air energy balance, 2 for the zone moisture

balance, and 3 for the wall energy balance. The expanded right hand sides of

the equations corresponding to equations (5-1 through 5-3) are presented

below. The variables which are designated by an upper case F with a subscript

are implemented in the ZONEFN subroutine in a special way. Rather than being

variables or constants, they are FORTRAN functions which return the value of a coefficient or energy ternw In the definition of terms following each equation

below, the FORTRAN function name used in the emulator model program will be

listed. For the zone air energy balance:

ZONEFN **sa + Fl(T0

- ^i) + ^(T^* T^) + FC g(T^g- T^) + F^g

F

cz

(5-8)

where: F

sa

cw

eg

ig

cz

: energy entering or removed from the zone by the air handling

unit (function SUPPLY).

: coefficient for infiltration of outdoor air into zone (function

INFILT)

: coefficient for convection between the exterior wall surface

and the zone air (function CONVWA)

: coefficient for convection between the surface of a window in

the exterior wall and the zone air (function CONVGL)

: energy entering the zone due to internal gains from lights,

equipment, and people (function INTERN).

: thermal mass of the zone air and furnishings (function MASSTH).

For the zone moisture balance:

40
ZONEFN (5-9) _ ^sv + + mg

'mz

where: F

sv

Im

mg

mz

: moisture added to or removed from the zone by the air handling

unit (function SUPPLW)

: coefficient for infiltration of moisture (or exfiltration) from

outdoor air into zone (function INFILW)

: moisture entering the zone due to internal moisture gain from

equipment and people (function MOIST).

: effective mass of zone air and furnishings which can absorb

moisture (function MASSA)

for the prediction of the rate of change of moisture in the zone, an

additional check is made to insure that the zone air does not become

oversaturated. This is accomplished by first using the function HUMRAT to

determine the humidity ratio at the zone air temperature, T^, and 100 percent

relative humidity. If the current indoor humidity ratio, W^, exceeds this

value and ZONEFN is a positive value, then ZONEFN is set to zero, indicating

that no further increase in zone moisture is possible. If is zero and

ZONEFN is negative, ZONEFN is set to zero indicating that zone moisture cannot

be decreased since the air in the zone is completely dry.

The right hand side function for the exterior wall mass heat balance is:

ZONEFN

WT8 o 1 ~ Tw>

+ ^ci^iw “ Tw)

cw

(5-10)

Where: F

co

F • ci

: coefficient for conduction and convection heat transfer between

the mass center of the wall and the outdoor air (function

CNDCTO).

: coefficient for conduction heat transfer between the mass

center of the wall and the inner surface of the wall (function

CNDCTI).

: the sol-air temperature (defined in section 3.3.4)

: thermal mass of the wall structure (function MASSW)

5.3.4 Surface Heat Balance Subroutine

The subroutine SURFAC is called by ZONE to solve the algebraic equations (5-4)

and (5-5) for the inner surface temperatures of the exterior wall and glass

window. The two equations after rearrangement to solve for the surface

temperatures are shown below. The variables which are designated by an upper

case F with a subscript are implemented in the SURFAC subroutine in a special

41
way. Rather than being variables or constants* they are FORTRAN functions

which return the value of a coefficient or energy term. In the definition of

terms following each equation below* the FORTRAN function name used in the

emulator model program will be listed. For the exterior wall inner surface

balance

T-1W

Where : FAw

Bw

Rw

Tr

? ts

fAw tw

* fBw Ti * fRw FTr + Fts (5-11)

FAw + FBw

+ FRw

coefficient for conduction heat transfer between inner surface

of wall and mass center of wall (function CNDCTI).

coefficient for convection heat transfer between inner surface

of wall and zone air (function CONVWA).

coefficient for radiation heat transfer between inner surface

of wall and zone radiant temperature (function WRADIA)

mean radiant temperature of zone (function TRN).

energy which passes through windows in the zone and is absorbed

at the wall surface (function WSOLAR)

For the window glass inner surface heat balance:

Where:

T-ig

Ag

Bg

?Rg

Tr

FAg To + FBg Ti + FRg FTr

FAg + FBg

+ FRg

(5-12)

coefficient for conduction and convection heat transfer between

the inner surface of the glass and the outdoor air (function

GCNDCT )

coefficient for convection heat transfer between the inner

surface of the glass and the zone air (function GCNVCT).

coefficient for radiation heat transfer between the inner

surface of the glass and the zone radiant temperature (function

GRADIA)

mean radiant temperature of zone (function TRN).

5.3.5 Fortran Coefficient Functions

The names of the special FORTRAN functions used to give coefficients or terms

in the equations (5-8) through (5-12) were listed in the definitions following

the equations. All of these functions have the same general form. The

functions can be designated as either time-dependent or t ime- independent.

Values of coefficients such as conduction coefficients may be considered timeindependent and will need to be calculated only once. In this case, the

function automatically determines if the function has already been evaluated

and if so*, returns the previously determined value.

42
Other coefficients such as a coefficient for infiltration heat transfer must

usually be considered time-dependent. In this case the coefficient vill be

calculated each time the function is called unless the function is called

twice for the same time step. In this case the function returns the previously

determined value.

The general form of the FORTRAN function for a coefficient function is shown

in table 5-1. The parameter TIMEDP is set to true for a time-dependent

function and false for a time-independent function. The specific equations to

evaluate the functions used in the emulator model program are listed in

section 5.4.

5.4 Equations for Heat and Mass Transfer Coefficients and Terms

The following sections contain the equations used to evaluate the terms and

coefficients to expand equations (5-8) through (5-12). The name of the

FORTRAN coefficient function (see section 5.3.5) will be given in parenthesis

following the title of the section. The equations listed would be placed in

the area marked Equations to evaluate function* in the generalized

coefficient function listing in table 5-2.

Table 5-2. General form of the FORTRAN function for a coefficient function

FUNCTION [name] (TIME)

LOGICAL TIMEDP. NOINIT

PARAMETER (TIMEDP = .TRUE.)

INTEGER OLD, TIME

SAVE NOINIT, [name], FUNC

DATA NOINIT/. TRUE./, OLD/ -1/

IF (NOINIT . OR .( TIME .NE . OLD . AND . TIMEDP ) ) THEN

FUNC = [ equations to evaluate function ]

OLD = TIME

NOINIT = .FALSE.

ENDIF

[name] = FUNC

RETURN

END

5.4.1 Sensible Heat Gain/Loss from Supply Air (SUPPLY)

The energy added to or removed from the zone by the supply air is determined

from the difference in enthalpies of return and supply air. In addition, the

energy which is added to or removed from the room by the local zone heating or cooling equipment is included. The equation used is:

43
F (5-13) sa

~ ^sal “sal + kSa2 msa2

” ^za mra+ ^£h

~

Where

^sal

“sal

hsa2

msa2

**za “ra

^5,h ^ JZ.C

: enthalpy of supply air

: mass flow rate of supply air

: enthalpy of hot deck supply air (optional)

: mass flow rate of hot deck supply air (optional)

: enthalpy of zone air

: mass flow rate of return air

: energy added to zone by local zone heating equipment

: energy removed from zone by local zone cooling equipment

This function is considered to be time-dependent.

5.4.2 Coefficient for Sensible Heat Gain/Loss due to Infiltration (INFILT)

Sensible gain/loss due to infiltration is assumed to equal a coefficient

multiplied by the difference in temperature between zone and outdoor air. The

coefficient is determined from:

**1 “ Cp.infl p infl ^i *air

where Cp

£nfi : specific heat of infiltrated air

p infl : density infiltrated air

Ia£ r

: air exchange rate (air changes/second)

: volume of zone space

(5-14)

The air exchange rate is considered to be a function of windspeed and indooroutdoor temperature difference [Kusuda and Saitoh (1980)* Achenbach and

Coblenz (1963)]. The equation for air exchange rate per hour is:

air

^•g.air

0.695 [0.15 + 0.013 Vw

+ 0.005 I T

q^ 1 ] (5-15)

where Vw : wind speed (mph)

I

g a £r

: standard air exchange rate (model parameter, see table 2-8)

= 1.5 for leaky building

1.0 for standard building

0.5 for modern type building

This function is considered to be time dependent.

5.4.3 Coefficient for Convective Heat Flow Rate from Walls (CONVWA)

The convective heat flow rate from the inner surface of the zone walls,

including floor and ceiling, is assumed to be equal to a coefficient

44
multiplied by the difference between the wall surface temperature and the zone

air temperature. The coefficient is given by:

F cw

~ b iw,c ^w

(5-16)

where h iu>c : convective heat transfer coefficient at the interior surface

(model parameter* see table 2-8)

: surface area of the wall (model parameter, see table 2-8).

This coefficient is assumed to be time- independent

5.4.4 Convective Heat Flow Rate from Glass Windows (CONVGL)

The convective heat flow rate from the inner glass surface of zone windows is

assumed to be equal to a coefficient multiplied by the difference between the

glass surface temperature and the zone air temperature. The coefficient is

given by:

lg.c A

g

where hig»c

A

g

convective heat transfer coefficient of inner glass

surface to zone air

the glass window area

(5-17)

This coefficient is assumed to be time- independent. The abridged zone model

effectively assumes that there are no windows in the zone, and therefore this

coefficient is not used.

5.4.5 Sensible Internal Gains (INTERN) and Latent Internal Gains (MOIST)

5. 4. 5.1 Internal Gain Schedule Table

The FORTRAN functions used to determine internal gains, INTERN and MOIST,

determine their values from a table which was created during the emulator

model initialization. The table contains internal gain as a function of timeof-day. Interpolation is not used with the table and therefore the internal

gain is not a continous function of time. Internal gains remain at a constant

level until changed by the next entry in the table at a particular time. The

internal gain schedule table is stored in five FORTRAN arrays. The arrays

contain time-of-day, sensible heat gain from equipment at the corresponding

time-of-day, sensible heat gain from people, moisture gain from equipment and

moisture gain from people. A function SCHED is used to determine a sensible or moisture gain from the internal gain tables. If moisture gain is requested (by

MOIST), SCHED returns the sum of the moisture gain from people and equipment.

45
If a sensible gain is requested (by INTERN)* SCHED returns the sum of the

sensible gains from people and equipment.

5. 4. 5. 2 Determining Values for the Internal Gain Schedule Table

The values in the internal gain schedule table are created during the

initialization process of the emulator model. An example input parameter file

was given in table 2-7. The following equations show how a row in the schedule

table is created from the information in a line of table 2-7.

5. 4.5. 2.1 Convective Heat Gain from People

The internal gain from a person is assumed to depend only on the activity in

which the person is engaged and not on the surrounding temperature. Table 5-3

is used as the basis for the sensible and moisture gains from people at

various activity levels. The number of occupants at each activity level at a given time is part of the information entered in the building use parameters

in table 2-7. The sensible heat generated by the occupants of a zone at a particular time is then assumed to be described by the following equation

[Building Environment, NBS (1977), p. A-21]:

Q people,

(1 n p»k W

p,s,k (5-18)

where r : ratio of radiation energy to total sensible heat of people

P (=0.4)

np^ : number of people who are doing the same (k-th) activity

Wp 8> k

: sensible heat generation by a person doing the k-th activity.

ft a

* t

: number of activities

Table 5-3. Occupant sensible and moisture gains

[Table 18 entitled 'Rates of Heat Gain from Occupants of Conditioned Spaces'

ASHRAE fundamentals (1981), p.26.25 and Threlkeld, page 363]

1. seated, very light work

2. seated, light work (typing)

3. standing, light work

4. light bench work

5. walking, 3 mph, light machine work

sensible gain moisture gain

(Btu / hr. ) (lb. water /

230. 0.20

255. 0.27

315. 0.34

345. 0.46

345. 0.74

46
5. 4. 5. 2. 2 Moisture Gain from People

Mpeople, lat

/

k=l

Wp,lat,k (5-19)

where Nact

r Vk

p.lat.k

: number of activity

: number of people who are doing the same (k-th)

: moisture generation by a person doing the k-th

activity

activity

5 .4. 5. 2. 3 Convective Heat Gain from Equipment and Lighting

The total sensible internal gain from equipment and lighting is determined

from the gains due to lighting added to the gain due to equipment. Parameters

for the lighting power, lighting type, and equipment sensible gain at various

times of the day are entered in the building use parameter file (see table 2-

7, for example). For building equipment, the sensible internal gain is

described by the following equation [ASHRAE (1981), pp. 26.26 - 26.29]:

^equip.c “ ^ ” r

e^ We,s

where r

g

: ratio of radiation energy to total equipment energy

(=0.2 to 0.8)

We>8 : total sensible energy produced by equipment

( 5 - 20 )

The sensible internal gain from lighting is assumed to be described by the

following equation:

Q light,c = (1 “ (5-21)

where : total electric input

r^ : ratio of radiation energy to total light energy

= 0.8 for incandescent lights [Building Environment, NBS (1977), p. A25]

= 0.5 for fluorescent lights

5. 4. 5. 2. 4 Moisture Gain from Equipment

Moisture gain from equipment is specified in the parameters for the building

use as shown in table 2-7. The equation for equipment moisture gain is:

[ASHRAE (1981), p. 26.26 - 26.28]

M

equip,lat “ We,lat

where W is the total moisture produced by equipment.

(5-22)
5.4.6 Thermal Mass of Zone Air and Furnishings (MASST&)

The thermal mass which stores energy and dampens fluctuations in zone air

temperature is contained in the air in the zone and the zone furnishings* such as desks* cabinets* and partitions. The shell parameter file* an example of

which is given in table 2-8* contains the volume of the zone and the thermal

capacity of the furnishings. The capacity of the zone air is given by:

C • = C • air p.i p.K V- i 1 (5-23)

where C

P* 1

p i

V i

specific heat of moist zone air

density of zone air

volume of zone air (volume of interior space) (parameter*

table 2-8)

The total thermal mass of the zone is then:

cz where Cfur :

Ce + C • fur air

thermal mass of furnishings (parameter* table 2-8)

(5-24)

5.4.7 Moisture Gain/Loss due to Infiltration (INFILW)

Moisture gain/loss due to infiltration is assumed to equal a coefficient

multiplied by the difference in humidity ratio between zone and outdoor air.

The coefficient is determined from:

FIm

= P inf 1 Vi *air (5-25)

where p infi : density of infiltrated air

I

fl £r

: air exchange rate (see equation 5-15)

: volume of zone (parameter, table 2-8)

The air exchange rate is determined by the same equation as in section 5.4.2.

The density of the infiltration air is determined by taking the average of the

densities of the zone air and the outdoor air. This function is considered to

be time-dependent.

5.4.8 Moisture Gain/Loss from Supply Air

The moisture added to or removed from the zone by the air handling unit is

determined from the difference in the humidity ratios of return and supply

air. The equation used is:

48
F (5-26)

sv

= W

sal m

sal

+ W

sa2 “sa2

” Wra mra where

"sal

“sal

W

sa2

“sa2

W

ra “ra

: humidity ratio of supply air

: mass flow rate of supply air

: humidity ratio of hot deck supply air (optional)

: mass flow rate of hot deck supply air (optional)

: humidity ratio of return air

: mass flow rate of return air

This function is considered to be time-dependent.

5.4.9 Mass of Zone Air and Furnishings that can Absorb Moisture (MASSA)

The mass in the zone which stores moisture and dampens fluctuations in zone humidity ratio • is contained in the air in the zone and the zone furnishings,

such as curtains, paper, and wood furnishings. The shell parameter file, an example of which is given in table 2-8, contains the volume of the zone and a factor which is used to produce an effective air mass for moisture absorption

by the combination of the air and the zone furnishings. This factor is

difficult to measure since the moisture absorption properties of office

materials are not only difficult to determine but depend on air flow and

existing moisture content. This function is considered time-independent. The

moisture absorbing capacity of the zone air and furnishings is described by

the following equation:

Fmz

= p

i Vi Smg (5-27)

where : density of zone air

: volume of zone air (volume of interior space)

R^g : a factor used to produce the effective air mass for moisture

absorption by the zone furnishings and air.

5.4.10 Overall Thermal Conductance between Wall and Outdoor Air (CNDCTO)

Heat transfer between the mass center of the wall and the outdoor air is a combination of conduction heat transfer within the wall and convection heat

transfer from the wall outer surface to the outdoor air. The overall

conductance is expressed by:

F co (5-28)

49
where LQW : thickness of wall between the mass center and the outer surface.

kQW : overall conductivity of the wall between the mass center and the

outer surface.

h

Q

: combined heat transfer coefficient (radiative + convective) for

heat transfer between the outer wall surface and the outside

air.

Aw : surface area of the exterior wall.

The four parameters in this equation are obtained from the shell parameter

file, exemplified by table 2-8. This is considered a time independent

function.

5.4.11 Thermal Conductance between the Wall and Indoor Air (CNDCTI)

Heat transfer between the mass center of the wall and the inner surface of the

wall is by conduction. The thermal conductance is expressed by:

Fci

where

*w

(5-29)

thickness of wall between the mass center and the inner surface,

overall conductivity of the wall between the mass center and the

inner surface.

surface area of the exterior wall.

The three parameters in this equation are obtained from the shell parameter

file, exemplified by table 2-8. This is considered a time independent

function.

5.4.12 Thermal Capacitance of the Wall

The mass in the exterior wall which stores energy and dampens fluctuations in

wall temperature is contained in the building materials of which the wall is

constructed. The wall is divided into an inner and outer section by the mass

center of the wall. The shell parameter file, an example of which is in table

2-8, contains the thickness, density, thermal conductivity and specific heat

of the inner and outer wall sections. This function is considered timeindependent. The thermal capacitance of the exterior wall is described by the

following equation:

cw

= (Cp.ow pow ^ow + Cp.iw P iw ^iw^ ^w (5-30)

50
where C p,ov* p.iv

P ov*

'ov»

p

IV

L

iv

average specific heat of outer and inner sections of

vail.

average density of outer and inner sections of vail,

thickness of outer and inner sections of vail,

area of the exterior vail.

5.4.13 Radiative Heat Transfer Coefficient (WRADIA)

Radiative heat transfer betveen the inner surface of the exterior vail and the

other surfaces in the zone is assumed to be equal to a radiative heat transfer

coefficient multiplied by the difference in the mean radiant temperature of

the zone and the vail surface temperature. For simplicity, this mode of heat

transfer vas not implemented in the abridged building emulator model program.

Determination of this mode of heat transfer involves the determination of viev

angle factors of surfaces vhich enclose the zone. [Walton (1980), Sovell and

Walton (1980)]. Using the viev angle, the radiant heat transfer coefficient

could then be determined [Walton (1980), p. 69]. Hovever, it is necessary to

knov the emissivities and the surface temperatures of all zone surfaces, as veil as the specific geometry. The abridged model effectively assumes that

the radiative coefficient is small or the difference betveen the surface and

mean radiant temperatures is negligible.

5.4.14 Solar Energy Entering the Zone Absorbed by Walls (WSOLAR)

If the zone contains a glass vindov, solar energy may enter the zone through

the vindov and be absorbed by the vails and furnishings. For simplicity, this

mode of heat transfer vas not implemented in the abridged building emulator

model program. Since the only zone vail vhich has inner surface temperature as

a state variable is the exterior vail, and this is the vail vhich is assumed

to have a glass vindov, no solar energy from the vindov vould be absorbed by

the exterior vail. Solar energy could be absorbed by the furnishings, hovever.

The abridged model effectively assumes that the zone vails are vithout

vindovs. The equation for solar energy absorption is included in this section

for completeness. If it is assumed that solar irradiation is not absorbed by

the glass vindov:

A0 r- S„ I g g c (5-31)

Av

: glass vindov area

: solar transmittance through glass [ASHRAE (1981), p. 27.38]

: total solar radiation on the glass vindov

: shading coefficient [ASHRAE (1981), p.27.38]

ts

A vhere g

r„

51
5.4.15

Thermal Conductance Through Glass to Outdoor Air (GCNDCT)

Heat transfer between the inner surface of the glass window and the outdoor

air is a combination of conduction heat transfer within the glass and

convection heat transfer from the glass outer surface to the outdoor air. This

mode of heat transfer is not included in the abridged emulator model program

since the abridged model assumes that there are no windows in the zone. The

overall conductance equation is included for completeness and is expressed by:

FAg

where

A

g

A, (5-32)

thickness of glass,

conductivity of the glass.

combined heat transfer coefficient (radiative + convective) for

heat transfer between the outer glass surface and the outside

air.

surface area of glass window.

5.4.16

Thermal Convection from Inner Glass Surface to Indoor Air (GCNVCT)

Heat transfer between the inner glass surface and the zone air is assumed to

be equal to a convective heat transfer coefficient multiplied by the

temperature difference between the inner glass surface and the zone air. The

heat transfer coefficient must be supplied in the building shell parameter

file. This mode of heat transfer is not included in the abridged version of

the building emulator model program since the abridged model assumes that

there are no windows in the zone.

5.4.17 Radiant Heat Transfer between Glass and Zone (GRADIA)

Radiative heat transfer between the inner surfaces of windows and the other

surfaces in the zone is assumed to be equal to a radiative heat transfer

coefficient multiplied by the difference in the mean radiant temperature of

the zone and the glass surface temperature. Since there are assumed to be no windows in the zone* this mode of heat transfer was not implemented in the

abridged building emulator model program. Determination of this heat transfer

is equivalent to the process described in section 5.4.13.

5.4.18 Mean Radiation Temperature (MRT)

The mean radiant temperature is used to determine radiation heat transfer

between the zone environment and wall or glass surfaces. The mean radiant

temperature is a function of wall and glass surface temperatures, the radiant

heat from people, the radiant heat from equipment, the properties of the

52
surfaces* and the geometry of surfaces [Walton (1983)* p. 69]. Radiation heat

transfer is not included in the abridged version of the emulator model program

and the equation for mean radiant temperature is not included here.

5.5 Property Equations for Hoist Air

The building emulator model program uses a number of functions to determine

properties of moist air. The equations used are not presented in this report,

but the names of the computer routines* the property determined, and a reference for the equations* is given in table 5-4.

5.6 Local Equipment Model

The local heating and cooling equipment is defined as the equipment located in

or adjacent to the zone which is designed to directly heat or cool only that

zone. The purpose of this equipment is to adjust the energy in the incoming

air from the air handling unit to match the heating/cooling requirements of

the zone. There are three different types of local equipment which may be

present. These are: 1. energy consuming local equipment such as reheat coils*

fan coils, heat pumps* perimeter radiation* or electric heaters; 2. air volume

reduction controls to adjust the volume of supply air to match the zone

requirements for VAV systems; and 3. air mixing controls to mix two air

streams of different temperature's to match the zone requirements for multizone

or dual-duct systems. Local equipment is usually under local control* which

means that the controller for the equipment is located adjacent to the zone

and the control sensor is located in the zone.

The subroutine LOCAL in the building emulator model program is used to

simulate local equipment and determine energy usage required to maintain the zone air temperature as close to a local setpoint as possible. This subroutine

is called by the ZONE subroutine for each zone model time step in the solution

of the zone model equations.

53
Table 5-4 Moist air property routines

CPAIR : specific heat of moist air as a function of humidity ratio [ASHRAE

fundamentals 1981* page 5.3]

RHOAIR : density of moist air as a function of humidity ratio and dry bulb

temperature [ASHRAE fundamentals 1981* page 5.3]

ENTHAL : enthalpy of moist air as a function of dry bulb temperature.

humidity ratio, and atmospheric pressure [ASHRAE fundamentals 1981.

page 5.4]

TEMP : dry bulb temperature of moist air as a function of enthalpy.

humidity ratio, and atmospheric pressure [ASHRAE fundamentals 1981*

page 5.4]

DEWPT : dev point temperature of moist air as a function of humidity ratio.

and atmospheric pressure [ASHRAE fundamentals 1981. page 5.4.

equation 40b]

HTJMRAT : Humidity ratio of moist air as a function of dry bulb temperature.

relative humidity, and atmospheric pressure [ASHRAE fundamentals

1981 . page 5.4]

RELHUM : Relative humidity of moist air as a function of dry bulb

temperature, humidity ratio, and atmospheric pressure [ASHRAE

fundamentals 1981. page 5.4]

PWSAT : Saturation pressure of water vapor as a function of absolute

temperature [ ASHRAE fundamentals 1981. page 5.2]

TDBSAT : dry bulb temperature of saturated air as a function of the enthalpy

of saturated air and the atmospheric pressure [ ASHRAE fundamentals

1981, page 6.5, table 2]

The zone model will require a value for the energy being added to or removed

from the zone due to the HVAC system. This energy, designated as q ent » which

enters or leaves the zone, consists of two parts, the energy removed or added

by the air from the air handling unit, q a £ r

» and the energy removed or added

by local energy consuming equipment, 9i 0 cal* The 9iocal ener &y ma y not be

present on certain systems. q a £r can be calculated by:

^air " “sal ^sal + msa2 ^sa2 ~ ^i mra ( 5 - 33 )

where m

8a are supply air mass flow rates (1 indicates the cold deck and 2 the

optional hot deck for a dual deck system), h

ga are the supply air enthalpies.

54
m

ra

is the return air mass flow rate, and h^ is the interior space enthalpy of

the zone.

The actions that the local equipment model simulates are those of the local

controller which will control the local equipment to maintain a zone air

temperature at a setpoint. If the actual zone air temperature deviates from

the setpoint, the local controller will attempt to return the temperature to

the setpoint by increasing or decreasing the local equipment output. A simple

proportional plus integral control (Pi-control) model will be used. For this

model, the adjustment the local controller attempts to make to q ent is

determined from:

d<*ent

= K

p

(e’ e

0

} + K

I (e dt) (5-34)

where e, error = T set p0 int ” T

i* e

0

t*ie Previous value of the error, Kj is an integral gain factor, K

p

is a proportional gain factor, and dt is the zone model time step. Once dqent is determined, the type of HVAC system determines

how q air and q^ oca ^ are to be determined.

If the system is dual deck, it is assumed that there is no local conditioning

equipment beyond the mixing box for the two deck air streams. Then:

d<Uir = dW (5 * 35)

^airCt+dt) ~ ^air(t) + d<*air (5-36)

The total mass flow rate of air from the air handling unit is assumed to be

constant at m

ga

. The supply air enthalpy after mixing must then be:

9 air a i sa + hi

msa

(5-37)

Since the flow rates of the two decks must sum to equal the total mass flow

rate, and the two deck air streams are assumed to mix ideally and have a known supply air temperature, humidity ratio, and enthalpy, the unknown mass

flow rates for the decks can be solved for:

msa2

m

sal

msa ^ ^sa

~ ^sal ^

^sa2 ” ^sal

msa

” msa2

(5-38)

(5-39)

also

If the

W

sa

“ msal ^sal + msa2 ^sa2

system is a single deck VAV, it is assumed that there is no

(5-40)

local

5 5
conditioning equipment beyond the VAV box for modulating the supply air flow

rate into the zone. This may not be true for some systems, for example VAV

reheat systems. Then:

d<lair = dW* (5’41)

^airCt+dt) " ^air(t) + d<Iair C5-42)

The enthalpy of the air from the air handling unit is assumed to be constant

at h

Sfl

. The mass flow rate, m

ga , can be varied. The supply air flow rate for

the new value of energy addition to the zone must then be:

msa (5-43)

The local zone VAV box will have a maximum and minimum air flow rating. If

these limits are exceeded, the air flow rate will be constrained to the

limits. If this is the case, the energy entering the room will have to be

recalculated based on the limited mass flow rate.

If the system is single deck constant volume, it is likely to have local

conditioning equipment. In this case:

d<

* local = d<*ent* (5‘ 44)

^localCt+dt) “ ^local(t) + dt*local (5-45)

If <li oca i exc eed8 the local capacity, it will remain at the local capacity. If

qiocal Become 8 negative, this represents a need for local cooling. If the

local equipment has no cooling capacity, such as a reheat coil, then qi oca i will be set to zero. If cooling capacity is present, then <l]_ oca ]_ can take on negative values up to the cooling capacity of the local equipment.

56
6. OCCUPANT COMFORT MODEL

6.1 Assumptions and Purpose of Comfort Model

The purpose of the occupant comfort model is to determine in a basic way if

the conditions in the zone are in accordance with ASHRAE standard 55-74 for

thermal comfort. This standard basically assumes that comfort is acceptable if

the zone dry bulb temperature and humidity ratio are within certain limits.

The test for comfort is only performed if the building zone is determined to

be in an occupied state. The zone is assumed to be occupied if the internal

gain schedule table, contained in the model function SCHED. predicts that the

sensible gain from people is non-zero for the current time-of-day. (see

section 5.4.5.

6.2 Comfort Envelope Equations

The equations to determine if the current state of zone air is within the

ASHRAE comfort limits are based on the simple two-dimensional geometry of the

psychometric chart. The upper and lower comfort limits for humidity ratio are constant and this allows a simple comparison of zone humidity ratio to the

limits. If the humidity ratio in the zone is greater than 0.012 or less than

0.0044 pounds of water per pound of dry air then the zone is automatically

outside of comfort limits.

The maximum and minimum comfort limits for dry bulb temperature are a linear

function of humidity ratio between the upper and lower humidity ratio limits

as expressed by:

Corner limit = 72 - 58 * 131 - 58 Vza • and (6-

1)

T

upper limit = 81 ' 73 - 384 - 74 wza

* < 6 -2>

where W za

is the zone humidity ratio in pounds of water per pound of dry air,

and T is the dry bulb temperature in degrees Fahrenheit. If the zone dry bulb

temperature is above or below these limits, then the zone is outside of

comfort limits.

The total time that the zone is outside of comfort limits during occupied

hours is determined by the comfort model subroutine. The maximum and minimum

values of zone dry bulb temperature and relative humidity are also recorded.

57
7. COMPILATION OF ENERGY USE

The predicted energy used by the emulated building during a test run of the

emulator is determined by the FORTRAN subroutine CMPILE. CMPILE uses the

values of the temperature* enthalpy* and air flow rate variables in the state

vector to calculate the energy used between the current and previous time

steps in several categories* and essentially integrates the energy use by

keeping cumulative totals which are placed in the state vector.

7.1 Categories of Energy Use

Energy used by a building is assumed to be categorized as follows:

1. energy used by the fans in the air handling unit (the power used by

the fans will be constant unless a variable air volume system is

used* in which case fan power will vary as a function of air flow).

2. energy used to heat or cool the air passing through the air handling

unit (this would not include cooling by outside air).

3. energy used by local space heating or cooling equipment (reheat

coils* heat pumps* fan coils* perimeter radiation* air

conditioners)

It is assumed that the calculated energy values are not in terms of any fuel

and do not include any plant efficiencies. The energies are either in energy

added to the air and extracted from steam or hot water for heating, or in

energy required to be mechanically extracted from the air by chilled water or refrigerant for cooling.

There are 4 primary sources of reduction in energy use which EMCS algorithms

may cause. Values for the energy reductions may be calculated from the above

types of energy use in a building. These are:

1. savings from not running electrical loads (fans* pumps).

2. savings from reducing the daily heating or cooling requirements of a building by reduction of operating hours or changing of local space

temperature setpoints.

3. savings by using the outside air for heating or cooling and avoiding

mechanical heating or cooling.

4. savings from reduction of reheating or recooling of air passing

though the air handler and building space for purposes of humidity

control and compensation for load variation between zones. These

causes of energy consumption are sometimes termed system effects,

since the amount of energy for this purpose will vary between

different HVAC systems, such as between reheat and multizone.

58
7.2 Equations for Energy Use

The first category of energy consumption is air handling unit fan energy

consumption. If the system is constant volume* the steady-state power

consumption of the fan can be used to calculate energy. If the system is

variable air volume, it is assumed that the fan power consumption will vary

with the volume flow rate through the unit. If the rated power consumption of

a fan is known at a rated flow rate* then the power consumption at another

flow rate can be determined based on fan lav number 3 from the ASHRAE

equipment handbook [1983* page 3.5* table 2]:

P — p fan rated ( ’a ' ’rated (7-1)

where: P

rate<j

• power consumption of fan at rated air flow rate.

Vrate d

1 rate<* fan volumetric air flow rate.

V„ : current volumetric flow rate,

a The fan energy between the current and last time step is determined by

multiplying Pf an

by the emulator model time step. The cumulative energy for

the fan is determined by adding the current time step energy to the sum of

values of energy for previous time steps.

The thermal energy used for heating and cooling the building should be divided

into the energy added to the building air for heating purposes and the energy

removed from the building air for cooling purposes. Therefore the energy

values will always have two parts* heating and cooling. The summation or subtraction of two energies will then consist of the separate summation or subtraction of the heating and cooling parts ( this is analogous to twodimensional vectors with x and y components ). For the following equations,

an energy rate is calculated. This rate is multiplied by the emulator model

time step to determine the energy used since the last time step. The

cumulative energy is obtained by adding this to the sum of energy values of

previous time steps.

The requirements energy is the energy required to bring the air leaving the

building zones through the return ducts and the ventilation air to the level

of the supply air enthalpy (including both air streams if dual deck).

The requirements energy rate (power) is:

qreq

~ m

sal ^ hsal~ hmc^

where the me subscript indicates the mixed air temperature in the air handling

unit with only minimum outside air or ventilation air (see section 4.2). If

q re „ < 0 this is cooling energy, if q^. > 0 this is heating energy. If the

system is dual deck there is an additional component of the requirements

energy rate for the. hot deck.

59
(7-3) ^req2 “ msa2 ^ ^sa2~ ^mc^

If q re g2 < 0 this is cooling energy* if q req 2

> 0 this is heating energy. The

hot deck component is added to the cold deck component of requirements energy*

keeping heating and cooling energies separated.

The actual energy is the energy required to bring the air passing through the

air handler to the level of the supply air enthalpy from the level of the

mixed air enthalpy. This is usually equal to the energy supplied or removed by

the heating and cooling coils. This actual energy will be less than the

requirements energy since part of the requirements energy can be met by an

economizer algorithm using outdoor air for cooling.

The actual energy rate (<l act i for single or cold deck and q act 2 f°r hot deck

if a dual deck system) is then:

^actl “ msal ^ ^sal ^ma ^ (74)

1act2 = Bsa2 ( h sa2

* h

ma ) (7 ' 5)

^actl an<* ^act2 are summe<*» keeping heating and cooling energies

separated.

The requirements and actual energy can be used to determine the economizer

energy. This is not usually considered an energy use, since the outdoor air

costs nothing to condition. However* for purposes of calculating energy

savings for economizer algorithms (savings category 3 above) it is useful to

consider this as energy consumption. The economizer energy rate is:

^econ “ ^req ~ ^act (7-6;

The actual energy can be divided into two parts which are the other two

categories of building energy use. The energy that would be required to bring

the mixed air to the actual condition of the air entering the building space

(after mixing for a dual deck system) can be termed the load energy. This will

equal the actual energy for all systems except dual deck. This represents the

energy requirements with system effects minimized (used for savings in

category 2, above). The load energy is:

^load “ ^ msal ^sal + m

sa2 ^sa2^ ^ m

sal

+ m

sa2^ ^ma

The final category of energy use is the energy used

recooling of air (used for savings in category 4, above).

^rhc " ^act ^load + ^local

This sum is of the heating and cooling components of the

(7-7)

for reheating and

This is defined as:

(7-8)

energies.

60
When a test run using the emulator is complete, the emulator model program

vill generate a summary report file containing energy usage for the emulated

building over the test period. This task is performed by subroutine REPORT.

The report will contain energy usage in the categories discussed above for

heating and cooling energy components. An example energy report is given in

table 7-1. The report file will also contain results on the amount of time

the zone was not within comfort limits and the maximum and minimum space

temperatures and relative humidities reached during occupied hours.

Table 7-1. Heating and cooling energy from emulator model

FAN ENERGY —*» +

-+ +

I HEATING | COOLING

ECONOMIZER ENERGY

LOAD ENERGY

REHEAT/RECOOL ENERGY

I +

I •+ •+

61
8. SUMMARY

The use of a building emulator is one of the approaches to testing of EMCS

systems without making actual connections to HVAC equipment or a building.

The NCEL/NBS emulator will be used for either factory or field acceptance

tests of EMCS algorithms. The building simulation program in the emulator is

the main part of this report. This program simulates the responses of a building to EMCS control actions.

The building emulator program was written in FORTRAN 77 and the source code is

appended in the APPENDIX. State variables in the emulator model are analog

sensor* energy* command* and comfort variables. The building emulator model

program contains the weather (WEATHR)* air handling unit (AHU)» zone (ZONE)*

and comfort (COMFRT) model. The use of energy during a simulation is

accumulated by the energy compilation subroutine (COMPILE). This energy use information is useful to evaluate the performance of an EMCS algorithm.

Prior to simulation* data files of initial values of state variables* building

shell parameters* HVAC system parameters, and weather data parameters must be

called. The air handling unit model is a steady-state model* while the zone model is a dynamic model. In the zone model* a system of first order ordinary

differential equations is solved by the explicit forward difference scheme.

The time step for the difference scheme is adjusted according to the changes

of state variables. By using this adaptive time step* which is uniformly

distributed in a given emulator simulation time period, savings in computation

time can be realized when the state variables approach steady state.

The air handling unit model simulates a single deck (or single duct)* constant

volume system without humidity control as shown in Figure 4-2. Damper

positions of the air handling unit are varied if economizer option is

selected. A single zone with exterior walls without window is modeled under

the assumptions described previously. For this abridged system, a local

equipment model is implemented. Heating/cooling energy delivered from the

local equipment to zone air is adjusted under the control law of the

proportional plus integral controller (Pi-controller). Capabilities of

simplified current models can be upgraded by providing necessary information

in the spaces already reserved for future extension.

62
9. REFERENCES

ASHRAE, 1981 Fundamentals Handbook. ASHRAE. 1981 .

Berglund, L., Control and Simulation of Thermal Comfort.T Workshop on HVAC

Control. Modeling, and Simulation. Atlanta. GA. Feb. 2-3. 1984.

Borresen, B.A., dermal Room Models for Control Analysis,* ASHRAE Trans.,

Vol . 87, Part 2. 1981.

Brokaw, R.S., ’Calculation of Flue Losses for High-Efficiency Furnaces and

Appliances,' ASHRAE Journal. Jan. 1979, pp. 49-51.

Building Environment Division of NBS, 'Technical Guidelines for Energy

Conservation,' AFCEC-TR-77-12 and NBSIR 77-1238, Nat'l Bureau of

Standards, June 1977.

Bullock, C.E., 'Dynamic Simulation Models for Commercial Air Conditioning and

Heat Pump Systems,' Workshop on HVAC Control, Modeling, and Simulation,

Atlanta. GA, Feb. 2-3, 1984.

Chi, J., 'Building Environmental Control System Analysis Computer Program - Volume 1 Component Mathematical Models,' Final Report submitted to NBS,

September 1982.

Fanger, P.O., Thermal Comfort . McGraw-Hill, 1970.

Ferziger, J.H., Numerical Methods for Engineering- AppJ.ig.a-tl.on, John Wiley,

1981.

Kao, J.Y., 'Strategies for Energy Conservation for a Large Office Building,'

NBSIR 83-2746, Nat'l Bureau of Standards, July 1983.

Kusuda, T., 'NBSLD, The Computer Program for Heating and Cooling Loads in

Buildings,' NBS BSS 69, Nat'l Bureau of Standards, 1978.

May, W.B., 'Time of Day Control and Duty Cycling Algorithms for Building

Management and Control Systems,' NBSIR 83-2713, Nat'l Bureau of

Standards, June 1983.

May, W.B., 'Control Algorithms for Building Management and Control Systems--

Hot Deck/Cold Deck/Supply Air Reset, Day/Night Setback, Ventilation

Purging, and Hot and Chilled Water Reset,' NBSIR 83-2713, Nat'l Bureau of

Standards, March 1984.

Modest, M.F. and Macken, N., 'Seasonal Performance Comparisons for an Apartment Building in New York State Climate Regions,' ASHRAE Trans.,

1977, pp. 881-892.

63
Park* C. ’Single-zone Computer Model for Residential Furnace Location

Analysis,' ASHRAE Trans., Vol. 87, Part 2, 1981.

Park, C., Kelly, G.E. and Kao, J.Y. 'Economizer Algorithms for Energy

Management and Control Systems,' NBSIR 84-2832, Nat'l Bureau of

Standards, Feb. 1984.

Park, C., 'An Optimum Start/Stop Control Algorithm for Heating and Cooling

Systems in Buildings,' NBSIR 83-2720, Nat'l Bureau of Standards, May

1983.

Park, C., 'Demand Limiting Algorithms for Energy Management and Control

Systems,' NBSIR 84-2826, Nat'l Bureau of Standards, February 1984.

Sowell, E.F. and Walton, G.N., 'Efficient Computation of Zone Loads,' ASHRAE

Trans., Vol. 86, Part 1, 1980.

Threlkeld, J.L., Thermal Environmental Engineering. 2nd ed., Prentice-Hall,

1970.

Walton, G.N., 'A New Algorithm for Radiant Interchange in Room Loads

Calculations,' ASHRAE Trans., Vol. 86, Part 2, 1980.

Walton, G.N., 'Thermal Analysis Research Program Reference Manual,' NBSIR 83-

2655, Nat'l Bureau of Standards, March 1983.

Wise, B.B., 'Energy Monitoring and Control System Software Testing Device,'

Ninth Energy Management & Controls Society Conference Proceedings, Nov.

1984.

64
APPENDIX

COMPUTER PROGRAM LISTINGS

65
12

3

4

5

6

7

8

9

10

1

12

13

14

13

16

17

18

19

20

21

22

23

24

25

26

27

28

29

30

31

32

33

34

33

36

37

38

39

40

41

42

43

44

43

46

47

48

49

30

31

32

33

34

33

36

37 o

oo

oo

o

oo

o

o

o

oo

o

o

o

o

o

o

o

o

oo

o

ooo

ooo

ooo

C Version 1.7 - FEBRUARY 03, 1983 - W.B. MAY .NATIONAL BUREAU OF STANDARDS

Ciiis3i9>3itiiiiii>i3333:i:333>iiiiiii3i33333:9::::ti>i::t:::::33:::::::

PROGRAM ALGTST

C3I333333333I3I3I3III3333II3333333I3133333333:33333333383I31333333333333

C This is the main program used to develop and test simulation routines

C for a simple building model to be used in an EMCS algorithm tester

C developed by the U.S. Navy Civil Engineering Laboratory. This test

C program takes the place of the EMCS tester software. Subroutine

C EMODEL is the emulator model. Subroutine UPDATE represents the I/O

C task, which determines if any commands have arrived from the EMCS and

C updates the analog and digital signals to the EMCS based on the

C simulation calculations. It is assumed that these routines will

C eventually run in parallel, UPDATE at a higher frequency. Parallel

C execution is approximated by serial calls, but with UPDATE being

C called more often.

C

C Determine start time

CALL SYNCCL

Update commands and outputs

000 CALL UPDATE

Call the emulator model

CALL EMODEL

Update clock.

CALL ECLOCK

GO TO 1000

END

38133333339>3l98l9l33333333339338333383333333333313333833tr:33:33:3

BLOCK DATA COMDAT

33333333333333333383333333333333388333333333333333333399333333333::

This module contains a list of all common block data and globals which

are referenced in the support routines for the EMODEL subroutine.

Additional common blocks used only with EMODEL are described in the

subroutine EMODEL or the subroutine PARAMA

common blocks are: ASENSR - contains all analog sensor readings in arrays. Sensor

readings are copied from the emulator model outputs after

each time step. Array index Indicates which of 5 possible

buildings the reading is from.

COMMND - contains all command digital inputs from the EMCS in arrays.

readings are directly available to the emulator model.

CPA - contains all control point adjustments from the EMCS in arrays.

DSENSR - contains all digital sensor readings in arrays.

EVENT - GLOBAL which is a logical variable Indicating an EMCS command.

PAST - contains day, hour, minute, and second that model subroutine

PFILES - contains the names of all parameter files, and the season. was last called, and the simulation major timestep.

STATE 1 - contains the state variables for the emulator model.

TIME - current time is broken into DAY, HOUR, MINUTE, SECOND. In this

test program, all timos ore zero at the start of the program,

the basic time step of the program is one second. ENDDAY , ENDHR

ENDMIN, and ENDSEC are the scheduled time to end the simulation.

66
38

39

60

61

62

63

64

63

66

67

66

69

70

71

72

73

74

73

76

77

78

79

80

81

82

83

84

83

86

87

88

89

90

91

92

93

94

93

96

97

98

99

100

101

102

103

104

103

106

107

108

109

110

1 1 1 112

113

1 14

113

oooooooooo

oooo

ooo

oooo

UNITS - contains the logical units used by the program for various

Inputs and outputs.

Temperature, Outside Air (F)

Relative Humidity, Outside Air (X)

Temperature, Space (F)

Relative Humidity, Space (X)

Temperature, Return Air (F)

Relative Humidity, Return Air (X)

Position, Reheat Valve (Btu/s)

Temperature, Supply Air (F)

Temperature, Mixed Air (F)

c TOUTA - EMCS ana l og input

c RHOUTA - EMCS ana l og input

c TZONEA - EMCS ana l og input

c RHZONA - EMCS analog input

c TRETA - EMCS ana l og input

c RHRETA - EMCS ana l og input

c ZHCDEM - EMCS ana l og input

c TSUPA - EMCS ana l og input

c

c TMIXA • EMCS ana l og input

REAL TOUTA , RHOUTA , TZONEA , RHZONA , TRETA

& RHRETA, ZHCOEM, TSUPA , TMIXA

COMMON / ASENSR / TOUTA, RHOUTA, TZONEA ( 3) , RHZONA ( 5 ), TRETA ( 5 )

& RHRETA ( 3 ) , ZHCDEM C 5 > , TSUPA ( 5 )

, TM I XA ( 5

DATA TOUTA, RHOUTA, TZONEA , RHZONA , TRETA, RHRETA, ZHCDEM/ 27* 0 . 0/

DATA TSUPA, TMIXA/ 10* 0.0/

ON - EMCS command, Supply Fan On/Off

ECON - EMCS command, O.A. Damper Off^Auto

VENT - EMCS command, Ventilation Damper Open/Close

SETBAK - EMCS command. Zone Thermostat Setback On/Off

LOGICAL ON, ECON, VENT, SETBAK

COMMON / COMMND / ON ( 3 ) , ECON ( 5 ) , VENT ( 5 ) , SETBAK C 5

DATA ON, ECON, VENT, SETBAK/20* . FALSE.

/

SUPPLY - EMCS CPA, Supply Air Temperature <F)

ZONE EMCS CPA, Zone Temperature (F)

REAL SUPPLY , ZONE

COMMON / CPA / SUPPLY ( 31 , ZONE (5)

DATA SUPPLY, ZONE/ 10*0. 0/

DPSSF - EMCS digital input, differential Pressure Switch, Supply Fan

LOGICAL DPSSF

COMMON / DSENSR / DPSSFC3)

DATA DPSSF/3*. FALSE./

EVENT - Logical variable which is true when a command from the EMCS

file has been received.

LOGICAL EVENT

GLOBAL EVENT

OLDOAY

OLDHR

OLDMIN

OLDSEC

STEP

Day number when EMODEL subroutine was last called.

Hour of day when EMODEL subroutine was last called.

Minute of hour when EMODEL subroutine was last called.

Second of minute when EMODEL subroutine was last called.

Major time step for simulation, model routines for AHU , weather, zone, comfort, and energy compilation are called

at this Interval. Not the same as the Euler timestep used

within the zone model.

67
116

117

118

119

120

121

122

123

124

129

126

127

128

129

130

131

132

133

134

133

136

137

138

139

140

141

142

143

144

145

146

147

148

149

130

131

132

133

134

133

156

137

138

139

160

161

162

163

164

163

166

167

168

169

170

171

172

173

OOOO

! NTEGER OLDDAY , OLDHR , OLDM I N , OLDSEC , STEP

COMMON / PAST / OLDDAY , OLDHR , OLDM I N , OLDSEC , STEP

DATA OLDDAY , OLDHR , OLDM I N , OLDSEC , STEP/4 * 0 , 30/

C

C HVFILE - Character filename of HVAC equipment parameter file.

C USFILE - Character filename of Building usage parameter file.

C CLFILE * Character filename of climate description parameter file.

C SHF I LE - Character filename of building shell parameter file.

C INFILE - Character filename of state vector initialization file.

C SEASON - Character name of weather season to use in test.

C

CHARACTER* 1 3 HVF I LE , USF I LE , CLF 1 LE , SHF I LE , I NF I LE , SEASON

COMMON / PFILES / HVF I LE , USF l LE , CLF I LE , SHF 1 LE , I NF I LE , SEASON

A - Real array containing state variables at the current time.

L - Logical array containing current logical state variables.

C

C

C

C

C

C

C

C

C

C

c

c

cc

c

c

c

c

LOGICAL L

REAL A

COMMON / STATE1 / At 1 00) ] LC 1 00)

DATA A/1 00*0. 0/, L/1 00* . FALSE.

/

DAY

HOUR

MINUTE

SECOND

ENDDAY

EMDHR

ENOMIN

ENDSEC

Current day number, first day of test is day 1. Current hour of day, from 0-23.

Current minute within the hour, 0-59.

Current second within the minute, 0-39.

Day number on which the test Is to stop.

Hour of the day on ENDDAY within which test is to stop.

Minute within ENDHR within which test is to stop.

Second within ENDMIN within which test is to stop.

I NTEGER DAY , HOUR , M I NUTE , SECOND , ENDDAY , ENDHR , ENDM I N , ENDSEC

COMMON /TIME / DAY. HOUR, Ml NUTE, SECOND, ENDDAY, ENDHR, ENDM I N, ENDSEC

CLU - Logical Unit for writes to console

FLU - Log 1 ca

Unit for read from EMCS command i nput file

I LU - Logical Un i t for read from console.

PLU - Log i ca

Un i t for read from parameter files.

I NTEGER CLU , FLU , I LU , PLU

COMMON / UNITS / CLU, FLU, I LU, PLU

DATA CLU/ 1 / , FLU/8/ , I LU/ 1 / , PLU/9/

END

BaaBaasaaaaaBaaaaaaBaaaaazaaaasaazBzaazaaaaaaaaaaaaaaaeaazxaaaaaaaasaas

SUBROUTINE SYNCCL

Caa<BaKaa * ax *saaaa8aaaa3aaaa3aaaaaaaaa333BS3aa3aaaa3eaaa3a3aasaaaaa3a8aa = 3 C This routine is used to determine the start time of the test. In a

C real-time system, this routine would set the real time clock, to the

C proper time.

C

C I nput var tables: none

C Output variables: none

C Procedures called: none

C Common blocks: UNITS, TIME

C

I NTEGER CLU , FLU , I LU , PLU

COMMON / UNITS / CLU, FLU, l LU, PLU

68
174

175

176

177

170

179

100

101

102

103

104

105

106

107

100

109

190

191

192

193

194

195

196

197

190

199

200

201

202

203

204

205

206

207

200

209

210

21 1 212

213

214

215

216

217

210

219

220

221

222

223

224

225

226

227

220

229

230

231

ooooo

ooooo

o

o

o

oooooooooo

o

o

oooo

oooo

I NTEGER DAY , HOUR . M I NUTE , SECOND , ENDDAY , ENDHR , ENDM I N , ENDSEC

COMMON / TIME / DAY , HOUR , M I NUTE , SECOND , ENCDAY , ENDHR , ENDM I N , ENDSEC

Local Variables:

CFFILE - Character name of file containing configuration information.

CHARACTER* 15 CFFILE

PARAMETER (CFFILE * 'BLDG. CNF')

Read the configuration file to determine the start time.

OPEN! PLU, FI LE*CFFILE, STATUS* 'OLD' ,ERR=9S00)

READ( PLU, * ) HOUR, MI NUTE, SECOND

CLOSE (PLU)

DAY * 1 RETURN

ERRORS

1000 WRI TEC 1,1)

1 FORMAT! IX, 'Cannot successfully open conf i gurat 1 on file')

STOP

END

BasssssssasssBasasasssassssssssassBzszsssssszssssssssssssxssssssssssss

SUBROUTINE CONPAR

8883838388X3838383833388833838388838333333X333333333333X3X333333333333

This routine is used to determine all of the emulator routine control

parameters including duration of test, parameter file

names and weather season.

Input variables: none Output variables: none Procedures called: none

Common blocks: UNI TS, TI ME, PAST, PFI LES

I NTEGER CLU , FLU , I LU , PLU

COMMON / UNITS / CLU, FLU, I LU, PLU

I NTEGER DAY , HOUR , M I NUTE , SECOND , ENDDAY . ENDHR , ENDM I N , ENDSEC

COMMON / TIME / DAY, HOUR, MI NUTE, SECOND, ENDDAY, ENDHR, ENDMI N, ENDSEC

I NTEGER OLDDAY , OLDHR , OLDM I N , OLDSEC , STEP

COMMON / PAST / OLDOAY, OLDHR, OLDM I N, OLDSEC, STEP

CHARACTER* 15 HVFI LE, USFI LE, CLFI LE, SHFI LE, I NFI LE, SEASON

COMMON / PFILES / HVFI LE. USFI LE, CLFI LE, SHFI LE, I NFI LE, SEASON

Internal Variables:

CFFILE - Character name of file containing conf 1 gurat i on Information.

CHARACTER* 15 CFFILE

PARAMETER ( CFF I LE = ' BLDG . CNF '

Read the configuration file to determine the test duration, HVAC type

file, building use file, weather data file, building type file,

inlt lal izat ion file, season, and simulation time step.

69
232

233

234

235

236

237

238

239

240

241

242

243

244

245

246

247

248

249

250

251

252

253

254

255

256

257

258

259

260

261

262

263

264

265

266

267

268

269

270

271

272

273

274

275

276

277

278

279

280

281

282

283

284

285

286

287

288

289

O

O

O

C

0PEN( PLU, FI LE=CFFI LEj STATUS* ' OLD ' , ERR=9800

REAO( PLU, *

READ ( PLU , * ) ENODAY , ENDHR . ENDM I N , ENOSEC

READ ( PLU, * ) HVFILE

READ (PLU,*) USFILE

READ (PLU, * ) CLFILE

READ (PLU,*) SHFILE

READ ( PLU, * ) INFILE

READ (PLU,*) SEASON

READ (PLU,*) STEP

CLOSE (PLU)

Compute end I no time

ENDOAY = ENDDAY DAY

ENDSEC = ENDSEC SECOND

I F ( ENDSEC . GT . 59 ) THEN

ENDSEC = 0

ENDM I N = ENDMIN 1 END IF

ENDMIN = ENDMIN MINUTE

I F ( ENDM I N . GT . 59 ) THEN

ENDMIN = 0

ENDHR * ENDHR 1 END IF

ENDHR = ENDHR HOUR

I F ( ENDHR . GT . 23 ) THEN

ENOHR * 0

ENDOAY a ENODAY 1 END IF

OLDDAY » DAY

OLDHR a HOUR

OLDMIN a MINUTE

OLDSEC = (-STEP)

RETURN

C-- ERRORS - - - 9800 WR I TE ( 1,1)

1 FORMAT! 1 X, ' Cannot successfully open configuration file*)

STOP

END

Ciaisiaiisiiiiiitai3ii33iii3i33a:xi«iE3C333:::33:3issi:ax23s:::3::::;:::

SUBROUTINE ECTROL( TIMEST, GO)

C* 3IIIII3lllllll3lll333II1333ai3:i:iI3333333:33::3l3i:illI3X3X33::: == -

C This routine Is used to control the execution timing of the emulator

C model routine. If GO is returned as TRUE, the emulator simulation

C should be executed. This occurs when 1. the time step duration has

C passed, and 2. a control event has been detected. TIMEST contains the

C time step that the model is to use in simulation.

C

C I nput var tables: none

C

C Output variables:

C TIMEST - The current value to be used for the major simulation time step.

C GO - A logical variable which indicates that the simulation

C subroutine should be executed If true.

LOGICAL GO

INTEGER TIMEST

70
290

291

292

293

294

299

296

297

298

299

300

301

302

303

304

303

306

307

308

309

310

31 1 312

313

314

313

316

317

318

319

320

321

322

323

324

323

326

327

328

329

330

331

332

333

334

333

336

337

338

339

340

341

342

343

344

345

346

347

oooo

oooooooo

o

oooooo

ooooo

ooooooooo

Procedures called: NUSTEP, TIMINC

Common blocks: none

Global Variables:

EVENT - Logical variable which is true when a command from the EMCS

file has been received.

LOGICAL EVENT

GLOBAL EVENT

Local Variables:

STEPDN - Logical variable, true if time for another simulation execution.

TIMEUP - Logical variable, true if test ending time has been reached

LOGICAL TIMEUP, STEPDN

Update the time and determine if the execution tlmestep has passed.

If so, STEPON will be true. Also determine if the ending time has

been reached. If so, TIMEUP will be true.

CALL T I MI NCI STEPDN, TIMEUP)

I F( TIMEUP) THEN

CALL REPORT

STOP ' Emu l at 1 on Comp late

END IF

I F ( EVENT . OR . STEPDN ) THEN

GO .TRUE.

STEPDN = .FALSE.

CALL NUSTEP ( T I MEST

ELSE

GO * .FALSE.

ENOIF

RETURN

END

1883333333388333X8383338833X3333838X333333333333333X3333833333333333333

SUBROUTINE NUSTEP ( T I MEST

1883338333383333383333383X333333333333333333333333333333333333333333333

This routine returns the amount of time in seconds since the last time

it was called.

Input variables: none

Output variables:

TIMEST * Time in seconds since routine was last called.

INTEGER TIMEST

Procedures called: SBTIME

Common blocks: TIME, PAST

I NTEGER DAY , HOUR , M I NUTE . SECOND , ENDDAY , ENDHR , ENDM 1 N , ENDSEC

COMMON / TIME / DAY, HOUR, MI NUTE, SECONC, ENDDAY , ENDHR , ENDM I N , ENDSEC

C

I NTEGER OLDDAY , OLDHR , OLDM I N , OLDSEC , STEP

COMMON / PAST / OLDDAY , OLDHR , OLDM I N , OLDSEC , STEP

71
346

349

330

331

352

333

354

353

336

337

356

359

360

361

362

363

364

365

366

367

368

369

370

371

372

373

374

373

376

377

378

379

380

381

382

383

384

383

386

387

388

389

390

391

392

393

394

393

396

397

398

399

400

401

402

403

404

403

oooo

ooooooooo

o

o

ooo Loca 1 Var tables: none — - - determ t ne time si nee 1 ast call

CALL SBT I ME ( DAY , HOUR , M I NUTE , SECOND , OLDDAY , OLDHR , OLDM I N , OLDSEC

& TIMEST)

....... reset last time called variables

OLDDAY = DAY

0LDHR=H0UR

OLDM! N=M I NUTE

OLDSEC= SECOND

RETURN

END

I3ii:>:i:33iiiiiiiii::3isc333»3iiiii:i:39::it:::3:93t>iin:3r:::::::::

SUBROUTINE ECLOCK

33933333333333I33ISI3333333333933I3I33333333I33333333333393333333333333

This routine Increments the time. In a real-time system, the time would

be obtained from a real-time clock..)

I nput var tables: none Output var i ab 1 es : none Procedures called: none Common blocks: TIME

I NTEGER DAY , HOUR , M I NUTE , SECOND , ENDDAY , ENDHR , ENDM I N , ENDSEC

COMMON / TIME / DAY, HOUR, MI NUTE, SECOND, ENDDAY , ENDHR, ENDM I N, ENDSEC

Local Variables: none

1 ncrement t i me SECOND * SECOND 1

I F ( SECOND . GE . 60 ) THEN

SECOND > 0

MINUTE = MINUTE 1

I F( MI NUTE. GE. 60) THEN

MINUTE = 0

HOUR = HOUR «• 1

I F ( HOUR . GE . 24 ) THEN

HOUR = 0

DAY * DAY 1 END IF END IF

END IF RETURN

ENO

C3339S3S333939I339I3333I3I33399933333333333333333393I339IC3I3333333333S:

c

c

c

c

cc

c

V*c

c

c

SUBROUTINE T I M I NC ( STEPDN , T I MEUP

399B393l339393933l333333933333339333t939i9S3339333393SI3lt3333333333333

This routine determines If the

desired duration of the execution has passed (TIMEUP=true) , and if the

simulation tlmestep has passed since the last execution of the

simulation subroutine (STEPDN=true )

I nput var tables : none Output variables:

STEPDN - Logical variable,

T I MEUP - Logical variable.

true

true

if

if

time for another

test ending time

simulation execut

has been reached

on

LOGICAL STEPDN, T I MEUP

72
406

407

406

409

410

41 1 412

413

414

413

416

417

418

419

420

421

422

423

424

423

426

427

426

429

430

431

432

433

434

433

436

437

438

439

440

441

442

443

444

443

446

447

448

449

430

431

432

433

434

433

436

437

438

439

460

461

462

463

ooooo

ooo

o

o

o

o

oo

o

o

oo

o

oo

o

o

oo

ooooo

oooo Procedures called: none

Common blocks: TIME, PAST

I NTEGER DAY, HOUR, MI NUTE, SECOND, ENDDAY, ENDHR, ENDMI N, ENDSEC

I NTEGER OLDOAY , OLDHR , OLDM I N , OLDSEC , STEP

COMMON / TIME / DAY, HOUR. Ml NUTE, SECOND, ENDDAY, ENOHR, ENDMI N, ENDSEC

COMMON / PAST / OLDOAY, OLDHR, OLDM I N, OLDSEC, STEP

Local Variables:

JTIME * Integer used to contain time differential in seconds.

INTEGER JTIME ................. check for end of program

CALL SBT I ME ( ENDDAY , ENDHR , ENOM I N , ENDSEC , DAY , HOUR , M I NUTE , SECOND

6 JTIME)

IF( JTIME. LT.O) THEN

T1MEUP = .TRUE.

ELSE

T1MEUP = .FALSE.

END IF ......... check, for end of time step ; CALL SBT I ME ( DAY , HOUR , M I NUTE . SECOND , OLDDAY , OLDHR , OLDM I N , OLDSEC

& JTIME)

IF( JTIME. GE. STEP) THEN

STEPON a .TRUE.

ELSE

STEPON a '.FALSE.

END IF

RETURN

END

iss:rx::iiii:i323SB:::33::3:i3st3e:K:s:2::::s:::r::2:sx3::::::::::::r:

SUBROUT 1 NE SBT I ME ( DAY , HOUR , M I N , SEC , MOAY , MHOUR , MM I N , MSEC , DSEC

This module subtracts two sets of time integers (xxx - Mxxx) and

puts the result in seconds Into DSEC.

I nput

DAY

HOUR

MIN

SEC

MDAY

MHOUR

MMIN

MSEC

var tables:

- Day number of first time from which second time Is subtracted.

- Hour of first time.

- Minute of first time.

- Second of first time.

- Day number of second time which is to be subtraced from first.

- Hour of second time.

- Minute of second time.

- Second of second time.

Output variables:

DSEC - differential in seconds between first and second times.

I NTEGER DAY , HOUR , M I N , SEC , MDAY , MHOUR , MM I N , MSEC , DSEC

Procedures called: none

Common blocks: none

Local Variables:

73
464

469

466

467

468

469

470

471

472

473

474

479

476

477

478

479

480

481

482

483

484

489

486

487

488

489

490

491

492

493

494

499

496

497

498

499

900

901

902

903

904

909

906

907

908

909

910

91 1 912

913

914

919

916

917

910

919

920

921

CDAY - Worn, i no copy of

CHOUR - Work. 1 ng copy of

CHIN - Work.

no copy of

CSEC - Work. 1 no copy of

DDAY - time differenti*

OHOUR - time differentii

OMIN - time different it

day from first time,

hour from first time,

minute from first time,

second from first time.

1 1 , day component

1 1 , hour component

il, minute component.

I NTEGER COAY , CHOUR , CM I N , CSEC , DDAY , DHOUR , DM I

CDAY * DAY

CHOUR s HOUR

CMIN » MIN

CSEC * SEC

DSEC » CSEC - MSEC

!F(DSEC. LT. OJTHEN

DSEC » DSEC 60

CMIN = CMIN - 1 END IF

DMIN a CMIN - MMIN

I F ( DM I N . LT . 0 ) THEN

DMIN * DMIN 60

CHOUR a CHOUR - 1 END I F

DHOUR a CHOUR - MHOUR

I F ( DHOUR . LT . 0 ) THEN

DHOUR * DHOUR 24

CDAY a CDAY -1 END I F

DDAY a CDAY - MDAY

DSEC a DSEC + (60«DMIN) (3600«0H0UR) (86400aDDAY)

RETURN

END

Cmiiassaaia i s r 3 a r x a 133333333331

SUBROUTINE UPDATE

: 9 a 9 s 9 a i iraaaaaaa:

This routine is intended to emulate a real time I/O task which is

connected to hardware I/O channels. Instead, it receives input from

a 'boundary value' file to define changes In sensor readings.

The following sensors and actuators are supported:

EMCS DIGITAL OUTPUTS VARIABLE NAME 0 1 Supply Fan On/Off ON off on

O.A. Damper Off /Auto ECON off auto

Ventilation Damper Open/Close VENT open close

Zone Thermostat Setback On/Off SETBAK off on EMCS ANALOG OUTPUTS

Control Po I nt Adjustment, Supply Air Temperature (F)

Control Point Adjustment, Zone Temperature (F)

VARIABLE NAME

SUPPLY

ZONE

EMCS ANALOG INPUTS

Temperature, Outside Air (F)

VARIABLE NAME

TOUTA

74
322

323

324

323

326

327

328

329

330

331

332

333

334

333

336

337

338

339

340

341

342

343

344

343

346

347

348

349

330

331

332

333

334

333

336

337

538

339

360

361

562

363

364

363

566

367

368

369

370

371

372

373

374

373

376

377

378

379

oo

o

o

C Temperature, Mixed Air (F) TMIXA

C Temperature, Supply Air (F) TSUPA

C Temperature, Return Air (F) TRETA

C Temperature, Space (F) TZONEA

C Relative Humidity, Outside Air (X) RHOUTA

C Relative Humidity, Return Air (X) RHRETA

C Relative Humidity, Space (X) RHZONA

C Position, Reheat Valve (Btu/s) ZHCDEM

C

C EMCS DIGITAL INPUTS

C Differential Pressure Switch, Supply Fan DPSSF

C

C I nput var tables: none

C Output variables: none

C Procedures called: RELHUM

C Common blocks: UNI TS, TIME, ASENSR, COMMND, DSENSR, CPA, STATE

C

INTEGER CLU, FLU, ILU, PLU

COMMON /UNITS / CLU, FLU, I LU, PLU

C

I NTEGER

COMMON /

C

REAL

&

COMMON /

&

C

LOGICAL

COMMON /

C

LOGICAL

COMMON /

C

REAL

COMMON /

C

REAL

LOGICAL

COMMON /

DAY, HOUR, MINUTE SECOND, ENDOAY, ENOHR, ENDM I N , ENDSEC

TIME / DAY, HOUR, MINUTE, SECOND, ENDDAY, ENDHR, ENDMIN, ENDSEC

TOUTA , RHOUTA , TZONEA , RHZONA , TRETA

RHRETA, ZHCDEM, TSUPA , TMIXA

ASENSR / TOUTA, RHOUTA, TZONEA (3), RHZONA ( 3 ), TRETA ( 3),

RHRETA C 3 ) , ZHCDEM ( 3 ) , TSUPA ( 5 )

, TM I XA i 5

ON, ECCN, VENT, SETBAK

COMMND / ON ( 3 ) , ECON ( 3 ) , VENT ( 3 ), SETBAK ( 3

DPSSF

DSENSR / DPSSF (3)

SUPPLY , ZONE

CPA / SUPPLY ( 5 ), ZONE ( 3

STATE/

A

L

A< 100) , L( 100)

Global Variables:

EVENT - Logical variable which is true when a command from the EMCS

file has been received.

LOGICAL EVENT

GLOBAL EVENT

C

C Local Variables:

C ANSWER - character variable to receive answer on debug print question

C DATE - day number for the next command event

C DBINIT - logical variable, true if debug print question has been asked

C DBP - logical variable, true if debug print is enabled.

C EX - exponential of ST.

C LEFT - seconds left before next command event will take place.

C NBLD - building number used as index in EMCS variable arrays.

C ST - negative reciprocal of TCSENS.

C TCSENS - time constant of EMCS analog sensors in seconds.

C TIME ~ time In hours, minutes, and seconds for the next command event

75
580

581

582

583

584

585

586

587

588

589

590

591

592

593

594

595

596

597

598

599

600

601

602

603

604

605

606

607

608

609

610

61

612

613

614

615

616

617

618

619

620

621

622

623

624

625

626

627

628

629

630

631

632

633

634

635

636

637

ooooo

oooooo

oooooooo

ooo

oooo WAIT - logical variable, True if command event has been read from

command event file and is being stored, Felse if no command

event is stored.

REAL EX, ST, TCSENS

LOGICAL WAIT, DBP, OBI NIT

INTEGER DATE, TIME (3), LEFT

CHARACTER* 1 ANSWER

PARAMETER ( NBLD= 1 , TCSENS=15.0, ST* ( - 1 . 0) /TCSENS)

DATA DBP /.FALSE./, DBI Nt T/ . FALSE. /, WA1 T/ . FALSE.

/

NAMELIST /DEBUG/ DAY, HOUR, MI NUTE, SECOND, TOUTA, RHOUTA, TZONEA,

&RHZ0NA , TRETA, RHRETA, ZHCOEM, TSUPA, TMI XA, DPSSF

Reed the boundary value file to see if EMCS has sent a command. THe

command will have a time associated with it which is the time the

command was given. The command data is not changed until the command

time. The variable WAIT, If true, indicates that the routine is waiting

for the command time to arrive before reading the boundary file again.

I F( .NOT. WAIT) THEN

500 READ! FLU, * , ERR=500, END* 900) DATE, TIME, ONI NBLD) , ECONC NBLD)

8

• VENT ( NBLD ) , SETBAK ( NBLD )

,, 8 SUPPLY (NBLD), ZONE (NBLD)

WAIT « .TRUE.

END IF

Determine if the time for the command obtained from the boundary file

has come. If so, then change the command data in common, and set the

event flag to indicate a command has occurred. LEFT is the time left

before the command becomes active.

100 I F ( WAI T) THEN

CALL SBT I ME ( DATE , T I ME ( 1 ), T I ME ( 2 ), TME ( 3 ) , DAY , HOUR , M I NUTE , SECOND

8 LEFT

I F ( LEFT . LE . 0 ) THEN

L( 1 )*0N(N3LD)

L ( 2 )

ECON ( NBLD

L ( 3 ) * VENT ( NBLD

L ( 4 ) *SETBAK ( NBLD

A ( 80 ) ^SUPPLY ( NBLD

A( 81 ) = ZONE (NBLD)

WAIT*. FALSE.

EVENT* . TRUE.

ELSE

EVENT*. FALSE.

END IF

ELSE

EVENT * .FALSE.

END IF Update the sensor data in common with the data from the state variables

calculated by the emulator model. To simulate sensors with finite time

constants, an exponential approach Is used.

76
630

639

640

641

642

643

644

649

646

647

640

649

690

691

692

693

694

699

696

697

690

699

660

661

662

663

664

669

666

667

660

669

670

671

672

673

674

679

676

677

670

679

660

661

662

603

664

609

666

607

606

609

690

691

692

693

694

699

OOOOOOOOO

DPSSF(NBLD) • L(10)

C

EX*EXP(ST)

TOUTA

TRETA(NBLD)

TZONEA ( NBLD ) TSUPA(NBLD) « TMIXA(NBLD) * ZHCDEM ( NBLD ) »

(1 . -EX)*A(01

(1 . -EX ) *A ( 20)

(1 . -EX) *A(06)

( 1

. -EX) *A( 32)

( 1

. -EX) *A(23)

A( 62)

EX* TOUTA

EX*TRETA< NBLD)

EX* TZONEA (NBLD)

EX* TSUPA ( NBLD

EX*TM1 XA( NBLD)

RHOUTA = (1 . -EX)*RELHUM(A(01 ) ,A(02) ,29.921 ) +EX*RHOUTA

RHRETAC NBLD) = ( 1

. -EX) *RELHUM( A( 20) , A ( 22) , 29 . 921 ) +EX*RHRETA ( NBLD

RHZONA ( NBLD ) = ( 1

. -EX) *RELHUM( A( 06) , A< 07) , 29 . 921 ) +EX*RHZONA ( NBLD

C

IRC .NOT. DBI NIT) THEN

WR I TE ( 1 , FMT*’ ( 1 X

, " DEBUG PRINTOUT?") ’)

READ( 1, FMT* ' (A 1 )

1 ) ANSWER

I F ( ANSWER . EQ . ' Y * ) DBP * .TRUE.

DBI NIT * .TRUE.

END IF

IF(DBP) WRITEt 1 , DEBUG)

RETURN

END

C Version 1.28 - FEBRUARY 13, 1989 - W.B. MAY .NATIONAL BUREAU OF STANDARDS

C***********************************************************************

SUBROUTINE EMODEL

C This Is the mein subroutine for the simulation routines for the EMCS

C algorithm tester. The simulation sections consist of the weather data

C section, the Air Handling Unit simulation, the zone simulation, the

C local heating and cooling equipment simulation, the comfort

C simulation, and the data compilation routines. NOTE: all common block.

C variables are described in subroutine PA RAMA.

C

C I nput var tables: none

C Output variables: none

C Procedures called: AHU, CMPI LE. COMFRT, CONPAR, ECTROL, LOCAL, PARAMA . WEATHR

C ZONE

C Common blocks: CL I MAT, STATE

C

REAL WOAVG , CW I ND , TDBAVG , TDBAMP , S I NORG , CSOLR

COMMON / CL I MAT / WOAVG, CW I ND, TDBAVG, TDBAMP, S I NORG, CSOLR

C

REAL A

LOGICAL L

COMMON / STATE1 / AC 1 00) , L( 1 00)

Loca l var tables

GO - logical variable; if true, indicates that EMODEL should be

executed

NOINIT - logical variable; if true, Indicates the the simulation

parameters have not been initialized.

TIMEST - The current value to be used for the major simulation time step

LOGICAL NOINIT, GO

INTEGER TIMEST

DATA NOINIT/. TRUE./

77
696

697

696

699

700

701

702

703

704

705

706

707

706

709

710

71

712

713

714

715

716

717

718

719

720

721

722

723

724

725

726

727

728

729

730

731

732

733

734

735

736

737

738

739

740

741

742

743

744

745

746

747

748

749

750

751

752

753

ooo

ooo

ooooo

ooooo

ooo

ooo

ooo Initialize all non - t 1 me-dependent parameters

I F ( NO I N I T ) THEN

CALL C0NPAR

CALL PARAMA

N0INIT * .FALSE.

CALL WEATHR

CALL ECTROL ( T I MEST , GO

ELSE

Determine if simulation routine is to be executed

CALL ECTROL (TI MEST, GO)

IF( .NOT. (GO) ) RETURN

Determine current weather parameters for this time step.

CALL WEATHR

Simulate the air handling unit equipment given return air temperature,

outdoor air conditions and the state of control signals. This is an algebraic/logical model assummlng a >ery small time constant.

CALL AHU

Simulate a building zone given the supply air conditions and

environmental conditions. This requires solution of differential

equat 1 ons

CALL ZONE! T I MEST)

Determine the comfort level in the building zone

CALL COMFRT ( T I MEST

Compile the data generated for later evaluation

END IF CALL CMP I LE( TI MEST)

RETURN

END

Cii>iiiii3iiiiiiiiiia33taii*ia*iti««tr33sss33t33:isxiisii«i«:333s=::::::

C’c

c

c

c

cc

c

c

c

c

cc

SUBROUTINE PARAMA

BaBBaBaaaaaBaBBBBBBSsaaaaaBaBBaaBBaaaxaaasasaassssBBaBBBBBassaasssssss

This routine is used to initialize the emulator building simulation

parameters. It is called once at the beginning of a one season test.

I nput var tables: none Output variables: none Procedures called: none

Common blocks: TIME, CL I MAT, HVAC, HVAC2, PFI LES, SHELL, STATE 1

, UNI TS, USE

Description of

TIME - contai

CL l MAT - contai

HVAC - contai

Common b l ocks

ns the current t i-me ns parameters which

ns parameters which

vector for the simulation,

describe climate around the

describe the HVAC system of

bu i l d i ng

the bu I Id i ng

73
754

755

756

757

756

759

760

761

762

763

764

765

766

767

768

769

770

771

772

773

774

775

776

777

778

779

780

781

782

783

784

785

786

787

788

789

790

791

792

793

794

795

796

797

798

799

800

801

802

803

804

805

806

807

808

809

810

81 1

C

C

C

C

C

C

C

C

C

C

C

C

C

cc

cc

c

c

cc

c

cc

c

cc

HVAC2

PFILES

SHELL

STATE1

UNITS

USE

contains character data describing The HVAC system of the building,

contains the character names of the files containing parameters,

contains parameters which describe tine building shell,

contains the state variables for the emulator model,

contains the logical unit assignments for program read and write,

contains parameters which describe building use schedules.

Common B 1 ock Var tables:

DAY

HOUR

MINUTE

SECOND

ENOOAY

ENDHR

ENDMIN

ENDSEC

Current day number , f i rst day of test i s day 1

Current hour of day, from 0-23.

Current minute within the hour, 0-59.

Current second within the minute, 0-59.

Day number on which the test is to stop.

Hour of the day on ENDDAY within- which test is to stop.

Minute within ENDHR within which test is to stop.

Second within ENDMIN within which test is to stop.

I NTEGER DAY , HOUR , M I NUTE , SECOND , ENDOAY , ENDHR . ENOM I N , ENDSEC

COMMON / TIME / DAY, HOUR, MI NUTE, SECOND, ENDDAY, ENDHR, ENOMI N, ENDSEC

WOAVG

CWIND

TDBAVG

TDBAMP

SINORG

CSOLR

Average outdoor humidity ratio in lb water/lb dry air.

Average windspeed, MPH

Average outdoor dry bulb temperature, F

Dirurnal amplitude of outdoor dry bulb temperature, F

Origin of sine wave describing dry bulb temperature variation in

seconds after midnight.

Peak, solar heat gain on horizontal.

REAL WOAVG , CW I ND , TDBAVG , TDBAMP , S I NORG , CSOLR

COMMON / CL I MAT / WOAVG, CW I ND, TDBAVG, TDBAMP, SI NORG, CSOLR

C

C DUDECK

C HUMCON

C RAFAN

C VAV

C VOLSAF

C CAPLOC

C PGA I NL

C TCTHSE

C QRAFAN

C FVAMIN

C QSAFAN

C HCAP1

C CCAP1

C DTSE7B

C POWSAF

C POWRAF

C GAINIL

C LOCEQT

C

logical variable which is true for a dual deck system,

logical variable which is true for a humidity control system,

logical variable which is true for a Return Air Fan in system,

logical variable which is true for a Variable Air Volume System.

Air handling unit fan rated supply air volume in ft* 3/sec.

Local equipment capacity for heating in Btu/sec.

Local equipment controller proportional gain in Btu/sec. F Local equipment controller sensor time constant in sec. Return Air Fan air heating rate in Btu/sec.

Minimum ventilation air, fraction of supply air mass Supply Air Fan air heating rate in Btu/sec.

Cold Deck Heating Coil Capacity in Btu/sec.

Cold Deck Cooling Coil Capacity in Btu/sec.

Local equipment controller setpoint setback in F power requirement of supply fan at rated volume in kW

power requirement of return fan at rated volume in kW

Local equipment controller integral gain in Btu/s2 F character description of Local equipment type

LOGICAL

COMMON / HVAC /

&

&

CHARACTER* 10

COMMON / HVAC2/

DUDECK , HUMCON , RAFAN , VAV

DUDECK , HUMCON , RAFAN , VAV , VOLSAF , CAPLOC , PGA I NL

TCTHSE, QRAFAN, FVAM I N , QSAFAN , HCAP1 , CCAP1

DTSETB

POWSAF , POWRAF , GAINIL

LOCEQT

LOCEQT

C

79
812

813

814

815

816

817

818

819

820

821

822

823

824

825

826

827

828

829

830

831

832

833

834

835

836

837

838

839

840

841

842

843

844

845

846

847

848

849

850

851

852

853

854

855

856

857

858

859

860

861

862

863

864

865

866

867

868

869

C

C

C

C

Cc

c

c

c

c

cc

c

c

c

c

cc

c

c

cc

c

c

c

cc

cc

cc

c

cc

c

c

cc

c

c

c

c

cc

c

c

cc

c

c

cc

HVFILE - Character filename of HVAC equipment parameter file.

USFILE - Character filename of Building usage parameter file.

CLFILE * Character filename of climate description parameter file.

SHFILE - Character filename of building shell parameter file.

INFILE * Character filename of state vector initialization file.

SEASON - Character name of weather season to use in test.

CHARACTER* 1 5 HVF I LE , USF I LE . CLF I LE , SHF I LE , I NF I LE , SEASON

COMMON / PFILES / HVFI LE, USF! LE, CLF1 LE, SHFI LE, I NFI LE, SEASON

AREAW

HI WALL

HOWALL

LENOW

LENIW

KOW

KIW

CFUR

VOLAIR

CPOW

RHOOW

CPIW

RHOIW

AIREXC

RMGFUR

Area of the building wall exposed to the outside in ft *2.

Inner wall film coefficient in Btu/s F ft *2.

Outer wall film coefficient in Btu/s F ft *2.

thickness of wall outside the mass center in ft.

thickness of wall inside the mass center in ft.

wall conductivity outside the mess center in Btu/s F ft.

wall conductivity Inside the mass center in Btu/s F ft.

thermal capacity of zone furnishings in Btu/F.

volume of air in zone in ft°3.

specific heat of wall outside the mass center in Btu/lb F. density of wall outside the mass center in lbm/ft*3.

specific heat of wall inside the mass center in Btu/lb F. density of wall inside the mass center in lbm/ft*3.

base air infiltration rate without wind and thermal in ACH

furnishings moisture absorption mass relative to air

REAL AREAW , H I WALL , HOWALL , LENOW , LEN I W , KOW , K I W , CFUR

& VOLA I R 4 CPOW , RHOOW .CPIW, RHO I W , A I REXC , RMGFUR

COMMON / SHELL / AREAW, HI WALL, HOWALL, LENOW, LEN I W, KOW, K I W, CFUR,

& VOLAIR, CPOW. RHOOW, CPIW, RHOIW, A 1 REXC, RMGFUR

contents of state vector

A( 1 ) : TEMPERATURE of outside air (F)

AC 2 ) : HUM t D I TY RATIO of outside air (lbw/lbda)

AC 3) : ENTHALPY of outside air CBtu/lb)

AC 4 ) : Windspeed CMPH)

AC 6): TEMPERATURE of zone air CF)

AC7): HUMIDITY RATIO of zone air C lbw/lbda)

AC 8): TEMPERATURE of zone wall interior CF1

AC 9): TEMPERATURE of Zone wall surface CF)

AC 10): TEMPERATURE of zone glass CF)

AC 20): TEMPERATURE of return air CF)

AC 21 ) : ENTHALPY of return air CBTU/lb)

AC 22) : HUMIDITY RATIO of return air C lbw/lbda)

AC 23) : TEMPERATURE of mixed air CF)

AC 24 ) : ENTHALPY of mixed air CBTU/lb)

AC 25) : HUMIDITY RATIO of mixed air C lbw/lbda)

AC26): TEMPERATURE of cooling coil discharge air CF) cold

AC27): ENTHALPY of cooling coil discharge air CBTU/lb) cold

AC 28): HUMIDITY RATIO of cooling coll discharge air C lbw/lbda) cold

AC 29): TEMPERATURE of heating coil discharge air CF) cold

AC 30): ENTHALPY of heating coll discharge air CBTU/lb) cold

AC31): HUMIDITY RATIO of heating coil discharge air C lbw/lbda) cold

AC 32): TEMPERATURE of supply air CF) cold

AC 33): ENTHALPY of supply air CBTU/lb) cold

AC 34): HUMIDITY RATIO of supply air C lbw/lbda) cold

deck

deck

deck

deck

deck

deck

deck

deck

deck

30
670

871

872

873

874

873

876

877

878

879

680

881

882

883

884

883

886

887

888

889

890

891

892

893

894

893

896

897

898

899

900

901

902

903

904

903

906

907

908

909

910

91 1 912

913

914

913

916

917

918

919

920

921

922

923

924

923

926

927

C

C

C

C

C

C

Cc

c

c

c

c

cc

c

c

c

c

cc

c

cc

cc

c

c

c

cc

c

c

c

c

N#c

cc

cc

c

cc

c

c

cc

c

c

cc

cc

AOS): ENTHALPY of return air plus ventilation air(BTU/lb)

A ( 3 1 ) : MASS FLOW RATE of return air Mb/s)

A<32): MASS FLOW RATE of supply air Mb/s)

A(33): MASS FLOW RATE of ventilation air Mb/s)

cold deck

A ( 60

A ( 61 ) A ( 62

A ( 63 ) A (64)

A ( 63)

A( 66)

A( 67)

AOS)

A( 69)

A( 70)

A( 71 ) A( 72)

AC 73)

AC 74

HEATING POWER to heating coil CBtu)

COOLING POWER from cooling coil (Btu)

HEATING POWER to local zone equipment CBtu/s)

COOLING POWER from local znne equipment CBtu/s)

HEATING ENERGY ( requ 1 rements ) CBtu)

COOLING ENERGY C requ i rements ) CBtu)

HEATING ENERGY (economizer) CBtu)

COOLING ENERGY (economizer) CBtu)

HEATING ENERGY (load) (Btu)

COOLING ENERGY (load) (Btu)

HEATING ENERGY (reheat) CBtu)

COOLING ENERGY C recool) (Btu)

ELECTRICAL ENERGY to fans (kWh)

HEATING ENERGY to local zone equipment (Btu)

COOLING ENERGY from local zone equipment (Btu)

cold deck

cold deck

AC 80)

AC 81 ) : AC 82) : AC 83) :

CONTROL POINT ADJUSTMENT supply air temperature (F)

CONTROL POINT ADJUSTMENT zone temperature (F)

CONTROL POINT ADJUSTMENT outside air damper position (F)

CONTROL POINT ADJUSTMENT supply air RH (X)

AC 90)

AC 91 ) AC 92)

AC 93)

A (94)

COMFORT

COMFORT

COMFORT

COMFORT

COMFORT

seconds outside of comfort range

maximum dry bulb temperature (F)

minimum dry bulb temperature (F)

relative humidity at maximum temperature

relative himidity at minimum temperature

(X)

(X)

LC 01 ) : LC02)

LC 03) : LC 04 )

COMMAND air handling unit on/ off

COMMAND economizer on/off

COMMAND ventilation on/off

COMMAND zone setback on/off

LC10): STATUS air handling unit on/off

L ( 1 1 ) : STATUS zone occupied/ unoccupied

REAL A

LOGICAL L

COMMON / STATE1 / AC 1 00) , LC 1 00)

CLU

FLU

ILU

PLU

Log leal Unit

Log i ca l Unit

Logical Unit

Log 1 ca l Unit

for writes to

for read from

for read from

for read from

console.

EMCS command

console,

parameter f

1 nput file. es

1 NTEGER CLU . FLU .ILU, PLU

COMMON /UNITS / CLU, FLU, I LU, PLU

LSCHED

NSCHED

TSCHED

SENS I

Index of current building usage state description in arrays.

Number of possible states in schedule arrays.

Start time of LSCHED 'th building usage state in hours past

m 1 dn i ght

Equipment sensible heat gain for LSCHED'th Dulldlng usage stato

01
928

929

930

931

932

933

934

933

936

937

936

939

940

941

942

943

944

943

946

947

948

949

950

931

952

933

934

953

936

937

958

939

960

961

962

963

964

963

966

967

968

969

970

971

972

973

974

973

976

977

978

979

980

981

982

983

984

983

oooooooooooooooooooooooooooooooooooo

oooooooo

In Btu/s.

LATENT - Equipment moisture gain for LSCHED'th building usage state

In Ibm/s.

PESENS - Sensible heat gain from people for LSCHED'th building usage state

In Btu/s.

PELATE - Moisture gain from people for LSCHED'th building usage state

In Ibm/s

PARAMETER ( NSCHED = 10)

I NTEGER LSCHED

REAL TSCHED , SENS I B , PESENS , LATENT, PELATE

COMMON / USE / LSCHED , TSCHED ( NSCHED) , SENSIB( NSCHED), LATENT (NSCHED)

& , PESENS ( NSCHED )

, PELATE ( NSCHED

Local

AINIT

BTUSPW

BULB

CHMAX

CTMAX

CTM1N

EQUIPL

EQUIPS

INDEX

LIGHT

LINIT

NACCAT

NACTI

OCCLAT

OCCSEN

RE

RL

RP

SEAS IN TOLAT

TOSEN

Var tables:

- Temporary value to hold state vector initializer read from

Initialization file.

- Conversion factor, Btu per Watt.

- Indicates type of lighting in describing building schedule.

- Time of maximum outdoor dry bulb temperature, hrs from midnight.

- Maximum daily outdoor temperature for a particular season in F.

- Minimum daily outdoor temperature for a particular season in F.

* Equipment moisture gain for a schedule period in lb water /hr.

- Electric power used by any equipment in zone for a schedule period

in Watts.

- Temporary value to hold state vector index read from initialization

file. 1 - Lighting power in the zone for a schedule period In watts.

- Temporary variable to hold logical state vector initializer read

from initialization file.

- Number of possible activity categories for occupants of zone.

- Array containing number of occuoants in NACCAT activity categories

for a schedule period.

- Occupant moisture gain per person for a schedule period in

lb water/hour.

- Occupant sensible heat gain per person for a schedule period in

Btu/hr

- ratio of radiated energy to total equipment energy for a schedule

per 1 od

- ratio of radiated energy to total lighting energy for a given

l i ght I ng type

- ratio of radiated energy tc total sensible heat of people.

- temporary variable to hold season name from climate file for

comparison with desired season in variable SEASON.

- total moisture gain from occupants for a schedule period in lb.

water /hour

- total sensible heat gain of occupants for a schedule period In Btu/hour

REAL EQU I PL , EQU I PS , RE , RL . RP , TOLAT , TOSEN

PARAMETER (NACCAT = 3, BTUSPW = 9. 4781 7E -4, RP = 0.4)

LOGICAL LINIT

INTEGER NACTI V( NACCAT) . INDEX

REAL AINIT, CHMAX , CTMAX , CTM I N , OCCSEN ( NACCAT )

, OCCLAT ( NACCAT » , L I GHT

CHARACTER*! BULB

CHARACTER* 5 SEAS IN C

DATA OCCSEN / 230.. 233., 313., 345., 343. /

82
966

987

968

989

990

991

992

993

994

999

996

997

998

999

1000

1001

1002

1003

1004

1009

1006

1007

1008

1009

1010

101 1 1012

1013

1014

1019

1016

1017

1018

1019

1020

1021

1022

1023

1024

1029

1026

1027

1028

1029

1030

1031

1032

1033

1034

1039

1036

1037

1038

1039

1040

1041

1042

1043

ooo

ooo

ooo

o

o

o

o

oo

o

o

o

o

o

o

DATA 0CCLAT / 0 20, 0.27, 0.34, 0.46, 0.74 /

Occupant sensible and mrlsture gains:

(from ASHRAE fundamentals, 1931, p 26.25 and Thretkeld, page 363)

sensible BTU/hr moisture LB/hr

1. seated, very light work 230. 0. 20

2. seated, light work, typing 255. 0. 27

3. standing, light work 315. 0.34

4 . 1 1 ght bench work 345. 0.46

5. walking, 3 mph, light machine work 345. 0. 74

Initialize State Variables

0PEN( PLU, FI LE= I NFI LE,STATUS= 'OLD 1 , ERR=994

DO 250 I = 1,100

READ ( PLU, *,END=280,ERR=260) INDEX, AINIT

A( INDEX) > AINIT

250 CONTINUE

260 DO 270 J a 1,100

REA0( PLU, *,END=280,ERR=993) INDEX, LI NIT

L(INDEX) a LINIT

270 CONTINUE

280 CLOSE (PLU)

Read Building Use parameters

OPEN ( PLU, FI LEaUSFI LE, STATUS= 'OLD ' , ERR=998)

LSCHED a o READ ( PLU , * , END* 500

00 400 I a 1.NSCHE0

TSCHED in hours, LIGHT and EQUIPS in W, EQUIPL in lb water/hr

READ ( PLU , * , END= 500 ) TSCHED ( I ) , LIGHT, BULB , EQU I PS, RE, EQUI PL,

& ( NACTI V(K )

, K= 1, NACCAT)

Determine lighting type, I = Incandescent, F = fluorescent

C

300

I F ( BULB . EQ . ’ 1' ) THEN

RL * 0.8

ELSE I F ( BULB . EQ . ' F

' ) THEN

RL a 0.5

ELSE

RL 0.5

ENOIF

LSCHED a LSCHED + 1 TOSEN a o.O

TOLAT = 0.0

DO 300 J = 1 , NACCAT

TOSEN = TOSEN + NACTIV(J) « OCCSEN(J)

TOLAT a TOLAT NACTI V(J) * OCCLAT(J)

CONTINUE

SENSIB(I) = ( (

1 -RL) «LIGHT (1 -RE) «EQUI PS)

LATENT (I ) = EQUIPL / 3600.

PESENS(I) = ( 1 - RP ) * TOSEN / 3600.

PELATE(I) = TOLAT / 3600.

« BTUSPW

33
1044

1045

1046

1047

1046

1049

1050

1051

1052

1053

1054

1055

1056

1057

1056

1059

1060

1061

1062

1063

1064

1065

1066

1067

1066

1069

1070

1071

1072

1073

1074

1075

1076

1077

1076

1079

1060

1081

1082

1083

108^ 1085

1086

1087

1088

1089

1090

1091

1092

1093

1094

1095

1096

1097

1098

1099

1 100

1 101

ooo

oooo

oooo

400 CONTINUE

Read Building She 1 1 Parameter

500 0PEN( PLU, F! LE=SHFI LE, STATUS* 1 OLD 1 , ERR=999)

READ( PLU, * ) AREAW

READ (PLU,*) HI WALL

READ (PLU,*) HCWALL

READ (PLU, * ) LENOW

READ ( PLU, * ) LENIW

READ ( PLU, * ) KOW

READ ( PLU, * ) KIW

READ (PLU,*) CFUR

READ ( PLU, * ) VCLAIR

READ (PLU,*) CPOW

READ (PLU,*) RHOOW

READ (PLU,*) CPIW

READ ( PLU , * ) RHOIW

READ (PLU,*) A1REXC

READ ( PLU, * ) RMGFUR

CLOSE (PLU)

Read HVAC system Parameters

OPEN ( PLU, FI LE=HVFILE, STATUS* •OLD' , ERR=997)

READ (PLU,*) VOLSAF

READ (PLU,*) POWSAF

READ (PLU,*) POWRAF

READ (PLU,*) LOCEOT

READ (PLU,*) CAPLOC

READ( PLU, * ) PGA I NL

READ ( PLU, * ) GAINIL

REACH PLU, * ) TCTHSE

READ (PLU,*) VAV

READ (PLU,*) RAFAN

READ (PLU,*) QRAFAN

READ( PLU, * ) FVAMIN

READ (PLU,*) DUDECK

READ( PLU, * ) QSAFAN

READ (PLU,*) HUMCON

READ (PLU,*) HCAP1

READ ( PLU, * ) CCAP1

READ (PLU,*) DTSETB

CLOSE (PLU)

Read C 1 1 mete Parameters

OPEN (PLU,

READ (PLU,

600 READ (PLU,

READ (PLU,

READ (PLU,

READ (PLU,

READ (PLU,

READ ( PLU,

READ (PLU,

FI LE=CLFI LE, STATUS* ’ OLD

*) SEASIN

* , END*995 ) SEASIN

* ) CTMAX

*) CTMIN

*) CHMAX

* ) WOAVG

* ) CSOLR

*) CWIND

ERR=996

84
1 102

1 103

1104

1 105

1106

1107

1106

1109

1110

1111

1112

1113

1114

1115

1116

1 117

1118

1119

1120

1121

1122

1123

1124

1125

1126

1127

1126

1129

1130

1131

1132

1133

1 134

1135

1136

1137

1138

1139

1140

1 141

1 142

1143

1144

1145

1 146

1147

1148

1149

1150

1151

1152

1153

1154

1155

1156

1157

1158

1159

C

c

c

I F ( SEAS I N ( 1 : 5 ) . NE . SEASON (1 : 5 ) ) GO TO 600

CLOSE (PLU)

TDBAVG C CTMAX CTMIN) /2.

TDBAMP * CTMAX - TDBAVG

S1NORG (CHMAX-30. )*3600.

RETURN

ERRORS

993 WRITECCLU, 7) l+J

7 FORMAT (IX,

STOP

' unable to

994 WRITECCLU. 6)

6 FORMAT (IX,

STOP

' unable to

995 WRITECCLU, 5)

5 FORMAT (IX,

STOP

' unable to

996 WRITECCLU, 4)

4 FORMAT C IX,

STOP

' unable to

997 WRITECCLU, 3)

3 FORMAT ( IX,

STOP

* unab l e to

998 WRITECCLU, 2)

2 FORMATC IX,

STOP

* unab lto

999 WRITECCLU, 1 ) 1 FORMATC IX,

STOP

END

'unable to

read INITIALIZATION file at line’, 16)

open INITIALIZATION file’)

locate requested season in CLIMATE file')

open CLIMATE parameter file')

open HVAC parameter file')

open USE parameter file')

open SHELL parameter file')

Iiii3s3s3iiiiiiimis33i8si:iis3itm3stsi3ssss:33sii3sisirs::33:::::::

SUBROUTINE WEATHR

C This routine is used to determine the current weather conditions. Wind

C and humidity ratio are considered constants. Outdoor temperature is

C assummed to vary sinusoidally with a period of 24 hours.

C

C AS of 2-15-85 there is no solar energy section.

C

C I nput var tables: none

C Output var 1 ab l es : none

C Procedures called: SINUSO

C Common blocks: STATE1 , CLIMAT, TIME

C

REAL A

LOGICAL L

COMMON / STATE1 / AC 1 00) , L( 1 00)

C

INTEGER DAY, HOUR, MINUTE, SECOND , ENDDAY , ENDHR , ENDM I N , ENDSEC

COMMON /TIME / DAY, HOUR, MI NUTE, SECOND, ENDDAY , ENDHR , ENOM I N , ENDSEC

C

COMMON / CLIMAT / WOAVG, CWI ND, TDBAVG, TDBAMP, S l N0R3, CSCLR

C

C Local Variables:

C

C PERIOD - Period of . s I nuso i da l temperature variation in seconds

C RTIME - Current time in seconds past midnight.

C TOUTA - Outside Air dry bulb temperature in F

35
1160

1161

1162

1163

1164

1169

1166

1167

1168

1169

1170

1171

1172

1173

1174

1173

1176

1177

1178

1179

1180

1181

1182

1183

1184

1189

1186

1187

1188

1189

1190

1191

1 192

1 193

1 194

1 193

1 196

1 197

1 198

1 199

1200

1201

1202

1203

1204

1203

1206

1207

1208

1209

1210

121 1 1212

1213

1214

1213

1216

1217

C VWIND - Wind velocity in MPH.

C WOUTA - Outside air humidity ratio in lb water / lb dry air

C

REAL PERIOD, RT I ME

PARAMETER ( PERI 00=86400 .

EQUIVALENCE ( A(1),T0U7A ),( A(2), WOUTA ) , ( A( 04 > , WI ND)

C

RTIME » HOUR* 3600 MINUTE*60 + SECOND

CALL S I NUSO ( TDBAVG , PER I OD , TDBAMP , RT I ME , S I NORG , TOUTA

WOUTA = WOAVG

WIND « CWIND

RETURN

END

SUBROUT I NE S I NUSO < TDBAVG , PER I OD , AMP , RT I ME , ROR I G , T

C3 31333 3ISI333333333I333 33 3333B3333333333 33.3 3333333333333333 = 3333333 = 33 = : C Returns a point on an approximate sinusoidal curve as a function of

C time. The curve shape is based on table data in SHAPE.

C ! nput var tables:

C TDBAVG - the average value or baseline of the curve.

C PERIOD - contains the period of the wave in seconds.

C AMP - the amplitude of the curve with reference to the baseline.

C RTIME - the current time in seconds since midnight.

C RORIG - the origin of the wave in seconds since midnite.

C Output variables:

C T - the current value of the wave at time RTIME.

C

C Procedures called: none

C Common blocks: none

C

C Local Variables:

C

C ICURP2

C ICURPE

C NPER

C NPIECE

C RANGE

C

C RCURPE

C RNPPS

C SHAPE

C

C

Index into SHAPE table for smallest time larger than current time.

Index into SHAPE table for largest time smaller than current time.

Number of whole curve periods passed by current time.

Number of entries In curve look-up table.

Difference between two table entries for times on either side

of current time.

Real number of curve pieces past origin for current time.

Real number of curve pieces per second.

Curve look-up table containing normalized curve values for one period of oscillation equally divided into NPIECE parts.

PARAMETER (NPIECE = 24)

REAL SHAPE (NPIECE)

DATA SHAPE/ +0. 0000, +0.2388, +0.3000, +0. 7071 , +0.8660, +0.9639,

* +1 . 0000, +0.9639, +0. 0660, +0.7071 , +0.3000, +0.2388,

& -0.0000, -0.2388, -0.3000, -0.7071 , -0. 8660, -0. 9659,

& -1.0000,-0.9639,-0.8660,-0.7071,-0 3000,-0.2508/

RNPPS = FLOAT (NPIECE) / PERIOD

RCURPE = RNPPS * (RTIME - RORIG)

ICURPE = A I NT (RCURPE)

NPER = ICURPE / NPIECE

ICURPE = ICURPE - NPER * (NPIECE)

RCURPE = RCURPE - NPER * (NPIECE)

36
1216

1219

1220

1221

1222

1223

1224

1223

1226

1227

1228

1229

1230

1 231

1232

1233

1234

1233

1236

1237

1236

1239

1240

1241

1242

1243

1244

1243

1246

1247

1248

1249

1250

1231

1252

1253

1254

1235

1236

1257

1250

1259

1260

1261

1262

1263

1264

1265

1266

1267

1268

1269

1270

1271

1272

1273

1274

1275

C

C

C

C

C

C

C

C

C

C

C

C

C

C

C

C

C

C

C

C

C

C

C

C

C

c

c

c

c

c

c

c

c

c

cc

cc

c

c

cc

ICURP2 = ICURPE 2

I F ( I CURP2 . 3T . NP I ECE ) I CURP2 = I CURP2 - NPIECE

RANGE » ( SHAPE ( I CURP2 ) - SHAPE ( I ClJRPE+1 ) ) T * ( ( RCURPE - ICURPE) * RANGE SHAPE ( I CURPE+ 1 ) ) * AMP + TDBAVG

RETURN

END

SUBROUTINE AHU

>3iiitii3isi333iiiiiiiiiiiiSBitiiii9f:iii>:3i«i:i>3i3iiti3i::33:s:::::3

This routine contains a steady state model of a generalized air handling

unit

Input variables: none Output variables: none Procedures called: CPAI R, DEWPT, ENTHAL, HUMRAT, RHOAI R, TDBSAT. TEMP

Common b lodes: HVAC, HVAC2, STATE

LOGICAL

COMMON / HVAC /

&

&

CHARACTER* 1

COMMON / HVAC2/

DUDECK, HIJMCON, RAFAN, VAV

DUDECK, HUMCON, RAFAN, VAV, VOLSAF, CAPLOC, PGA I NL,

TCTHSE, QRAFAN, FVAMI N, QSAFAN, HCAP1 , CCAP1 , DTSETB,

POWSAF , POWRAF, GA I N I

LOCEQT

LOCEQT

REAL A

LOGICAL L

COMMON / STATE1 / A( 1 00) , LC 1 00)

Local Variables:

AHSTAT

DELTA * DTRF

DTSF

ECONO

HCA1

HHA1

HMA

HMC

HMO

HOUTA

HRA

HSA1

MRA

MSA1

PRESS

QC00L1

QHEAT1

TCA1

TDPMA

THAI

TMA

TMC

TMO

TOUTA

TRA

TSAI

TSET

TZA

STATUS air handling unit on/off

Difference between supply air setpoint and fan entering

Temperature differential across the return fan.

Temperature differential across the supply fan.

COMMAND economizer on/off

ENTHALPY of cooling coil discharge air (BTU/lb)

ENTHALPY of heating coil discharge air (BTU/lb)

ENTHALPY of mixed air (BTU/lb)

ENTHALPY of return air plus ventilation air (BTU/lb)

Enthalpy of mixed air If outside air dampers are open

ENTHALPY of outside air (Btu/lb)

ENTHALPY of return air (BTU/lb)

ENTHALPY of supply air (BTU/lb)

MASS FLOW RATE of return air ( Ib/s)

MASS FLOW RATE of supply air ( Ib/s)

Atmospheric pressure (assummed constant) in inches

COOLING POWER from cooling coll (Btu)

HEATING POWER to heating coil (Btu)

TEMPERATURE of cooling coil discharge air (F)

Dewpoint temperature of mixed air (F)

TEMPERATURE of heating coil discharge air (F)

TEMPERATURE of mixed air (F)

Temperature of return air plus ventilation air (F)

Temperature of mixed air if outside air dampers are TEMPERATURE of outside air (F)

TEMPERATURE of return air (F)

TEMPERATURE of supply air (F)

CONTROL POINT ADJUSTMENT supply air temperature (F)

TEMPERATURE of zone air (F)

lir T.

col

CO 1 d deck,d deck

(Btu/lb)

cold deck

cold deck

of Hg

cold deck

cold deck

cold deck

cold deck

i open (F)

co 1 d deck

37
1276

1277

1278

1279

1280

1281

1282

1283

1284

1289

1286

1287

1288

1289

1290

1291

1292

1293

1294

1293

1296

1297

1298

1299

1300

1301

1302

1303

1304

1309

1306

1307

1308

1309

1310

131 1 1312

1313

1314

1313

1316

1317

1318

1319

1320

1321

1322

1323

1324

1323

1326

1327

1328

1329

1 330

1331

1332

1 333

OOOO

O

O

O

O

OOOO

C UN I TON

C VENOFF

C WCA1

C WHA1

C WMA

C WMC

C WMO

C WOUTA

C WRA

C WSA1

C WZA

C

COMMAND air handling unit on/off

COMMAND ventilation on/off

HUMIDITY

HUMIDITY

HUMIDITY

Humidity

Hum 1 d 1 ty

HUMIDITY

HUMIDITY

HUMIDITY

HUMIDITY

RATIO of

RATIO of

RATIO of

ratio of

rat I o of

RATIO of

RATIO of

RATIO of

RATIO of

cooling coil discharge air (Ibw/lbda) cold deck,

heating coll discharge air (Ibw/lbda) cold deck

mixed air (Ibw/lbda)

return air plus ventilation air (Ibw/lbda)

mixed air If outside air dampers are open (Ibw/lbda

outside air (Ibw/lbda)

return air (Ibw/lbda)

supply air (Ibw/lbda) cold deck zone air (Ibw/lbda)

C

C

REAL MRA, MSA

EQUIVALENCE ( A( 01 ), TOUTA ( A (02) .WOUTA ). ( A( 20) , TRA ) EQUIVALENCE ( A ( 22 )

, WRA ( A( 80) . TSET ), ( A ( 24 ) . HMA ) EQUIVALENCE ( A( 06) , TZA >, ( A(07) , WZA >, ( A( 33) . HSA1 ) EQUIVALENCE ( A( 32) , TSAI ), ( A( 34J , WSA1 ( A( 52) , MSA1 ) EQUIVALENCE ( A ( 21 )

, HRA ). ( A( 51 ) .MRA )

( A (25) .WMA ) EQUIVALENCE ( A ( 26 )

, TCA1 >, ( A (28) . WCA1 ). ( A( 27) , HCA1 ) EQUIVALENCE ( A( 29) , THAI ), ( A( 31 )

. WHA1 >. ( A ( 30 )

. HHA1 ) EQUIVALENCE ( A ( 03 )

, HOUTA ). ( A ( 23 ), TMA )

( AOS) . HMC ) EQUIVALENCE ( A( 60 ) , QHEAT1 >, ( A (61 ) .QC00L1 ) LOG 1 CAL UN 1 TON , VENOFF . ECONO , AHSTAT

EQUIVALENCE ( L ( 01 ), UNI TON >. ( L( 03) .VENOFF ). ( L( 02) . ECONO ) EQUIVALENCE ( L( 10) .AHSTAT )

i PARAMETER (PRESS * 29.921)

Air handling unit is operating - determine air flow rates

I F( UNI TON) THEN

AHSTAT = . TRUE

HOUTA = ENTHAL ( TOUTA , WOUTA, PRESS)

IF(VAV) THEN

STOP 'NO VAV IMPLEMENTED"

ELSE

MSA1 = VOLSAF * RHOA I R ( WSA 1, TSET

MRA * MSA 1 END IF

Determine heating of air by return fan

IF(RAFAN) THEN

DTRF * QRAFAN / (MRA CPA I R( WZA) ) TRA » TZA DTRF

ELSE

TRA a TZA

END IF WRA * WZA

HRA * ENTHAL (TRA, WRA, PRESS)

88
1334 C

-

1 333 C Determine mixed air conditions

1336 C

1337 C state of air for outside air dampers closed

1338 C minimum outside air dampers closed----

1339 IF(VENOFF) THEN

1340 HMC a HRA

1341 WMC = WRA

1342 TMC * TRA

1343 C minimum outside air dampers open

1 344 ELSE

1343 HMC = (FVAMIN * HOUTA) + ( 1 - FV AM IN)* HRA

1346 WMC » (FVAMIN « WOUTA) ( 1 -FVAMl N) *WRA

1347 TMC * TEMP ( HMC , WMC , PRESS

1346 ENDIF

1349 C determine mixed air temperature----

1330 IF( .NOT. DUDECK) THEN

1331 C single deck.

1332 C supply fan temperature rise

1353 DTSF = QSAFAN / (MSA1 * CPA I R( WRA))

1334 DELTA * TSET - DTSF

1333 I F ( . NOT . HUMCON ) THEN

1336 C economizer enab led-determ i ne damper position

1357 I F( ECONO) THEN

1358 HMO = HOUTA

1339 TMO * TOUTA

1360 WMO » WOUTA

1361 C ----------outside air damper should be closed

•1362 IFCTMC.LE. DELTA) THEN

1363 TMA * TMC

1364 WMA * WMC

1365 HMA » HMC

1366 C ---- outside air damper should be open

1367 ELSE I F(TMC.GT. DELTA. AND. TMO. GT. DELTA) THEN

1366 TMA = TMO

1369 WMA « WMO

1370 HMA = HMO

1371 C outside air damper should be part open

1372 ELSE I F( TMC. GT. DELTA. AND. TMO. LE. DELTA) THEN

1373 WMA = (WMC - WM0)/(TMC - TMO) *( DELTA -TMO ) WMO

1374 TMA a DELTA

1373 HMA a ENTHAL(TMA, WMA, PRESS)

1376 ENDIF

1 377 ELSE

1378 HMA HMC

1379 WMA « WMC

1360 TMA « TMC

1381 ENDIF

1382 C*************** ******************************** ********C

1383 C damper position control or report could be Installed C

1 384 C*******************************************************C

1385 C - - - 1386 C Determine heating or cooling required

1387 C

1388 C no heating or cooling required

1389 I F ( TMA . EQ . DELTA ) THEN

1390 TCA1 a TMA

1391 V/CA1 * WMA

89
1392

1393

1394

1395

1396

1397

1398

1399

1400

1401

1402

1403

1404

1405

1406

1407

1406

1409

1410

141 1 1412

1413

1414

1415

1416

1417

1418

1419

1420

1421

1422

1423

1424

1425

1426

1427

1428

1429

1430

1431

1432

1433

1434

1435

1436

1437

1438

1439

1440

1441

1442

1443

1444

1445

1446

1447

1448

1449

C

C

C

HCA1 a HMA

THAI * TMA

WHA1 WMA

HHA1 HMA

QHEAT1 * 0.0

QC00L1 a 0.0

---heating required

ELSE I F ( TMA . LT . DELTA ) THEN

THAI a DELTA

WHA1 a WMA

HHA1 a ENTHAL ( THA 1 ,WHA1 .PRESS)

QHEAT1 a MSA1 * ( HHA1 - HMA)

l F ( QHEAT 1 . GT . HCAP 1 ) THEN

HHA 1 a HCAP1/MSA1 +HMA

THAI a TEMP ( HHA 1 ,WHA1. PRESS)

QHEAT1 a HCAP1

END IF HCA1 a HHA1

TCA1 THAI

WCA1 a WHA1

QCOCL l = o .

cooling required

ELSE I F( TMA. GT. DELTA) THEN

THAI » TMA

WHA1 a WMA

HHA1 a HMA

QHEAT1 a o.O

TCA1 a DELTA

TDPMA a DEWPT( WMA, PRESS)

I F(TCA1 . QT. TDPMA) THEN

HCA1 » WCA1 a ELSE

WCA1 a HCA1 a END l F QC00L1 =

l F ( QCOCL

HCA1 a QCOCL. 1

I F ( HCA

TCA1

WCA1

ELSE

TCA1

WCA1

END IF END IF END IF

ENTHAL ( TCA 1

, WMA , PRESS

WMA

HUMRAT ( TCA1

ENTHAL ( TCA

100.0, PRESS)

WCA 1 , PRESS)

MSA 1 a ( HMA - HCA 1 ) .GT.CCAP1) THEN

HMA - CCAP1/MSA1

= CCAPI

. GT . ENTHAL ( TDPMA , WMA , PRESS )

a TEMP ( HCA 1

, WMA, PRESS)

a WMA

THEN

TDBSAT ( HCA1

HUMRAT ( TCA

PRESS)

100.0, PRESS)

determine supply air temperatures

WSA 1 a WCA

TSAI a TCA1 DTSF

HSA 1 a ENTHAL ( TSA 1

, WSA 1

, PRESS)

ELSE

STOP 'HUMIDITY CONTROL NOT IMPLEMENTED'

END IF

• ELSE

STOP 'DUAL DECK NOT IMPLEMENTED*

90
1430

1431

1432

1433

1434

1433

1436

1437

1438

1439

1460

1461

1462

1463

1464

1463

1466

1467

1468

1469

1470

1471

1472

1473

1474

1473

1476

1477

1478

1479

1480

1481

1482

1483

1484

1485

1486

1487

1488

1489

1490

1491

1492

1493

1494

1495

1496

1497

1498

1499

1500

1301

1302

1303

1304

1305

1306

1507

oooooooooooo

o

o

o

oo

o

o

oo

o

o

o

o

oo

oo

o

oooo

END IF

Air Handling Unit Off

ELSE

AHSTAT * .FALSE.

MSA1 * 0.

MRA MSA1

TSAI » TZA

TRA » TZA

TMA TZA

WSA1 a WZA

WRA > WZA

WMA a WZA

END IF

RETURN

END

SUBROUTINE LOCAL (TIMEST)

This routine contains a simple model of the local heating and/or

cooling equipment in the zone. It is assummed that the change in output

of the local equipment is proportional to the change in error between

the local controller setpoint and the room air temperature and also is

proportional to the Integral of the error ( error *t imestep) . Local

equipment may be of several types as configured with the variables in common blocks HVAC and HVAC2.

Input variables:

TIMEST - The current value to be used for the simulation time step.

Output variables: none Procedures called: ENTHAL

Common blocks: HVAC, HVAC2, STATE

LOGICAL

COMMON / HVAC /

&

&

CHARACTER* 1

COMMON / HVAC2/

DUDECK, HUMCON, RAFAN, VAV

DUDECK, HUMCON, RAFAN , VAV, VOLSAF, CAPLOC, PGA I NL

TCTHSE, QRAFAN, FVAMIN, QSAFAN, HCAP1 , CCAP1

DTSETB

POWSAF , POWRAF , GA I N I

LOCEQT

LOCEQT

REAL A

LOGICAL L

COMMON / STATE1 / A(100),L(100)

Local Variables:

DQAIR

DQENT

DQLOC

ELOCC

ELOCH

ERROR

Required change in heat contained in supply air entering zone to

reduce zone temperature error. (Btu/s)

Required change in total supplied heat entering zone to reduce

zone temperature error. Equals DQAIR DQLOC .( Btu/s

Required change In auxi lllary heat added to supply air by local zone equipment to reduce zone temperature error (Btu/s)

COOLING ENERGY from local zone equipment (Btu)

HEATING ENERGY to local zone equipment (Btu)

difference between local zone setpoint and zone temperature at

91
1508

1509

1510

151 1 1512

1513

1514

1515

1516

1517

1518

1519

1520

1521

1522

1523

1524

1525

1526

1527

1528

1529

1530

1531

1532

1533

1534

1535

1536

1537

1538

1539

1540

1541

1542

1543

1544

1545

1546

1547

1548

1549

1550

1551

1552

1553

1554

1555

1556

1557

1558

1559

1 560

1561

1562

1563

1564

1565

O

O

current time.

ERR0R0 - difference between local zone setpoint and zone temperature at

prev i ous t i mestep

EX - Transition ratio for local controller temperature sensor

HSA1 - ENTHALPY of supply air (BTU/lb) cold deck

HSA2 - Enthalpy of hot deck (if existing) supply air (Btu/lb)

HZA - Enthalpy of zone air (Btu/lb)

MRA - MASS FLOW RATE of return air (Ib/s)

MSA1 - MASS FLOW RATE of supply air (Ib/s) cold deck

MSA2 - Mass flow rate of hot deck (if existing) supply air (Ib/s)

PRESS * Atmospheric pressure (assummad constant) in inches of Hg.

QAIR - Heat energy entering zone space through supply air.

QL0CH - HEATING POWER to local zone equipment (Btu/s)

SETBAK - COMMANO zone setback on/off

TCTHSE - Time constant of temperature sensor in local controller (s).

TREF - Current local zone setpoint temperature (F).

TSETLO - CONTROL POINT ADJUSTMENT zone temperature (F)

TTHSE - Temperature measured by temperature senscr in local controller (s)

TZA - TEMPERATURE of zone air (F)

WZA - HUMIDITY RATIO of zone air (Ibw/lbda)

REAL MRA, MSA1 , MSA2

LOGICAL SETBAK

EQUIVALENCE

EQUIVALENCE

EQUIVALENCE

EQUIVALENCE

A( 81 ) A( 06)

TSETLO

TZA

A( 51 ) , MRA

A( 73) , ELOCH

) , ( L ( 04 )

, SETBAK ) ), ( A( 07) , WZA ) ),( A ( 62) , QLOCH )

) , ( A< 74 ) , ELOCC )

A ( 52

A ( 33)

MSA 1 HSA1

PARAMETER (PRESS * 29.921)

DATA TTHSE/75 . / , ERRGRO/O .

/

I F( SETBAK) THEN

TREF » TSETLO - DTSETB

ELSE

TREF a TSETLO

END IF EX = EXP( -T I MEST/ TCTHSE)

TTHSE « ( 1 -EX ) * TZA EX TTHSE

ERROR * TREF - TTHSE

DQENT a PGA I NL « ( ERROR -ERR0R0) GAINIL a ERROR * T I MEST

ERRORO a ERROR

HZA a ENTHAH TZA, WZA, PRESS)

THE FOLLOWING VALUES ARE FOR DUAL DECK (NOT IMPLEMENTED)

MSA2 0.0

HSA2 - 0.0

QAIR ( HSA1 * MSA1 ) ( HSA2 * MSA2 ) - (HZA a MRA)

--dual deck----------

I F ( LOCEQT . EQ. ’DUAL DECK 1 ) THEN

STOP ' DUAL DECK NOT IMPLEMENTED'

C--VAV

ELSE IF(VAV) THEN

DQAIR = DQENT

QAIR a QAIR DQAIR

MSA1 = QAIR / ( HSA1 - HZA)

PAUSE 'ZONE VAV BOX AIR FLOW LIMITS NOT IMPLEMENTED'

C

- -constant vo l ume reheat

92
1966

1S67

1368

1369

1370

1371

1372

1373

1374

1373

1376

1377

1378

1379

1380

1381

1382

1383

1384

1383

1386

1387

1388

1389

1390

1391

1392

1393

1394

1393

1596

1397

1398

1399

1600

1601

1602

1603

1604

1605

1606

1607

1608

1609

1610

161 1 1612

1613

1614

1613

1616

1617

1618

1619

1620

1621

1622

1623

onooooooooooooooooooo

o

ooo

oooo

ELSE IF(L0CEQT.EQ. ‘REHEAT* ) THEN

DQLOC » DQENT

QLOCH « QLOCH DQLOC

I F ( QLOCH . GT . CAPLOC ) QLOCH a CAPLOC

I F ( QLOCH . LE . 0 . ) QLOCH a 0.

ELOCH « ELOCH QLOCH * TIMEST

ELSE

STOP ‘LOCAL EQUIPMENT TYPE NOT IMPLEMENTED’

END IF

RETURN

END

SUBROUTINE ZONE (TIMEST)

CII3I33SI3I3IIII3333I33333I33II33333333333333333333I33S3S3II33333333333:

C ZONE solves for the values of the variables (in array Y) describing a

C building zone. The current values in Y are used as the starting point.

C New values of Y at one time step in the future are calculated and

C replace the old values of Y. C

C I nput var tables:

C TIMEST - The current value to be used for the major simulation time step.

C

INTEGER TIMEST

Output variables: none

Pr ocedures called: EULER , SURFAC , ZONE , ZONEFN , LOCAL

EXTERNAL ZONEFN

Common blocks: STATE1,TIME

REAL A

LOGICAL L

COMMON / STATE1 / A( 1 00) , L( 1 00)

INTEGER DAY, HOUR, MINUTE, SECOND, ENDDAY, ENDHR, ENDMIN, ENDSEC

COMMON / TIME / DAY, HOUR, MI NUTE, SECOND, ENDDAY , ENDHR, ENDM I N, ENDSEC

Local Variables:

ADSTEP

DELMAX

El SMI

El STEP

EUSTEP

FMAX

I TIME

Logical variable which if true, enables the use of the algorithm

to internally determine a minimum Euler timestep. If false, a fixed Euler timestep is used.

The maximum change in state variable allowed for any of the

differential equations describing the zone within one Euler

timestep in degrees F. Theoretically, each differential equation

should have its own DELMAX. A simplification of this is used.

Minimum Euler timestep which can be used.

Current Euler timestep to be used In solving the differential

equations, determined from dt = dYmax / CdY/dt), where dYmax Is DELMAX

Euler timestep to be used if ADSTEP is false.

Magnitude of the largest dY/dt. for all of the first order

differential equations used to describe

the zone.

A time value which is a combination of the current minutes and

seconds of the time of day. Used to detormlno if t I me - dependent

functions should be recalculated.

93
1624

1625

1626

1627

1628

1629

1630

1631

1632

1633

1634

1635

1636

1637

1638

1639

1640

1641

1642

1643

1644

1645

1646

1647

1648

1649

1650

1651

1652

1653

1654

1655

1656

1657

1658

1659

1660

1661

1662

1663

1664

1665

1666

1667

1668

1669

1670

1671

1672

1673

1674

1675

1676

1677

1678

1679

1 680

1681

ooo

ooo

oo

ooo

oooooooooooooooooo

NDEQNS - Number of variables in the zone model state vector which are defined by differential equations. These variables will be the

first NDEQNS variables in the vector. The remaining variables

will be defined by algebraic equations.

In the state vector

NITER - Number of iterations •equlred to solve the zone differential

equations for one major simulation timestep.

NSTVAR - Total number of variables in the zone model state vector.

STEP - Floating point equivalent of the integer major timestep In the

simulation.

Y - Zone model state vector which is a subset of the simulation

CF)

PARAMETER (ADSTEP = . TRUE. , EUSTEP = 0.1)

PARAMETER (NSTVAR = 5, NDEQNS = 3, DELMAX = 0.10, EISMIN = 0.001)

INTEGER I TIME

REAL Y( NSTVAR)

EQUIVALENCE ( A(6),Y(1) )

state vector. Elements are defined as : Y( 1 )

zone dry bulb temperature (F)

Y ( 2)

zone hum i d i ty rat i o Y(3) : zone 1 nternal wall temperature at mass

Y ( 4 ) : zone l nter 1 or wall temperature (F)

Y(5) : zone I nter I or glass temperature (F)

LOGICAL ADSTEP

-determine Euler -time stepI Ft ADSTEP) THEN

FMAX =0.0

DO 900 I * 1, NDEQNS

FMAX AMAX1 ( ABS< Z0NEFN( I , SECOND, NSTVAR, Y )

900

I F ( FMAX . LE

CONTINUE

1 . E-8) FMAX * 1 . E-8

STEP = El STEP

El STEP

El STEP

NITER

El STEP

ELSE

NITER = El STEP

END IF

1 oop

FLOAT ( T I MEST

= DELMAX / FMAX

= AMI N1 ( El STEP,

= AMAX1 ( El STEP,

= AINT( STEP / El STEP

= STEP / FLOAT (NITER)

STEP ) EISMIN

)

T I MEST ,

' EUSTEP

EUSTEP

through Internal IterationsFMAX )

I TIME * SECOND MINUTE«60

DO 1000 I * 1, NITER

solve NDEQNS simultaneous linear differential equations

CALL EULER (NDEQNS, I T1 ME. Y , E I STEP , ZONEFN, NSTVAR

solve algebraic equations for wall and glass surfaces

CALL SURFAC ( NDEQNS , NSTVAR . Y ,

I T 1 ME

C - C Simulate the local zone space heating and cooling equipment. This Is C an algebraic approximation.

94
1682 C

1683 CALL LOCAL ( E I STEP

1684 1000 CONTINUE

1685 C

1 686 RETURN

1 687 END

1 688 (Iiiiiit3iitiitiiiiis33i>i»3i3ts>:iiiin>:::3s:::>siiiiici:::33:::::::

1689 SUBROUTINE EULER l NDEQNS , TIME, Y. TIMEST, FNC, NSTVAR)

1 690 C> *** 313a<>laix|>ls * aa:axi:3sx * l,l:13x3i:ia::s::s::ic>3axi,>t:::a::: = 3::3:

1691 C This subroutine advances the system of NDEQNS ODE'S described by:

1692 C dY(l)/dt » FNCC i

, t, Y ( 1 ), Y< 2) .... Y ( NDEQNS) by one step of size

1693 C TIMEST using the explicit Euler method. Y contains the dependent

1694 C variables and is overwritten with the new values.

1695 C

1696 C Input variables:

1697 C NDEQNS - Number of differential equations to be advanced.

1698 C TIME - A time value which must be distinct from the time value on 1699 C the previous call to EULER.

1700 C Y - contains the current values of the dependent variables.

1701 C TIMEST - The Euler tlmestep size.

1702 C FNC - The address of the function describing dY/dt.

1703 C NSTVAR - dimension of the Y array. Not all elements may be advanced.

1704 C

1705 REAL Y( NSTVAR), TIMEST

1706 INTEGER TIME

1707 C

1708 C Output variables: none

1709 C Y - contains the new value of the dependent variables.

1710 C

1711 C Procedures called: FNC

1712 C

1713 EXTERNAL FNC

1714 C

1715 C Common b l ock.s : none

1716 C

1717 C- advance dependent variables

1718 C

1719 DO 900 J * 1 , NDEQNS

1720 Y ( J ) * Y ( J ) + TIMEST * FNC(J, TIME, NSTVAR, Y)

1721 900 CONTINUE

1722 C

1723 C-- the independent variable, time is advanced externally

1724 C

1 725 RETURN

1 726 END

1 727 Cxxaxxx3XxaaxaxxaxX3X3333B,x:,xaxxxX3xa33X3X3xaxx3X3XXXXXxxxa=3xxa3:x:: =3

1728 SUBROUTINE SURFACC NDEQNS, NSTVAR, Y, TIME)

1 729 Cxaxa33aa3xaaaxxaaxaa33xa3x3xxaaa3axaxxaxxx33xxx3xxaxaxxxaixx:xxlx 3ic:::

1730 C This subroutine is used to calculate the next value of the zone state

1731 C variables which can be calculated from an algebraic equation.

1 732 C

1733 C Input variables:

1734 C NDEQNS - Number of variables in the zone model state vector which are

1735 C defined by differential equations. These variables will be the

1736 C first NDEQNS variables in the vector. The remaining variables

1737 C will be defined by algebraic equations.

1738 C in the state vector

1739 C NSTVAR - Total number of variables in the zone model state vector

95
1740

1741

1742

1743

1744

1745

1746

1747

1746

1749

1750

1751

1752

1753

1754

1755

1756

1757

1758

1759

1760

1761

1762

1763

1764

1765

1766

1767

1768

1769

1770

1771

1772

1773

1774

1 775

1776

1777

1778

1779

1780

1781

1782

1783

1784

1785

1786

1787

1788

1789

1790

1791

1792

1793

1794

1793

1796

1797

o

o

o

ooo

o

oo

o

o

oo

oo

o

o

o

o

ooooooooo

ooooo

Y - contains the current values of the dependent variables.

TIME - A time value which is a combination of the current minutes and

seconds of the time of day. Used to determine if t i me -dependent

functions should be recalculated.

REAL Y(NSTVAR)

INTEGER TIME

Output variables:

Y - contains the new value of the dependent variables.

Y(4): zone Interior wall temperature (F)

Y(5): zone Interior glass temperature (F)

Procedures called: CNDCT I , CONVWA , GCNDCT , GCNVCT . GRAD l

Common blocks: STATE

REAL A

LOGICAL L

COMMON / STATE1 / A( 1 00) , L( 1 00)

Local

A41

A42

A43

A51

A52*

A53

TITOUTA

TRN

TW

var 1 ab l es

- conduction component of heat transfer between interior of wall and zone air.

- convection component of heat transfer between interior of wall and zone air.

- radiation component of heat transfer between interior of wall and zone air.

- conduction component of heat transfer between interior of wall and

outside air.

- convection component of heat transfer between Interior of wall and

outside air.

- radiation component of heat transfer between interior of wall and

outside air.

- zone dry bulb temperature (F)

outside air dry bulb temperature (F)

Mean radiant temperature in the zone (F)

- zone Internal wall temperature at mass center CF)

EQUIVALENCE ( Ad). TOUTA ) C

TI Y( 1

TW » Y( 3)

TRN * MRT (TIME)

C----------wal

interior surface temperature thermal energy balance

A41 CNDCT I ( TI ME)

A42 CONVWA C T I ME

A43 » WRADI A( TI ME)

Y ( 4

« ( A41 * TW + A42 « TI A43 « TRN WS0LAR( TIME) ) /

ft ( A4 1 A42 «• A43 ) C-- glass interior surface temperature thermal energy balance

A31 » GCNDCT ( TIME)

A32 = GCNVCT ( TIME)

A33 » GRAD I A( TI ME)

I F ( A31 . NE. 0. AND. A52. NE . 0 . AND . A33 . NE . 0

ftY ( 5

a ( AS 1 a TOUTA A52 « TI A53 « TRN ) /

8 ( A3 1 A32 A33 ) C

RETURN

96
1790

1799

1000

1001

1002

1003

1004

1005

1006

1007

1000

1009

1010

101 1 1012

1013

1014

1015

1016

1017

1010

1019

1020

1021

1022

1023

1024

1025

1026

1027

1926

1029

1030

1031

1032

1033

1034

1035

1036

1037

1030

1039

1940

1041

1942

1043

1044

1945

1046

1047

1940

1049

1050

1051

1052

1853

1854

1855

ooo

oooooooooooo

ooo

oooooooo

o

oo

o

o

oo

o

o

o

o

o

o

o

o

END

FUNCTION Z0NEFN( I, TIME, NSTVAR, Y)

133il33ll3l333l33ir313l33333133833333333333333333333333333333333333::

This function returns the right hand side of the linear ordinary

differential equation for the building zone model. This should be a linear function of time, the set of dependent variables, boundary

conditions, and parameters.

I nput var

ITIME

NSTVAR - Y

tables:

Index of differential equation for which function is being

determ i ned

A time value which is a comb i nat

seconds of the time of day. Used

functions should be recalculated

Total number of variables in the

Zone model state vector which is

state vector.

on of the current minutes and

to determine if time-dependent

zone model state vector,

a subset of the simulation

INTEGER TIME

REAL Y(NSTVAR)

Output variables:

ZONEFN - right hand side of linear ODE.

Procedures called: CNDCTI , CNDCTO, CONVGL, CONVWA, HUMRAT, I NFI LT, I NF I LW,

I NTERN, MOIST, MASSA, MASSTH, MASSW, SOLAR, SUPPLW,

SUPPLY

REAL CNDCTI , CNDCTO, CONVGL, CONVWA, HUMRAT, INFILT, INFILW, INTERN,

& MOIST, MASSA, MASSTH, MASSW, SOLAR, SUPPLW, SUPPLY

Common blocks: STATE

REAL A

LOG I CAL L

COMMON / STATE1 / A( 1 00) , L( 1 00)

Local

MAXW

PRESS

Tl

TIG

TIW

TOUTA

TW

WI

WOUTA

Var tables:

- Humidity ratio in the space at current temperature and 1 00X RH.

- Atmospheric pressure (assummed constant) in inches of Hg.

- zone dry bulb temperature (F)

- zone interior glass temperature (F)

- zone interior wall temperature (F)

- TEMPERATURE of outside air (F)

- zone internal wall temperature at mass center CF)

- zone hum i d i ty rat 1

- HUMIDITY RATIO of outside air (Ibw/lbda)

PARAMETER (PRESS = 29.921)

REAL MAXW

' EQUIVALENCE ( A(1), TOUTA ),( A(2), WOUTA ) Zone Air Energy Balance

I F( I . EQ. 1 ) THEN

TI * Y ( 1 ) WI * Y ( 2)

97
1856

1857

1858

1859

1860

1861

1862

1863

1864

1865

1866

1867

1868

1869

1870

1871

1872

1873

1874

1875

1876

1877

1878

1879

1880

1881

1882

1883

1884

1885

1886

1887

1888

1889

1890

1891

1892

1893

1894

1895

1896

1897

1898

1899

1900

1901

1902

1903

1904

1905

1906

1907

1908

1909

1910

191 1 1912

1913

oooo

oooooooooo

o

oo

ooo

ooo

ooo

ooo

&

&

&

&

TIW * YC4

TIG * Y ( 5)

ZONEFN * (SUPPLY (TIME)

I NFI LT( TI ME, TI , WI ) * (TOUTA - TI ) + CONVWA ( T I ME ) * (TIW - TI ) CONVGL ( T I ME ) * (TIG - TI ) + INTERN (TIME) )/ MASSTH(TIMEj WI , TI

Zone Air Moisture Ba 1 ance

ELSE I F ( I . EQ . 2 ) THEN

WI » Y ( 2)

TI * Y(1 ) ZONEFN * ( INFILW(TIME,TI ) * (WOUTA - WI ) & SUPPLW ( T I ME ) & MCI ST( TI ME) ) / MASSA( TI ME, TI

If zona air reaches saturation, prevent further Increase

MAXW » HUMRAT(TI , 100, 0, PRESS)

IF(WI .GT.MAXW, AND. ZONEFN. GT. 0. 0) ZONEFN = 0.0

IF(WI .LE. 0.0. AND. ZONEFN.LT. 0.0) ZONEFN a 0.0

Wall Thermal Energy Balance

ELSE I F ( I . EQ . 3 ) THEN

TW > Y ( 3)

TIW Y ( 4

NOTE: TOUTA in the equation below should be the so lair temperature.

* (TOUTA - TW)

» (TIW - TW)

/ MASSWT I ME

ELSE

ZONEFN a 0. END IF

ZONEFN ( CN0CT0( TIME)

CNDCTI (TIME)

SOLAR ( TIME)

RETURN

END

FUNCTION SUPPLY (TIME)

laaaaBBaaaaaaBaaaBBaasaaasBBaasaaaaaxaaaassaaaarxsaBaaBaaBsasaasassssar

Returns the amount of energy entering or being removed from the zone

by the HVAC system. This Includes energy from the elr handling unit

end the local conditioning equipment.

UNITS: returns Btu/sec.

I nput ver tables:

TIME - Used to determine if 1 1 me -dependent functions should be recalculated.

INTEGER TIME

Output variables:

SUPPLY - Energy entering or being removed from the zone by the HVAC system.

98
1914

1919

1916

1917

1916

1919

1920

1921

1922

1923

1924

1923

1926

1927

1928

1929

1930

1931

1932

1933

1934

1933

1936

1937

1938

1939

1940

1941

1942

1943

1944

1943

1946

1947

1948

1949

1930

1931

1932

1933

1934

1933

1956

1937

1938

1959

1960

1961

1962

1963

1964

1965

1966

1967

1968

1969

1970

1971

oo

o

o

o

o

o

o

o

o

o

o

o

o

o

o

o

o

o

o

oo

oooo

REAL SUPPLY

Procedures called: ENTHAL

Common blocks: STATE

REAL A

LOGICAL L

COMMON / STATE1 / A( 1 00) , L( 1 00)

Local

FUNC

HSA

USA2

HZA

MRA

MSA1

MSA2

NOINIT

OLD

PRESS

QLOCH

QLOCC

TIMEDP

TZA

WZA

Variables:

- Holding variable for function return value.

- enthalpy of the supply air in BTU/lbm

- enthalpy of the hot deck, (if existing) supply air Btu/lbm

- enthalpy of zona air in Btu/lbm

- mass flow rate of return air (Ibtn/s)

* mass flow rate of supply air in LBm/sec.

- mass flow rate of hot deck, (if existing) supply air in LBm/sec.

- logical variable, true if function has never been calculated.

- time when function was last calculated.

- Atmospheric pressure (assummed constant) in Inches of Hg.

- energy added to the zone by local equipment during the time step.

- energy taken from the zone by local equipment during the time step,

- logical parameter , true i f funct ion is time dependent

- temperature of zone air in degrees F. humidity ratio of zone air in Ibm water/ Ibm dry air.

LOGICAL TIMEDP, NOINIT

PARAMETER (TIMEDP > .TRUE., PRESS * 29.921)

INTEGER OLD

DATA NOINIT/. TRUE./, 0LD/-1/

REAL MRA , MSA 1 , MSA2

EQUIVALENCE ( A(06),TZA ),( A(07),WZA ),( A(33),HSA1 ) EQUIVALENCE ( A(31),MRA ),( A(62), QLOCH ),( A(63),QL0CC ) EQUIVALENCE ( A(52),MSA1 )

I F ( NOI NIT. OR . (TIME. ME. OLD. AND. TIMEDP) ) THEN

THE FOLLOWING VALUES ARE FOR DUAL DECK (NOT IMPLEMENTED)

MSA2 s 0.0

HSA2 >0.0

C

HZA « ENTHAL (TZA, WZA, PRESS)

FUNC » ( HSA1 *MSA1 ) (HSA2«MSA2) - ( HZA*MRA ) QLOCH - QLOCC

C

OLD » TIME

NOINIT * .FALSE.

END IF

SUPPLY = FUNC

RETURN

END

CassssssasaassaaaaaasssaaasaaaaaaaaBasassaasssasassaaaacsaaxsrsxssxarxrz

FUNCTION I NF I LT ( T I ME , T I

, W I

C33>333333SS33S33Sa38333a33333S333383333333333333333aaaaaaa*3333333S3SX3

C Returns the Infiltration coefficient for energy gain/ loss due to

C outside air inf i Itrat ion.

C

C UNITS: returns Btu/sec F.

99
1972

1973

1974

1975

1976

1977

1978

1979

1980

1981

1982

1983

1984

1985

1986

1987

1988

1989

1990

1991

1992

1993

1994

1995

1996

1997

1998

1999

2000

2001

2002

2003

2004

2005

2006

2007

2008

2009

2010

201 1 2012

2013

2014

2015

2016

2017

2018

2019

2020

2021

2022

2023

2024

2025

2026

2027

2028

2029

oo

oooooo

o

ooooooooooo

o

o

oooo

oooo

c

C I nput var tables:

C TIME - Used to determine if time -dependent functions should be recalculated

C TI - Dry bulb temperature of zone air in F. C Wl - Humidity ratio of zone air in Ibm water/lbm dry air.

C

INTEGER TIME

Output var i ab 1 es

INFILT - the infiltration coefficient for energy fla in/ loss.

REAL INFILT

Procedures called: AXRATE , CPA I R , RHOA I

Common blocks: STATE1 , SHELL

REAL A

LOGICAL L

COMMON / STATE1 / A( 1 00) , L( 1 00)

EQUIVALENCE ( A( 1 ) , TOUTA ),( A(2),W0UTA ) REAL LENOW,LENIW,KOW,KIW

COMMON / SHELL / AREAW, Hi WALL, HOWALL, LENOW, LENI W, KOW , K I W , CFUR,

& VOLAIR, CPOW, RHCOW, CPIW, RHOI W, A I REXC, RMGFUR

Local Variables:

CP I A is specific heat of Infiltration air (average of inside and outside)

FUNC - Holding variable for function return value.

NO I NIT - logical variable, true if function has never been calculated.

OLD - time when function was last calculated.

RHOI A - density of Infiltration air (average of inside and outside)

TIMEDP - logical parameter, true if function is time dependent.

TOUTA - the outdoor air temperature in degrees F. WOUTA - the outdoor air humidity ratio in Ibm water/lbm dry air.

INTEGER OLD

LOGICAL TIMEDP, N01 NIT

PARAMETER (TIMEDP » .TRUE.)

DATA NO I NIT/. TRUE./, OLD/- 1/

I F ( NO I NI T . OR . ( T I ME . NE . OLD . AND . TIMEDP) ) THEN

RHOA I R is a function returning the density of moist air

AXRATE is a function for the air exchange rate in eir changes per second

VOLAIR is the volume of the zone air

RHOI A = ( RHOA I R( WOUTA, TOUTA) RHOA I R ( WI , T I )

/ 2. CPIA a ( CPA I R( WOUTA ) CPAIR(WI)) / 2. FUNC = CPIA * RHOI A VOLAIR « AXRATE ( TI ME , T I

OLD * TIME

NOINIT = .FALSE.

END IF INFILT * FUNC

RETURN

100
2030

2031

2032

2033

2034

2035

2036

2037

2036

2039

2040

2041

2042

2043

2044

2045

2046

2047

2048

2049

2050

2051

2052

2053

2054

2055

2056

2057

2058

2059

2060

2061

2062

2063

2064

2065

2066

2067

2068

2069

2070

2071

2072

2073

2074

2075

2076

2077

2078

2079

2080

2081

2082

2083

2084

2085

2086

2087

oooooooo

o

ooooooooo

o

o

oooo

oooo

oooooooooo

END

FUNCTION AXRATE(TIME, TI

Returns the air exchange rate for a building zone. UNITS: returns air exchange rate In air changes per second

TI and TO must be in degrees F

Input variables:

TIME - Used to determine if t 1 me -dependent functions should be recalculated.

TI - Dry bulb temperature of zone air in F. INTEGER TIME

Output variables:

INFILT - the infiltration coefficient for energy gain/ loss.

REAL AXRATE

Procedures called: none Common blocks: STATE 1

, SHELL

REAL A

LOG I CAL L

COMMON / STATE1 / A( 1 00) , L( 1 00)

EQUIVALENCE ( All), TOUTA* ) , ( A(04),WIND ) REAL LENOW, LENIW, KOW, KI

COMMON / SHELL / AREAW, HI WALL. HOWALL, LENOW, LENI W, KOW, KI W, CFUR,

& VOLA I R , CPOW , RHOOW , CPI W, RHO I W , A I REXC , RMGFUR

Local Variables:

FUNC - Holding variable for function return value.

NOINIT - logical variable, true if function has never been calculated.

OLD - time when function was last calculated.

TIMEDP - logical parameter, true if function is time dependent.

TOUTA - the outdoor air temperature in degrees F. VWIND - the outdoor wind velocity in miles per hour.

INTEGER OLD

LOGICAL TIMEDP, NOINIT

PARAMETER (TIMEDP = .TRUE.)

DATA NOI NIT/. TRUE./, OLD/- 1/

I F( NOI NIT. OR. (TIME. NE. OLD. AND. TIMEDP) ) THEN

A l REXC is the standard air exchange rate without the effects of wind

and temperature differential in air changes per hour:

= 1.5 for leaky building

3 1.0 for standard building

= 0.5 for modern type building

AXRATE=A I REXC/0 . 695* (0. 15+0. 01 3* VW I ND+0 . 005* ABS ( TOUTA- T I )

FUNC 3 AXRATE / 3600.

101
2088

2089

2090

2091

2092

2093

2094

2095

2096

2097

2098

2099

2100

2101

2102

2103

2104

2105

2106

2107

2108

2109

2110

21 1 1 2112

21 13

21 14

2115

21 16

2117

2118

21 19

2120

2121

2122

2123

2124

2125

2126

2127

2128

2129

2130

2131

2132

2133

2134

2135

2136

2137

2138

2139

2140

2141

2142

2143

2144

2145

noon

o

ooooooo

oooo

oooo

OLD = TIME

N0INIT = .FALSE.

END IF

AXRATE = FUNC

RETURN

END limn

FUNCTION CONVWA ( T I ME

C« l>f3t3l331tll>lllll33>llc:3Sllllllll3333l3333333338333:iI333:33:33:::3

C Returns the convection coefficient for convection heat transfer between

C the zone air and the wall (not including glass) surfaces enclosing the

C zone.

C

C UNITS: returns Btu/s F. C

C I nput var tables

C TIME * Used to determine if t 1 me -dependent functions should be recalculated

C

INTEGER TIME

Output var i ab l es

CONVWA - the convection coefficient between zone air and wall.

REAL CONVWA

Procedures called: none Common blocks: SHELL

REAL LENOW , LEN I W , KOW , K I

COMMON / SHELL / AREAW , H I WALL , HOWALL , LENOW , LEN I W , KOW , K I W , CFUR

& VOLA I R , CPOW, RHOOW, CPI W, RHO I U, A I REXC , RMGFUR

Local Variables:

FUNC - Holding variable for function return value.

NO I NIT - logical variable, true if function has never been calculated.

OLD - time when function was last calculated.

TIMEDP - logical parameter, true if function Is time dependent.

INTEGER OLD

LOGICAL TIMEDP, NO! NIT

PARAMETER (TIMEDP = .FALSE.)

DATA NO! NIT/. TRUE. /.OLD/- 1/

I F( NOI NIT. OR . (TIME. NE. OLD. AND. TIMEDP) ) THEN

HI WALL is the interior wall surface film coefficient

AREAW is the surface area of the well exposed to the outside

FUNC * HI WALL » AREAW

C - - - - OLD * TIME

NO! NIT « .FALSE.

END IF CONVWA = FUNC

RETURN

END

C aaas==a3 -saaassassas3aaaar=asaassaaazza3=s==s3=ss=saaaasaas=ss3s=zr==:==

FUNCTION C0NVGL( TIME)

102
2146

2147

2148

2149

2150

2151

2152

2153

2154

2155

2156

2157

2158

2159

2160

2161

2162

2163

2164

2165

2166

2167

2168

2169

2170

2171

2172

2173

2174

2175

2176

2177

2178

2179

2180

2181

2182

2183

2184

2185

2186

2187

2188

2189

2190

2191

2192

2193

2194

2195

2196

2197

2198

2199

2200

2201

2202

2203

oooooooo

oooooooooo

o

oooo

oooooooooooo

o

ooo

C Returns the convection coefficient for convection heat transfer between

C the zone air and the glass (windows) surfaces enclosing the zone.

C

C UNITS: returns Btu/s F. C

REAL CONVGL

THIS FUNCTION NOT IMPLEMENTED (NO GLASS IN WALL).

CONVGL

RETURN

END

0.0

FUNCTION I NTERN ( T I ME

Returns the sensible internal gain in the zone from lighting, people and ecu i pment

UNITS: returns Btu/s.

Input variables:

TIME - Used to determine if t i me -dependent functions should be recalculated.

Output variables:

INTERN * sensible internal gain.

REAL I NTERN

Procedures called: SCHED

Common b l ock.s : none

INTERN = SCHED(TIME, 'S’

RETURN

END

333313333333333333133333333333333333333333333333339333331333333333333:

FUNCTION SCHED ( T I ME , SORL

I131331333ll93l333313333313333333:i333939:3333333311133333i:::3i::::33:

This routine is used to determine the current values of sensible or latent Internal gains in the zone. The values of the gains are determined

from a schedule of occupancy , lighting and equipment usage.

Input variables:

SORL - Character used to select the return of latent gain if ‘ L‘ and

sensible gain’ if ‘S'.

TIME - Used to determine if t line -dependent functions should be recalculated

INTEGER TIME

CHARACTER* 1 SORL

Output variables: SCHED

SCHED - Contains either sensible internal gain In Btu/sec. or latent

moisture gain In Ibm water /sec.

Procedures called: none Common blocks: STATE 1

, TI ME, USE

PARAMETER (NSCHED = 10)

103
2204

2205

2206

2207

2206

2209

2210

221

2212

2213

2214

2215

2216

2217

2218

2219

2220

2221

2222

2223

2224

2225

2226

2227

2228

2229

2230

2231

2232

2233

2234

2235

2236

2237

2238

2239

2240

2241

2242

2243

2244

2245

2246

2247

2248

2249

2250

2251

2252

2253

2254

2255

2256

2257

2258

2259

2260

2261

oooo

oooooo

ooo

o

o

ooooooooooo

REAL

COMMON /

&

C

I NTEGER

COMMON /

C

REAL

LOGICAL

COMMON /

TSCHED, SENSIB, LATENT

USE / LSCHED , TSCHEO ( NSCHED )

, SENS I B ( NSCHEO ) , LATENT ( NSCHED

, PESENS ( NSCHED )

, PELATE ( NSCHED

DAY. HOUR, MINUTE, SECONO, ENDDAY, ENDHR, ENDMIN, ENDSEC

TIME / DAY .HOUR, MINUTE, SECOND, ENDDAY .ENDHR .ENDMIN, ENDSEC

A

L

STATE1 / A< 100) , L( 100)

Loca 1 Var tables:

DHOUR - current hour of the day In decimal hours.

ISCHED - Index Into schedule tables for sensible and latent loads for

current time of day.

NOINIT - logical variable, true if function has never been calculated.

OCCUPY - STATUS zone occupied (. true. ) /unoccup led (.false.).

OLD - time when function was last calculated.

TIMEDP - logical parameter, true if function is time dependent.

INTEGER OLD

LOGICAL TIMEDP, NOINIT, OCCUPY

PARAMETER (TIMEDP * .TRUE.)

DATA NOI NIT/. TRUE./, OLD/- 1/

EQUIVALENCE ( L( 11), OCCUPY )

I F( NOI NIT OR. (TIME. NE. OLD. AND. TIMEDP) ) THEN

Find row of schedule table for current time of day

DHOUR * HOUR MINUTE/60. SEC0ND/3600.

DO 1 000 I « 1 , LSCHED

I F ( TSCHED ( I ) . GT . DHOUR ) GO TO 2000

1000 CONTINUE

ISCHED *1-1

OLD * TIME

NOINIT * .FALSE.

END IF

Determine sensible or latent load from tables.

The following assumes that the sensible and moisture gains of people

ere constant under ell conditions. This Is not strictly correct.

I F ( SORL . EQ. 'L* ) THEN

SCHED » LATENT ( ISCHED) PELATE ( I SCHED

ELSE I F ( SORL . EQ . ‘S') THEN

SCHED * SENSIB ( ISCHED) PESENS ( I SCHED

determine If building Is occupied

I F ( PESENS ( I SCHED ) . GT . 0 ) THEN

OCCUPY = .TRUE.

ELSE

OCCUPY » .FALSE.

104
2262

2263

2264

2269

2266

2267

2266

2269

2270

2271

2272

2273

2274

2279

2276

2277

2278

2279

2280

2281

2282

2283

2284

2289

2286

2287

2288

2289

2290

2291

2292

2293

2294

2299

2296

2297

2298

2299

2300

2301

2302

2303

2304

2309

2306

2307

2308

2309

2310

231 1 2312

2313

2314

2319

2316

2317

2318

2319

oooooo

o

oooooooo

oooo

oooo

END IF

ELSE

SCHED > 0.

END IF

RETURN

END

FUNCTION MASSTH ( T I ME , W I , T I ) Cl3I99l3933l8SI33333i333333333313S3i3333S333333333Blil33ll3C333333333333

C Returns the thermal mass of the zone air and furnishings (not including

C walls).

C

C UNITS: returns Btu/F.

C

C Input variables:

C TIME - Used to determine if t i me -dependent functions should be recalculated.

C TI - Dry bulb temperature of zone air in F. C WI - Humidity ratio of zone air in Ibm water /Ibm dry air.

C

INTEGER TIME

Output var tables:

MASSTH - the thermal mass of the zone air and furnishings in Btu/F.

REAL MASSTH

Procedures called: CPAIR,RHOAIR

Common b 1 ock.s : SHELL

REAL LENOW , LEN I

, KOW , K I W

COMMON / SHELL / AREAW, HI WALL, HOWALL, LENOW , LEN I W , KOW , K I W, CFUR

6 VOLA I R , CPOW , RHOOW . CP I W , RHO I W , A I REXC , RMGFUR

Local Variables:

CAIR - thermal mass of the zone air.

FUNC - Holding variable for function return value.

NOINIT - logical variable, true if function has never been calculated.

OLD - time when function was last calculated.

TIMEDP - logical parameter, true if function is time dependent.

INTEGER OLD

LOGICAL TIMEDP, NOINIT

PARAMETER (TIMEDP = .TRUE.)

DATA NOINIT/. TRUE. /,0LD/-1/

I F( NOINIT. OR. (TIME. NE. OLD. AND. TIMEDP) ) THEN

VOLA I R is the volume of the zone air

RHOAIR is a function returning the density of moist air

CPA I R is a function returning the specific heat of moist air

CFUR is the thermal mass of the zone furnishings

CAIR » CPAIR(WI) « RHOA I R ( W I , T I ) * VOLAIR

FUNC = CFUR + CAIR

C - - OLD » TIME

NOINIT = .FALSE.

END IF

105
2320

2321

2322

2323

2324

2329

2326

2327

2328

2329

2330

2331

2332

2333

2334

2339

2336

2337

2338

2339

2340

2341

2342

2343

2344

2349

2346

2347

2348

2349

2390

2391

2392

2393

2394

2399

2396

2397

2398

2359

2360

2361

2362

2363

2364

2369

2366

2367

2368

2369

2370

2371

2372

2373

2374

2379

2376

2377

onooon

o

ooooooooo

o

oooo

ooooooo

MASSTH = FUNC

RETURN

END

C>ii»3»i>iiiiciii33:33iiiiiiii3mit33i3i:3:3 3::B3iinitii::3 3 3:::::::

FUNCTION I NF I LW ( T I ME , T l , W I

I3I33333I3313333I333333S33:33333333333333:3333S33I33333I333333:3::333:

C Returns the infiltration coefficient for moisture gain/ loss by outside

C air infiltration.

C

C UNITS: returns Ibm dry air/second.

C

C I nput var tables:

C TIME - Used to determine if 1 1 me -dependent functions should be recalculated.

C TI - Dry bulb temperature of zone air In F. C WI - Humidity ratio of zone air in Ibm water/lbm dry air.

C

INTEGER TIME

Output variables:

INFILW - the mass flow rate of dry air from outside to inside which

carries water vapor with it. The actual change in moisture will

be assumed to be: INFILW times the i ndoor -outdoor humidity ratio

different ial

REAL INFILW

Procedures called: AXRATE , RHOA i

Common blocks: SHELL, STATE

REAL A

LOGICAL L

COMMON / STATE1 / A( 1 00) , LC 1 00)

REAL LENOW , LEN I W , KOW , KIW

COMMON / SHELL / AREAW, H I WALL , HOWALL , LENOW , LEN I W , KOW , K I W , CFUR

& VOLA I R , CPOW , RHOOW , CP I W , RHO I W , A I REXC , RMGFUR

Local Variables:

FUNC - Holding variable for function return value.

NOINIT - logical variable, true if function has never been calculated.

OLD - time when function was last calculated.

RHO I A - density of DRY infiltration air (average of inside and outside)

TIMEDP - logical parameter, true if function is time dependent.

TOUTA - the outdoor air temperature in degrees F. INTEGER OLD

LOGICAL TIMEDP, NOINIT

PARAMETER (TIMEDP = .TRUE.)

EQUIVALENCE ( A(1), TOUTA ) DATA NOINIT/. TRUE. /,0LD/-1/

! F( NOINIT. OR. ( T I ME . NE . OLD . AND . T I MEDP ) ) THEN

RHOA I R Is AXRATE is

VOLA I R is

a function returning the density of moist air

a function for the air exchange rate in air changes per second

tho volume of the zone air.

106
2378

2379

2380

2381

2382

2383

2384

2383

2386

2387

2388

2389

2390

2391

2392

2393

2394

2393

2396

2397

2398

2399

2400

2401

2402

2403

2404

2405

2406

2407

2408

2409

2410

2411

2412

2413

2414

2413

2416

2417

2418

2419

2420

2421

2422

2423

2424

2423

2426

2427

2428

2429

2430

2431

2432

2433

2434

2433

oooo

oooo

ooooooooo

RH0IA » ( RHGA I R ( 0 . 0 , TOUTA ) + RH0AIR(0. 0, TI ) ) / 2.

FUNC « RHOIA’ * VOLAIR * AXRATE ( T 1 ME , T I

OLD « TIME

NO! NIT » .FALSE.

END IF

INFILW » FUNC

RETURN

END

FUNCTION SUPPLWC TI ME)

Returns the net amount of water entering the zone from the air handling

unit. If negative, water Is being extracted from the zone. UNITS: Ibm water/ second.

Input variables:

TIME - Used to determine If 1 1 me -dependent functions should be recalculated.

INTEGCR TIME

Output variables:

SUPPLW - water entering or being removed from the zone by the HVAC system.

REAL SUPPLW

Procedures called: none Common blocks: STATE!

C

C

C

C

C

C

C

C

C

C

C

C

C

C

REAL A

LOGICAL L

COMMON / STATE! / A( 1 00) , L( 1 00)

Local

FUNC

MRA

MSA1

MSA2

NOINIT

OLD

PRESS

TIMEDP

WRA

WSA1

Variables:

- Holding variable for function return value.

- mass flow rate of return air (Ibm/s)

- mass flow rate of supply air in LBm/sec.

- mass flow rate of hot deck, (if existing) supply air in LBm/sec.

- logical variable, true if function has never been calculated.

- time when function was last calculated.

- Atmospheric pressure (assummed constant) in Inches of Hg.

- logical parameter, true if function is time dependent.

- humidity ratio of return air in Ibm water/ Ibm dry elr.

- humidity ratio of supply air in Ibm water/ Ibm dry elr.

INTEGER OLD

LOGICAL TIMEDP, NOINIT

PARAMETER (TIMEDP = .TRUE.)

DATA NOINIT/. TRUE. /.0LD/-1/

REAL MRA , MSA 1, MSA2 , WRA , WSA 1, WSA2

EQUIVALENCE ( A(22),WRA ),( A(31),MRA ) EQUIVALENCE ( AC34),WSA1 ),( A(32),MSA1 ) C

I F( NOINIT. OR. (TIME. NE. OLD. AND. TIMEDP) ) THEN

C

107
2436

2437

2438

2439

2440

2441

2442

2443

2444

2449

2446

2447

2448

2449

2450

2431

2452

2453

2454

2455

2456.

2437

2458

2459

2460

2461

2462

2463

2464

2465

2466

2467

2468

2469

2470

2471

2472

2473

2474

2475

2476

2477

2478

2479

2480

2481

2482

2483

2484

2485

2486

2487

2488

2489

2490

2491

2492

2493

oooo

oooo

oooooooooo

o

oooo

ooooooooooo

o

o

ooo

C THE FOLLOWING VALUES ARE FOR DUAL DECK (NOT IMPLEMENTED)

MSA2 0.0

WSA2 >0.0

NOTE: this is approximate, since MRA and MSA1 are in Ib/sec of moist air.

FUNC « ( WSA1 a MSA1 ) + ( WSA2 a MSA2) - (WRA a MRA)

OLD TIME

NO! N! T a .FALSE.

END I F

SUPPLW a FUNC

RETURN

END

FUNCTION MOISTCTIME)

Returns the amount of water generated in the zone by people and equipment.

UNITS: Ibm water/ second.

Input variables:

TIME - Used to determine if t 1 me -dependent functions should be recalculated.

Output variables:

MOIST - moisture gain from people and „equ I pment . REAL MOIST

Procedures called: SCHED

Common b lodes: none MOIST SCHEDITIME, 'L'

RETURN

END

33332323838388333383338323333X38833X3333323333333398382332233322332233

FUNCTION MASSA ( T I ME , T I

999339999aa9lt39III3a3l9I33i9999399l9ai3193339333339I39IIi:3333::393::

Returns the mass of the zone air which can absorb and release water

vapor, (mass of dry air)

UNITS: returns Ibm of dry air.

I nput var tables:

TIME - Used to determine if t I me -dependent functions should be recalculated.

TI - Dry bulb temperature of zone air in F. INTEGER TIME

Output variables:

MASSA - the infiltration coefficient *o- energy gain/ loss.

REAL MASSA

Procedures called: RHOAIR

Common blocks: SHELL

REAL LENOW , LEN I W , KOW , K1W

108
2494

2499

2496

2497

2490

2499

2900

2901

2902

2903

2904

2909

2906

2907

2900

2909

2910

291 1 2912

2913

2914

2919

2916

2917

2910

2919

2920

2921

2922

2923

2924

2929

2926

2927

2926

2929

2930

2931

2932

2933

2934

2939

2936

2937

2930

2939

2940

2941

2942

2943

2944

2949

2946

2947

2946

2949

2990

2991

COMMON / SHELL / AREAW , H I WALL , HOWALL , LENOW , LEN I W , KOW ( KIW, CFUR

& VOLA I R , CPOW , RHOOW ,CP1W, RHO I W , A I REXC , RMGFUR

C

C Local Variables:

C FUNC - Holding variable for function return value.

C NOIN1T - logical variable, true If function has never been calculated.

C OLD - time when function was last calculated.

C TIMEDP - logical parameter, true If function is time dependent.

C WI - zone humidity ratio, set to zero for dry air.

C

INTEGER OLD

LOGICAL TIMEDP, NOI NIT

PARAMETER (TIMEDP « .TRUE.)

DATA NOI NIT/. TRUE./, OLD/- 1/

C

I\ NOI NIT. OR. (TIME. NE. OLD. AND. TIMEDP) ) THEN

C

- - C VOLA I R Is the volume of the zone air.

C RHOAIR is a function returning the density of moist air.

C RMGFUR is a factor used to produce the effective air mass for moisture

C absorption by the zone furnishings and air.

C

WI > 0.0

FUNC » RHOA IR(WljTI) * VCLAIR * RMGFUR

C

OLD a TIME

NOI NIT a .FALSE.

END IF

MASSA a FUNC

RETURN

END

FUNCTION CNDCTO(TIME)

C3ll3333l:ll3l333333933:3333:3333:8S3t33333:::3::33>333333II3333::::::3:

C Returns the coefficient for thermal conductance between the mass center

C of the wall and the environment outside of the building.

C

C UNITS: returns Btu/second F. C

C I nput var tables:

C TIME - Used to determine if t i me -dependent functions should be recalculated

C

INTEGER TIME

C

C Output var tables:

C CNDCTO - conductance coefficient between outside air and wall mass center

C

REAL CNDCTO

C

C Procedures called: none

C Common blocks: SHELL

C

REAL LENOW, LENIW, KOW, KIW

COMMON / SHELL / AREAW , HI WALL , HOWALL, LENOW , LEN I W, KOW , K I W , CFUR

0 VOLA I R , CPOW , RHOOW , CP I W , RHO I W , A I REXC , RMGFUR

C

C Loca l Var tables:

C FUNC - Holding variable for function return value.

109
2552

2553

2554

2555

2556

2557

2556

2559

2560

2561

2562

2563

2564

2565

2566

2567

2568

2569

2570

2571

2572

2573

2574

2575

2576

2577

2578

2579

2580

2581

2582

2583

2584

2505

2586

2507

2500

2509

2590

2591

2592

2593

2594

2595

2596

2597

2598

2599

2600

2601

2602

2603

2604

2605

2606

2607

2608

2609

oooo

oooo

ooooooooo

o

o

ooooooo

C N0INIT - logical variable, true If function has never been calculated.

C OLD - time when function was last calculated.

C T1ME0P - logical parameter, true if function is time dependent.

C

INTEGER OLD

LOGICAL TIMEDP.NOINIT

PARAMETER ( T I MEDP * . FALSE .

DATA NOINIT/. TRUE./, OLD/- 1/

C

I FC NOINIT. OR. (TIME. NE. OLD. AND. TI MEDP) ) THEN

HOWALL is the surface film coefficlont for the outside wall

LEN0W is the thickness of the wall outside of the mass center

KOW is the thermal conductivity of the wall outside of the mass center

AREAW is the area of the wall exposed to the outside

FUNC * 1 / (1 /HOWALL LENOW/KOW) « AREAW

OLD » TIME

NOINIT = .FALSE.

END IF CNDCTO » FUNC

RETURN

END

iasaii:aiiiiiaixsi3sxais:i:ii3iiiii3:siis:3:33aizissssxic:::33::::::z

' FUNCTION CNDCTI (TIME)

8383333333833333X33333333338883383833333333333333338388883333333333333

Returns the coefficient for thermal conductance between the mass center

of the wall and the Inner wall surface surrounding the zone. UNITS: returns Btu/second F.

I nput var tables •' TIME - Used to determine if 1 1 me -dependent functions should be recalculated.

INTEGER TIME

Output variables:

CNDCTI - conductance coefficient between zone air and wall mass center.

REAL CNDCTI

Procedures called: none Common blocks: SHELL

REAL LENOW , LEN I W , KOW , K I

COMMON / SHELL / AREAW, H I WALL , HOWALL , LENOW . LEN I W , KOW , K l W, CFUR

4 VOLA I R , CPOW , RHOOW . CP I W , RHO I W , A I REXC , RMGFUR

C

C Local Variables:

C FUNC - Holding variable for function return value.

C NOINIT - logical variable, true If function has never been calculated.

C OLD - time when function was last calculated.

C T I MEDP - logical parameter, true if function Is time dependent.

C

INTEGER OLD

LOGICAL TIMEDP.NOINIT

110
2610

261 1 2612

2613

2614

2619

2616

2617

2616

2619

2620

2621

2622

2623

2624

2623

2626

2627

2628

2629

2630

2631

2632

2633

2634

2639

2636

2637

2638

2639

2640

2641

2642

2643

2644

2643

2646

2647

2648

2649

2690

2631

2632

2653

2634

2633

2656

2637

2638

2639

2660

2661

2662

2663

2664

2665

2666

2667 oo

oooo

oooo

oooooooo

o

oooo

ooooo

o

o

oooooo

PARAMETER (TIMEDP = .FALSE.)

0ATA N0IN1T/. TRUE. /.OLD/- 1/

C

I F( NOI NI T . OR. ( TIME. NE. OLD. AND. TIMEDP) ) THEN

KIW is the thermal conductivity of the wall ins Ida of the mass canter

LENIW is the thickness of the wall inside the mass center

AREAW is the surface area of wall exposed to the outside

FUNC » KIW / LENIW * AREAW

OLD > TIME

NOINIT * .FALSE.

END IF

CNDCTI a FUNC

RETURN

END aasaaaaaaaaaaassasasaasssaassxssaaasaxsaasssaassasaxsasassssaasssssxss

FUNCTION SOLAR l T I ME

3aa3a3a33aiSBa33aa3333aSS33ai33I3l3X33333333333333333SS33333333333S=3:

Returns the amount of solar heat gain through the windows in the zone. UNITS: returns Btu/second.

REAL SOLAR

THIS FUNCTION NOT IMPLEMENTED (NO GLASS IN WALL).

SOLAR >0.0

RETURN

END

iiiaz:223ixsss3i3Bi3:2si:3xi2C3is:xx22t3s:::i3:23Eii3ii3i3::::2:::::::

FUNCTION MASSWCTIME)

333333X38X333333383333333333X333333X33333333333333333333333=333=======

Returns the thermal capacitance of the wall surrounding the zone. UNITS: returns Btu/F

I nput var tables:

TIME - Used to determine if time -dependent functions should be recalculated

INTEGER TIME

Output variables:

MASSW - therma l capac i tance of wall.

REAL MASSW

Procedures called: none

Common blocks: SHELL

REAL LENOW, LENIW, KOW, KIW

COMMON / SHELL / AREAW. HI WALL. HOWALL, LENOW , LENI W, KOW , K I W

CFUR , & VOLAi R, CPOW, RHOOW , CPI W. RHO I W , A I REXC , RMGFUR

Local Variables:

111
2668

2669

2670

2671

2672

2673

2674

2675

2676

2677

2678

2679

2680

2681

2682

2683

2684

2685

2686

2687

2688

2689

2690

2691

2692

2693

2694

2695

2696

2697

2698

2699

2700

2701

2702

2703

2704

2705

2706

2707

2708

2709

2710

271 1 2712

2713

2714

2715

2716

2717

2718

2719

2720

2721

2722

2723

2724

2725

ooo

oooo

o

ooo

oooo

FUNC - Holding variable for function return value.

N0INIT - logical variable, true if function has never been calculated.

OLD - time when function was last calculated.

TIMEDP - logical parameter, true if function is time dependent.

INTEGER OLD

LOGICAL TIMEDP, NOI NIT

PARAMETER (TIMEDP = .FALSE.)

OATA NOI NIT/. TRUE./, OLD/- 1/

I F ( NO I N I T . OR . (TIME. NE. OLD. AND TIMEDP) ) THEN

c-

c CPOW is the i c RHOOW is the

c LENOW is the

c CPIW is the i c RHOIW is the

c LENIW is the

c AREAW is the

c FUNC » (I

1 inside of mess center,

'all inside of mass center,

ill exposed to the outside

a (CPOW « RHOOW * LENOW CPIW » RHOIW * LENIW) « AREAW

OLD » TIME

NOI NIT * .FALSE.

END IF

MASSW FUNC

RETURN

END

aaBassBsasBBBasssaaas:

FUNCTION WRADIA(TIME)

laaasaaaaaaaa:

Returns the coefficient for radiation heat transfer between the wall and

any other walls or furnishings (at room mean radiant temperature )

REAL WRADIA

THIS FUNCTION NOT IMPLEMENTED (no radiation assumed).

WRADIA > 0.0

RETURN

END aaaaaaaaaaaaaaaasaaaaaaaaaaaaaaaaazassaaaasassasraaaaaaaaaaaasaaas

FUNCTION WSOLAR(TIME)

BaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaseaasaaaaasaaaaaaaaaassaBasaaaaaaasaaaa

Returns the amount of solar energy entering the zone which is absorbed

by the wall surface.

REAL WSOLAR

THIS FUNCTION NOT IMPLEMENTED (no windows assumed).

WSOLAR a 0.0

RETURN

END

Qiai:a:3233Jnii3ji3i:::jjs3iTi33i333r:«S22:3232:2:c«33« = r.ii = ::a = ::::::: FUNCTION GCNDCT ( TIME)

Cl*SSSa3a9aaaSISai53ISS2SSB3ISS2BS29SSS9SS933:SSSSSSSSaSS2ISB3S53S:S5S=::

C Returns the coefficient for conductance through the glass In the zone

112
2726

2727

2728

2729

2730

2731

2732

2733

2734

2735

2736

2737

2738

2739

2740

2741

2742

2743

2744

2745

2746

2747

2748

2749

2750

2751

2752

2753

2754

2755

2756

2757

2758

2759

2760

2761

2762

2763

2764

2765

2766

2767

2768

2769

2770

2771

2772

2773

2774

2775

2776

2777

2778

2779

2780

2781

2782

2783

oooo

ooo

o

oooo

oooo

o

oooo

oooo

o

oooo

oo from the Inner wpll surface and the environment outside the building.

REAL GCNDCT

THIS FUNCTION NOT IMPLEMENTED (no windows assumed).

OCNDCT 0.0

RETURN

END

aaaaaaaaaasBasssaBassaaasassasssaaassasaszassasassaasasasssszsssssssss

FUNCTION GCNVCTCTIME)

aaB3I3a3B8aaS333aaaS33a383333asaaS3333333333333333sasa383S333333333333

Returns the coefficient for convection between the inner glass surface

and the air In the zone. REAL GCNVCT

THIS FUNCTION NOT IMPLEMENTED (no windows assumed).

GCNVCT 0.0

RETURN

END

8883388838883838388338838338833888888383833333333888883833333333333333

FUNCTION GRAD I A( TIME)

8888838388888888888333338833883388883833333333333833383333333833333333

Returns the coefficient for radiation heat transfer between the glass

and any other surfaces or furnishings in the zone. REAL GRAD I

THIS FUNCTION NOT IMPLEMENTED (no windows assumed).

GRADIA « 0.0

RETURN

END

SS3S3SSSaZXSa3B8SSS3SSSSSS2aSSaaS3SSSSSSS3S33S3SSBae3SSSX=3=S===S=:==

FUNCTION MRT (TIME)

88888388838838333888333833333333338333333333333333388333883333333333333

Returns the value of the mean radiant temperature in the zone. REAL MRT

THIS FUNCTION NOT IMPLEMENTED (no radiation transfer assumed).

MRT = 0.0

RETURN

END

^338333333333333333383333333333333333333333333333333333833 83333 3333333 S3

SUBROUTINE COMFRT(TIMEST)

033333333333333333333333333333333333333333333333333333333333333333333333

C This routine is used to determine if a zone space Is within the limits

C of comfort. Comfort is defined in accordance with ASHRAE standard

C 55-74 as being acceptable when the dry bulb temperature and humidity

C ratio are within certain limits. The check, of comfort is only performed

113
2784

2785

2786

2787

2788

2789

2790

2791

2792

2793

2794

2795

2796

2797

2798

2799

2800

2801

2802

2803

2804

2805

2806

2807

2808

2809

2810

281 1 2812

2813

2814

2815

2816

2817

2818

2819

2820

2821

2822

2823

2824

2825

2826

2827

2828

2829

2830

2831

2832

2833

2834

2835

2836

2837

2838

2839

2840

2841

oooo

oooo

o

o

oo

o

o

oo

o

oo

oo

o

o

ooooo

C if the building zone is occupied, as determined by the logical variable

C OCCUPY in the state vector.

C

C I nput var tables:

C TIMEST - The current value to be used for the major simulation time step.

C

INTEGER TIMEST

Output variables: none Procedures called: RELHUM

Common blocks: STATE

REAL A

LOGICAL L

COMMON / STATE1 / A( 1 00) , L< 1 00)

Local

OCCUPY

OUT

PRESS

RHAMAX

RHAMIN

SECOUT

TCLEFT

TCRITE

T2A

T2MAX

T2MIN

WZA

Var tables:

- STATUS zone occupied C . true . ) /unoccup 1 ed (.false.).

- Logical variable which is true when comfort Is outside of limits.

- Atmospheric pressure (assummed constant) in Inches of Hg.

- COMFORT relative humidity at maximum temperature (occupied) (X)

- COMFORT relative humidity at minimum temperature (occupied) (X)

- COMFORT seconds outside of comfort range

- Minimum temperature comfort limit, function of humidity ratio.

- Maximum temperature comfort limit, function of humidity ratio.

- temperature of zone air in degrees F.

- COMFORT maximum dry bulb temperature during occupied period (F)

- COMFORT minimum dry bulb temperature during occupied period (F).

- humidity ratio of zone air in Ibm water/ Ibm dry air.

PARAMETER (PRESS 29.921)

LOGICAL OCCUPY, OUT

EQUIVALENCE ( L( 1 1 ). OCCUPY ) , ( A(06),TZA ) , ( A(07),WZA ) EQUIVALENCE ( A(90), SECOUT ) , ( A(91),TZMAX ) , ( A(92),TZMIN ) EQUIVALENCE ( A (93), RHAMAX ) ,

( A (94), RHAMIN )

Test for conditions outside of comfort limits

I F( OCCUPY) THEN

TCLEFT » 72.58 - 131. 58* (WZA)

TCRITE * 81.73 - 394. 74* (WZA)

I F( WZA. GT. 0.01 2. OR.WZA.LT. 0.0044) THEN

OUT * .TRUE.

ELSE I F(TZA.LT. TCLEFT. OR. TZA.GT. TCRITE) THEN

OUT * .TRUE.

ELSE

OUT « .FALSE.

END IF IF(OUT) SECOUT » SECOUT TIMEST

Test for max and min space conditions

I F ( TZA . GT . TZMAX ) THEN

TZMAX a TZA

RHAMAX = RELHUM ( TZA, WZA, PRESS)

ELSE IF(TZA.LT.TZMlN)THEN

114
2842

2843

2844

2843

2846

2847

2848

2849

2830

2831

2832

2833

2834

2833

2856

2837

2838

2839

2860

2861

2862

2863

2864

2863

2866

2867

2868

2869

2870

2871

2872

2873

2874

2873

2876

2877

2878

2879

2880

2881

2882

2883

2884

2885

2886

2887

2888

2889

2890

2891

2892

2893

2894

2893

2896

2897

2898

2899

O

O

O

O

O

O

TZMIN TZA

RHAMIN » RELHUM ( TZA , WZA , PRESS

ENOIF

END IF

Test for max rate of temperature change

NOT IMPLEMENTED

RETURN

END

SUBROUTINE CMP I LE ( T I MEST

C3II333333339339S333933333333333399333333333333333333t333333::333:3::333

C This routine is used to compile cumulative totals of energy use. C Four types of energy use are compiled. These are:

C 1 . Fan Energy Use

C 2. Economizer Energy = energy provided by economizer.

C 3. Load Energy = energy required to heat or cool air to space entering

C condition.

C 4. Reheat /recool Energy = energy used to reheat and recool air for space

C temperature control.

C Total Thermal Energy Use is also compiled. All thermal energy values have

C a heat i ng component and a coo l i ng component

C

C Input variables:

C T I MEST - The current value to be used for the major simulation time step.

C

INTEGER T I MEST

C

C Output variables: none

C Procedures called: RHOA I R , SHOTS

C Common blocks: HVAC, STATE

LOG I CAL DUDECK , HUMCON, RAFAN , VAV

COMMON / HVAC / DUDECK, HUMCON, RAFAN, VAV, VOLSAF, CAPLOC, PGA I NL,

& TCTHSE, QRAFAN, FVAMI N, QSAFAN, HCAP1 ,CCAP1 , DTSETB,

& POWSAF , POWRAF , GA I N I

REAL

LOGICAL

COMMON / STATE1 / A ( 1 00 ) , L ( 1 00

Local Variables

EACTC

EACTH

EECOC

EECOH

EFANS

ELOADC

ELOADH

ELOCC

ELOCH

EREC

EREH

EREQC

EREQH

actual energy for cooling (Btu)

actual energy for heating (Btu)

COOLING ENERGY (economizer) (Btu)

HEATING ENERGY (economizer) (Btu)

ELECTRICAL ENERGY to fans (kWh)

COOLING ENERGY (load) (Btu)

HEATING ENERGY (load) (Btu)

COOLING ENERGY from local zone equipment (Btu)

HEATING ENERGY to local zone equipment (Btu)

(recool) (Btu)

(reheat) (Btu)

(requirements) (Btu)

( requ 1 rements ) (Btu)

COOLING ENERGY

HEATING ENERGY

COOLING ENERGY

HEATING ENERGY

115
2900

2901

2902

2903

2904

2905

2906

2907

2908

2909

2910

291 1 2912

2913

2914

2913

2916

2917

2918

2919

2920

2921

2922

2923

2924

2925

2926

2927

2928

2929

2930

2931

2932

2933

2934

2935

2936

2937

2938

2939

2940

2941

2942

2943

2944

2943

2946

2947

2948

2949

2930

2931

2932

2933

2934

2955

2936

2937 oo

ooooooo

c

C

C

C

C

C

C

C

C

C

C

c

cc

c

cc

cc

cc

c

c

cc

HMA

HMC

HSA1

HSA2

MRA

MSA 1 MSA2

POWR

POWS

QACT1

QACT2

QLOAO

QREQ1

QREQ2

RVOLR

RVOLS

TRA

TSAI

TSA2

VRA

VSA

WRA

WSA1

WSA2

ENTHALPY of mixed air (BTU/lb)

ENTHALPY of return air plus ventilation airCBTU/lb)

ENTHALPY of supply air (BTU/lb) cold deck

ENTHALPY of hot deck (if existing) supply air in (Btu/lbm) mass flow rate of return air (Ibm/s) mass flow rate of supply air in LBm/sec. mass flow rate of hot deck (if existing) supply air in LBm/sec.

electric power to return air fan (kW)

electric power to supply air fan (kW)

actual energy for cold deck (or single deck) (Btu).

actual energy for hot deck (if existing) (Btu).

load energy (Btu), may be negative for energy extraction (cooling)

requirements energy for cold deck (or single deck) (Btu).

requirements energy for hot deck (if existing) (Btu).

ratio of current return volume flow rate to rated volume flow rate

ratio of current supply volume flow rate to rated volume flow rate

TEMPERATURE of return air (F)

TEMPERATURE of supply air (F) cold deck

Temperature of hot deck (if existing) supply air (F)

Volume flow rate of return air in ft* 3/sec.

Volume flow rate of supply air in ft* 3/sec.

HUMIDITY RATIO of return air (lbw/lbda)

HUMIDITY RATIO of supply air (lbw/lbda) cold deck

Humidity ratio of hot deck (if existing) supply air (lbw/lbda)

REAL MSA1 . MSA2. MRA

EQUIVALENCE

EQUIVALENCE

EQUIVALENCE

EQUIVALENCE

EQUIVALENCE

EQUIVALENCE

EQUIVALENCE

A( 52) , MSA1

A( 64 ) , EREQH

A ( 66 ) , EECOH

A ( 69 ) , ELOAOC

A( 71 ), EREC

A( 32) , TSAI

A (31 )

, MRA

) , ( A( 33) , HSA1

) , ( A( 65) , EREQC

), ( A( 67) , EECOC

) , ( A( 73) , ELOCH

),( A( 72) , EFANS

) , ( A ( 20) , TRA

) , ( A(74) , ELOCC

, ( *A ( 24 ) , HMA )

, ( AOS) , HMC )

, ( A( 68) , ELOADH -)

, ( A (70) , EREH )

, ( A ( 34 )

, WSA 1 )

, ( A(22) , WRA )

Fan Energy

These calculations are based on fan law number 3 from ASHRAE equipment

handbook 1983, page 3.3, table 2. MSA2 * 0.

WSA2 > 0.

TSA2 0.

I F ( MSA 1 . GT .0.0) THEN

VSA MSA1 / RHOA I R ( WSA 1 , TSA 1

I F ( MSA2 . GT .0.0) VSA * VSA MSA2 / RHOA I R ( WSA2 , TSA2

RVOLS » VSA / VOLSAF

POWS POWSAF * RVOLS * RVOLS * RVOLS

EFANS * EFANS POWS TIMEST / 3600.

IF(RAFAN) THEN

VRA = MRA / RHOA I R ( WRA , TRA

RVOLR = VRA / VOLSAF

POWR = POWRAF * RVOLR « RVOLR « RVOLR

EFANS » EFANS POWR « TIMEST / 3600.

END IF END IF

116
2996

2999

2960

2961

2962

2963

2964

2969

2966

2967

2968

2969

2970

2971

2972

2973

2974

2979

2976

2977

2978

2979

2980

2981

2982

2983

2984

2989

2986

2987

2988

2989

2990

2991

2992

2993

2994

2999

2996

2997

2998

2999

3000

3001

3002

3003

3004

3009

3006

3007

3008

3009

3010

301 1 3012

3013

3014

3019 no

oooo

ooooo

oooo

oooo

o

oo Requ i rements Energy

MSA2 >0.0

HSA2 >0.0

QREQ1 > MSA1 * (HSA1 - HMC) * TJMEST

QRE02 » MSA2 * ( HSA2 - HMC) « TIMEST

I F(QREQ1 . 0T. 0. ) THEN

EREQH > EREQH QREQ1

ELSE

EREQC > EREQC QREQ1

END IF

IFCQREQ2.9T.0. ) THEN

EREQH > EREQH * QREQ2

ELSE

EREQC > EREQC QREQ2

END IF

Actua 1 Energy

QACT1 > MSA1 « ( HSA1 - HMA) * TIMEST

QACT2 > MSA2 * ( HSA2 - HMA) * TIMEST

I F(QACT1 . GT. 0. ) THEN

EACTH EACTH * QACT1

ELSE

EACTC » EACTC QACT1

ENDIF

I F ( QACT2 . GT . 0 . ) THEN

EACTH « EACTH QACT2

ELSE

EACTC « EACTC QACT2

ENDIF

Econom i zer Energy

EECOH = EREQH - EACTH

EECOC = EREQC - EACTC

Load Energy

QLOAD ( (MSA1 *HSA1 MSA2*HSA2) - (MSA1 +MSA2) *HMA ) > TIMEST

I F ( QLOAD . GT . 0 . ) THEN

ELOADH > ELOAOH QLOAD

ELSE

ELOADC > ELOADC QLOAD

ENDIF

Reheat/Recool Energy

EREH = EACTH - ELOADH ELOCH

EREC = EACTC - ELOADC ELOCC

X

117
3016

3017

3016

3019

3020

3021

3022

3023

3024

3025

3026

3027

3028

3029

3030

3031

3032

3033

3034

3033

3036

3037

3038

3039

3040

3041

3042

3043

3044

3043

3046

3047

3048

3049

3030

3031

3032

3033

3034

3033

3036

3037

3038

3039

3060

3061

3062

3063

3064

3063

3066

3067

3068

3069

3070

3071

3072

3073

CALL SHOTS

RETURN

END

Ciiii:iiiiiiiiiiiiiii:33iii3iiiiiiiiii>ii>i333ii:iii9ii3tin:3i3::::::::

SUBROUTINE REPORT

C>II>3S11BI1IX1 9I3333333I38333I33933I3333333333333III33133I33333333S3333

C This rout

no Is used to output The fine l energy usege numbers aT The

C conclusion of a tost.

C The output is current senT To The console buT should ovenTually be

C Transm iTTed To a file for furTher analysis by an addTlonal program.

C

C I npuT var tables: none

C Output variables: none

C Procedures called: none Common blocks: PFI LES, STATE 1

, UNI TS

INTEGER CLU, FLU, ILU, PLU

COMMON / UNITS / CLU, FLU, I LU. PLU

CHARACTER* 1 3 HVF I LE , USF I LE , CLF I LE , SHF 1 LE , I NF I LE , SEASON

COMMON / PFI LES / HVFI LE, USFI LE, CLFI LE, SHFI LE. I NFI LE, SEASON

REAL A

LOGICAL L

COMMON / STATE1 / A( 1 00) . L( 1 00)

Loca l Var tables:

EECOC

EECOH

EFANS

ELOADC

ELOADH

EREC

EREH

EREQC

EREQH

RHAMAX

RHAMIN

SECOUT

T2MAX

TZMIN

COOLING

HEATING

ELECTRI

COOLING

HEAT I NG

COOLING

HEATING

COOLING

HEATING

COMFORT

COMFORT

COMFORT

COMFORT

COMFORT

ENERGY (economizer) (BTu)

ENERGY (economizer) (Btu)

CAL ENERGY To fans (kWh)

ENERGY (load) ( BTu

ENERGY (load) ( Btu

ENERGY (recool) (Btu)

ENERGY (reheat) ( Btu)

ENERGY (requirements) (BTu)

ENERGY (requirements) (Btu)

relaTlve humidity at maximum Temperature (occupied)

relative humidity at minimum Temperature (occupied)

seconds outside of comfort range

meximum dry bulb temperature during occupied period

minimum dry bulb temperature during occupied period

(X)

(X)

(F) . (F) .

C

C

EQUIVALENCE ( A (64), EREQH ) EQUIVALENCE ( A (66), EECOH ) EQUIVALENCE ( A (69), ELOADC ) EQUIVALENCE ( A (90). SECOUT ) EQUIVALENCE ( A (93), RHAMAX )

( A ( 63 ) , EREQC ) ,

( A(67), EECOC ),(

( A ( 70 ) , EREH ) ,

, ( A ( 9 1 )

, TZMAX )

,

( A ( 94 ), RHAM I N )

A (72), EFANS ) A (68), ELOADH ) A( 71), EREC )

, ( A(92), TZMIN )

WRITE (CLU, 7) HVF I LE

WRITE (CLU, 8) USF I LE

WR I TE ( CLU , 9 ) CLFI LE

WR I TE ( CLU ,10) SHF I LE

WR I TE ( CLU ,11) INFILE

WR I TE ( CLU ,12) SEASON

WRITE ( CLU, 1

118
3074

3073

3076

3077

3078

3079

3000

3001

3002

3003

3004

3003

3006

3007

3000

3009

3090

3091

3092

3093

3094

3093

3096

3097

3096

3099

3100

3101

3102

3103

3104

3103

3106

3107

3100

3109

3110

31 1

31 12

3113

3114

3113

3116

3117

3110

31 19

3120

3121

3122

3123

3124

3123

3126

3127

3120

3129

3130

3131

C

C

c«

WRITE ( CLU, 2

WRITE ( CLU, 1

*1 WRITE (CLU.3)

WRITE ( CLU , 1

WRITE (CLU, 4)

WRITE (CLU, 1)

WRITE (CLU, 3)

WRITE (CLU,1)

WRITE (CLU, 6)

EFANS

EECOH, -EECOC

ELOAOH, -ELOADC

EREH, -EREC

WRITE (CLU, 13) SEC0UT/60.

WRITE (CLU, 14) TZMAX , RHAMAX

WRITE (CLU, 15) TZMI N, RHAMI

RETURN

1 FORMAT ( / , IX, 18X, 4! '- '),• + ' ,2(16! ' -

' ),'') ) 2 FORMAT( IX, ’FAN ENERGY I

1 .012. 5, ' KWH I

'

3 FORMAT (IX, 22X ,

' I HEATING I COOLING I')

4 FORMAT! IX, 'ECONOMIZER ENERGY I ',2(012.5, ' BTU I

' )

3 FORMAT! IX, 'LOAD ENERGY I*. 2(012. 5, ' BTU I

' )

6 FORMAT! IX, ‘REHEAT/RECOOL ENERGY |•,2(G12.3, , BTU I ' )

7 FORMAT! IX, 'HVAC PARAMETERS FROM: ', A1 5)

0 FORMAT! IX, 'USE PARAMETERS FR0M:',A15)

9 FORMAT (IX, ' CL I MATE PARAMETERS FR0M:',A13)

10 FORMAT! IX, 'SHELL PARAMETERS FR0M:',A15)

11 FORMATdX, ' INITIALIZING DATA FR0M:',A15)

12 FORMAT! IX, 'WEATHER SEASON USED:',A13)

13 FORMAT!/, IX, 'MINUTES OUTSIDE COMFORT ENVELOPE = '.F10.2)

14 FORMAT! IX, 'MAXIMUM ZONE TEMPERATURE * ',F10.2, ' (

' , F3 . 1

,

' XRH )'

13 FORMAT! IX, 'MINIMUM ZONE TEMPERATURE = '.F10.2, ' (

' , F5 . 1

,

' XRH )’

END

SUBROUTINE SHOTS

Cisii:3:i:ii33i3:ii33:::3i:3i3ti3t:i:i:i>33:::::::3:i3:3i333:::3::::::::

C This is a special subroutine used for development of the emulator model

C which will eventually be eliminated. Its purpose is to allow snapshots

C of the state vector to be taken and stored in a file for analysis.

C Instructions on which state vector variables are to be shot and what

C Interval the snapshots are to be taken at, are read from the file

C assigned to logical unit 7. The snapshots are stored in the file

C assigned to logical unit 10.

C

I NTEGER DAY , HOUR , M I NUTE , SECOND , ENDDAY , ENDHR , ENDM I N , ENOSEC

COMMON / TIME / DAY, HOUR, MI NUTE, SECOND, ENDDAY, ENDHR , ENDM I N, ENDSEC

LOGICAL L

COMMON / STATE1 / A! 1 00) , L( 1 00)

C

LOGICAL NO! NIT

INTEGER I VECP! 21

DATA NOINIT /.TRUE./

I F ( NO I N I T ) THEN

READ(7, *,ERR=900,END=900) REP 1 NV, NCHAN . ( IVECP<

), 1=1 , NCHAN

I F ( NCHAN . GT . 21

& STOP 'TOO MANY STATE VECTOR POSITIONS IN SNAPSHOT 1 GO TO 1000

900 WRITE! 1, FMT= 1 ( IX, "REPORT INTERVAL IN T I MESTEPS? '* )

' ) READ! 1

, * ) REP I NV

119
3132

3133

3134

3133

3136

3137

3136

3139

3140

3141

3142

3143

3144

3143

3146

3147

3148

3149

3150

3131

3132

3153

3134

3133

3156

3157

3136

3139

3160

3161

3162

3163

3164

3163

3166

3167

3168

3169

3170

3171

3172

3173

3174

3173

3176

3177

3176

3179

3180

3181

3182

3183

3184

31 83

3186

3187

3188

3189

ooooonoooooooo

o

o

oooooooo

WRITE! 1

, FMT= *( IX, “NUMBER OF VECTOR POSITIONS TO REPORT?")*)

READ! 1

, * ) NCHAN

WRITE! 1 ,FMT= '( IX. “POSITION INDICES?") *) READ! 1 j * ) ! I VECP (l) 4 I*1, NCHAN

1000 REPLEF * 1 NOINIT a .FALSE.

END IF

REPLEF = REPLEF - 1 IF! REPLEF. LE.O) THEN

RTIME » FLOAT! DAY- 1 >*24. HOUR + & FLOAT ! M I NUTE ) /60 . + FLOAT! SECOND)/ 3600.

WRITE! 10, 1 ) RTIME, !A! IVECP! I ) ) , 1*1, NCHAN)

1 FORMAT! 1X,Q1 1 . 5, 6G1 1 .3: ,/, ! 1X.7G1 1 .3)

REPLEF * REP1NV

ENOIF

RETURN

END

C Version 1.3 - FEBRUARY 4, 1983 - W.B. MAY , NATIONAL BUREAU OF STANDARDS

C<itxxssi<sxi3iixiift:333si:iiiiiiiiit:iasss3:3333isxiixiEic:::3::::::::

FUNCTION CPA I R!W)

II3iai333IIlSIIlIXI333II333l3II3ltl>33Ila33:3333ll3I33BIII3::3I3:3:33:

This function returns the specific heet of moist sir ss a function of

the humidity ratio.

UNITS: specific heat Is returned in Btu/lbm F

REFERENCE: ASHRAE Handbook, of Fundamentals 1981, page 3.3

CPA I R 0.24 0.444 * W

RETURN

END

aaanaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaanaaaaaaaaaaaaaaaaaaaarraararaa

FUNCTION RHOAI R! W, TF) aaaaaaaaaaaaaaaaaaaaaaaaazaaaaaaaaxasaaaaaaaaaaaaaaaaaasaasaaaaaaarraaa

This function returns the density of moist air as a function of

the humidity ratio and the dry bulb temperature.

UNITS: density Is returned In lbm/ft*3

temperature TF must be In degrees F REFERENCE: ASHRAE Handbook of Fundamentals 1981, page 3.3

PRESSR is atmospheric pressure in lbf/ft*2 (assumed standard atmosphere)

RAIR Is the gas constant for dry air In ft-lbf/lbm R

TA Is the absolute temperature of the air In degrees R FTOR is the conversion factor from F to R degrees

PARAMETER (PRESSR * 2116.2, RAIR * 33.332, FTOR = 439.67)

C

TA * TF FTOR

RHOAI R a PRESSR / (RAIR « TA * (1 - 1.6078 * W) ) RETURN

END

C>iBsaBsissB***is«i3aaas33S33sia3B3ii&3si3aass33:3iias323aisasss33Ss:3

FUNCTION ENTHAL ( TDB f WPPDA , PA

CniBsaaasaaaaaaaaaaaaaasaaaaaaaaeaaaaaaaaaassaassaaaaaaaaaasssassss::

C from ASHRAE 1981 Fundamentals handbook, Page 3.4

120
3190

3191

3192

3193

31 94

3199

3196

3197

3198

3199

3200

3201

3202

3203

3204

3209

3206

3207

3208

3209

3210

321 1 3212

3213

3214

3219

3216

3217

3218

3219

3220

3221

3222

3223

3224

3229

3226

3227

3228

3229

3230

3231

3232

3233

3234

3239

3236

3237

3238

3239

3240

3241

3242

3243

3244

3249

3246

3247

C

C TDB average dry bulb temperature (n degrees F

C PA atmospheric pressure in inches of mercury

C WPPOA « humidity ratio in pounds of water per pound of dry air

C ENTHAL » enthalpy of mixture in BTU per pound of dry air

C

ENTHAL » 0.240 * TDB WPPDA « (1061. + 0.444*TDB)

RETURN

END

(it

FUNCTION TEMPI ENTHAL, WPPOA, PA)

C Based on equation for enthalpy In ASHRAE 1981 Fundamentals

C handbook.. Page 9.4, solved for TEMP.

C

C TEMP » average dry bulb temperature in degrees F

C PA a atmospheric pressure in inches of mercury

C WPPDA = humidity ratio in pounds of water per pound of dry air

C ENTHAL a enthalpy of mixture in BTU per pound of dry air

C

TEMP a (ENTHAL - 1061.*WPP0A) / (0.240 0.444* WPPOA)

RETURN

END

.Ciiii3i:ti9ii9ii3ii(i3i3i33:9ii33i3ia33i:3i33i:3:i3i3i:33iti::3i::::::

FUNCTION DEWPT( WPPDA, PA)

CltII3IBB33IIIIII3l3i:3:iIt3l3l3lll3,l3l33l33133:33I33:iItl3:::3:::::3

C from ASHRAE 1981 Fundamentals handbook, Page 5.4, equation 40b

C equation Is valid in the range of TDB from 32 to 150 F. C

C PA atmospheric pressure in Inches of mercury

C WPPDA a hum i d i ty rat i o in pounds of water per pound of dry a i

C DEWPT a dewpoint temperature of mixture in degrees F. C PW = water vapor partial pressure for ;he mixture in inches Hg.

C

PW a PA*WPPDA / (0.62198 WPPDA)

ALPHA a ALOQ(PW)

DEWPT a 79.047 30. 5790* ALPHA + 1 . 8893*ALPHA*ALPHA

RETURN

END

Ciiii>3i>il3iiiaii3ii!:3i3ii:3i33>3iii3::i3:::ia::3ti::aiiii::::::::::

FUNCTION HUMRAT ( TDB , RH , PA

3331

C from ASHRAE 1981 Fundamentals handbook. Page 5.4

C

C TDB » average dry bulb temperature in degrees F

C PA > atmospheric pressure In inches of mercury

C RH > relative humidity in percent

C HUMRAT a humidity ratio in pounds of water per pound of dry air

C

PASCAL a PA a 3.38638E+3

TDBK a ( ( (TDB-32)/1 . 8)+273. 15)

PWAT a RH a PWSAT(TDBK) / 100.

HUMRAT a 0.62198 * ( PWAT / (PASCAL - PWAT) ) RETURN

END

C33133333333113I3333I333331333S33III3333:3I333:3::::i33:iii::::>:::::3

FUNCTION RELHUM(TDB, WAIR, PA)

Cllll331333l3:t33133i:33S13333a3:3ll3I3331I33aii:i33IIII3I l l:ill:13:il

121
3246

3249

3250

3251

3252

3253

3254

3255

3256

3257

3256

3259

3260

3261

3262

3263

3264

3265

3266

3267

3266

3269

3270

3271

3272

3273

3274

3275

3276

3277

3276

3279

3280

3281

3282

3283

3284

3285

3286

3287

3288

3289

3290

3291

3292

3293

3294

3295

3296

3297

3298

3299

3300

3301

3302

3303

3304

3305

C from ASHRAE 1981 Fundamentals handbook, Page 5.4

C

C TDB » average dry bulb temperature in degrees F

C PA » atmospheric pressure in inches of mercury

C RELHUM a relative humidity in percent

C WAIR a humidity ratio in pounds of water per pound of dry air

C

PASCAL a PA a 3.38638E+3

TDBK a (( (TDB-32)/1 . 81+273. 15)

PWAT a PASCAL * WAIR / (.62198 WAIR)

RELHUM a PWAT / PWSAT ( TDBK ) a 100.

RETURN

ENO

FUNCTION PWSAT ( T

Ciiti3xzti:iiiii<ixii33:ii3aiti3axii3i3<3:3:3333:s:is3::x33ts::3::::::

C from ASHRAE 1981 Fundamentals handbook, Page 5.2

C

C TDB = dry bulb temperature In degrees K. C PWSAT a saturation pressure of water vapor in pascals

C

DATA Cl / -5674 . 5359 /, C2/6. 3925247 /, C3/-0 . 9677843E-2 /

DATA C4/ 0. 62215701E-6/,C5/0. 20747825E-8/, C6/-0 . 9484024E- 1 2/

DATA C7/ 4.1635019 /, C8/-5800 . 2206 /, C9/ 1.3914993 /

DATA CIO/ -0.04860239 /, Cl 1 /0. 41 764768E-4/, Cl 2/-0. 14452093E-7/

DATA Cl 3/6. 5459673/

T2»T*T

T3*T2aT

T4»T3aT

IF(T.0T, 273. 15)00 TO 100

PWSLN « Cl /T+C2+C3*T+C4»T2+C5*T3+C6aT4+C7*AL0Q( T)

GO TO 200

100 PWSLN » C8/T+C9+C1 0*T+C1 1 *T2+C 12»T3+C1 3« ALO0 ( T)

200 PWSAT= EXP (PWSLN)

RETURN

END

Ciisss=sssitiisisiiiB33SSBsaa9sisi8«asass3sa333s:33ass:iiiss::333:=r:=

FUNCTION TDBSAT ( HSAT, PA

CaXII3XlS3XiaXXI3S3Xa33aai333IIIXIXXI3333X333333333Z33>X3ll33333333333

C This function uses a table look-up using data from ASHRAE 1981

C Fundamentals handbook. Page 6.5, table 2. C function Is valid in the range of TDB from 32 to 100 F. C

C TDBSAT dry bulb temperature in dogrees F at saturated conditions

C PA atmospheric pressure in Inches of mercury assummed to be 29.921

C HSAT « enthalpy of saturated mixture in BTU per pound of dry air

C

REAL

DATA

&

&

&

&

&

&5

6

DATA

H( 69 ), T ( 69

H / 11.760,12.170,12.587,13.010,13.441,13.878,14.322,14.773,

15.233, 15.700, 16.175, 16.660, 17.152, 17.653, 18.164,18 685,

19.215, 19.756, 20.306, 20. 868, 21.441, 22.025,22.621 , 23.229,

23.850, 24.484,25. 131 , 25.792, 26.467,27. 157,27. 862, 28. 582,

29. 318, 30. 071 , 30. 840, 31 . 626, 32.431, 33.254, 34 . 097, 34.959,

35.841, 36.743, 37.668, 38.615,39. 583,40. 576,41 . 592,42. 633,

43.701,44. 794, 45. 91 3, 47. 062, 48. 238,49 445,50. 681,51. 949,

53.250,54 . 582, 55.951,57. 355, 58. 794,60.271 , 61 . 787, 63 343,

64 . 940,66.578,68.260, 69.988. 71.761/

T / 32.0, 33.0, 34.0,. 35.0, 36.0, 37.0, 38.0, 39.0,

i nhg

122
3306

3307

3300

3309

3310

331 1 3312

3313

3314

3313

3316

3317

3310

3319

3320

3321

3322

3323

3324

3323

3326

3327

3329

3329

3330

3331

3332

3333

3334

3333

3336

3337

3336

3339

3340

3341

3342

3343

3344

3343

3346

3347

3346

3349

3330

3331

3332

3333

3354

3335

3336

3337

3330

3339

3360

oooo

oooo

oooo

40.0. 41 . 0. 42.0, 43.0, 44.0, 45.0, 46.0, 47.0,

49.0. 49.0, 30.0, 51 .0, 52. 0, 33.0, 54.0, 55.0,

56.0. 37.0, 56.0, 39.0, 60.0, 61 . 0, 62.0, 63.0,

64.0. 63.0, 66.0, 67.0, 68.0, 69.0, 70.0, 71 . 0,

72.0, 73.0, 74.0, 73.0, 76.0, 77.0, 78.0, 79.0,

00.0. 01 .0, 82.0, 63. 0, 84.0, 85.0, 86.0, 87.0,

00.0. 09.0, 90.0, 91 .0, 92.0, 93.0, 94.0, 95.0,

96.0. 97.0, 98.0, 99.0, 100.0/

TDBSAT

RETURN

END

TERPOL ( T , H , HSAT , 69

FUNCTION TERPOL (Y, X, XI , LTAB)

C This Is a general function to Interpolate a value from a table.

C XI Is the independent value < X Is a list of Independent values

C Y Is a list of dependent values. LTAB Is the length of the table.

C

REAL Y ( * ), X( *

Locate position In table

DO 1000 I > 1 , LTAB

IF( CXI .0T.X( 1 ) .AND. XI . LT. X( I +1 ) ) . OR.

& (XI .LT.XC I ) . AND. XI .OT.XC 1+1 ) ) ) 00 TO 2000

1000 CONTINUE

independent value supplied is not In table

IF(XCLTAB) .GT.XC1 )) THEN

IFCXI . LT.XC 1 ) ) THEN

TERPOL = Y ( 1 ) ELSE IFCXI .OT.XCLTAB) ) THEN

TERPOL = YCLTAB)

END IF

ELSE IFCXC1 ) .OT.XCLTAB) ) THEN

IFCXI 0T.XC1 )) THEN

TERPOL » Y ( 1

ELSE IFCXI .LT.XCLTAB)) THEN

TERPOL » YCLTAB)

END IF

ELSE

TERPOL * 0.0

END IF

RETURN

Interpolate

2000 TERPOL a (XI-XCD) / (X( I +1 ) -X( I ) ) a (Y( I +1 ) -Y( I ) ) Y(!)

RETURN

END

123
NBS-114A (rev. 2-ao) U.S. OEPT. OF COMM. 1. PUBLICATION OR 2. Performing Organ. Report No. 3. BIBLIOGRAPHIC DATA

REPORT NO.

SHEET (See instructions

)

NBSIR 85-3291 JANUARY 1986

4. TITLE AND SUBTITLE

BUILDING EMULATION COMPUTER PROGRAM FOR TESTING OF ENERGY MANAGEMENT AND

CONTROL SYSTEM ALGORITHMS

5. AUTHOR(S)

William B. May, Jr. and Cheol Park

6. PERFORMING ORGANIZATION (If joint or other than NBS, see instructions)

national bureau of standards

DEPARTMENT OF COMMERCE

WASHINGTON, D.C. 20234

7. Contracv'Grant No.

I. Type of Report & Period Covered

9.

SPONSORING ORGANIZATION NAME AND COMPLETE ADDRESS ( Street. City. State, ZIP)

Naval Civil Engineering Laboratory

U.S. Department of Defense

Port Hueneme, California 93043

Office of Building and Community Systems

U.S. Department of Energy

1000 Independence Avenue, SW

Washington, DC 20585

10.

SUPPLEMENTARY NOTES

| | Document describes a computer program; SF-I8S, FlPS Software Summary, is attached.

11.

ABSTRACT (A 200-word or less factual summary of most significant information. If document includes a significant

bibliography or literature survey, mention it here)

A building emulator can be used to test energy management and control systems (EMCS)

The emulator uses a computer program to simulate the responses of a building

including the equipment, building space, and building envelope to EMCS commands.

Building model software for the emulator has been developed at the National Bureau

of Standards (NBS) in an effort to assist the United States Naval Civil Engineering

Laboratory (NCEL) , which is developing a sophisticated building emulator.

The concept of the building emulator and the building emulator computer program

are described in this report. The program includes the weather, the air handling

unit, the zone, and the comfort model. In addition, the energy compilation routine

is also included. The models presented here are simplified models. With these

abridged models, a single zone building with exterior walls and a single deck air

handling unit are simulated. A complete FORTRAN source code of the building

emulator computer program is appended.

12.

KEY WORDS (Six to twelve entries; alphabetical order; capitalize only proper names; and separate key words by semicolon s) air handling unit simulation model; building emulator; building space zone model;

computer simulation; EMCS algorithms; energy management and control systems (EMCS);

HVAC emulator; local equipment simulation model; testing EMCS algorithms.

13. AVAILABILITY

1 Uni imi ted

| | For Official Distribution. Do Not Release to NTIS

i [ | Order From Superintendent of Documents, U.S. Government Printing Office, Washington, D.C.

20402.

IX 1 Order From National Technical Information Service (NTIS), Springfield, VA. 22161

14. NO. OF

PRINTED PAGES

132

15. Price

$ 16.95

USCOMM-OC 9043-PS0