Environrrwnt at Protection
Environmental Research
Laboratory
Corvallis OR 97330
EPA-600/3-78-074
July 1978
Research and Development
6EFW
User Guide for the
Hydrodynamical-
Numerical Model
-------
REPORTING SERSES
Research reports of the Office of Research and Development, U.S. Environmental
Protection Agency, have been grouped into nine series These nine broad cate-
gories were established to facilitate further development and application of en-
vironmental technology. Elimination of traditional grouping was consciously
planned to foster technology transfer and a maximum interface in related fields.
The nine series are
1, Environmental Health Effects Research
2. Environmental Protection Technology
3 Ecological Research
4 Environmental Monitoring
5 Socioeconomic Environmental Studies
6, Scientific and Technical Assessment Reports (STAR)
7. Interagency Energy-Environment Research and Development
8. "Special" Reports
9 Miscellaneous Reports
This report has been assigned to the ECOLOGICAL RESEARCH series This series
describes research on the effects of pollution on humans, plant and animal spe-
cies, and materials. Problems are assessed for their long- and short-term influ-
ences. Investigations include formation, transport, and pathway studies to deter-
mine the fate of pollutants and their effects This work provides the technical basis
for setting standards to minimize undesirable changes m living organisms in the
aquatic, terrestrial, and atmospheric environments.
I his document is available to the public through the National I echnical Informa-
tion Service, Springfield, Virginia 22161.
-------
EPA-600/3-78-074
July 1978
USER GUIDE FOR THE
ENHANCED HYDRODYNAMICAL-NUMERICAL MODEL
by
A. D. Stroud
and
R. A. Bauer
Compass Systems, Inc.
San Diego, California 92109
Contract Number 68-03-2225
Project Officer
R. J. Callaway
Marine and Freshwater Ecology Branch
Corvallis Environmental Research Laboratory
Corvallis, Oregon 97330
Corvallis Environmental Research Laboratory
Office of Research and Development
U. S. Environmental Protection Agency
Corvallis, Oregon 97330
-------
DISCLAIMER
This report has been reviewed by the Corvallis Environmental Research
Laboratory, U. S. Environmental Protection Agency, and approved for pub-
lication. Approval does not signify that the contents necessarily reflect
the views and policies of the U. S. Environmental Protection Agency, nor
does mention of trade names or commercial products constitute endorsement
or recommendation for use.
11
-------
FOREWORD
Effective regulatory and enforcement actions by the Environmental Protection
Agency would be virtually impossible without sound scientific data on pollu-
tants and their impact on environmental stability and human health. Respon-
sibility for building this data base has been assigned to EPA's Office of
Research and Development and its 15 major field installations, one of which
is the Corvallis Environmental Research Laboratory (CERL).
The primary mission of the Corvallis Laboratory is research on the effects
of environmental pollutants on terrestrial, freshwater, and marine ecosystems;
the behavior, effects and control of pollutants in lake systems; and the de-
velopment of predictive models on the movement of pollutants in the biosphere.
A. F. Bartsch
Director, CERL
111
-------
PREFACE
This Users Guide to the Optimized Hydrodynamical-Numerical (HN) Model
system reported by Bauer and Stroud (1) details the information needed to use
the model on CDC equipment. The reported efforts focused on a set of spe-
cific objectives and as a result the Users Guide and program library that
evolved exclude the display program and other programs that are logically
related to the HN model but that were not needed.
HN model work at Fleet Numerical Weather Central (FNWC), Monterey, was
begun in 1966 by Dr. T. Laevastu and was continued at various levels of effort
into this project. In 1971 the model work was shifted to the Naval Environ-
mental Prediction Research Facility (NEPRF), when an Oceanographic Department
was established under Dr. Laevastu. During the early years of the use of the
HN model at FNWC and NEPRF, each new application resulted in the creation of
a new program deck and often each new application also required features in
the model that were unique to the problem under investigation. As a result
by 1973 the number of separate models and special features became unmanageable
in terms of both man time and computing resources and a major effort was
undertaken to create a single program with selectable options that^ould be
used on any geographic configuration. This effort resulted in the Optimized
Model. Unfortunately the timing of this effort caused the programs related to
the HN model to be distributed between the NEPRF CDC 3100, where quick
turn around was possible, the Lawrence Berkeley Laboratory (LBL) CDC 7600, where
cost effective calculations could be done, and FNWC's CDC 1604, where existing
routines could be used without program conversion. The major elements of the
iv
-------
optimized system referred to in Bauer (2) as Phases I, II and III were sent to
Mr. Callaway at Environmental Protection Agency (EPA) Corvallis in 1975, where
Phase I, a CDC 3100 program, is run on Oregon State's CDC 3300 computer and
Phase II is run on the Bonneville Power Authority (BPA) CDC 6500. The original
CDC 3100 program was converted to run on BPA's 6500 and extended to contour
water heights by EPA. Since this program has concentrated on efforts of specif-
ic calculations rather than the description of a particular area, Phase III was
not used and as a result is not included in this documentation.
Other major related programs that have not been included are the Varian
contour plot, velocity vector plot program written for the CDC 3100 by J.
Harding at NEPRF (3), the wave program written by Bauer (4) for the CDC 6500
that with minor modification could be used to generate the HN bathymetric
data from contour charts, and finally the AT model sequence (Bauer et. al.
[5]) that can be used to produce three hourly wind fields using the FNWC
hemispheric fields to bound the local area grid solution forecasts using
Bengtsson's limited area 3 parameter model.
The capabilities of the HN model have been considerably enhanced by the
effort to develop and test specific modifications to the basic optimized model
and can be considered operational. The goal of having a single documented
model with many optional computations that can be used in any area has been
realized. Additional effort is required to convert the related programs and
merge them into the program library; however, with the disestablishment of
the NEPRF Oceanographic Department and the reorganization of the EPA Corvallis
Laboratory that occurred during the project, it is not clear what facilities
and agencies would be involved in future efforts. What is clear is that the
level of understanding and development of this model evolved only after many
man years of trial and error, so that we strongly recommmend that anyone who
-------
has never used the model arrange to have help in setting up his first model
runs. We also strongly recommend against converting the program to other
computer hardware unless there are long term commitments to use the model,
because the conversion would automatically lose the flexibility and structure
of the present program library and because the conversion and debugging
effort required would certainly be many times the cost of producing the re-
quired results with the present set of programs.
-------
ABSTRACT
This guide provides the documentation required for use of the Enhanced
Hydrodynamical-Numerical Model on operational problems. The enhanced model
is a multilayer Hansen type model extended to handle near-shore processes
by including:
•Non-linear term extension to facilitate small-mesh studies of
near-shore, including river inflow dynamics;
•Layer disappearance extension to enable appropriate procedures in
tidal flat and marshy regions, as well as some down/upwelling cases;
•Thermal advection enhancement for treatment of thermal pollution
cases by method of moments coupled with heat budget procedures;
•Monte Carlo diffusion enhancement to deal with dispersion via
statistical methods and comparison to the method of moments experi-
ments .
The guide includes a description of the model system, source code main-
tenance procedures, notes on implementation, functional descriptions of rou-
tines and option sets, listings and example runs for a test data set.
This report was submitted in fulfillment of contract No. 68-03-2225 by
Compass Systems, Inc. under the sponsorship of the U. S. Environmental Pro-
tection Agency. This work covers the period from June 1975 to March 1977
and was completed March 31, 1977.
VII
-------
CONTENTS
Foreword iii
Preface iv
Abstract vii
Figures x
Tables x
Acknowledgment xi
1. Introduction 1
2. Recommendations and Conclusions 2
3. The Enhanced HN Model Program System 4
FORTRAN program structure 4
The library maintenance system 7
Program library structure and usage 9
4. HN Model Implementation Notes 17
Implementation process 17
Phase II timing notes 20
Notes on experimental/proposed features/methods. . .22
References 43
Appendices
A. Phase I Program Notes 45
Description of Phase I program 45
Listing for Phase I sample run 52
B. Phase II Program Notes 68
Description of program LAYER3S 68
Listing for Phase II sample run 77
C. Non-Linear Modifications 93
CONVCT and CONVCU series 93
CONPRT, CONPRU and CONPRV series 94
Listing for non-linear correction sets 95
D. Layer Disappearance Modifications 98
FLAT, FLATA and FLATB series 98
FLATIO and FLATIP series 101
Listing for layer disappearance correction sets. . 103
E. Monte Carlo Modification
MONTE 110
Listing for Monte Carlo correction sets 112
F. Phase IV Program Series 116
Thermal advection program 116
Listing for Phase IV sample run. 122
IX
-------
FIGURES
Number Pa9e
1 The HN model system of programs 5
2 Intermediate data set format 6
3 HN Phase I program structure 10
4 HN Phase II basic program structure 12
5 HN Phase II enhancement impacts 13
6 Phase IV program structure • 15
TABLES
Number
1 Master Control Card Set 25
2 HN Phase I Control/Data Cards 28
3 UPDATE Directives Examples 31
4 CDC UPDATE Control Card Parameter Examples 31
5 HN Library: Phase I Modules 32
6 HN Library: Phase II Modules 33
7 HN Library: Phase IV Modules 36
8 HN Library Usage Job Control Deck Setup for CDC6500 38
9 HN Phase I Control/Data Card Set 39
10 HN Master Control Card Set (All Phases) 40
11 HN Library Usage By CDC UPDATE Directives Set 41
12 HN Phase II Timing Performance Chart 42
D-l Tidal Flat Code Conventions 100
-------
ACKNOWLEDGMENT
We wish to express our appreciation to R. Callaway, Project Officer,
for providing the proper environment to perform theoretical studies of the
model and wish to note the thoughtful user oriented feedback from C.
Koblinski, also of the Corvallis laboratory.
We wish to thank Dr. T. Laevastu for his guidance and to acknowledge all
of his innovative enhancements to the original Hansen model, which formed
the basis for the optimized model.
We wish to thank the Naval Environmental Prediction Research Facility
for allowing use of their facilities; S. Larson, K. Rabe and J. Harding for
their many constructive suggestions; and the NEPRF computer support personnel
for their assistance.
We wish to acknowledge the use of the ERDA computing facilities at
Lawrence Berkeley Laboratory, University of California at Berkeley, California.
Finally, we wish to express our thanks to S. Haines for preparing the
illustrations and tables and to D. Bauer for typing the report.
XI
-------
SECTION 1
INTRODUCTION
The Hansen-type optimized multi-layer Hydrodynamical-Numerical (HN)
model is a numerical difference solution of the vertically integrated form of
the hydrodynamic equations. Analytic details are presented in the enhanced
model described by Bauer and Stroud (1). Included in this document are de-
scriptions of the Data Initialization (Phase I), Hydrodynamic Computational
(Phase II) and Auxiliary Computational (Phase IV) HN programs; the program
library structure, usage and maintenance; and notes on implementation of the
model. Program notes, listings, sample runs and module listings are included
as appendices.
-------
SECTION 2
RECOMMENDATIONS AND CONCLUSIONS
The HN Model program designed to permit options, boundary specification,
grid size and number of layers to be specified by control card's is operational,
The model consists of four' separate programs, Phases I, II, III and IV. All
are now coded to run on a CDC 6000-7000 series computer. A single program '••'•'
library documented in this Users Guide includes Phases I, II and IV. Phase I
produces the initial, geographic data set. Phase II, the main computational
routine, solves the hydrodynamic equation''set whiclvcan include as'options
the non-linear convective term and the treatment of"areas where a layer dis-
appears (tidal flat). As an option it also advects a set of points using-
Monte Carlo techniques. Phase IV uses the HN velocity fields as input to
the Pedersen-Prahm advection method and decays the heat content of the ad-
vected volume.
Both Phase I and II have been tested on both synthetic data developed to
test boundary and layer disappearance cases and on real data cases. Results
on all cases attempted have been satisfactory. Results from Phase IV reflect
the assumption of.pure advection with no diffusion. In a real data test the
plume developed by Phase IV was more compact and hotter downstream from the
source than was observed by IR overflights.
From the experience gained using the HN model it is clear that valuable
information can be obtained about the character of an area and that model
-------
results can aid both in the design of field experiments and in interpretation
of field data. With the present model initial results for most areas can be
obtained with one or two days of effort. While the thermal advection part
of the model did not verify as well as hoped because of the assumption of no
diffusion, this single example of the use of the HN derived velocity field
should not discourage others who need to know spatial and temporal variation
in the flow fields to design and interpret measurement or to derive input
for other models.
-------
SECTION 3
THE ENHANCED HN MODEL PROGRAM SYSTEM
The enhanced HN model is composed of a series of FORTRAN programs main-
tained on a program library which is accessed through a library maintenance
system. The structure of the program system, of the library system and of
the resulting program library are presented. More detailed documentation of
the individual programs, sample runs and enhancements are included as appen-
dices.
FORTRAN PROGRAM STRUCTURE
The structure of the enhanced HN model as shown in Figure 1 is composed
of four separate programs linked by the master control card and intermediate
data sets. The four phases are: Phase I~Data Initialization; Phase II-Hydro-
dynamic Computations; Phase Ill-Data Display and Phase IV-Auxiliary Computa-
tions. The design of the sequence of programs allows the user to verify, save
and reuse key intermediate results and minimizes the computer resources re-
quired during the computational phases.
The Phase I data initialization program is designed to check the master
control card set defined in TABLE 1 and to create the intermediate data set
with type 1 and type 2 time 0 records as defined by Figure 2 based on the water
depths at U and Vpoints and the Phase I data set, defined in TABLE 2.
The Phase II hydrodynamic computation program is usually the most ex-
pensive program to run in the sequence and, since it must be run out for at
least a half tide cycle before the results become usable, one of the major
goals in designing the system was to allow the Phase II program to be stopped
-------
Master Control
Card Set
Phase I
Data Initialization
Intermediate
Data Set
Phase II
Hydrodynamic
Computations
Phase IV
Auxiliary
Computat ion
Phase III
Data Display
CDC 3100
6500
Phase I
Data Set
CDC
6500
7600
Restart Loop
Intermediate\
Data ]
Set /
CDC
6500
7600
CDC
3100
6500 t
Phase III
Data Set
Conversion of Phase III from CDC 3100 to CDC 6500 was per-
formed at EPA Ccrvallis. The CDC 3100 program is document-
ed in (2) .
Figure 1. The HN model system of programs.
-------
ARRAY
CHARACTER
BYTE
HEADER
DATA FIELD
PHYSICAL RECORD
LOGICAL RECORD
GROUPED RECORDS
A 2 dimensional matrix used in the HN model system with
rows and columns indexed from the upper left hand
corner
6 bit CDC display code
24 bits integer value
First 5 bytes of a physical record
Byte containing an element of an array as a scaled ±. in-
teger value
Header +1 to 1245 data fields
1 or more physical records containing an entire array
stored in standard FORTRAN array element order
2, 3 or 6 logical records having the same headers where
the order of the logical records in the group is the
only key to the contents of the logical record
Records are written and read with CDC FORTRAN Buffer statements
in binary mode.
When the data set is on tape the tapes are unlabelled, 556 BPI, 7 track
in stranger format 1 file per model.
Header format (first 5 bytes of every physical record.) Headers on the tape
are in ascending sort order. All header bytes are + integers.
TIME
LAYER
TYPE
ROW-COL
MODEL NO
TIME
LAYER
TYPE
1
3
4-7
8
9
10
20
Time in seconds
Layer number
Symbolic HTZ, HTU and HTV arrays (must be first record
in set. Placed on the set by Phase I.
Computed £, U and V arrays. (Time 0 fields are placed on
the set by Phase I. Records with time > 0 are placed
on the set by Phase II.)
Used by FLATIP to save FLAT modified £, U and V arrays
Reserved for other groups of triple logical records
Double logical record group used to input variable winds
from the AT model
Single logical record used by FLATIO to save char-
acter array for display, length=(MAXDIM+9)/10+2
Single logical records used by CONPRV to save convective
terms (one separated set per IO step)
Hex grouped logical record used by Phase IV on a separate
data set to save results and provide restart capability
ROW-COL
MODEL NO
(Number of rows)*10000B+(number of columns)
Number used to identify data set and checked against the
MODEL number given in the master control card set
Figure 2. Intermediate data set format.
-------
and restarted easily. This feature allows the initialization period results
to be produced once and used to restart the model with different features
activated in each restart. This is accomplished by having the Phase II pro-
gram read the intermediate data set type 1 record to verify the correct data
set is being read and to set up the topography, then to continue reading the
intermediate data set until the type 2 record with a time matching the
start time specified in the master control card deck is found.
The intermediate data set in conjunction with the master control card
deck allows both Phases III and IV to produce products based on the HN compu-
tations after the hydrodynamic computational data rather than as an integral
part of it, simplifying the logic and minimizing the memory and memory swap-
ping required to obtain results.
The Phase III display program is not included in this document since it
is highly dependent on the hardware available to the user. The Phase III pro-
gram available to the authors uses the CDC 3100 and CalComp plot routines in
a highly overlaid structure required to fit within the 26,000 24 bit words of
the 3100 memory. This Phase III program has not been changed since it was re-
ported by Bauer (2). The CDC 6500 version reprogrammed and extended by EPA
at Corvallis would serve as a better program base for future efforts.
The auxiliary computations, Phase IV, use the velocity fields in the
type 2 record from the intermediate data set to drive the duplex thermal ad-
vection model. The reasons for separating the auxiliary computations are
that the advection processes can be run on longer time steps than the hydro-
dynamical model without producing noticeable differences in the results and
the thermal advection model has its own relatively large memory requirements.
THE LIBRARY MAINTENANCE SYSTEM
This section describes the general structure of the CDC library
7
-------
maintenance system, UPDATE (6), as used for the HN programs. The HN program
library is usually maintained on a magnetic tape or in permanent disk file
accessed by UPDATE. UPDATE supports the following modules: decks, common
decks and correction sets. Each module consists of a specific set of card
images and is uniquely named. Each card image is identified by reference
to the module name and a sequence number assigned to the card. As used with
the HN program library, a deck usually consists of one or more routines. A
common deck is a subset of card images that appears several times in the
program but only once in the library. A simple reference to the common deck
at the required location causes UPDATE to replace the reference card with
the common deck cards referenced. This technique greatly reduces the effort
required to maintain FORTRAN (7) COMMON blocks and other repeated sections of
FORTRAN code.
Correction sets include library maintenance directives and associated
new or replacement card images which modify existing modules. Correction
features available include text insertion before, after or in place of speci-
fied card images, deactivation or reactivation of card images or prior cor-
rection sets and permanent structural modifications (renaming, renumbering,
adding or purging modules). Most correction actions are fully recoverable
since deactivated modules or card images are not removed from the library un-
less specifically directed during new library generation. Card images on the
library are maintained in a compressed format which allows retention of cur~
rent and historical status for each image and module.
After an initial library creation, use of the system merely requires that
access be provided to the library and additional desired correction sets.
Provided the appropriate access, the library system tailors the library based
on the correction directives to form an expanded format subset of the library
-------
suitable for input to the CDC FORTRAN EXTENDED compiler, usually on file
COMPILE. Use of the library system avoids the much less efficient process
of maintaining and using actual card decks with large systems of programs.
A more detailed and precise documentation of the CDC UPDATE library system
features may be found in (6).
In order that correction sets may be presented and discussed, TABLE 3
presents the most frequently used CDC UPDATE directives, and TABLE 4 pre-
sents some frequently used parameters required for control of the CDC UP-
DATE program.
PROGRAM LIBRARY STRUCTURE AND USAGE
The HN program library is maintained on magnetic tape in a form compat-
ible with UPDATE for CDC 6000/7000 et. al. series computer systems. The HN
library currently consists of several FORTRAN programs, sub-programs and
variations as UPDATE modules. Some of the modules retained in the present
library version are the result of prior evolution rather than current devel-
opment effort. Such modules, e.g. VARWIND, VARCORL, have been retained as
they provide valuable options for future HN efforts. Other potential modules
with similar or greater value have not been specially added due to separate
evolution and lack of specific requirement. As the variety of modules present
on the HN library involves overlapping deck requirements or multiple deck im-
pact, the library has been maintained with "simplest" Phase II usage in mind.
Phase I
The library modules associated with Phase I have been presented in TABLE
5, including two decks and various correction sets. Figure 3 presents the re-
sultant program structure of Phase I as a block diagram with references to the
individual card image identification associated with each block. Appendix A
-------
PHASE I
MAIN PROGRAM
LAYER3S A
INITIALIZATION
MASTER CONTROL
PROCESS B
BATHY -DATA
OLDT,NEWT,
SETU,SETV,SETZ
LAND/SEA
BOUNDARIES
GRDZ
AUXILIARY
PROCESSES
NEWP,OLDP,
VOLM
OUTPUT
PROCESSING
TAPE
SUBROUTINES
TIDCUR
BLOCK KEY
A
B
C
D
E
F
G
H
I
J
K
DECK
HNI
HNI
HNI
HNI
HNI
HNI
HNI
HNI
UTIL
UTIL
UTIL
MODULE ID LIMITS
HNI.2,HNI.54
HNI.55,HNI.76 '
HNI.85,HNI.91
HNI.77,HNI.84
HNI.122,HNI.143
HNI.92,HNI.121
HNI.237,HNI.292
HNI.194,HNI.236
HNI.144,HNI.193
HNG5.28,UTIL.66
HNG8.127,UTIL.98
PR3.6,UTIL.63
Figure 3. HN Phase I program structure.
10
-------
provides documentation of the Phase I program and associated subprograms in-
cluding a FORTRAN source program listing with complete card image identifica-
tion as provided by UPDATE on file COMPILE and printed output examples.
Phase II
The library modules associated with Phase II have been presented in
TABLE 6, including two common decks, three decks and many correction sets
showing the logical development of the computational model and various en-
hancements. Note that the deck UTIL is associated with both Phases I and II.
Figure 4 presents the resulting program structure of the basic (simplest, en-
hanced) Phase II as a block diagram with references to the individual card
image identification associated with each block. Appendix B provides docu-
mentation of the basic Phase II program and associated subprograms including
a FORTRAN source program listing as provided by UPDATE on file COMPILE.
Phase II Enhancements
Included in TABLE 6 are the modules associated with Phase II enhancements
and the module series associated with suppressing those enhancements until
specific usage is required. Correction set names beginning with the letter A_
have been used to identify modules that exist merely to suppress features that
are not a part of the basic Phase II program. If these suppression modules
are inactivated through the UPDATE YANK directive or if specific sub-
directives of such a module are inactivated by means of the UPDATE DELETE
directive, the desired features becomes a part of the resulting FORTRAN source
deck image on the COMPILE file.
Figure 5 presents the impact of each of the enhancement module groups on
the Phase II program in a block diagram similar to Figures 3 and 4. The
hexagonal blocks identify logical module groups that should not be separated
for independent use. The one exception to this rule is that CONPRT may be
11
-------
PHASE II
MAIN PROGRAM
SUBPROGRAMS
LAYERS A
INITIALIZATION
C COMPUTATION:
TIDE WIND SORC
PRESCRIPTION
Hu,Hv
COMPUTATION
U,V
COMPUTATION:
FLOW PRESCRIP-
TION
OUTPUT
PROCESSING
D
SETHUV
BAR
XMIT
CXMIT
UVCAL
UVSTAR
PRTMAX
INPOUT
H
K
M
IIOERR
BARUV
I
BLOCK KEY
A
B
C
D
E
F
G
H
I
J
K
L
M
DECK
HNG
HNG
HNG
HNG
HNG
HNG
BARSTAR
UTIL
BARSTAR
HNG
BARSTAR
UTIL
UTIL
MODULE ID LIMITS
HNG.2,HNG.78
HNG8.26,HNG.133
HNG.134,BC8.6
BARMOD.1,HNG8.48
HNG.191,HNG.246
HNG.249,HNG8.71
HNG8.105,BARSTAR.36
HNG5.28,BC8.16
BARMOD.8,BARSTAR.50
HNG8.72,HNG8.94
HNG8.95,BARSTAR.2 3
HNG8.127,UTIL.98
PR3.6,UTIL.63
Figure 4. HN Phase IV basic program structure.
12
-------
PHASE II
MAIN PROGRAM
LAYERS A
INITIALIZATION
B
C COMP:
TIDE,WIND,SORC
Hu,Hv
C
'.
U,V:
FLOW
D
•
OUTPUT
E
AFFECTED
BLOCKS
AFFECTED
MODIFICATIONS SUBROUTINES
SUBROUTINES
CONPRT
CONPRU
CONPRV
SETHUV
F
XMIT
CXMIT
H
BARUV
(CONUV)
PRTMAX
INPOUT
M
IOERR
DISPLE N
(DSTART)
(DISPRT)
(DISPIN)
UVCAL
J
WAYTUV
O
BLOCK KEY
A-B
C
D-E
F
H-M
N
O
DECK
HNG
HNG
HNG
HNG
HNG.UTIL
BARSTAR
BARSTAR
MODULE ID LIMITS
(SAME AS FIGURE 4)
FLAT.7,MONTE.61(or,FLAT.51)
(SAME AS FIGURE 4)
FLAT.52,HNG8.71
(SAME AS FIGURE 4)
FLATIO.4,FLATIO.66
MONTE.64,MONTE.7 6
Figure 5. HN Phase II enhancement impacts.
13
-------
used with CONVCT and CONVCU without CONPRU and CONPRV. Obviously the CONPRT
module group should never be used without the CONVCT module group. Similarly
the FLATIO group should never be used without the FLAT group.
Any combination of Features 1, 2 and 4 as identified in TABLE 6 may be
used with or without their I/O options.
Appendices C, D and E provide library module listings of the Features
1, 2 and 4, respectively, with associated documentation.
Phase IV
The library modules associated with Phase IV have been presented in
TABLE 7, including five common decks, two decks and various correction sets.
Figure 6 presents the resultant program structure of Phase IV as a block di-
agram with references to the UPDATE card image identification. Appendix F
provides documentation of the Phase IV program and associated subprograms in-
cluding a FORTRAN source program listing as provided by UPDATE on file COM-
PILE. Note that routines with the same names as those in deck UTIL are in-
cluded in deck AD. Though these routines are somewhat different Phase IV
should ultimately evolve toward use of deck UTIL. That it does not do so
currently is due to the separate initial library evolution of the thermal
advection procedures.
Usage
TABLE 8 presents a basic CDC 6500 SCOPE (8) job control deck setup
for use of Phases I, II and IV including a list of the various logical unit
designations required by the job stream. TABLES 9 and 10 present the HN con-
trol/data card decks that are required for the particular HN phase desired in
the deck setup of TABLE 8. TABLE 11 presents various UPDATE directive
sets to achieve particular PIN objectives. Both the YANK directive and the
DELETE directive have been employed to remove the suppression directives
14
-------
PHASE IV
MAIN PROGRAM
ADVECT A
INITIALIZATION
B
THERMAL
TREATMENT
DUPLEX
ADVECTION
D
SOURCES;
TEMPERATURE
COMPUTATION;
WIND EFFECT
OUTPUT
PROCESSING
SUBPROGRAMS
PRTMAX
THERMAL
BDVECT
XMIT
INPOUT
IIOERR
BLOCK KEY
A
B
C
D
E
F
G
H
I
J
DECK
AD
AD
AD
AD
AD
AD
TH
AD
AD
AD
MODULE ID LIMITS
A2.1,AD.39
A3.5,A3.11
A2.47,A2.52
A2.53,A2.78
A2.79,A1.27
A1.28,A1.61
T1.9,TH.324
A2.97,A3.26
A1.62,A1.65
A1.81,A1.144
Figure 6. Phase IV program structure.
15
-------
associated with correction modules whose names begin with the letter A,
Note that to achieve a basic Phase II run, no YANK or DELETE directive is
required. The "Q" option on the UPDATE control card must be used as shown in
TABLE 8 in conjunction with the indicated COMPILE directive to access only
the required decks. The correction set idents used in TABLE 11 examples are
arbitrary.
16
-------
SECTION 4
HN MODEL IMPLEMENTATION NOTES
This section discusses the process of implementing the HN model for a
proposed application including the estimation of cost for the hydrodynamic
computational phase. Also notes on various experimental or proposed fea-
tures, methods and results are presented.
IMPLEMENTATION PROCESS
Throughout the implementation process consideration must be given to
resources available for the application, including funds, computer resources,
data resources and personnel. Success of the application will depend pri-
marily on the quality and quantity of these resources.
Preliminary Processes
Site selection is an important step that interacts with available data
resources. Features of interest for the study will play an important role in
determining the number of layers as well as determining the extent of the area
to be covered and the grid mesh length to be used. Interest should interact
with practicality here to ensure economic implementation. Bathymetric fea-
tures, tidal character, storm surge and river inflow all need to be considered
in determining how best to cover the area of interest.
Orientation of the grid with respect to land is not material in itself;
however, economical use of computer memory would point toward aligning at
least one side of the grid with the coastline. Inclusion of river estuaries
is preferable to attempting to prescribe river inflow across the estuary
mouth. Some additional notes on river inflow methods may be found in the
17
-------
proposed methods section. Note that specific area selection determines
the high tide maximum depth, and selection of the grid mesh size determines
the maximum practical time step through the Courant-Fredricks-Lewin (CFL)
stability criterion:
J,
VT > (2gHmax) .
Once the time step is known an estimate of cost of Phase II on the CDC 7600
at Berkeley can be made by referring to TABLE 12 and the notes on Phase II
timing.
The actual grid preparation is straightforward once the approximate
area has been outlined on the appropriate chart and the grid size has been
selected. Merely rule the z grid carefully on the chart; then read depth
values from the chart row by row for the hu and hv fields respectively. The
first hu value is taken one-half grid point to the right of the z(l,l) point
and thus along each row for ME values. When a u or v point overlies land,
a zero value is used. During this digitizing process it may be important to
shift or rotate features slightly in order to produce realistic cross sec-
tions in the depth fields. For this reason hand digitization is the most
efficient method of preparing the depth fields; however, programs do exist
which could be used to automate this process (4) .
Phase I
When the hu and hv depth decks are prepared along with the necessary
Phase I control cards (TABLES 8, 9 and 11) and the data are gathered to
prepare the master control card set (TABLES 1 and 10), Phase I can be run to
verify the depths and control card data.
Section 3 defines the data requirement for the various HN model phases.
TABLE 9 shows the data card sequence for these Phase I control cards and the
18
-------
depth data. TABLE 8 contains the computer job control card deck setup for
Phase I and TABLE 11 provides the required library directive set to access
the Phase I program on the HN program library.
Phase I provides various printed output options to assist the implemen-
tation. All Phase I control and data cards are verified as to format and
printed. A printed plot is presented for each TIDE card present with an in-
dication of appropriate offset to initialize properly the tidal input to
Phase II. A printed symbolic array of the land-sea boundaries and set of er-
ror messages indicate required correction to the depth fields in order to de-
fine valid land-sea boundaries. Correction can then be made in subsequent runs
by either using change cards or by repunching the depth cards. The CFL de-
rived time step is printed for both the maximum depth specified on the LAYS card
and the actual maximum depth in the field as verification of that computed dur-
ing the grid preparation phase. The VOLM card allows computation of regional
or section volumes. The TAPE card, in addition to causing the depth values to
be written on the intermediate data set, provides various print options for the
symbolic height and depth fields and reinitializes the Phase I program so
additional Phase I control/data card sets can be processed in the run.
Phase II
The intermediate data set prepared by Phase I and the master control
card set are used to execute Phase II. TABLE 8 shows a Phase II computer
job control deck setup. TABLE 10 presents the HN Phase II control card set
which provides access to and verification of the intermediate data tape during
Phase II initialization. TABLE 11 provides the essentials for determining the
appropriate library directives set to access the Phase II program on the HN
program library.
19
-------
Phase III
Use of HN Phase III on any computer system other than the CDC 3100 would
require a conversion effort to the local plotting capabilities as well as to
the local computer system. The documentation (2) of data requirements should
prove adequate for most applications.
Phase IV
In its present form Phase IV is a thermal advection model using the
method of moments based on Pederson and Prahm (9). With some additional ef-
fort to remove the thermal calculations the model could be used for the
advection of other properties. TABLE 11 presents the required library direc-
tives set required to access the thermal advection package on the HN program
library. TABLE 8 provides the computer job control deck setup for the thermal
advection, and TABLE 10 presents the required HN control cards used to con-
trol and access the intermediate data tape for thermal advection implementa-
tion. The variable ITS on the TIME card should be set to the earliest start
time of any source.
PHASE II TIMING NOTES
TABLE 12 presents timing information in the form of a "performance index."
The HN performance index, !„„/ is a comparative measure of job performance for
HN
HN Phase II. The index is basically the length of time in CDC 7600 central
processing unit (CPU) seconds required to perform all type computations for a
single grid node and computational step. The total number of CDC 7600 CPU
seconds required to accomplish a particular Phase II objective can be esti-
mated by using a composite index to describe the job. For example, a compos-
ite index for a non-linear HN version including tidal flat capability and tidal
flat I/O without saving any other data, including convective data, would be
-4
formed by using references 25, 23, 16 from the table (1.56+.37+.16=2.09(X10 ')).
20
-------
The index is then multiplied by the total number of sea points for all layers
and the number of computational steps for the job to compute the CDC 7600 CPU
second estimate.
Each of the indices for references B though 10 of TABLE 12 was computed
directly from the accounting statistics for HN Phase II jobs as described.
References 0-10 included in excess of 1000 computational steps for 260 sea
points. References B and D included 567 sea points computed for in excess of
8000 computational steps. The number of sea points for a particular area may
be computed as the summation of the symbolic height array from Phase I, as:
NO
Ps = £ HTZ(I).
1=1
The number of computational steps may be computed readily from information on
the Phase II TIME control card for the job, as:
Sc = (ITE-ITS)/ITINC.
Currently estimated CDC 7600 job cost at LBL Berkeley for about 1000 steps,
-4
1000 sea points, an index of 2.0X10 and approximate cost per CPU second equal
to 23ft would be:
COST = .23P SI = $46.00.
s c HN
The cost estimate will generally be conservative in the sense that it will
overestimate the cost for long jobs but may underestimate for some very short
jobs. A recent job with 3125 steps, 716 sea points and an estimated composite
-4
index of 2.25X10 resulting in a 505 CPU second estimate actually used 334 CPU
seconds with resulting actual index of 1.50. The index discrepancy resulted
primarily from the low tidal flat involvement relative to the index base.
A set of such indices for the particular computer system used are readily
derived based on several actual jobs. Several of the component indices are
21
-------
extremely area specific, e.g., the tidal flat computations should vary exten-
sively from area to area based on actual tidal flat involvement. Less than
25 % of the sea points in the tidal flat grid ever become involved in the
tidal flat computations directly and even then for less than 50% of the model
time. Note that the tidal flat I/O depends markedly on whether other data
are saved during the job (MODOUT on TIME card ^ °°.) This is due primarily to
the requirement to reread the entire output data file in order to print the
tidal flat change steps.
NOTES ON EXPERIMENTAL/PROPOSED FEATURES/METHODS
The following set of notes is intended to encourage additional thought
and work on the HN model and its application to various areas. The existence
of a successful model need not imply that the model is static. Much work needs
to be done in areas that have been little explored, if at all, numerically or
analytically.
Friction Term Modification
Experiments were performed which showed that replacing the smoothed with
uhsmoothed velocity elements of the friction term as indicated by the form of
the differential equations in (1) was insignificant for tests performed on the
tidal flat grid. After three hours of model time the maximum differences be-
tween comparable runs with the smoothed terms in the friction computation and
Mo
with old velocities in the friction computation were 8.X10 cm for water
_2
heights and 2.8X10 cm/sec for velocity components.
Removal of a constraint restricting the magnitude of the friction computa-
tion to be no larger in absolute value than the magnitude of the smoothed veloci-
ty component from the prior step also proved insignificant. A difference between
comparable runs with and without the constraint produced zero differences to the
22
-------
precision of the data saved on tape (10 cm, 10 cm/sec) after three hours
of model time for the tidal flat area.
Actual computations for the various enhancement tests retained both the
smoothed velocity term and the friction constraint to avoid introducing minor
differences between most recent and earlier runs. Since both modifications
were shown insignificant numerically, no major impact need be expected from
inclusion of these changes in library correction set HNG9 for future usage.
Further tests could be performed by using YANK feature on HNG9 set.
Notes on Modeling River Inflow
Extensive testing of river inflow specification has not been accom-
plished. If the river mouth lies on an open boundary without tidal speci-
fication, the boundary acts as an aperture to an infinite capacity channel.
The problem is alleviated if actual river mouth heights can be specified. If
the tidal values are not available, perhaps the most practical solution is to
include a portion of the estuary as part of the grid specification leaving no
opening for river inflow (a dead end channel) and then to specify reasonable
values for volume on a FLOW card at the end-point of the channel (the last
z point.) The channel should be least one Z cell in length to achieve the
desired flow effect into the adjacent water mass. If the estuary volume and
friction forces are modeled, tidal flows at the mouth should be reasonable
even if the estuary is realigned in the grid to prevent an excessively large
grid size or land-sea ratio.
Notes on Introducing Diffusion into Advectiye__Schemes
Vertical diffusion could be readily incorporated into the existing model
of moments by determining a volume increment to be applied at each step to
all affected cells. Extremely arbitrary experiments were run testing the
23
-------
feasibility by using an increment that was proportional to the absolute value
of the base ten logarithm of the total volume within the cell. The factor
chosen was too large for a first guess and the volumes increased too rapidly.
The factor probably should be based on the time step and the affected area with-
in each cell.
Horizontal diffusion effects would require a complete rewrite of the
method of moments advection scheme if they were to be incorporated therein.
24
-------
TABLE 1. MASTER CONTROL CARD SET (page 1 of 3)
Card Type
Program Variable
Description
NAME Card
1-4
5-80
GRID Card
1-4
5-7
8-10
11-20
21-30
31-40
41-50
51-60
TIME Card
1-4
5-10
11-20
21-30
31-40
41-50
51-60
61-70
FACT Card
1-4
5-10
11-20
21-30
31-40
NE
ME
DL
ROTANG
C
TK
R
MODOUT
ITO
ITS
ITE
ITINC
IPO
IPMOD
Blank
TIDSPD(l)
TIDSPD(2)
TIDSPD(3)
Format (20A4)
NAME
Comments or titles (any number of cards
may be used to identify the output)
Format (A4,213,7E10.3)
GRID
Number of rows
Number of columns
Distance between grid points (cm)
Counter-clockwise angles between north
and the +Y axis of the computational
grid (deg)
Wind drag coefficient
Mid latitude of grid (deg)
Bottom friction coefficient
Format (A4,I6,7I10)
TIME
Number of seconds between saved results
Time of first output in seconds
Start time for computations in seconds.
If this time is greater than zero,
the intermediate tape is assumed to
contain U,V and Z arrays for neces-
sary levels at the time specified.
End time for computations in seconds
Time step in seconds
Time of first printout
Number of seconds between printouts
Format (A4,6X,7E10.3)
FACT (when present must precede TIDE
cards in Phase I)
Speed of the tidal constituents
(Match order of component phase angle
and amplitude on TIDE card)
Kl = 15.0410686 (deg/hr)
01 = 13.9430356
M2 = 28.9841042
S2 = 30.0000000
N2 = 28.4397295
25
-------
TABLE 1 (Page 2 of 3)
Card Type
TIDE Card
1-4
5-7
8-10
11-13
14-16
17-22
23-29
30-36
37-43
44-50
51-56
57-62
63-68
69-74
75-77
78-80
WIND Card
1-4
5-7
8-10
11-13
14-16
17-20
21-30
31-40
41-50
51-60
LAYS Card
1-4
5-10
11-20
21-30
31-40
Program Variable
ITIDE(1,J)
ITIDE(2,J)
ITIDE(3,J)
ITIDE(4,J)
ITIDE(5,J)
TIDE(1,J)
TIDE(2,J)
TIDE (3, J)
TIDE (4, J)
TIDE (5, J)
TIDE(6,J)
TIDE (7, J)
TIDE (8, J)
TIDE (9, J)
TIDE (10, J)
IWIND(1,J)
IWIND(2,J)
IWIND(3,J)
IWIND(4,J)
Blank
IWIND(5,J)
IWIND(6,J)
IWIND(7,J)
IWIND(8,J)
HL(1)
HL(2)
HL(3)
pescription
Format (A4 , 413 , 16 , 4F7 . 2 , 4F6 . 2 , 2F 3 . 2 )
TIDE (O^J<4)
Minimum row (N) index for tidal input
Maximum row (N) index for tidal input
Minimum column
Maximum column
(M)
(M)
index for tidal input
index for tidal input
Time offset between time of computations
and time 0 tide
Tidal constituents
(deg)
Tidal constiuents
phase angle (sec)
1 phase angle
2
3
4
1 amplitude (cm)
2
3
4
Weight factor for second layer tide
heights (no
if 0)
tide set in second layer
Weight factor for third layer tide
heights (no
if 0)
tide set in third layer
Format (A4 ,413,4X4110)
WIND (OSJ<4) .
Minimum row (N) index for wind field
Maximum row (N) index for wind field
Minimum column
Maximum column
(M)
(M)
index for wind field
index for wind field
Wind start time (sec)
Wind end time (sec)
Wind direction
blows)
(degs true from which wind
Wind speed (m/sec)
measurement
Note: This is the only
in the input given in
meters rather than cm.
Format (A4,I6,3F10
LAYS
0,8F5.0)
Number of layers
Maximum depth
of layer
1
2 if zero the layer is
3 ignored in the compu-
+
•a 4" i /-\v» e
26
-------
TABLE 1 (Page 3 of 3)
Card Type
41-45
46-50
51-55
56-60
61-65
66-70
71-75
SORC Card
1-4
5-7
8-10
11-13
14-16
17-18
19-20
21-30
31-40
41-50
51-60
FLOW Card
1-4
5-40
41-50
51-60
61-70
71-80
COM? Card
1-4
5-10
Program Variable
ALP HA (1)
ALPHA (2)
ALPHA (3)
RHO(l)
RHO(2)
RHO(3)
HMIN
ISORC(1,J)
ISORC(2,J)
ISORC(3,J)
ISORC(4,J)
ISORC(5,J)
ISORC(6,J)
ISORC(7,J)
ISORC(8,J)
SORC(1,J)
SORC (2, J)
IFLOW( ,J)
FLOW(1,J)
FLOW (2, J)
FLOW (3, J)
FLOW (4, J)
MODEL
Description
Smoothing parameters for layer 1
2
3
Density for layer 1
2
3
Minimum allowable thickness of layer
in Phase II
Format (A4,4I3 ,2I2,2I10,4E10.3)
SORC (O^J^IO)
Minimum row (N) index for source
Maximum row (N) index for source
Minimum column (M) index for source
Maximum column (M) index for source
Minimum layer (LM) index for source
Maximum layer (LP) index for source
Start time for source
End time for source
Concentration introduced per time step
Volume of water introduced per time step
per level (cm^)
Format (A4 ,413, 212,2110, 4E10. 3)
FLOW (O^J^4)
Same as SORC card
Direction of flow (deg true)
Velocity increment (cm/sec/step)
Second layer factor
Third layer factor
Format (A4,I6)
COMP This card is used to start the
computations in the Phase II pro-
gram. It is the last card read in
the control card deck by the pro-
gram Phase II.
Model number used in record header on
intermedate data set
27
-------
TABLE 2. HN PHASE I CONTROL/DATA CARDS (page 1 of 3)
Card Type
NEWT Card
Program Variable Description
1-4
5-10
11-20
OLDT Card
FAC
SEALEV
1-4
5-10
11-20
PAC
SEALEV
New Format HTU and HTV Tables
1-6
IU or IV
Old Format HTZ, HTU and HTV Tables
1-6 IZ, IU or IV
SETU Card
1-4
Format (A4,F6.3,F10.0)
NEWT Read in new format HTU and HTV
tables that follow card
Conversion factor used to change depth
units to cm. If blank assumed = 1.
Sea level adjustment to be added to
water depths in cm. (No adjustment
or conversion occurs until TAPE
card processing).
Format (A4,F6.3,F10.0)
OLDT Read in old format HTZ, HTU and
HTV tables that follow card
Conversion factor used to change depth
units to cm. If blank assumed = 1.
Sea level adjustment to be added to
water depths in cm. (No adjustment
or conversion occurs until TAPE card
processing.)
Format (1216)
Follows NEWT card. Each row of the ar-
ray is read with a separate READ
statement. The array element (1,1)
is punched on the first field of the
first card and the last element in
the row (1,ME) is the last value on
the (ME-1)/12-KL card. The elements
(2,1)... (NE,1) all begin new card
sets. The IU array is read first
followed by the IV array.
Format (1216)
Follows OLDT card. Each array is read
using logic that duplicates an array
READ using a double indexed implied
DO loop with the column index vary-
ing faster. Only the last card in
each array se^t can contain unused
fields. HTZ array is read first
followed by the HTU and HTV arrays.
Format (A4,4I3,4X,I10)
SETU Set the value ICC in HTU within
the limits
28
-------
TABLE 2 (page 2 of 3)
5-7
8-10
11-13
14-16
17-20
21-30
SETV Card
SETZ Card
VOLM Card
.1-4
5-20
21-30
GRDZ Card
1-4
5-10
11-20
OLDP Card
1-4
5-16
NEWP Card
1-4
5-16
TAPE Card
ICC
Minimum row (N) dimension
Maximum?, row i.N) dimension; if zero set
to min value
Minimum column(M) dimension
Maximum column (M) dimension; if zero
set to min value
Blank
Value to be placed in army (same units
as HTU field prior ho conversion)
Same as SETU. Modify HTV array.
Same as SETV Modify HTZ array..
Format (A4,413,4X,110}
VOLM
Same as SETU card
Reset flag
1 = volume to zero before computing
vo 1 urne
0 = add computed volume to previous
volume
Format (A4,6X110)
GRDZ Set HTZ and -1 values in HTU and
HTV field; test boundary configura-
tion; provide print option
Blank
Flag to bypass combined symbolic HTZ,
HTU and HTV printout
0 = print
1 = bypass printout
Format (A4,4I3)
OLDP Punch old format HTZ, HTU and
HTV cards
Same as SETU card, except if left blank
limits set to 1, ME or NE as ap-
propriate
Format (A4,4I3)
NEWP Punch new format HTU and HTV
table
Same as OLDP card
FjiTf^.t (M, 16 ,110)
29
-------
TABLE 2 (page 3 of 3)
1-4
5-10
11-20
ICC
ICD
TAPE
1) Combine u,v depth print
2) Convert and adjust sea level as
required
3) Print HTZ, HTU and HTV
4) Set HTZ, HTU and HTV values less
than 0 to 0 and scale by 10
5) Write HTZ, HTU and HTV fields on
tape
6) Set Z, U and V fields to 2 for all
cells over water, -0 for all cells
over land. Write on tape for lay-
er 1.
7) Repeat 6 for layer 2 and then 3 if
HL(2) and HL(3) are non-zero
8) Write double EOF on tape, back-
space 1 EOF
9) Return to reset for next data set
Skip parameter
>0 Print HTZ, HTU and HTV but do not
write tape
0 Print and write tape
>0 Request U/V grid print of depth
field
No print
30
-------
TABLE 3. UPDATE DIRECTIVES EXAMPLES
DIRECTIVE
*DECK A
*COMDECK B
*IDENT C
*DELETE A.3,A.7
*BEFORE A.8
*INSERT A. 8
*RESTORE A.5,A.7
*CALL B
*COPY A,A.9,A.11
*YANK C.D
*ADDFILE
*COMPILE A
DESCRIPTION
The first card of deck A
The first card of common deck B
The first card of correction set C
Deactivates 3rd through 7th cards of deck A and substi-
tutes following nori-update text after card A.7
Inserts following text immediately prior to card A,8
Inserts following text immediately after card A.8
Reactivates cards A. 5 through A.7 and inserts following
text after card A.7
Replaced by active cards in common deck B
Copies cards A.9 through A.11 of deck A into current
correction set
Deactivate correction sets C through D
Insert following decks and common decks at the end of the
last deck on new program library
Specifies decks to be prepared for compilation on file
COMPILE
TABLE 4. CDC UPDATE CONTROL CARD PARAMETER EXAMPLES
CONTROL
PARAMETER
C
D
F
I
K
L
N
O
P
Q
FUNCTIONAL, DESCRIPTION
The
Allows specification of compilation file name.
default name is COMPILE.
Allows 80 columns to be used for data sets, otherwise
72 columns
Full upate (used here only for creation of new library).
All decks are processed.
Allows specification of an alternate input file name.
The default name is INPUT.
Provides decks to compile file in sequence specified on
the compile directive
List options for printed output
New program library . The default name is NEWPL
Allows specification of an alternate printed output file
name. The default name is OUTPUT.
Old program library. The default name is OLDPL
Quick update. Processes only decks listed on compile
directives.
31
-------
TABLE 5. HN LIBRARY: PHASE I MODULES
TYPE/
NAME
DECK/
UTIL
HNI
CORREC-
TION SET/
HNJ
BKYPRT
AHNIJ1
HNK
HNL
PHASE1
STATUS
A
A
I
I
A
A
A
E
DEPEN-
DENCY
PR3
BC8
BKY
FACILITY
INCOMPAT-
IBILITY
Phase II
Other
sites
FEATURE
DESCRIPTION COMMENTS
Includes routines : INPOUT ,
IOERR , XMIT , PRTMAX
Includes routines LAYER3S,
GRID , GEOPOS , TIDCUR
Phase II to Phase I conver-
sion of deck UTIL and 3100-
6500 conversion of deck HNI,
Inactive due to incompati-
bility.
Required for correct over-
printing at BKY for land/
sea boundary print feature
in routine GEOPOS .
Inactivates HNJ, BKYPRT
through YANK feature to
make library Phase II com-
patible and 6500 site com-
patible .
Corrects indexing error in
VOLM control card support-
ing code.
Corrects page break suppres-
sion so it works on 6500 and
7600.
Reactivates HNJ for non-BKY
site use of Phase I code.
Status Code:
A-Active on the. library
I-Inactive on the library
E-Extra-library correction set
32
-------
TABLE 6. HN LIBRARY: PHASE II MODULES (Page 1 of 3)
TYPE/
NAME
COMDECK/
CONTROL
BLANK
DECK/
HNG
BARS TAR
UTIL
CORREC-
TION SET/
HNG4
HNG 5
HNG6
HNG7
HNG 8
VARWIND
VARWINE
VARCORL
HNG 81
STATUS
A
A
A
A
A
A
A
A
I
A
I
I
I
A
DEPEN-
DENCY
HNG4
HNG 5
HNG6
HNG6
HNG8
VAR-
WIND
HNG8
HNG8
INCOMPAT-
IBILITY
HNG8
HNG7
UTILA,
PR3,
CONPRU
HNG 9
FEATURE
B
B
B
B
B
B
B
B
B
W
W
C
B
DESCRIPTION COMMENTS
COMMON/CONTROL/ (includes 500
positions for local routine
use)
COMMON// (the major array
variables for all layers)
Includes routines LAYERS,
SETHUV,UVCAL
Includes routines UVSTAR,
BAR,BARUV
Includes routines INPOUT,
IOERR , XMIT , CXMIT , PRTMAX
HNG4 through HNG6 are some
of the earliest corrections
to the initial creation of
the HN library. Earlier
correction sets have been
irrecoverably incorporated
by resequencing decks.
Minor corrections to the
original optimized code
superceded by HNG8 .
Major structural changes in
optimized code toward com-
plete layer generality. In-
troduced UVCAL, affects all
decks.
Allows wind array input as
Type 8 data records from
tape logical unit 2. Re-
quires review before use.
Revision of VARWIND for com-
patibility with more recent
enhancements. Corrects ro-
tation error.
Allows variable coriolis
parameter for arbitarily ro-
tated grid. Requires review
prior to use.
Corrects rotation error in
basic WIND feature, produces
z centered velocity print,
and Hu/Hv correction for
layers 2+.
33
-------
TABLE 6 (Page 2 of 3)
BC8
BARMOD
PR3
UTILA
UTILB
LIB2
HNG9
VALDEZ
SINGLE
CONTRM
BLANKA
CONVCT
CONVCU
CONPRT
A
A
A
A
A
A
A
I
A
A
A
I
I
I
HNG8
HNH8
HNG8
UTILA
CONTRM
BLANKA
HNG5
SINGLE
BLANKA
SINGLE
BARMOD
CONVCT
CONVCT
HNJ
HNJ
FLATIO
VARWIND
VARCORL
B
B
B
B
B
B
B
V
B
B
B
1
1
1
Uniform flow boundary condi-
tion to match that on first
row and column
Changes boundary conditions on
2nd order computations (BAR,
BARUV.) Land/sea amplitude
damping: ?l,U/Vt
Precision feature: sets preci-
sion of data on output tape
-3 -1
to 10 instead of 10 cm
(cm/sec)
Data type processing revision
to INPOUT supporting CONPRU
modification to IOERR for 6500
compatibility
Data type processing revision
to INPOUT supporting CONPRV
and VARWIND
Corrects minor coding errors
Corrects thickness check for
bottommost layer; removes
friction constraint; reverts
from bar term usage to old
velocity in friction term;
removes non-linear interpola-
tion scheme from star term
Changes single layer version
for array storage
Removes arrays for 2nd and 3rd
layer; leaves single layer
version. For array storage:
MAXDIM=800
Changes single layer version
such that CONTROL arrays di-
mensioned 300
Changes single layer version
for array storage: MAXDIM=300
Feature 1 : non-linear, term
computations. Forms and uses
ENTRY CONUV to routine BARUV
for the procedure
Corrects prior version boun-
dary condition to consistent
boundary treatment normal vs.
tangential
Feature 1 OUTPUT option: non-
linear terms printed in array
form using PRTMAX (2 arrays)
34
-------
TABLE 6 (Page 3 of 3)
CONPRU
CONPRV
FLAT
PLATA
FLATIO
FLATIP
MONTE
AHNG71
AVARC1
AVARW1
AVARW2
ACONV1
ACONV2
AFLATl
AMONTl
AVALDEZ
HNII
(series)
I
I
I
I
I
I
I
A
A
A
A
A
A
A
A
E
UTILA
CONPRT
UTILB
CONPRU
HNG8 . 1
BC8
FLAT
FLAT
FLATIO
UTILA
1
1
2
2
2
2
4
BC
BC
BW
BW
Bl
Bl
B2
B4
BV
Modify print option replacing
printed arrays with data rec-
ords TYPES on LU1
Modify data save option replac-
ing TYPES with TYPE10 data rec-
ords for VARWIND compatibility
Feature 2 : layer disappearance
enhancement follows thickness
comp. , allows tidal flat
phenomena
Corrects detail feature of re-
computation flag field to as-
sure proper tidal flat handl-
ing
Feature 2 output option: tidal
flat flag fields saved on LUl;
printed at normal end of job,
ends EOF LUl.
Modifies FLATIO to simplify
INPOUT processing in conjunc-
tion with UTILA
Feature 4 : Monte Carlo enhance-
ment follows thickness comp.
(FLAT comp. if present) abbre-
viated printout of "plume"
HNG7 The correction sets
VARCORL AHNG71 through AVALDEZ
VARWIND inactivate the correc-
VARWINE tion set(s) indicated
CONPRT immediately to the
CONPRU left through the YANK
CONVCT feature to maintain
CONVCU the HN library orient-
CONPRV ed toward the basic
FLAT Phase II run. Revers-
FLATA al of directives by
FLATIO YANK or DELETE effec-
FLATIP tively reactivates the
MONTE desired features
VALDEZ
Reactivates appropriate fea-
tures and feature options for
particular Phase II objective
FEATURE CODE:
1 Non-linear enhancement B BASIC HN PHASE II
2 Layer disappearance enhancement C VARIABLE CORIOLIS
4 Monte Carlo extension V MAXDIM=1200
W VARIABLE WIND
35
-------
TABLE 7. HN LIBRARY PHASE IV MODULES (page 1 of 2)
TYPE/
NAME
COMDECK/
A
B
C
ABC
THERM
DECK/
AD
TH
CORREC-
TION SET/
Al
Tl
COM
A2
T2
A3
T3
STATUS
A
A
A
A
A
A
A
A
A
DEPEN-
DENCY
Al
Tl
A2
T2
INCOMPAT-
IBILITY
FEATURE
3
3
3
3
3
3
3
3
3
DESCRIPTION COMMENTS
COMMON/A/ (control
linkage)
COMMON// (advect array
linkage)
COMMON/1/ (utility
linkage)
COMMON/ ABC/ (thermal
array linkage)
COMMON/THERM/ (ther-
mal control)
Includes routines AD-
VECT , BDVECT , INPOUT ,
IOERR , XMIT , PRTMAX
Includes routines
THERMAL , SETFLD , WATER
ANOMALY, CALC
Introduced routines
PRTMAX, XMIT and PRT-
MAX, provided 3100-
6500 conversion
Condensed TH routines
to single routine
THERMAL
Restructured common
decks B, ABC, THERM
Restructured linkage
to thermal and utili-
ty routines for du-
plex advect ion limit-
ed printout
Adjusted computation-
al boundaries , limited
the printed output
Restart/data save,
interpolated velocity,
v direction correc-
tions and general im-
provement of coding
structure
General improvement of
coding structure, com-
ments and corrected
anomaly treatment
36
-------
TABLE 7 (page 2 of 2)
A4
T4
HNIVTA
A2
T3
Replaces wind stress
with true wind in m/
sec units
Replaces dynamic
moisture exchange
computation with de-
termination from
relative humidity
Dummy ident
Status Code: A-Active module on library
E-Extra-library correction set
I-Inactive module on library
Feature Code: 3-Thermal advection
37
-------
TABLE 8. HN LIBRARY USAGE JOB CONTROL DECK SETUP FOR CDC6500
[JOBNAME],TU2.[USER NAME]
VSN(OLDPL=#l)t
REQUEST,OLDPL,MT,HI,NORING.
UPDATE,Q,L=124. (L=124 PROVIDES ESSENTIAL PRINTED OUTPUT)
UNLOAD,OLDPL.
FTN,L=0,OPT=2,I=COMPILE.
COMMENT. REPLACE L=0 WITH L,R=2 FOR REFERENCE LISTING.
RETURN,COMPILE.
VSN(TAPEl=#2,TAPE2=#3)
REQUEST,TAPEl,MT,HI,RING.*
REQUEST,TAPE2,MT,HI,NORING.*
LGO.
RETURN,TAPEl,TAPE2.
EXIT.
7
8 (binary end of record)
9
[CDC UPDATE DIRECTIVES SET]
7
8 (binary end of record card)
9
[HN PHASE CONTROL CARD SET]
6
7
8 (binary end of information card)
9
t Logical Units
2 Input for wind fields Phase II/output for advection fields
Phase IV (as required)
1 Intermediate file 556 read/write tape
LGO Load-and-Go unit
INPUT System control cards, update directive sets and master
control/data card decks
OUTPUT Control card list and horizontal printed output
COMPILE Corrected FORTRAN source code
OLDPL Program library tape
* Tape numbers should be reversed for Phase IV jobs. (Tape2 re-
quest optional based on indicated usage.)
38
-------
TABLE 9. HN PHASE I CONTROL/DATA CARD SET
IF USED REQUIRES TIME,FACT,PRECEDING)
NAME (OPTIONAL)
GRID (REQUIRED)
LAYS (REQUIRED)
TIME (OPTIONAL)
FACT (OPTIONAL)
TIDE (OPTIONAL,
SORC (OPTIONAL)
COMP (REQUIRED)
OLDT NEWT
{symbolic Z data} old {u depths}
or new
{U depths} format data {V depths}
cards re-
{V depths} quired
SETU (OPTIONAL)
SETV (OPTIONAL)
GRDZ (OPTIONAL WITH OLD FORMAT-REQUIRED WITH NEW FORMAT)
SETZ (OPTIONAL FOR OLD FORMAT-NOT REQUIRED FOR NEW FORMAT)
VOLM (OPTIONAL)
TAPE (REQUIRED AS THE LAST CARD OF THE SET)
Note: Multiple sets may be processed by Phase I. If multiple sets
are written on the intermediate data set, they are separated
by an end-of-file mark. User is expected to preposition in-
termediate data set file when stacking data sets with mul-
tiple Phase I runs.
39
-------
TABLE 10. HN MASTER CONTROL CARD SET (ALL PHASES)
NAME
GRID
LAYS
TIMEt
FACT
TIDE
SORCt
WIND
FLOW
COMP*
(OPTIONAL)
(REQUIRED)
(REQUIRED)
(REQUIRED)
(OPTIONAL)
(OPTIONAL)
(OPTIONAL)
(OPTIONAL)
(OPTIONAL)
(REQUIRED AS
LAST CARD OF SET)
t For thermal advection applications the following control card
usage variations apply:
TIME Card-The computational time step (ITINC) is replaced
by the thermal advection time step which must be
a multiple of the data save interval of the Phase
II run which generated the intermediate data set
in use.
SORC Card-The quantity of heat introduced per second is en-
tered in place of the concentration per step; the
volume per second is entered in place of volume
per step. These deviations facilitate run com-
parisons with various steps.
WIND Card-The direction and wind speed must be floating
point quantities.
* There is no required order for any control card except COMP
for Phase II or thermal advection applications.
40
-------
TABLE 11. HN LIBRARY USAGE BY CDC UPDATE DIRECTIVES SET
OBJECTIVE
Phase I
(at other than BKY site)
Phase II (Basic)t
Phase IV
Phase II with Tidal Flat
and Tidal Flat I/O
Phase II with non-linear
terms and non-linear I/O
Phase II with Tide Flat
non-linear and monte carlo
capabilities without FLAT
or CONVECTIVE I/O
Phase II basic modified to
1200 array positions (for
VALDEZ)
Phase II basic modified to
800 array positions (for
PRUDHOE)
Phase I (at BKY site)
REQUIRED UPDATE DIRECTIVE SET
*IDENT PHASE1
*DELETE AHNIJ1.1
*COMPILE HNI,UTIL
*IDENT HNII
*COMPILE HNG.UTIL
*IDENTI HNIVTA
*COMPILE AD,TH
*IDENT HNIIFO
*YANK AFLATl
*COMPILE HNG.UTIL
*IDENT HNIICO
*YANK ACONV1,ACONV2
*COMPILE HNG.UTIL
*IDENT HNIIFCM
*DELETE AFLATl.1
*DELETE ACONV1.1
*DELETE AMONTl.l (or *YANK AMONT1
*COMPILE HNG.UTIL
*IDENT PETROX
*YANK AVALDEZ
*COMPILE HNG.UTIL
*IDENT ICE
*YANK CONTRM,BLANKA
*COMPILE HNG.UTIL
*HNPHASEI
*YANK AHNIJl
*COMPILE HNI,UTIL
The Basic HN Phase II is the three layer capability version
with arrays restricted to single layer applications using < 300
grid points per array-
41
-------
TABLE 12. HN PHASE II TIMING PERFORMANCE CHART
JOB/
REF
B
D
O
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
BASIC
HN
X
X
X
X
X
X
X
X
X
X
X
X
X
X
X
X
X
CONVCT
X
X
X
X
X
X
X
-
-
-
X
X
X
X
X
X
CONPRT
X
X
X
-
-
-
-
-
-
-
X
X
X
X
X
FLAT
X
X
X
X
X
X
X
X
X
X
X
X
X
FLATIO
_
X
X
X
-
-
-
-
X
X
-
X
X
X
X
MONTE
POINTS
1000
1000
1000
1000
1000
1000
1000
1000
1000
1000
4000
1000
1000
MODOUT
240
240
120
120
120
120
120
CO
oo
oo
120
120
00
120
120
'HN104
1.18
1.50
3.27
3.27
3.19
2.94
2.97
2.93
2.49
2.65
2.81
3.03
5.35
.95
.08
.24
.03
.44
.16
.16
.24
.38
.08
.32
1.15
.37
1.54
1.56
2.11
REF*/COMMENT
*Performance
indices be-
low the
dashed line
are esti-
mates based
on the ref-
erences
given.
1,10
1,2 (CONPRT
Index=0,if
MODOUT= «.)
2,3,4
(MODOUT=120)
3,4,5 Basic
output
5,6
6,7
(MODOUT= «=)
7,8 (BASIC+
FLATIO)
1,9
(MODOUT=120)
2,8
3,4,9
B,D (older
^ersion
CONVCT)
B,14
22,11,6
22,23
22,15,19
22,11
42
-------
REFERENCES
1. BAUER, R. A. and A. D. Stroud. Enhanced Hydrodynamical-Numerical Model for
Near-Shore Processes. In press. Environmental Protection Agency-Corvallis
Environmental Research Laboratory, Corvallis, Oregon.
2. Bauer, R. A. Description of the Optimized EPRF Multi-Layer Hydrodynamical-
Numerical Model, ENVPREDRSCHFAC Technical Paper No 15-74, Environmental
Prediction Research Facility, Monterey, California, 1974. 66 pp.
3. Harding, J. Unpublished program notes for HNCNTR(CDC3100). Naval Environ-
mental Prediction Research Facility, Monterey, California, 1976.
4. Bauer, Roger A. Digitization Procedure for Wave Refraction Bathymetry.
Final report on Contract N00228-76-C-3121. Naval Environmental Prediction
Research Facility, Monterey, California. 42 pp.
5. Bauer, R. , S. Larson, T. Laevastu and A. Stroud. Linkage of the Bengtsson
Limited Area Forecast Model and the Optimized Hydrodynamical-Numerical
Model of the W. Hansen Type. Northwest Fisheries Center Processed Report,
July 1976. Northwest Fisheries Center, National Marine Fisheries Service,
Seattle, Washington.
/R)
6. Control Data Corporation. Update Reference Manual for the Control Data^
Cyber 170 Series Computer Systems, Cyber 70 Series Computer Systems, 6000
Series Computer Systems and 7600 Computer Systems, Publication No
60342500 Rev G, 1975.
7- Control Data Corporation. FORTRAN Extended Version 4 Reference Manual,
Publication No 60305600 Rev H, 1975.
8. Control Data Corporation. Scope 3.4 Reference Manual, Publication No
43
-------
60307200 Rev L, 1975.
9. Pedersen, L. B. and L. P. Prahm. A Method for Numerical Solution of
the Advection Equation 3C/8t=-V«VC, Meteorological Institute of Den-
mark, as submitted to TELLUS, August, 1973. 36 pp.
44
-------
APPENDIX A
PHASE I PROGRAM NOTES
DESCRIPTION OF PHASE I PROGRAM
Program LAYER3S (Blocks A-E:HNI.2,HNI.143)
This program was modified for CDC 6500 FORTRAN EXTENDED compiler from
the original CDC 3100 program and retains much of the integer arithmetic
logic used to reduce core requirements on the 3100. The program includes the
main program LAYER3S and the routines XMIT, GRID, INPOUT, IOERR, TIDCUR,
GEOPOS and PRTMAX. The main program is divided into a number of sections
that process master and Phase I control cards. Although sections are acti-
vated independently, the effect varies depending on the order of the control
cards. To describe the program the sections that process individual cards
are grouped into the blocks shown in Figure 3. The sequence numbers refer
directly to the UPDATE sequence numbers in the program listing at the end of
this appendix to the right of the FORTRAN code.
This program uses the following logical units:
LU Description
INPUT Card input file
OUTPUT Printed output file
PUNCH Punched card output file
1 Intermediate data set file
2 Scratch file
3 = INPUT
61 = OUTPUT
45
-------
Initialization (Block A:HNI.2,54)
This block provides the linkage to subroutines for standard control para-
meters in COMMON/DATA/, which also includes space saving scratch storage (IBUF)
for subroutines. The major array variables are linked to other routines by
blank COMMON, which does not require dimensional modification unless grids of
more than 1000 positions are to be processed. Local LAYER3S arrays are speci-
fied in the DIMENSION statement and the DATA statement is used to store maxi-
mum array dimensions and constants. Additional constants are initialized in
the first executable statements which are executed at the beginning of the pro-
gram and after a TAPE card is processed.
All labelled control cards are read from file INPUT and printed immediate-
ly on OUTPUT in "A" format. The four character labels are matched with an in-
dex by table look up, and the index is used to branch to the appropriate sec-
tion for further processing. The index key, card name and statement number for
routing are included as a comment table in the program.
The remainder of the block contains a section for the decoding and re-
printing of each Master control card to provide verification of the punching
and scaling of each field. When present TIDE cards must be preceded by both
TIME and FACT cards. After each TIDE card is interpreted by DECODE, TIDCUR is
called to produce a printer plot tidal curve (see description of routine TIDCUR.)
Unidentifiable cards, GRID specifications exceeding MAXDIM, or too many TIDE,
WIND, SORC or FLOW cards cause early termination.
Process Bathymetric Data (Block B:HNI.55,76; HNI.85,91)
After Phase I control cards have been identified in Block A the following
depth control, edit and data deck cards are processed in this logical block:
NEWT, OLDT, SETU, SETV, SETZ cards and OLDT and NEWT decks. The NEWT and OLDT
cards cause branching to the sections that read in the indicated type of depth
46
-------
card decks. After the tables are read the SETU, SETV cards are used to edit
the fields in the hu and hv fields in original units. The SETZ cards can be
used to edit symbolic hz if old style cards are to be punched after GRDZ pro-
cessing. Usually hu and hv editing take place before a call to GRDZ which pro-
duces an hz field from the hu and hv arrays. GRDZ must be called after NEWT
tables are read. No conversion or level adjustment occurs until the TAPE card
is processed. Fewer than the required number of depth table cards following the
OLDT or NEWT card will abort the run. Extra cards will cause an abort as an
unidentified label card. Normal processing returns contol to Block A to read
the next control card.
Land-Sea Boundaries (Block C: HNI.77,84)
GRDZ card is identified in Block A, routed to Block E for DECODE inter-
pretation and then to Block C for generation of the symbolic hz field by calling
GRID. On return the CFL limited time step is computed and printed for the max-
imum depth specified on the LAYS card and in the hu or hv field. See routine
GRID for details of symbolic height field development and implied land-sea
boundary verification procedures. Within GRID, GEOPOS is called to print a com-
bined geographic display of the three fields. Normal processing returns con-
trol to Block A to read the next control card.
Auxiliary Processes (Block D:HNI.122,143)
After identification of the cards VOLM, OLDP and NEWP in Block A and in-
terpretation in Block B, control is transferred for auxiliary processing to
Block D. Both OLDP and NEWP allow subsets of the hu, hv arrays to be punched:
OLDP produces an hz, hu and hv depth table deck; NEWP produces only the hu
and hv decks. OLDP can only be used for single layer model fields and normally
requires editing with SETZ cards to achieve the proper old style boundary
setting. VOLM computes the cm volume of a specified region to assist in
47
-------
prescribing external boundaries and transport computations. Normal processing
returns control to Block A to read the next control card.
Output Processing (Block E:HNI.92,121)
Both GRDZ and TAPE cards are interpreted in this block by the DECODE func-
tion. The GRDZ card is routed to Block C for further processing. Logical unit
2 is employed in this block to avoid destroying the bathymetric data during
processing of the TAPE card. If specified the routine GEOPOS is accessed
through entry point POGO to print a combined hu, hv depth field in input units
on the OUTPUT file. Conversion of the depth field to cm and adjustment of depth
reference to mean sea level (MSL) is then accomplished and the symbolic height
and hu and hv fields are printed on file OUTPUT by routine PRTMAX.
If tape output is specified all arrays are scaled by a factor of 10; neg-
ative boundary flags are removed; and routine INPOUT is called to buffer the
fields out to logical unit 1 as type 1 records. The £, U and V initialization
fields are prepared as integers scaled by 10 and routine INPOUT is again used
to buffer these fields out to LU1 as time 0 type 2 records.
Normal processing results in return of control tq Block A to reinitialize
data set dependent constants for the next data set after writing two end of
file marks and backspacing over the second on LUl. Since the bathymetric data
still exist for the current data set, if no bathy data are read in a subsequent
set, the same bathy data are available for reuse.
Routine TIDCUR (Block F;HNI.237,292)
The formal parameters represent actual parameters transferred from
LAYER3S as follows: ITS, ITE, IOFF and MODOUT are as prescribed on the TIME
card; SE provides the angular velocities from the FACT card; and PD provides the
amplitude and phase information from the TIDE card. The scratch file, LU2, is
used to save the major arrays and the arrays are used to store the computed
48
-------
tidal heights. Height computations occur for each MODOUT step starting at
(IOFF-1)/MODOUT for (ITE-ITS-1)/MODOUT steps. The computed heights are printed
and the variable height scales for the printer plot computed. The number of
tidal constituents, input tidal parameters, tidal height scale and title are
printed followed by the height scale and a line for each MODOUT step with
approximate minutes scale and height indicated. A binary search technique is used
to locate and print a reasonable zero time offset value for early Phase II con-
vergence. Then the major arrays are restored and return is effected to LAYER3S.
Routine GRID (Block H:HNI.144,193)
Subroutine GRID computes the symbolic grid, hz, and checks the combined
hz, hu and hv arrays for consistency in defining land and sea area boundaries.
A boundary passes vertically through U points and horizontally through V points
and must complete a division between groups of land and sea z grid points.
The formal parameters represent actual parameters transferred from
LAYER3S as follows: LH provides access to the maximum depth values for each
layer from the LAYS card adjusted to the same units as the bathymetric input
under OLDT or NEWT, and ICC provides the key to call or not call routine GEOPOS
for symbolic print and is used to return the maximum depth value to the main
program.
The symbolic height field is set to -0 for the grid size. For each layer
the following steps occur. The minimum layer depth is defined. All hu or hv
points less than the minimum layer depth are set to -0. Those hu or hv equal
to -0 and next to positive depth values in the same row for hu or the same
column for hv are reset to -1. This last step is restricted so that it may not
define a channel one cell wide. During the definitions phase of the symbolic
heights, additional -0 hu or hv points are reset to -1 when channel conditions
are detected. A channel exists when a pair of hu values bracketing a z point
49
-------
equal -0 or -1 and the corresponding z-bracketing hv points are positive,
or vice versa. A -1 is set in the symbolic height field if the layer equals
1 and all of the four adjacent hu or hv values are 0 or -1. Each symbolic
height field position is set to the layer number if one adjacent hu or hv
value is positive. For the first layer the hu and hv fields are then saved
on the scratch file, LU2. After all layers are processed the symbolic hz
field shows the land-sea boundaries and the horizontal extent of each layer.
Since it is possible that erroneous land-sea boundary configurations may
be generated, the boundary is tested for consistency. The only allowed con-
ditions are that there be 0, 2 or 4 boundary points (-1) in the set of four
values between four adjacent z points. A printed message occurs for each
condition where 1 or 3 boundary points (-1's) are included around a z grid
square. Corrections of such conditions may be made on subsequent Phase I
runs using SETU and SETV cards. The testing is done for each layer prior to
print of the combined symbolic z u v grid including the "land-sea" boundary
for that layer. If the print is selected, routine GEOPOS is called to per-
form the function.
Once all layers have been processed and the symbolic height field is in
final form, the hu and hv values for the first layer are read back from the
scratch file. The maximum depth in the hu and hv fields is saved in ICC, and
transfer is effected back to LAYER3S.
Routine GEOPOS; ENTRY POGO (Block G;HNI.194,236)
The formal parameter II provides GEOPOS the layer number for print head-
ing from GRID. The hu, hv and symbolic hz are printed in a quasi-geographic
relationship symbolically. Positive depth values at u and v points are replaced
with +1; all -0 values print as 0; all -1 values in the hu or hv fields print
as a vertical.or horizontal line of asterisks, respectively; and -1 values in
50
-------
the symbolic height field print as -1. Overprinting is employed to rep-
resent properly intersecting land-sea boundary segments as lines of asterisks.
Refer to the sample output for Phase I following the program listing. The
overprinting on the CDC 7600 at BKY required special treatment. A library
module has been .included to remain compatible with the BKY site (BKYPRT) as
required. GEOPOS returns control to GRID.
Entry POGO is called from LAYER3S for a single print on file OUTPUT of
actual depth values for both hu and hv in a quasi-geographic structure in
the same units as input under OLDT or NEWT control. A variable FORMAT (MCL)
is used to handle the variable line lengths required throughout the process.
Return from POGO is effected to LAYER3S.
Routine XMIT (Block I;HNG5.28J,UTIL.66)
The routine XMIT is used to transfer values from a single position or
array into another array. An array can be set to a constant, one array can
be transferred into another or a row or column can be shifted in an array.
The first parameter is the transfer from address; the second is the transfer
to address. The third parameter indicates the number of values to be trans-
ferred and the fourth and fifth parameters control the respective increments
of the transfer from and to addresses.
Routine PRTMAX (Block J;HNG8.127JJTIL.98)
Routine PRTMAX provides a single integer array print capability with
titling. The array to be printed is specified as the first parameter. The
second parameter is unused in this version but provides compatibility with
Phases II and IV. The third parameter allows control over spacing of the
printed array rows ("1HA" for single spacing). The title information is
specified in the fourth parameter and may include up to 80 characters of in<-
formation. The fifth parameter specifies the layer of the array to be printed.
51
-------
The sixth specifies the number of characters included in parameter four.
The seventh parameter specifies the time associated with the array to be
printed.
Arrays are printed in segments of NE row subsets of up to 15 columns.
All segments except perhaps the last will print 15 column subsets. Arrays
of 15 or fewer columns are printed in a single segment. Even ordinal seg-
ments will be printed on the same page as odd ordinal segments when single
spacing has been selected and there are fewer than 25 rows as a paper saving
feature. Array values are printed as 8 position integers. Return is ef-
fected to the main program on completion.
Routines INPOUT: IOERR (Block K:PR3.6,UTIL.63)
Routine INPOUT for Phase I provides the means of placing scaled, com-
pressed data on the intermediate data tape file, LU1. An indicative print-
out occurs on entry to the routine showing the time, layer, record type, grid
dimensions, model number and output record count effective with the last
output operation on the intermediate data file. The format parameters specify
respectively the first array to be written on the data tape, the layer, the
record type requested, the time and the scale factor which has been applied
prior to entry. The data are output in accordance with the data type specifi-
cations provided in Figure 2 as binary buffered output. Routine IOERR pro-
vides error processing when the buffer operation does not function normally.
LISTING FOR PHASE I SAMPLE RUN
The list is a reproduction of a PHASE I run. The first page of this
listing is generated by the UPDATE processing of the library file using the
list option L=124, as shown in the SCOPE control card example given in TABLE 8.
The FORTRAN compiler listing of the Phase I program shows the UPDATE module
52
-------
identification and card numbers on the right edge. The run results begin
after the PRTMAX listing and show an example run called TEST BOX MODEL. Con-
trol cards read by the program appear in the listing with the card labels at
the end rather than the beginning of the line.
The reproduced listing shows the tide curve on two separate pages
when in fact it is produced as continuous print lines spanning as many pages
as are needed to depict the entire run time interval and the resolution of
the number of seconds between saved results that are specified on the TIME
card. Following the tide listing is the suggested offset for the TIDE card
that would allow the model to be started with an increasing tide of near zero
amplitude.
The NEWT table cards are followed by a GRDZ symbolic HTZ listing with a
single land point separated from the rest of the area by lines of *** charac-
ters. The next page gives the CFL time step and volume calculation examples.
Water depths at U and V points are specified in a geographic relationship
with the Z points indicated by + signs. The individual HTZ, HTU and HTV arrays
are printed before the INPOUT lines printed when the type 1 and type 2 records
are written to the intermediate data set.
END LAYER3S is printed when an end of file is deleted on input or an ab-
normal condition causes early termination.
53
-------
IDENT
PHASE!
BKY 76/i6 UPDATE 1.2/110-1H 09 MAR 77 01.15.ZS
PAGE
*****
*****
*****
YANKSSS
YANKSSS
YANKSSS
YANKSSS
YANKSSS
YANKSSS
YANKSSS
YANKSSS
YANKSSS
YANKSSS
YANK***
»IDENT PHASE1
»DELETE AHN1J1 .1
_ »COnPILE HNl/JTIL
MODIFICATIONS / DIRECTIVES
»YANK VARU1NE
•YANK CONPRV
»YANK HNGT
•YANK VARCORL
•YANK VARWIND
•YANK VALDEZ
•YANK CONVCT.LONVCU
•YANK CONPRT.CONPRU
•YANK FLAT, PLATA
•YANK FLATIO,FLAT1P
• YANK P10NTE
•YANK BKYPRT
AVARUI2
ACONV2
AHNG71
AVARC1
AVARU1
AVALDEZ
ACONV1
ACONV1
AFLAT1
AFLAT1
AIMONT1
AHMJ1
\
1
1
1
1
1
2
1
2
1
2
-------
PROGRAM LAYER3S 7600--7600 OPTr2 FTN 1.6+139/010 09 MAR 77 01.15.264 BKY PAGE
PROGRAM LAVER35( INPUT,OUTPUT,PUNCH,TAPE 1,TAPE 2,TAPE 3=INPUT,TAPE61 = HNI .2
I.OU.TP.UTJ s coMi-ioN izt loog 1,1 u( iopo),iv(iooo) HN..L_3
COMMON/DATA/WE,HE,NO,MAX01M.MAXLAY,MODEL, IBUF( 1250 > HNI .1
DIMENSION NAME( 20 ), IC( 8 >, 1TIDE( 5 ), T 1DE( 10 ), 1UIIND( 8),ISORC( 8 ), HNJ.9
JL 5 0 R CX_2 ) ,_N AME Z (_3_>, NAMEUJ j ),NAP!EV(.3) ,1 N( 1_2_LJLLQ_SP D It) . L HJJ ) HN I J>_
DATANE,nE,DL,BOTANGJC)ftLftT.R,AUS » DDL=DL*OL» . 25 HNI. 27
30 1001 FORMATCMX2I3,7E10.3"> HNI.28
G = COS(2.*ALAT»RAD) $ G = ((G»5-9E-6-.0026373 )*G + 1.)»980.616 HNI.29
P R I N T_ J_0 0 tJJ E, Pi E,DL,ROTA N G_, C . ALAT . R, AUS. G HN I_=_30_
NO=NE*ME S IF(NO.GT.MAXDIM)900,1 HNI.31
5 DECOOEl 80,1005,1C ) PIODOUT, I TO, ITS, I TE , ITINC , IPO, IPMOD HNI.32
35 LPJ)5_J_ORri ATJ3 X. 16 j 7_110 ) HNJU 3_
PRINT" 10bT,"nbD"OU"T,ITO, ITS, ITE, ITINC, IPO, IPMOD $ GO TO 1 HN1.3H
6 DECODEC80,1006,1C )TIDSPD,TIDLAG,SLOPEX HNI.35
L016_F 0 R M_AJT (J_Q X_,_7_E_10 . 3J . HN I_._3.6_
PR'lNT lo6"6,TID"SP"b,TIDLAG, SLOPEX * GO TO 1 " HNI.37
tO 10 ITNO=ITNO + 1 4 IF(1TNO.GT.MAXTIDE ) GO TO 900 HNI.38
DE C p D EJ_8 0 ^ 10_10_,1C )_ITIDE ,TIDE S PRINT 1010 . 1T1DE,T IDE H N !_. 3 9_
1010 FORnAt(MX1I3,i6,MF7.2,')F"6.2,2F3.1 ) HNI.10
CALL TIDCUR( ITS,ITE,MODOUT,TIDSPD,TIDE,ITIDEC5 )) * GO TO 1 HNI.11
J,6_JWfilO=_IUN_0+l ?_J_F( lUlNO.GT.riAXUIND ) GO TO 900 HNI .12
15 DECODEC 80, 1016,10") luTlND " $ PRINT 1016,IUIND $ GO TO 1 HN"l.13
1016 FORMATLPH3f RH01 .RH02f BH03f HNI .15
,Y,Z HNI.16
1018 FORMAT<1XI6,3F10.0,8F5.3) HNI.17
.1° PR I NT 1018,1 D.H1L.H21.H3L, ALPHA ,ALPH2.ALPH3.RH01.RH02,RH03fY,Z HN l_J\i
GO TO 1 HNI.19
20 ISNO=ISNO+1 * 1F(ISNO.GT.MAXSORC) 900.22 HNI.50
L0.20_FORMAJ( !XJLU»?J 2^2U_Q_t1E10.._3 ) HNI_.5_1
21 IFNO-lFNO+1 * iF('lFNO.GT.riAX'FLOU) GO TO 900 HNI.52~
55 22 DECODE(80,1020,IC) ISORC.SORC 4 PRINT 1020,ISORC,SORC * GO TO 1 HNI.53
31_J)E_CpPJE( SO.Jip^.S^lC ) MODEL » PRINT 1005 .MODEL 4 GO TO 1 H_NK5.1_
C READ PHASE I CONTROL CARDS »»»••»»».»»»»»»»»»»»»»»»»»»»»»»»»»»»»»HNI.55
-------
PROGRAM LAYER3S 7600-*7600 OPT=2
FTN 1.6+139/010
09 MAR 77 01.15 .26$ BKY PAGE
60
65
80
Ui
CTi
51 DECODE(80,1051,IC)FAC,SEALEV
1051 FORMftTC.1 X.F6 .3.7F10.3)
READ IN NEW TYPE U AND V TABLES
1F( I.EQ.12)GO TO 175 S DO 60 N = 1,NE S MP = C P1E-1 )»NE+N
READ 105S .( IU(I ).lrN,MP,NE)
1059 FOflnATC 1216 )
60 PRINT 1060,N,(IU(I),I=N.MP,NE) * PRINT 1060
J.060 FoRnaT_i\xjj./Lfem. i_6j> -
HNI.56
_HN-LJ?_7_
HNI .58
HNI. 59
HNI.60
DO 120 N = 1,NE t> MP = (ME-1 )»NE+N
120 PRINT 1060,N,( IV(I >,I=N.MP.NE )
.CORRECT TABLES
KNI.61
HNI.62
_H(LLi L.
$ READ 1059,( IV(I ),1 = N,MP,NE) HNI.61
S PRINT 1060 * GO TO 1 HNI.65
HNI.66
70
75
160 DECODEl 80, 1016, 1C )NM ,NP,MM,MP, ICC $ 1 = 1-11
GO TO (166,166,161,100,360,300),!
161 LM=0 $ GO TO 170
166 Lfi=MAXDIn»I
170 NH=MAXO(NP,NM) S MP=HAXO< MP , MM )
IF( MM.LE .0 .OR .MP .GT.ME.OR .NM.LE.O .OR.NP .GT.NE > PRINT 1170, 1C
1170 FORI",ftT( 11K**»ERROR*»»20A1 )
DO 172 M=MM,MP * IP=LM+( M-l >»NE S IM=1P+NM * IP = IP+NP
DO 172 I-1M.1P
HNI .67
HNI .68
HNI .69
HNI .70
HNI .71
HNI .72
HNI. 73
HNI. 71
HNI .75
172 1ZCI) = ICC * GO TO 1 "" " ~ HNI.76
C SET 2 U AND V ARRAYS BOUNDARY VALUES HNI.77
173 IF( FAC.EO.O.) FAC = 1 . $ LH=H1L/FAC $ LH(2)=H2L/FAC $ LH( 3 )=H3L/FAC HNI ,J_6_
IF( H3L.EQ.O. ) P1A~XLAY = 2 S IF( H2L .EQ . 0 . ) MAXLAY=1 HNI.79
IF(ID.EG.O) ID=MAXLAY S IF( ID-NE .MAXLAY ) GO TO 900 KiMI.80
85 171
1171
C
175
90
176
1176
C
95 180
100 181
105
162
C
110
181
CFL = LH(ID) S CALL GHlDf LH J CC )
DO 171 J = l,2 « CEL=DL/SQRT(CFL*FAC»2.»G ) S PRINT 1171,CEL,CFL
CFL=ICC $ GO TO 1
FQRMAT(15H MAX TIME STEP=F9 . 2 , 1 1H MAX DEPTH=F9.0)
READ OLD TYPE TABLES
DO 176 L = l,3 S LM=(L-1 )*MAXDIM $ LP = ( ME-1 )»NE+LM « IND=13
DO 176 I = IM,1P,NE $ IND = IND + 1 $ IF( IND .LE . 12 ) GO TO 176
IND=1 $ READ 1059, IN S PRINT 1176, IN
1Z( I ) = IN( IND) $ GO TO 1
FORMAT( 6X1216 )
WRITE OUTPUT TAPE
DECODEf 80. 1005, 1C ) ICC^ICD
IFCI.EQ.13) GO TO 173 $ REWIND 2 S WRITE (2) IZ,IU,IV
IF(ICD.GT.O) CALL POGO(l) * DO 181 1 = 1, NO
IF(IUU).GT.O) IU(I) = AMINO(IU(I).LH( ID ) )*FAC + SEALEV
IF( 1V( I ).GT.O) 1V( I ) = AMINO( IV( I ),LH( ID ) )*FAC + SEALEV
CONTINUE $ CALL PPTMAXC IZ,1HNONE, 1H .NANEZ, 1 , 30, 0 )
CALL PRTdAX< lUjlHNONE.lH . NAMEU , 1 , 29 . 0 )
CALL PRTnAXf IV,1HNONE,1H ,NAP1E« , 1, 29 ,0 )
IF( ICC.LT.O) GO TO 250
DO 182 J=l jNO
IZ( J ) = IZ( J )»10 $ IF( IZ( J).LE.O >IZ( J) = -0
IU( J ) = IU( J )»10 S IF( IU( J ).LE.O)IU( J ) = -0
IV( J ) = IV( J )*10 $ IFC IV( J ).LE.O)1V( J ) = -0
CONTINUE » CALL INPOUT( IZ, 1 , 1 ,0, 10 )
OUTPUT INITIAL COMPUTATIONAL FIELDS FOR Z U ftND V
DO 195 1 = 1. MAXLAY S IF( LH( I ) ,EO .0 ) GO TO 200
L = 0 S IFC I .EO.l ) GO TO 181
REWIND 2 S READ(2)IZ,IU,IV $ L=LH( 1-1 )
DO 190 J = 1,NO $ IFdZUl.LT.l) GO TO 185 * IZ(J) = 2 * GO-TO 186
HNI .81
HNI. 82
HNI .83
HNI .81
HNI .85
HNI. 86
HNI. 87
HNI .88
HNI .89
HNI. 90
HNI .91
HNI. 92
HNI. 93
HNI. 91
HNI. 95
HNI. 96
HNI. 97
HNI. 98
HNI .99
HNI. 100
HNI. 101
HNI .102
HNI. 103
HNI. 101
HNI .105
HNJ.19
HNI. 107
HNI. 108
HNI. 109
HNI .110
HNI. Ill
185 IZ(J)=-0
HNI. 112
-------
PROGHAB LAYER3S 7600-*7600 OPT=£
FTN 1.6+139/OtO
09 PIAR TT 01.15.264 BKY PAGE
115 186
187
188
189
190
120 195
200
250
1250
C
125 300
310
320
130 C
360
135 370
390
too
t05
ito
It05
900
1900
its
IF(IUCJ).LE.L) GO TO 187 4 IU( J )=2 4 GO TO 188
IUI J )=-0
1F( IV( JKLE.L) GO TO 189 4 IV(J) = 2 4 GO TO 190
IVC J >=:-0
CONTINUE 4 CALL INPOUTl IZ . 1^0 JO )
CONTINUE
ENDFILE1 4 ENDFILE1 4 BACKSPACEl
REWIND 2 4 READC2) 1Z.IU.IV 4 PRINT 1250 4 GO TO 111
FORMAT; IHI >
PUNCH NEU U AND V ARRAY CARDS
MM=MAXO; PIPIJ ) 4 P1P = PI1NO< PIP ,P1E ) 4 NPI=PIAXO( NPI,_L> 4 NP=PIINOC NP, NE )
DO 310 N = NP1,NP 4 PINP\=; PIPI-1 )»NE+N 4 PIPIP = ( PIP-1 )»NE + N
PUNCH 1059, C IU( n,I=PlPIPI,MPIP,NE)
DO 320 N = NPI,NP 4 PIP1P1=( PIPl-1 )»NE+N 4 PIPIP = ( PIP-1 >»NE + N
PUNCH 1059, ( 1V( I ),I = PIPin,PIPIP,NE) 4 GO TO 1
PUNCH OLD TYPE TABLES
DO 390 L=l_,3 4 LPI = ; L-l )»PIAXDIPI 4 INO = 0
DO 370 N=NM,NP 4 IP1=LP1+N 4 DO 370 PI=MPI,PIP
I = ln + (Pl-l )»NE 4 IND=IND+1 4 IF( IND .LE . 12 ) GO TO 370
PUNCH 1059. IN 4 IND=1
IN( IND )=IZ( I )
PUNCH 1059, ( INC I), 1 = 1, IND) 4 GOTO 1
IF(ICC.NE.O) SUP1 = 0. 4 DO t05 PliPinjlP 4 P1PIPI=( PIAXO( 1 .M-l )-l )«NE
NN = ( Pl-l )»NE4DOt05N = NPl,NP 4 NNm=MAXO( N-l , 1 ) 4 PIP1N=NN+N 4 PlnlPIN=Plt1PI-
SUM^i/UPI+FLOATCPlAXO; IU( ilflN ), 0 )+PIAXO( IV( MMN ), 0 HRAXOt IU( PIPIPIN ), 0 ) +
» nAXOC IV(NNPI).O)) 4 SUPI=SUP1*DDL*FAC
PRINT 1105, SUP1 4 GO TO 1
FORMAT; 3oxi6HVOLunE CUBIC cn=Eii.t>
PRINT 1900
FORPIAT; i2H END LAYERSSI
END
HNI. 113
HNI .lit
HNI. 115
HNI .116
HNJ.20
HNI . 118
HNI . 119
HNI .120
HNI .121
HNI. 122
HNI . 123
HNI . 12t
HNI . 125
HNI .126
HNI .127
HNI .128
HNI .129
HNI .130
HNI . 131
HNI .132
HNI. 133
HNI.Ut
HNI .135
HNl!l37
HNI .138
HNI .139
HNI. 140
HNI .142
HNI . It3
CARD NR. SEVERITY DETAILS
I
6
6
80
_SYPIBOJ._IC_
UJfiTER DE
WATER DE
LH
DIAGNOSIS OF PROBLEPI
_H DLL E R I_TH_C Q NSJ A N T .GT. 10 CHARflCTERS. EXC ES S_C.HA RAE_TER_S__I_N_IT IALIZ E 0_1 NT 0_SU C C E E 0 I N G_U 0 R D S .
HOLLERITH CONSTANT .GT. 10 CHARACTERS, EXC'ESS CHARACTERS INITIALIZED l"NTO SlTCCEEDTN'G WORDS'."
HOLLERITH CONSTANT .GT. 10 CHARACTERS, EXCESS CHARACTERS INITIALIZED INTO SUCCEEDING WORDS.
ARRAY NAP1E OPERAND NOT SUBSCRIPTED. FIRST ELEPIENT WILL BE USED.
-------
SUBROUTINE GRID 7600-7600 OPT=2 FTN 1.6+139/010 09 mAR T7 01.1$.26$ BKY PAGE
1 SUBROUTINE GRID(LH,ICC)
COPimON 1Z( 1000 ).IU( 1000 ),IV( 1000)
COPlRON/DATA/NE,rcE,NO,nAXDII'l/l'IAXLAY,inODEL,IBUF( 1250)
D1P1ENS10N LHC 1 )
5 CALL XflIT< -O.IZ.NO.0.1 )
REWIND 2
C SET HTU BOUNDARIES
DO 800 11 = 1 .P1AXLAY
LHM=0 $ IF(II.NE.l) LHP1=LH( 11-1 >
10 DO 110 N=1,NE * MP = (nE-l )»NE+N $ IFG=0
DO 100 I=N.MP.NE S IF< IU( I > .LE .LHM >IU< I >=-0$ IF< IFG .Ed .0 >GOT080
IF( 1U< 1 ).GT.O ) GO TO 100 $ IU( I ) = -! * IFG=0 S GO TO 100
80 IF( IU( 1 KLE.O ) GO TO 100 S IFG=-1 S IF(I.EQ.N) GO TO 100
1NE=I-NE $ IU( INE)=-1
15 100 CONTINUE
110 CONTINUE
C SET HTV BOUNDARIES
DO 150 n=l,PlE $ IP = (n-l)»NE $ ln=IP+l * IP = IP+NE $ IFG = 0
DO 110 1=IM,IP* IF( !V< I ).LE.LHP1)IV( I ) = -0 s IF( IFG -EQ .0 )GOT01 30
20 IF( IV( D.GT.O) GO TO 110 $ IV( I ) = -! S IFG=0 $ GO TO 110
130 IF( IV( I).LE.O) GO TO 110 $ IFG=-1 S 1F( I .NE . Id ) IV( 1-1 )=-!
110 CONTINUE
150 CONTINUE
DO 200 M=l,nE S IP = (ft-l)»NE
25 IWP1=IP-NE $ IF(PIIW.LT.O) MPim=IP
C SET HTZ BOUNDARIES
DO 200 N=I,NE s I=IP+N s NM=MAXO( N-i , i )+IP s mn=nnn+N
IF< IU IZ(I)=-1 $ GO TO 200
HNI.
HNI.
HNI.
HNI.
HNI.
HNI.
HNI.
HNI.
HNI.
HNI .
HNI .
HNI.
HNI .
HNI .
HNI .
HNI.
HNI .
HNI.
HNI.
HNI .
HNI .
HNI.
HNI .
HNI.
HNI.
HNI .
HNI.
HNI .
HNI.
HNI.
HNI.
HNI .
Ill
115
116
117
118
150
151
152
153
151
155
156
157
158
159
160
161
162
163
161
165
166
167
168
169
170
171
172
,173
,171
,175
167 IU(nn)=IU( I ) = -! * IF( II .EO.l .AND.I . NE .FIN >I Z< Hit ) = -! * GOTO 190HNI.176
170 1F( IV(NFI).GT.O.OR.IV( D.GT.O ) GO TO 190 HNI.177
IV(Nf1)=lV( I ) = -! ? 1F( II.EO.l.AND.I .ME.Nil) 1Z(NM) = -1
190 1Z( I )=II HNI.179
200 CONTINUE S IF(II.EQ.l) URITE(2) IU,IV HNI.180
_£ B Oil N D HY _E_RR OR. TES TING H NI .18_1_
DO 300 fn = l,ME * MP = C|ilINO(Pl + l,nE)-l )»NE S CIPI=< m-1 )*NE HNI.182
10 DO 300 N=1,NE * NCNT=0 $ NN=Pin+N HNI.183
F.< IUJ.NN )^EQ_,jilI NCNT=1 $ IF( 1V(NN).EQ.-1 )NCNT=NCNT*1 HNI.181
$ IF(lU(NN).EQ.-l) NCNT=NCNT+1 HNI.185
IV(NN).EQ.-1 ) NCNT=NCNT+1 HNI.186
IFJ N C N T_. E 0.1 .OR . N C NT.EO .3) PR IJJT_ 1_299_,_N .M,NCNT H N Kl 8J_
15 1299 FORPSATf 29H »*»ERRO~R LANDSEA »»*N,P1,NCNT3I1 ) HNI.1~88
300 CONTINUE S IF(ICC.EQ.O) CALL GEOPOS(II) HNI.189
800 CONTINUE HNI. 190
900 REWIND 2 $ READ (2) IU,W * LHP1=0 $ DO 910 1 = 1, NO HNI.191
LHi1=M6XO( IU( I ), IV< I >,LHn> HNI. 192
50 910 CONTINUED lCC=LHn » RETURN ? END HNI.193
-------
SUBROUTINE GEOPOS 7600-7600 OPT = 2 T™ 1.6+139/010 09 MAR 77 01.15.26$ BKY PAGE
1 SUBROUTINE GEOPOS(II) HNI.m
1Z( 1000),IU( 1000),IV( 10001 = H-NJ_._195
>£>
15 ), IBUK 15 ),IBU2( 15 HNI. 196
1 ),IBU3( 15 >,IBU1t 1190) ' ' HNI.197
DIMENSION nCL< 3 ) _ HN.LJ.?A.
~DATft"~(nCL = 29HTl5,HlH + 15(Ifc,lH+)//2X15I7/)),(LU=61) HNI .199
, ,
PRINT 1000 $ NOP = (ME-1)/15 + 1 $ MF = -H HNJ.21
1000 FORMftTUHO) , . ., . _ ... „_ . . Hr!Jj-4-94-
~DQ~i"00~U=T7NOP S MF = MF + 15 4 plL=niNO( PlF + lt ,ME ) * URITE( LU, 1001 ) IIHNI.20~2
10 1001 FORMftTC1H112X6HLAYER= 13) _ . HNI.203
UIRITE(LU inn;M, IBUZ( K>,K=l, l >
FORI-lflTC I1,15( 13,2X11 ))
URITE( LU,1016 X IBUFt K ),K = 1 , I)
UIRITEI LU.1016 )( IBUF( K >,K = 1, I )
WRITE(LU,1050)( IBUK K ), IBU3C K ), IBUK K ), K = l , I )
FORMftT( IH + SXlSfRTjJljRS) )
CONTINUE '* GO TO 900
ENTRY POGO PROVIDES U/V STRUCTURED OUTPUT FOR DEPTH VALUES, ...
CNTRV pnr,n s NOP=( ME+11 )/15»15-11 S PRINT 1000 $ K=15
DO 200 nF = i,NOP,15 * NN=(1F + H $ IF( NN .LE .PIE )GO TO 110
K=ME-P1F + 1 $ NN=C1E * GO TO 300
WRITE(LU 1110) II (I I=MF NN ) * NN=( NN-1 )*NE 4 P1M=( MF-1 )»NE
FORMAT( 1H1//1X31HUIATER DEPTH AT U/V POINTS LEVELI2// 15 I7/ )
DO 150 J = I,NE $ n=ntn+j * N=NN+J
UR!TE(LU.I
CONTINUE * IF(K.EQ.IS) GO TO 900 $ K=15
STATEMENTS LINE 300 ARE EXTENDED RANGE OF DO 200
ENCOOEC 30,1 300, P1CL ) nCL( 1 ), K ,nCL< 2 ),(nCL( 3 ) , K ,C1CL( 3 )
FORP1AT( A9,I2,R9,A3,I2,R5)
IF( K.NE.15 ) GO TO 110
FORHATC 1HR)
PRINT 1900 * RETURN * END
HNI .
HNI .
HNI.
HNI .
HNI .
HNI .
HNI .
HNI .
HNI .
HNI.
HNI .
HNI .
HNI .
HNI .
HNI.
HNI.
HNJ.
HNI .
HNI.
HNI.
HNI.
HNI.
HNI.
HNI.
HNI.
HNI.
HNI .
HNI.
HNI.
208
209
210
211
212
213
211
215
216
217
218
219
220
221
222
223
22
225
226
227
228
229
230
231
232
233
231
235
236
CARD NR. SEVERITY DETAILS DIAGNOSIS OF PROBLEP1
HOLLERITH CONSTANT .GT. 10 CHARACTERS, EXCESS CHARACTERS INITIALIZED INTO SUCCEEDING blORDS.
13 I 30 NON-INNER LOOP BEGINNING AT THIS CARD IS ENTERED FROfl OUTSIDE ITS RANGE.
-------
SUBROUTINE TIDCUR 7600-7600 OP1-2 FTN 1 .6+1 39/010 09 F1AR 77 01.15.26$ BKY PAGE
1 SUBROUTINE TIDCUR( ITS, ITE,P10DOUT, SE,PD, IOFF ) 4 COPIPION Z( 3000 )
COMMON /DATA/ I A< 6).A<1).XS(1),JF(10)
1,IC( 1232)
DIMENSION SE( 1 ),PD( I ),LINE( 13)
5 EQUIVALENCE ( JEJAXTIDE )
DATA (IB-1H ) ( ID = 10H > (LINE=13(1H ) ), ( NAXTIDE=1 ),
»( RPD=. 01715329252), ( JF = 97H< A1,R9) (2A1,R8) ( A2, Al ,R7 )( A3, A 1, R6
»< A1 .AljR5 >< AS^Al .R1)(A6.A1.R3)(A7.A1.R2)(A8.A1.R1XA9.A1»
REWIND 2 4 WRITE ( 2) Z
10 T=NS=nODOUT 4 IE = C ITE-ITS-1 1/NS+l 4 Xn = 0 . 4 Nl = ( IOFF-1 )/NS-l
DO 25 J = 1.JE $ K = J + JE 4 A( J ) = PD( K )
25 CONTINUE 4 TN=T/3600.
DO 30 1 = 1, IE 4 Z(I) = 0. 4 T=FLOAT( I+NI )»TC1
DO 28 J=1.JE 4 Z( I ) = COS( ( ( SE( J )*T-PD< J ) )»RPD ) )»A( J ) + Z( I )
HNI. 237
HNI .238
HNJ.23
HNI .210
HNI .211
HNI .212
1HNI.213
Mil .211
HNI .215
HNI .216
HNI .217
HNI .218
HNI. 219
HNI .250
15 28 CONTINUE 4 XS=ABS(Z(I>) 4 IF(XS.GT.XM) Xfl=XS HNI.251
30 CONTINUE 4 Tn=Tn»60. 4 DO 36 J=1,IE,12 4 K=J+11 i T=J 4 J J=TH»THN I . 25 2
36 PRINT 1036 .JJ .( Z( I ). I = J.K )4S=Xn/50. + .84IF( S.GE.2. )GO TO 704SSS=10 .HNI . 25 3
1036 FORMATi 112, ( 8X12F9.3/ ) )
61 SS = SSSr.SSS/10. 4 LS=1
20 63 IF< S.GE.SS )GO TO 694 IF( LS .EQ . 3 >GO TO 6 14SS=SS/2 . 4LS=LS+14GO TO 63
69 S=SS 4 GO TO 75
70 S=J=S
75 XS(1)=S»50.S XS( 3) = XS(1 )».5 4 XS(2) = -XS(3) 4 XS(1)=-XS(1)
PRINT 1075, JE,( SEC J ),PD( J ),A(J ),J = 1,JE),S 4S=1./S
25 1075 FORP1AT< 1HQ/1 12, 1 1X6HUK 0/H )5X1HK( 0 )1X5HH( CP1 )//( 20X3F9 . 3 ) )
PRINT 1077. XS,< ID,J=1, 13)
1077 FORMAT! 65X6HHEIGHT/2( 13XF9 . 3, 3X U8X 1H03X2C 16XF9 . 3 1/6X12A10, A6 )
IH=1H. 4 DO 80 J=13.113,25 4 GO TO 200
30 80 CONTINUE 4 DO 90 1 = 1, IE 4 Jr.p|INO( MAXK Z( I >»S+63 .5 , 1 . ), 126 >
IH=1HX 4 GO TO 200
85 IF(MOD( Ij3).NE.O) GO TO 88 4 T=I-1 4 JJ = T»Tfl
PRINT 1086, JJ, LINE 4 GO TO 89
1086 FORC1AT( I6,13A10)
35 88 PRINT lOSS^LINE
1088 FORflAT( 6X13A10)
89 IH = IB 4 IF(tnOD( J + l 2, 25 ) .EO . 0 ) IH = ID 4 GO TO 200
90 CONTINUE
I = J = 11818 4 DO 130 K = l,18 4 XH=0 . 4 T=FLOAT( I )/3600 .
10 DO 100 L = 1,JE 4 XP1=COS( ( SE(L )*T-PD(L))*RPD)»A(L) + XM
100 CONTINUE 4 IF(J.NE.l) GO TO 101 4 PRINT 1100. XM, I 4 GO TO 105
101 J=J/2
1100 FORP1AT(8H HEIGHT=E9 . 3, 7H OFFSETI6 )
105 IF(XM-.2) 110,130,120
15 110 I=I+J 4 GO TO 130
120 I=I-J
130 CONTINUE 4 REWIND 2 4 PRINT 1130 4 READ ( 2 ) Z 4 GO TO 900
1130 FORflATC 1HH)
C LINES 200+ ARE EXTENDED RANGE OF DO 80, 90
50 200 JU=( J-l )/10 + l 4 JP=P10D( J-l, 101 + 1 4 IF(JP.EQ.IO) 60 TO 220
IF( JP.EO.l ) GO TO 210
ENCODEC 10,JF( JP),L1NE( JU » LINEC JU ), IH,LINE( JU ) 4 GO TO 230
210 ENCODE( 10.JFC JP),LINE( JU ) ) IH,LINE(JU) 4 GO TO 230
220 ENCODE< 10,JF( JP ),LINE( JU ) ) LINE(JU),IH
55 230 IF( IH.EQ.ID.OR.IH.EQ.1H ) GO TO 90 4 1F( 1H .EO . 1HX ) 45,80
900 RETURN 4 END
HNI.251
HNI. 255
HNI. 256
HNI .257
HNI. 258
HNI. 259
HNI. 260
HNL.l
HNJ.21
HNL.2
HNJ.26
HNJ.27
HNJ.28
HNI. 267
HNI .268
HNI. 269
HNI. 270
HNI. 271
HNI. 272
HNJ.29
HNJ.30
HNI. 276
HNI. 277
HNI. 278
HNI. 279
HNL.3
HNI. 281
HNI. 282
HNI. 283
HNL.l
HNL.5
HNI. 285
HNI. 286
HNI. 287
HNI. 288
HNI. 289
HNI. 290
HNI. 291
HNI. 29 2
-------
SUBROUTINE INPOUT 7600-7600 OPT=2
FTN t.6+t39/010 09 MAR 77 01.15.26* BKV PASE
I
5
1
C
1000
10
C
15 100
110
120
20 130
110
150
170
25
175
180
190
30 195
500
600
SUBROUTINE !NPOUT( ft, LAYER, I TV , IT.FACOUT )
COMMON/DATA/ NEfIHE,NO ItAXDIfVMAXLAY, MODEL . IU< 1250 )
DIMENSION A( 1 )
INTEGEP A,FACOUT
OflTA( I REC=0),( PI 21 = 7777 7777B ) , ( Pll 2=77776 ), . OR . SH IFTC ME , 21 ). OR, MODEL
KP=MftXO( 1,MTN6( 3, 10-ITY ) )
00 600 K=1,KP $ KK=(K-1 )»MAXDIM
OUTPUT
NU=2 * IU7600 OPT = 2 FTN 1.6+139/010
SUBROUTINE IOERR< L , IREC ) S DIMENSION LAB(I)
DATA(LAB=10HEOF ERR INr10HPAR ERR IN,10HEOT ERROUT . 10HPAR ERR_OUT_)
,(K = 0)
PRINT 1000, LAB(L), IREC * IFCL-2) 1,3,2
FORNflT( A11.I10)
BACKSPACE 1 5 ENDFILE 1 * BACKSPACE 1 * 1FCL.EQ.3) 60 TO 1
ENOFILE 1 * BACKSPACE 1
IREC-IREC-1 $ K=K+1 $ IF(K.LT.20) RETURN
STOP
END
09 MAR 77 01.15.26$ BKY PAGE 1
UTIL. 55
UTIL, 56 ...
UTILA.5
UTILA.6
UTILA.7
UTILA.8
UTILA.9
UTILA.10
UTIL. 62
UTIL. 63
-------
SUBROUTINE XfllT
7600-7600 OPT=2
FTN M.6+H39/010
09 MAR 77 01.15.264 BKY PAGE
1 SUBROUTINE XM1T( IA, IB ,N, JA, JB ) ^DIMENSION IA( 1) . 1B( 1 >*l=J = l
DO 5 NN=1.N * 1BU)=IAU) $ J=J+JB * I=I+JA
5 CONTINUE $ RETURN
END
SUBROUTINE PRTMAX 7600-7600 OPT=2 FTN 1-. 6 +139/010
1 SUBROUTINE PRTflAX< S.D . I S,NAME1 ,LE,NC, IT )
DIMENSION S( 1),D( 1)
INTEGER S
5 COnnON/DATA/ NE. PIE. fJOrHAXDinfPIAXLAVTPIODELrIW( 1250)
DATA (LU=61 ),(nAXNAnE=8),(NCPW=10)
NU=(NC-1 )/NCPU)+l
NOP = (P1E-1 1/15 + 1 S nF=-11
CALL XniT(NAnEl,IW,NU,l,l) $CALL XHIT( 1H , IUK NUI + 1 J.MAXNAnE-NU.O,
10 DO 100 IK = 1,NOP * HF=P1F+15 S PlL=niNO( MF + lt ,ME )
IF(NE.LT.25.AND.MOD( IK .2) .EQ .0 . AND .IS.EQ. 1H ) GO TO 30
20 WR1TE(LU,1000)( IU( I ), 1 = 1 ,P1AXNAME ) ,LE, IT
1000 FORFIATC 1H 120X8A10, 1 3, 7X5HTlnE=I 10 )
30 URITEfLU.1001 )( I, I=rtF,flL)
15 1001 FORMAT< //6X15I8 )
Pl(nF=CplF-l )*NE * PIL=(P1L-1 )»NE $ IF( D( 1 ) .NE.tHNONE ) GO TO 60
DO 50 J=1.NE $ JK=J+nMF * NN=ML+J
50 WRITECLU^OOOJIS^.CSCKJ.KsJK.NN.NE)
2000 FORCIAK Al, It, 1518)
20 60 CONTINUE
100 CONTINUE $ RETURN
END
HNG5.2B
HNG5.29
HNG5.30
UTIL.66
09 MAR 77 01.15.26* BKY PAGE 1
HNG8.127
UTIL.68
UTIL.69
HNJ.t
KMJ.5
HNJ.6
UTIL.73
UTIL.76
1 1HNG6.13
UTIL.77
UTIL.78
HHG3.128
HNG5.35
UTIL.81
UTIL.82
UTIL.83
UTIL.81
UTIL.85
HNJ.7
NNJ.8
UTIL.97
UTIL.98
-------
LOAD NAP.
LINK - BKV 6000/7000 8.1
09 MflR 77 01.15.30
PAGE
CT>
UJ
FL REQUIRED
FL REQUIRED
INITIAL IRAN
TO LOAD 11100
TO RUN 34300
SFER TO LAVER3S
7726
BLOCK ASSIGNMENTS.
BLOCK
/DftTA/
LAYER3S
GRID
GEOPOS
TIDCUR
1NPOUT
IOERR
XM1T
PRTMAX
/STP.END/
/FCL.C./
/OLINLNT/
/QJOBSF.L7 .
/OASAFLS/
/.QCOP1P./
/Q8.IO./
FTN1LIB
BACKSP=
BUFIO=
BUFOUT=
C10=
conio=
DECODt=
ENCODE:
ENDFIL=
EOF
ER[>1S4
FECP1SK =
FLTIN=
FLTOUT=
FMTAP=
FORSYS=
/QOVERL7/
FORUTL=
6ETFIT=
GOTOE8=
_ _INCOP1 =
/IO.BUF./
INPB =
INPC =
KODER=
KRAKER=
OUTB =
ADDRESS
100
2150
11775
12352
13056
13710
11310
11220
11210
11517
11520
11511
11512
11511
11515
11516
15012
15012
15071
15210
15315
15330
15111
15526
15672
15762
16031
17522
17563
17711
20256
20650
21631
21632
21650
21713
21727
22211
22110
23023
23271
23760
21115
LENGTH
2350
7325
355
501
632
230
60
20
257
1
21
1
2
1
1
211
0
62
111
55
13
61
115
111
70
52
1166
11
156
315
372
761
1
16
13
11
262
227
363
216
167
135
277
FILE
LGO
LGO
LGO
LGO
LGO
LGO
LGQ
LGO
FTN1LIB
FTN1LIB
FTN1LIB
FTN1LIB
FTN1LIB
FTN1LIB
FTN1LIB
FTN1L1B
FTN1LIB
FTN1LIB
FTN1LIB
FTN1LIB
FTN1LIB
FTN1L1B
FTN1LIB
FTN1LIB
FTN1LIB
FTN1LIB
FTN1LIB
FTN1LIB
FTN1LIB
FTN1LIB
FTN1LIB
FTN1LIB
FTN1LIB
21711
201
FTN1L1B
-------
LOflD NAP.
LINK - BKY 6000/7000 8.1
09 NAB 77 01.15.30
PAGE
CT\
BLOCK
OUTC =
RDL =
RDO=
RDU=
REUIND=
SINCOS=
SORT
SVSA1D:
SY5=
SVS=1ST
UNIT
UZERn..
UITH =
UITL =
UITU =
II
TEST BOX MODEL
5 53333i.am
5 5 .333E+05 0
1200 7777777
1200 7777777
1 77777
1 77777.
15 nuifiAotn
.150E+02
5515
5515 0
19
259
199
739
979
1219
1
-25 .000
39 !
99 '.
ADDRESS
25120
25355
25 til
25t37
25516
25621
25715
25763
2576H
26020
26102
261&7
26,171
26221
26250
26336
0.
0
0
0.
.94303562
.139E+02
117.
117.00 1
3.890
1 .t28
-2.821
7.307
-.135
-9.0tt
U( 0/H )
15.011
13.9t3
28.981
30.000
.500
LENGTH
235
37
23
107
53
71
16
1
31
62
65
2
30
27
66
5670
3.2E-6
.320E-05
86100
86100
0.
8 981 101 2
.290E+02
121. 359
21.00 359.
1.566
551
-2.217
7.677
-1 .508
-8.580
K(0)
117.000
121.000
359.000
52.000
-
FILE
FTN1LIB
FTNtLIB
FTfilLIB
FTN1L1B
FTN1L1B
FTNtLIB
FTN1LIB
FTNtLIB
FTW1L1R
FTN1L1B
FTN1LIB
FTM1L1B
FTN1LIB
FTN1L1B
FTN1LIB
NAME
12 003 GRID
.120E+02 .300E-02 0. .980E+03
100 0 300 T1P1E
100 0 300
.99 1.022 LAYS
.9900 0000 0001 0220.0000.0000.0000.000
30 FACT
.300E+02 0. 0.
052. 1.0 2.0 3.0 1.0 TIDE
no •;? nn i nn ?.oo 3.00 i. 000. 00.0
"5.078 " V.tit 5.570 5.516 5.318 1.987 1.178
- 311 -1 128 -1 872 -2-511 -3.033 -3.110 -3.629
-1.179 -.629 .309 1.306 2.332 3.353 1.338
7.860 7.816 7.629 7.211 6.598 5.802 1.812
-2 862 -1 161 -5 379 -6.176 -7.127 -8.207 -8.796
-7.927 -7.103 -6-133 -5.017 -3.875 -2.651 -1.109
men)
1.000
2.000
3.000
t.ooo
HEIGHT
12.500 0 l*.»uu
X
X
X
X
X
X
X
3.813 3.101 2.289
-3.681 -3.562 -3.273
5.255 6.071 6.766
3.710 2.522 1.220
-9.180 -9.350 -9.301
-.181 .989 2.079
25 .000
-------
159
219
279
339
.X
X.
399
-------
HEI6HT= .200E+00 OFFSET 33tfct
1
2
3
*
5
1
2
3
t
5
•)•)
It
30
30
30
35
to
35
35
35
to
to
20
20
25
30
35
25
25
30
35
35
10
10
20
25
35
15
20
25
30
30
0
5
15
25
35
5
10
20
30
30
COflP
NEWT
0
5
15
25
35
0
10
20
30
30
GRDZ
LAYER= 1
i
2
3
t
5
MAX
MAX
1 ' i
1
1 1
1
1 1
1
1 1
1
1 1
1
TINE
TIME
1 1
1 1
1
1 1
1
1 1
1
1 1
1
1 1
1
STEP =
STEP =
1 t
111 -10
1 1 *»*»»»
11111
1 1 1
11111
I 1 1
11111
1 1 1
11111
1
1
1
1
1 1 1
2.70 MAX DE?TH= T77T7 .
119.03 MAX DEPTH= tO.
1
VOLM
VOLUME CUBIC CM= .1028E+12
2 $
0
5 5
1
VOLUME CUBIC CM= .1097E+12
1
VOLP1
TAPE
-------
WATER DEPTH AT U/V POINTS LEVEL 1
12345
1 * 30 +
35 25
2 + 30 +
35 25
3 + 30 +
35 30
1 + 35 +
tO 35
5 » tO+
10 35
CTi
1
1 1
2 1
3 I
t 1
5 1
1
1 30
2 30
3 30
t 35
5 tO
1
1 35
2 35
3 35
1 to
5 10
INPOUT
INPQUT.
20+ 10
15
20+ 10
+ -1+ 0+
5 -1
+ 5+ 5 +
20 10 10
25+ 20+ 15+ 15 +
25
30+ 25
30
35+ 35
30
20 20
+ 25+ 25 +
30 30
+ 35+ 35 +
30 30
SYMBOLIC Z WATER HEIGHT TABLE
2
1
1
1
1
1
WATER
Z
20
20
25
30
35
WATER
Z
25
25
30
35
35
0
0
3 t 5
1 1 -1
111
1 1 1
1 1 1
1 1 1
DEPTH AT U POINTS (CM)
3 t 5
10 -1 0
10 5 5
20 15 15
25 25 25
35 35 35
DEPTH AT V POINTS (CM)
3 t 5
15 5 -1
20 10 10
25 20 20
30 30 30
30 30 30
1 1 5 5 99 0
1 2 5 5 99 3
1 TIME= 0
1 TIME= 0
1 TIP1E= 0
LRVF.R3S
-------
APPENDIX B
PHASE II PROGRAM NOTES
DESCRIPTION OF PROGRAM LAYER3S
This program was originally written for the CDC 6500 FORTRAN EXTENDED
compiler (7). The program includes the main program LAYERS and the routines
SETHUV, BAR, XMIT, CXMIT, BARUV, UVCAL, UVSTAR, PRTMAX, INPOUT and IOERR.
The program performs numerical computations simulating hydrodynamic processes
by difference solution of the basic hydrodynamic equations described in (1) .
To describe the program arbitrary logical blocks have been used which
were defined in Figure 4. The sequence numbers refer directly to the UP-
DATE sequence numbers in the program listings at the end of this appendix to
the right of the FORTRAN code.
This program uses the following logical units:
LU Description
INPUT Card input file
OUTPUT Printed output file
1 Intermediate data tape file
3 =OUTPUT
Initialization (Block A;HNG.2,78)
This block provides the linkage to subroutines for standard control para-
meters in COMMON/CONTROL/. To minimize the number of parameters passed in
subroutine calling sequences the order of the arrays in the control common is
essential as is the MAXDIM size of each. The dimension of the NNM, NNP arrays
must be greater than or equal to NE. Additionally the 500 positions following
68
-------
A6 on card sequence HNG8.1 are presumed destroyed through scratch usage by
subroutines, so all permanent usage variables precede variable NM on the
next card. COMMON// also has order dependent structure: the seven main compu-
tational arrays must be in the sequence shown and size MAXDIM £ NO for proper
program function. A program comment shows the required dimension of the dummy
array variable ZLAYK when multiple layers are desired. ZLAYK(7*MAXDIM*(K-l))=
Z2(MAXDIM),...,HGV2(MAXDIM),Z (MAXDIM),...,HGV (MAXDIM) where Z2 represents
the second layer Z array, etc. Zi, etc., will be used in subsequent computa-
tional descriptions to imply any and all layers of the indicated computation.
Three single layer basic array sizes have been provided on the library with
the primary (active) size being 300 positions (refer to TABLE 6). The rela-
tive position of the XK, YK arrays and MAXDIM size is also important.
After initialization of run constants and variables, each master control
card is read and printed in "A" format, identified by name in the index table
NAME and routed to subsections of the initialization block for DECODE inter-
pretation and preliminary processing. The following conditions cause immedi-
ate termination of the run during master control card processing: no COMP
card (end of INPUT file); unidentified control card; or too many TIDE, WIND,
SORC or FLOW cards.
After normal processing of each master control card except COMP control
is returned to read the next control card. Following the COMP card routine
INPOUT is used to obtain the type 1 data records from LU1 for the HTZ, HTU and
HTV arrays. Boundary condition 4 defined in (1) is applied to the depth
fields via CXMIT call and then additional program constants are defined in-
cluding the uppermost and lowermost sea positions for each column of the sym-
bolic array for each layer. These column sea limit values are used to
69
-------
minimize the range of computation DO loops. The layer dependent initial or
restart type 2 Z, U and V fields are read using INPOUT, constants are comput-
ed and the thicknesses HGU, HGV are computed in SETHUV. The program termin-
ates if records matching the header computed from the INPOUT call are not
found on LUl.
£ Computation and Prescriptions (Block B:HNG8.26,133)
After initiating the time loop, routine BAR is called to compute £ into
array ZMi; then ?i is computed into array Zi. Boundary condition 5 is im-
plicit in the £i, condition 3 is applied during the computation for U and V
values outside the open boundary and condition 4 is applied to Hu and Hv
values outside the open boundary. Conditions 3, 4 and 5 as applied consti-
tute the dynamic form of prescribing £ under condition 2. Any available ex-
ternal tidal information replaces the dynamically computed values to give the
alternate form of condition 2 prescription. Volume sources are then pre-
scribed as modifications to positions of Zi as specified on SORC cards. The
first time or any time thereafter specified for a wind to start or end the
wind support section defines the wind stress effect fields in XK and YK. Any
subsequent wind card must start at or subsequent to termination of the prior
wind specification or at the same time as the prior wind specification. If
a wind is specified to start prior to the end of a prior WIND termination
specification, but after the prior wind initiation, the wind specified will
not begin until the prior wind is terminated. Currently up to four constant
wind subfields may be applied simultaneously, or in sequence. Modification to
extend the maximum number of WIND, SORC, FLOW and TIDE cards can be accom-
plished by simply changing the array dimensions and the MAX variable asso-
ciated with it in the DATA statement.
70
-------
Hu, Hv Computation (Block C:HNG.134,BC8.6)
Routine SETHUV is called to perform the computation of both Hu and Hv
fields into arrays HGUi and HGVi, respectively. Routine CXMIT is employed
to apply boundary condition 4 at the lower boundary of Hvi fields and at the
right boundary of Hui fields.
U, V Computation, Prescription (Block D;BARMOD.1,HNG8.48)
* —
Routine UVSTAR is called to compute Vi into array ZSTi. Ui is calcu-
lated into the Zmi array by calling routine BARUV. The routine UVCAL call com-
putes Ui into the ZMi array, whereupon routine CXMIT is used to establish
boundary condition 3 on the right boundary of the ZMi array. The ZSTi array is
*
then used to store the Ui values computed by call to routine UVSTAR. The Ui
values from array ZMi are transferred into array Ui by calling routine XMIT.
Routine BARUV is then employed to calculate Vi into the ZMi array. Vi values
are computed into array Vi by calling UVCAL, whereupon boundary condition 3 on
the lower boundary of Vi is effected by call to CXMIT.
Having now computed the Ui and Vi values, incremental FLOW card specifi-
cations are added to the fields. Actual velocity values could be specified
through simple modification, though there has been little success in specifying
velocity in this or the incremental manner due to the dearth of external infor-
mation and inadequate theoretical tools.
Output Processing (Block E:HNG.191,246)
In accordance with output specification on the TIME card, routine INPOUT
is called to save Zi, Ui and Vi arrays on LUl. Similarly the TIME card spec-
ifies the printed output interval for the Z and vector velocity fields. Rou-
tine PRTMAX is called for each layer to print these fields on file OUTPUT. At
the end of this block determination is made whether the computation is over
based on the end time on the TIME card. Normal completion of.the program
71
-------
is indicated by END LAYERS in the job DAYFILE. If the job ends abnormally
either STOP or some other message will appear in the dayfile. If the message
is not indicative as to the cause, check the input parameters and dimensions.
Is NO greater than MAXDIM? Is MAXDIM equal to the dimension of all major
arrays? Is NE greater than the dimension of NNP and NNM? Is ZLAYK sufficient
to hold multiple layer arrays?
Routine SETHUV (Block F:HNG.249,HNG8.71)
Routine SETHUV requires no formal parameters as all required variables ap-
pear in either /CONTROL/ or //common. The routine depends on the sequence and
size of HTU, HTV and HGU, HGV in the single layer case and on the entire struc-
ture of //common for multilayer cases. HTV and HGV are never referenced by
name in this routine, only by index modification to name references to HTU and
HGU respectively, similarly in two layer cases Z , HGU are not referenced ex-
£ £
cept through Z, HGU by index modification. (Z and HGU are the first and sixth
set of MAXDIM memory positions in array ZLAYK for the 2 layer case.) One call
causes computation of all layers of HGU, followed by computation of all layers
of HGV. No analytic boundary conditions are imposed during the thickness compu-
tation per se. The Hv computation is v-symmetric; the Hu computation is u-sym-
metric. Computations on the last row of Hv and the last column of Hu are re-
placed externally as boundary condition 5 solely for dynamic £ computation.
Routine BAR (Block:HNG.105,BARSTAR.36)
Routine BAR requires no formal parameters and refers only to HTZ and first
layer arrays Z and ZM: reference to succeeding layers when required is by index
modification from the first layer array. For each position of array Zi the
four nearest Zi points are used with the central value to compute a smooth-
ed value which is returned in array ZMi. Positions over land (HTZ=0) are not
defined in BAR since they are preset to zero. Further, no use of the ZMi array
72
-------
changes the zero values over land: if the U, U or V value of the same index is
defined, the ZMi point is defined, since HTZ^O). Boundary condition 5 is ef-
fected by index substitution prior to the computation for both land-sea and
open boundaries. Smoothing of the central value is interpreted as eddy dif-
fusivity effect over the last time step. The computation is symmetric in space
around the central point except at boundaries where computation is quasi-
symmetric. Independent smoothing parameters may be specified at each level on
the LAYS card.
Routines CXMIT, XMIT (Block H:HNG5.28,BC8.16)
Routine CXMIT provides a special controlled transfer. Routine XMIT is
the same as that described in Phase I program notes. Routine CXMIT is similar
to XMIT in function with the extended feature that transfer as described for
XMIT occurs from A to B if and only if C is greater than D. The last parameter
of CXMIT provides an index adjustment to the array C so that a control array
(such as HTU) may be referenced to assure positive depth before transfer is
effected. The parameters N, JA and JB function precisely as do the parameters
of the same name in routine XMIT.
Routine BARUV (Block I:BARMOD.8,BARSTAR.5Q)
Routine BARUV is similar to routine BAR in function since the routine com-
putes smoothed u or v velocity component fields to achieve eddy viscosity ef-
fect for the last time step. Three parameters are required: the source array,
the control array and a key to indicate whether u or v type computations are
being performed.
For each position of array Bi the four nearest Bi points are used with
the central value to compute a smoothed value which is returned in array ZMi.
Positions over land (Hi=0) are not defined in BARUV for the same reasons given
in routine BAR. In actuality some values defined by the BAR routine are reset
73
-------
to zero by BARUV to insure lower v land-sea boundary points and right hand u
land-sea boundary points remain zero. At normal land-sea boundaries, boun-
dary condition 1 is effected. At tangential open or closed boundary points
condition 6 is effected. The computation is symmetric in space except at
open or closed boundaries where it is quasi-symmetric. Independent smoothing
parameters may be specified at each level on the LAYS card.
Routine UVCAL (Block J:HNG8.72,94)
Routine UVCAL computes the u or v velocity component values for the next
time step using two input parameters: the source array (UV) and the key (Al)
to indicate u or v computations. The UVi array is used to store either the Ui
or Vi computational values calculated from the ZMi, ZSTi and the HGUi (HGVi)
arrays. References to HGUi also access required HGVi values using index mod-
ification that relies on the array structure of blank common. Wind stress ef-
fect terms are applied to the topmost layer of thickness greater than zero.
In the present version the friction term is constrained to be less than or
equal to the smoothed velocity component. The smoothed velocity component is
included in the friction term computation where a prior time step velocity com-
ponent is expected. These numerical deviations showed no significant effects
on experimental runs.
Boundary conditions 3 and 6 at open boundaries are implicit in the ZMi and
ZSTi arrays. At closed boundaries condition 6 is applied by index modification
and condition 1 is applied by default since zero velocity values persist at
land-sea boundaries. Although values are computed in the last column of the Ui
calculation and in the last row of the Vi computation, they have no impact on
the calculation for this or future steps as these values are replaced externally
prior to use in the £ computation. The computation is quasi-symmetric in space
near boundaries and symmetric otherwise, after the external boundary adjustment.
74
-------
Routine UVSTAR (Block K:HNG8.95,BARSTAR.23)
* *
The routine computes the "star" terms Vi and Ui by linearly interpolating
each velocity component to the computational grid points of the normal compon-
ent. The star terms are used in the coriolis and frictional effect terms of
the velocity component computations. Routine UVSTAR requires three parameters
designated UV, H and N. The first provides access to the variable to be in-
terpolated, the second is no longer required and the third provides a key
used to differentiate between v or u interpolation. The interpolated values
from the UVi array are stored in the ZSTi array. All zero values at land-
sea boundaries are used in the computation making the computation symmetric
in the general case and at land-sea boundaries. The routine is quasi-symmet-
ric at open boundaries where condition 3 is applied by index substitution.
Routine PRTMAX (Block L:HNG8.127,UTIL.98)
Routine PRTMAX prints a simple floating point array with a title row and
column identification similar to PRTMAX of Appendix A. It also prints a
velocity field as speed and geographic direction values separated by a comma.
Seven parameters are required: S, an array of floating point values or of v
velocity components; D, either an indicator that no array is present or the
u velocity component; IS, the carriage control for printer vertical spacing;
NAMEl, the variable length title in display code; LE, the level associated
with the array printed in the title line; NC, the number of display code
characters NAMEl; and IT, the time printed on the title line.
The array is printed in segments of 15 columns by NE rows, except the
last segment which prints only to the ME column. When NE £ 25, ME £ 15 and
single spacing is selected, two segments of up to 15 columns are printed on a
single page.
75
-------
Routines INPOUT; IOERR (Block M;PR3.6,UTIL.63)
Routine INPOUT provides access to the intermediate data set (LUl) . It
uses the last 500 positions of control common for an input-output buffer and
depends on the blank common structure. INPOUT uses five parameters: A, the ar-
ray or arrays to be processed; LAYER, the layer associated with the arrays
(positive layer number indicates output; negative layer number indicates in-
put) : ITY, type of data; IT, the time associated with the arrays; and FACOUT,
a factor which for most data types is used to scale the output data and to re-
verse the output scaling on input. Output data are packed in a compressed for-
mat with a 24 bit byte assigned to each data field, and input data are unpack-
ed expanding the data field into a full 60 bit word for most data types. Out-
put arrays are partitioned into 1245-byte or less sections, resulting in a
maximum physical records size of 500 words. Each physical record has a two
word header constructed or verified from input parameters and control common
elements.
Routine IOERR provides error processing when normal access to logical unit
one was not accomplished by routine INPOUT. The routine requires two para-
meters: L, an index which indicates the particular problem encountered; and
IREC, the record number where the error occured. The routine was designed
primarily to handle physical tape problems. Parity conditions should never
occur when the program reads and writes to disk. The routine prints a single
line to indicate the problem and the record number. Input end of file condi-
tion causes program termination after the printed line. A total of twenty in-
put or output parity errors is allowed prior to termination. Output end of
tape condition causes a backspace over the last record output, an end of file
to be written, another backspace and termination.
76
-------
LISTING FOR PHASE II SAMPLE RUN
The list is a reproduction of a Phase II run. The first page of the
listing was produced by the UPDATE processing the library file with list
option L=124 as shown in the TABLE 8 example of SCOPE control cards. The
FORTRAN compiler listing of Phase II shows the UPDATE module identification
and card number to the right. The run results are presented immediately
following the listing of routine PRTMAX starting withthe NAME card. The
NAME through COMP lines are control card images. They are followed by the
INPOUT lines printed by the initialization calls to access the intermediate
data set file. The following pages show the output of the time loop execu-
tion. Water heights are printed as specified on the time card with veloci-
ties following. The INPOUT line reflects the call to save data on the inter-
mediate data set.
77
-------
CURDS ENCOUNTERED IN INPUT
BKY 76/66 UPDATE 1.2/110-1H 09 MAR 77 01.15.30
PAGE 1
-J
oo
»».»» .COMPILE HNG. UTIL
MODIFICATIONS / DIRECTIVES
YANKS*$
YANXSSS
YANK$$$
YANKS44
YANK$$*
YANKSSS
YANKSSS
YANK4S*
YANKSS*
YftNKSS*
HNG
HNG
HNG
HNG
HNG
HNG
BARSTAR
BARSTAR
BARSTAR
BARSTAR
BARSTAR
BARSTAR
UTIL
UTIL
• YANK VARUIINE
•YANK CONPRV
•YANK HNG7
•YANK VARCORL
•YANK VARU1ND
•YANK VALDEZ
•YANK CONVCT.CONVCU
•YANK CONPRT,CONPRU
•YANK FLAT, PLATA
•YANK FLATIO.FLATIP
•YANK MONTE
•YANK HNJ
•YANK BKYPRT
•CALL CONTRL
•CALL BLANK
•CALL CONTRL
•CALL BLANK
•CALL CONTRL
•CALL BLANK
•CALL CONTRL
•CALL BLANK
•CALL CONTRL
•CALL BLANK
•CALL CONTRL
•CALL BLANK
•CALL CONTRL
•CALL CONTRL
AVARUI2
ACONV2
AHNGT1
AVARC1
AVARU1
AVALDEZ
ACQNV1
ACONV1
AFLAT1
AFLAT1
AMONT1
AHNIJ1
AHNIJ1
HNG
HNG
HNG
HNG
HNG8
HNG8
BARSTAR
HNG8
BARSTAR
HNG8
BARSTAR
HNG8
UTIL
UTIL
1
1
1
1
1
1
1
2
1
2
1
1
2
3
7
250
252
7t
76
3
97
25
107
38
117
3
70
-------
PROGRAM LAYER3 7600->7600 OPT=2
FTN t.6+t39/OtO
09 P1AR 77 01.15.3U BKY PAGE
1 PROGRAM LAYER3( INPUT, OUTPUT, TAPE1, TAPE3=OUTPUT )
COMMON /CONTROL/HTZC 300 ) . HTU( 300 ) . HTV( 300 ).NNP( 20 >.NNM< 20 ).NE.
lME,NO,MAXDIM,MftXLAY,MODEL,ROT
A,ZZ(3),HLO( 1 ),HL(3),ALPHA( 3),BETA( 3),LEV( 3),A5RHO( 3 ) , A5RH( 3 ) , A3,
5 2.NFl,!in.MP,IPlM.II .11 . I 2 . IM r IP .TK . XU, YU .MAXNAWE .MAXSORC .MAXFLOW .
3MAXT1DE LM,LP,ICi( 9),NAME( 10) IUI(t63)
C DIMENSION ZLAYK(MAXDIM»7*(K-1 ) ) (WHERE K - ACTUAL NO. LAYERS).
COMMON ;/ Z( 300), U< 300), V( 300), ZM( 300).ZST( 300)
A ,HGU( 300),HGV( 300)
10 K ,ZLAYK( 1 )
U .XK( 300)fYK( 300)
l.COR
DIMENSION IUIND(8,t ),rfIDE(5,t ),TIDE( 10,t ),IC(8)
1 .TIOSPD! t >.IFLOU( 8,30 >,FLOW( t,30 IjISORCCS.lO ).SORC( 2.10),
15 2RHOO( 1 ),RHO( 3)
EQUIVALENCE! A3,R),.EO.NAmE(lM GOTO 3
2 CONTINUE 4 GO TO 900
3 GO TO ( l,t, 5, 6, 10, 16, 18, 20, 21, 31 ),I
35 t DECODE! 76, 1004. 1C )NE.f!ErDL .ROTANG,C P ALAT,R r AUS. CMPDLftT
lOOt FORfiftT( 2I3,7E10.3)
ROT=810.-ROTANG $ TK=ALAT*RPD $ NO=NE*ME
TK = COS( 2.»TK ) $G=((TK*5.9E-6-.0026373)*TK + l . )»980.616 * GO TO 1
5 DECODE! 76, 1005, 1C )MODOUT, ITO, ITS, ITE , ITINC, IPO, IPMQD
tO IF(MODOUT.EQ.O) MODOUT=600 $ IF( 1PP10D .EQ .0 )IPMOD = 3600 $ GO TO 1
1005 FORMAT! 16^7110)
6 DECODE(76, 1006, 1C )TIDSPD,TIDLAG, SLOPEX
1006 FORMAT(6X,7E10.0)
DO 7 K = 1.ITCN * TIDSPD(K)=TIDSPD(K)»RPD/3600.
t5 7 CONTINUE * GO TO 1
10 ITNO=ITNO+1 $ IF( ITNO.GT.HAXTIDE) 60 TO 900
OECODF! 76. 10 10. 1C )( ITIDE( I, UNO), 1 = 1.5 ).(TIOE( 1 r ITNO ) r 1 = 1 . 10 )
1010 FORMAT(tI3,I6,tF7.2,tF6.2,2F3.1 >
DO 11 K=1,ITCN * TIDE(K,1TNO)=TIOE(K,ITWO)»RPD
50 11 CONTINUE * SO TO 1
16 IUNO=IUNO+1 « SF( lUNO.GT.flAXWIND ) GO TO 900
DECODE! 76, 1016, ICM IW1ND( I , IUNO >, 1 = 1,8 > * GO TO 1
1016 FORriAT(tI3,tX6I10)
18 DECODEC76, 1018, 1C )HL, ALPHA, RHO « GO TO 1
55 1018 FORmAT(6X3F10.0.8F5.3)
20 ISNO=ISNO+1 * IF< ISNO.GT.PIAXSORC ) 60 TO 900
HNG.2
CONTRM.l
HNG6.1
A6HNG8.1
HNG0.3
HNGB.t
BLANKA.l
BLArJKfl.2
BLANKS. 3
BLANKA.t
Bl.ftNKA.5
HNG8.2
HNG.8
HNGS.5
HMG8.6
HNG8.7
HNG8.8
i) HNG5.9
) HNG8.9
HMG8. 10
HiMG.lt
HNG.15
HMG81.2
t:;;G3l.3
HW68.11
HN6.18
HNG.19
BLANKA.6
HNG.20
HFJG.21
HNG.22
HNG.23
HNG.2t
HNGt. 3
HNG8.12
HNG.27
HNG8.13
HNG5.12
HNG6.5
HNG6.6
HNG.31
HNG6.7
HNG6.8
HNGt .6
HNGt. 7
HNGt. 8
HNGt .9
HNGt. 10
HNGt .11
HNGt .12
HNGt .13
HNGt.lt
HNGt. 15
HNG8.lt
HNGt. 17
HNG.t3
DECODE! 76,1020.SCXISORCCI,ISNO),1 = 1,8),(SOHC(I,ISNO),1 = 1,2)
HNG.tt
-------
PROGRAm LAYER3 7600-7600 OPT = 2
FTN M .6+N39/OMO
09 PIAR 77 01.15.316 BKV PAGE
CO
o
niNSrnC=MlNO(niNSORC, I50RC( 7,ISNO) )* GO TO 1
21 1FNO=IFNO+1 S IF< IFNO.GT.PIAXFLOU) GO TO 900
60 DECODE! 76 , 1 020 , 1C M I F LOU! 1 , 1FNO > , I = 1 , 8 > , 1 t- LOW! 1 , 1FNO > , I = 1 , t>
1020 FORIUT
TK = ( ROT-FLOW! I IFNO))»RPO
FLOUC 1, IFNO) = C05( TK >«FLOU( 2, 1FNO)
FLOUI(Z,IFNO> = SIN!TK)»FLOU(2,lFfJO>
65 MR'FI.O'.l=nlNO( nINFLOW_, 1FLOW! T^IFNO) )5 GO TO 1
31 DECODE! 76, 1005, 1C 1P100EL t CALL 1NPOUT! HTZ , - 1 , 1 , 0 , 1 0 . )
CALL C ~t PIIT< H TU( N 0-NE ), H TU( NO >,HTU! NO ),0. ,NE,-1 ,-!,-!)
CALL C'.'rilT! HTV! NO-1 ).HTV(NO) HT V< NO ) . 0 . . PIE . -NE -NE -NE )
A1=ITINC » A3=R»A1 S AM=A1/DL S A5=G*A1 S A6-A3».5 $C 3=C»A 1» 10000 .
70 Dl-1 ./( DL»DL ) * CFAC-10. S IF(ITS.NE.O) CFAC=1000.
DO 33 I=ln, IP
IF( HTZ! I > .GT.O. .OR.HTU! I ) .GT.O . .OH .HTV( I ) .GT.O. >GO TO 32
lF(K!if'Hn).EO.I ) NNn $ IF(IP.LT.IM) £0 TO 5_1
90 IMfl=(pi-l )»NE + 1 $ DO 50 I = Id,IP S IF! HTZ! I ) .LE .0 . ) GO TO 50 $ ZZ = 0
MCI=I-NE s iF(nn.LE.O) nn-i s NM=MAXO( 1-1, inni $ JL=HTZ(I)
tl L-LEV1JL) i II=L+I 6 I1=L+MM S I2=L+Nn
ZZ = (HGU( II )«U( II )-HGU( II )»U! 11 HHGV( 12)*V( 12)-HGV( II )»V( II ) ) + ZZ
t5 Z( in = Z,1( II >-ZZ»At 4 IF(JL.EO.l) GO TO 50 $ JL=JL-1 6 GO TO 1 1
95 50 CCWTKx'UE
HNG.H5
HNG.M6
HNG8.15
HNG6.11
HNG5 .15
HNG5 .16
HNG. 5 2
PR3.1
BC8.1
BC8.2
HNG8. 16
PR3.2
K(OG .55 _
HNG .56
HNG. 5 7
HNG. 5 8
HNG. 59
HNG. 60
HNG8.18
PR3.3
HNG8.20
HNG8.21
HNG8.22
HNG8.23
HNG8.2H
PR3.M
HNG8.25
fHNG.77
HNG.7*
HNG8.26
HNG8.27
.HNG8.28
HNG8.29
HNG8.3.0
HMG8.31
HNG8.32
HNG. 96
51 CONTINUE HNG. 97
IF! ITNO.EO.O) GO TO 70 S DO 60 J=1,ITNO S TK^FLOAT! IT IDE( 5 , J >+lT )HNG .98
ZZ = 0. 5 DO 55 K = 1,ITCIM S L=K + ITCN HNGM.18
ZZ = COS(TIDSPD(X)»TK-TIDE(K,J))»TIDEIL,J) + ZZ
100 55 CONTINUE S IFCZZ.EQ.O.) ZZ=ZEROP S ZZ3=ZZ»TIOE( 10, J )
ZZ2 = ZZ»TIDE< 9, J) 5 Pin=ITIDE< 3. J ) S MP= I TIDE( M, J )
DO 60 n-nn,np * ip= $ iP=ip+niDE(2,j)
DO 60 I = I(H,1P S DO 60 JL = l,flAXLAY $ L=LEV( JL )
K = HT^ I ) $ IF(K.GE.JL.AND.ZZ( JL ).NE.O. ) Z( I+L )=ZZ( JL )
HNGH.19
HNG81 .H
HNG. 101
HNG. 105
HNG81 .5
HNG81 .6
105
110
60 CONTINUE HNG.109
70 IF( H.LT.fllNSORC) GO TO VO $ DO 79 J = 1,15NO HNG.110
lF(SORC(2.J).EQ.O..OR.IT.LT.lSOPC(T.J).OP.lT.GT.150BC(BfJ))GOT079 HN6.111__
nn=isoRC(3,J) s np=isoRC(t,j> HNG.112
DO 78 n=nn,np * IP=»NE * in=ip+i50RC( i.J) t IP=IP+1SOHC(2,J) HNG.113
-.n=150BC(j'j) S LP=ISORC(6,J) » DO 77 l = in JP » DO 77 L=LP.LP ' HN68.33,
IL = L£V(L) + I » 1FCHTK I ) .GT.FLOATC L-l ) )Z( JL ) = Z< JL )+SOBC( 2, J )»01 HNG8.3t
77 CONTINUE
_-._
79 CONTINUE
HNG.118
_HNG.119
HNG.120
-------
PROGRAM LAYER3 7600-7600 OPT=2
FTN 1.6+139/010
09 MAR 77 01.15.3U BKV PAGE
115 90
91
100
120
120
125
130
135
130 110
IF( IT.NE.ITU) GO TO 110 S IF( IUNO.EQ .0 ) GO TO 91
IFt IT.EO.lUINDt 5 .IUNO)) GO TO 120
CALL XPIITt 0. ,XK,nAXDin»2 ,0,1) $ ITU = 0
1UNO=IUNO+1 4 1F< IWNO.GT.mAXUIND) GO TO 110
IFl IT.GE.IUlNDf 6.IUNO)) GO TO 100
ITW = 1WIND(5,IUNO) 4 1F( ITU .ST. IT) 60 TO 110
YU=t -ROTANG-90.-FLOAT( IUIND( T, I UNO ) ) )»RPD
TK=IWIND< 8. IUNO) S TK=TK»TK»C3 4 XU=COSt YW )»TK 4 YU=S1N(Y
nn=iuiND( 3,iuiNO) 4 PIP=IUIND(I,IUNO)
DO 130 n=nn,pip * ip=tn-i >»NE 4 IM=IP + IUIND( I.IUNO)
IP = 1P + IUIND( 2.IUNO) * DO 130 UIPI.IP * XK( I ) = XU
IF( lUNO.GE.nAXUIND) GO TO 135
IF( JT.GE.IUINDC5 .IUNO + 1 )) GO TO 100
ITUI = lUIINDt6,lUNO>
CALL SETHUV
DO 111 JL = 1J1AXLAV 4 K=LEVUL) + NO
HNG.121
HNG.122
HNG.123
HNG.121
HNGt .17
HNG6.18
HNG81 .7
'U )»TK HUG. 127
HNG.128
HNG.129
HNG .1 30
HNG.131
HNG6.20
HMG6.21
HNG. 133
HNG.131
BC8.3
135
00
110
115
150
155
160
165
111 CONTINUE_
CALL CXiniTCHGUC K-NE),HGU(K ),HGU+1 HNG8.36
HNG8.37
CALL BARUV(V,HG«,0) 4 CALL UVCAL( V ,-1.) BAR(n00.2
DO 302 JL = 1,H1AXLAY 4 K = LEV( JL UNO BC8.10
CALL CXMIT(V (K-l ).V (K ).H6V(K ) .0..HE.-NE.-NE.-NE ) BC8.1
202 CONTINUE
CALL UV5TAR(U,HTU,1 ) $ 00 250 L=1,MAXLAY
250 CALL XnlTlZm JLKUl JLKNOjKl >
302 CONTINUE
IF(IT.LT.niNFLOU) GO TO 500 4 DO 190 J=1,IFNO
IF( IT.LT.IFLOUI7.1 ) .OR.IT.6T.IFLOU( 8.J » 60 TO 190
BfTT.
H.!>: 3
K;.G8
Pin=IFLOU< 3,J ) S I-1P::IFLOUK1,J ) 4 L(1-IFLOU( 5 , J ) S LP=IFLOUI( 6 , J ) HNG3
TK = 1. * 00 180 L=LP1,LP * IF(L.NE.l) TK=FLOU( L + l, J ) HNG8
_DG 180 n=nn.MP,t ip-(m-i>»NE s in-ip + iFLOuc i.j ) $ IP-IP-HFLOUK 2. J )HNGS
DO 180 I = lm,"lP i JL=LEV(L)+1 HNG8
IFt HTU( I l.GT.HLO(L)) Uf JL )t=Uf JL HFLOUt 1, J )»TK
!F<_HTVt ] KGT.HLOt D) V( JL )-V( JL HFLOlJt 2. J )»TK
180 CONTINUE
190 CONTINUE
_500 1 fJin0pilLJI000UTJKNE.O.OR. IT.LT.1TO) GO TO 600
DO 5i"6 L = 1,MAXLAY 4 JL=LEVtL ) + l
510 CALL INPOUTt Z( JL ), L , 2, IT, CFAC )
HNG8
_BTJS.8
HNG8
Hi.GS
_KW3.
HNG8
PR3.
HUG 8
700 1F( i10D( IT, IPTIOD J.NE.O.OR.IT.LT.IPO) GO TO 800 " HNGT
DO 705 L=l,tnAXLAY « JL=LEV< L HI HNS8
CALL PRT nAX( Z( JL ),1HNONE,1H ,31HUATER HEIGHT DEVIATION ( CH ) LAVERK K G 8
'
12
.39
._1_0_
.11
.12
^13_
.11
.15
,._1_6_
.17
.18
ill_
.19
5
•51_
2312
.52
.53
705 CALL PRTMAX( V( JL>,U( JL),1H ,18HCURRENT SPEED AND DIRECTION ( CPI/SECHNG8 .55
_ L !!LUt )_ LAVE R .L.18.1T) _ H NGf.. 5_6_
600 CONTINUE HNG.216
900 CONTINUE HNG.217
_ END _ H PIG. 218
-------
SUBROUTINE SETHUV T600-«7600 OPT=2
FTN 1.6+139/010
09 BAR 7T 01.15.3U BUY PAGE
oo
to
1 SUBROUTINE SETHUV
COMMON /CONTROL /HTZ( 300),HTU< 300KHTVI 300>,NNP( 20),NNM< 20),NET
1ME, NO, MAXD1M,MAXLAY, MODEL, ROT
A,ZZ( 3>,HLO( 1>,HL(3),ALPHA< 3 ) ,BETA< 3 ) ,LEV( 3 ) , A5RHOI 3 >, A5RHC 3 >, A3,
5 2,IM.IP.NP,NK.KL.JL,ZZ2.IPP.IUCt92)
C DIMENSION ZLAYK(MAXDIM»T»(K-1» (WHERE K - ACTUAL NO. r^YERS).
COMMON // Z< 300) U( 300),V( 300),ZM< 300),ZST( 300)
A .HGU< 300>.HGV< 3&0)
K ,ZLAYK( 1 )
10 U ,XK< 300), YK( 300)
l.COR
IPP=NO * NK=NE S KL=0 4 DO 120 KK=1,2
IF(KK.EQ.l) GO TO 10 S NK=1 * KL=MAXDiri
10 DO 110 KUl.ME $ Im=NNM(in) 4 IP=NNP(P1) * IF(IP.LT.in) GO TO 110
15 IF(XK.EQ.Z) IPP=M»NE
DO 100 I=IM,1P $ ZZ=HTU(I+KL) * IF( ZZ.LE.O. > GO TO 100
NP-I+NK * IF(NP.GT.IPP)NP = I $ JL=P1AXLAY
60 IFC ZZ.GE.HLO( JD) GO TO 80 $ JL=JL-1 $ GO TO 60
80 I2=LEV(JL) * 11=12+1 * I2=I2+NP S I3=I1+KL $ JL=JL-1
20 ZZ2=HL( JL)-< Z( 11 ) + Z( I2))< .5
HGU( 13)=ZZ-ZZ2 S IF(JL.EO.O) 60 TO 100 * ZZ=ZZ2 $ GO TO 80
100 CONTINUE
110 CONTINUE
120 CONTINUE $ RETURN
2$ END
HNG.219
CONTRd.l
HNG6.1
A6HNG8.1
HNG8.57
BLANKA.l
BLANKA.2
BLANKA.3
BLANKA.l
BLANKA.5
HNG8.2
HNG8.58
HNG8.59
HNG8.60
HNG8.61
HNG8.62
HNG8.63
HNG8.61
HNG8.65
HNG8.66
HNG8.67
HNG8.68
HNG8.69
HNG8.70
HNG8.71
CARD NR. SEVERITY DETAILS
DIAGNOSIS OF PROBLEM
I
_zz_
LIBSCRIPTED. FIRST ELEMENT UIILL BE USED.
16 I ZZ ARRAY NAdE OPERAND NOT SUBSCRIPTED, FIRST ELEMENT WILL BE USED.
18 I ZZ ARRAY NAME OPERAND NOT SUBSCRIPTED, FIRST ELEMENT WILL BE USED.
21 I ZZ ARRAY NADE OPERAND NOT SUBSCRIPTED^ FIRST ELEMENT WILL BE USED.
21
ZZ
ARRAY NAME OPERAND NOT SUBSCRIPTED, FIRST ELEMENT WILL BE USED.
-------
SUBROUTINE UVCAL
7600-T600 OPT=2
FTN t .6+t39/010
09 ClAR TT 01.15.3U BKY PAGE
00
(jj
SUBROUTINE UVCAL(UV,A1>
DIMENSION UV(1)
COMMON /CONTROL/HTZ( 300),HTU(300),HTV(300 ),NNP( 20 ),NNW 20 >,NE,
HMG8.T2
HNG8.73
CONTRtt.l
HNG6.1
_AJ.I2J_3J_tHL 0( 1 )rHL(3),ALPHA(3),BETA( 3),LEV( 3) A5RHO( 3).A5RH( 3 ) , A3 . A6HTJG8 . 1 _
2,IPP,KK,Nk,J~L,lK,JP,AA,Il,I2,I3,ZZ2,lOuKt8fi) HNG8.75
DIMENSION ZLAYK (WHERE K - ACTUAL NO. LAYERS). BLANKA.l
COnFiOM // Z( 3 00 1. U( 300 ).V( 300 ).ZCU 300).ZST( 300 ) _ BL A N XA_..2_
A ,HGU( 300),HGV( 300)
10 K ,ZLAYK( 1)
U .XK( 300) YK( 300)
l.COR
IF< ftl .LT.O. )GO TO
15 20 DO 101 n=l,ME $ I
IF(KK.NE.O) 1PP=N
DO 100 I=IH.IP S
10 * IFP=NO * KK=0 <
IK=I+KK $ IF( HTU( IK )
> NK=NE
S 1FCIP
.LE.O. )
$ SO TO 20
.LT.ld) GO TO 101
GO TO 100
BLANKA.3
BLANKA.t
BLflJKA.5
HNGB.2
HHGB.77
K.'.G8.78
KKiO.79
HNG3.80
HNG8.81
DO 50 L=1.MAXLAY * JL=riAXLAY-L + l * IF( HGUt LEV( JL ) + IK ) .GT . 0 . )GOTOfcOHNG8 . 82
50 CONTINUE HNG8.83
HNG8.81
DO 80 L=1,JL S I2=LEV(L) « I3=I2+IK 4 11=12+1
ZZ2=A1*(Z< 12 + JP)-Z( 11 )) * IF(L.EQ.JL) AA = A3
65 IF( HGU( I3KGT.O. ) GO TO 70 $ UV< 1 1 1 = 0 . * GO TO 80
70 UVl 11 ) = Zn< 11 )-SQRTC U( 13 )*U( I 3 ) + ZST( 1 1 )»ZST( 11 ))»AA»U( I3)/HGU( 13)
25 1 -A5RH(L)»ZZ( 1 >-A5RHO< L >+ZZ2+ZST< 11 )»A1»COR
IF(K.EQ.O) GO TO 80 $ UV( 11 ) = UV( 1 1 ) + XK( IK )/HGU( 1 3 ) $ K=0
80 ZZ=ZZ2
100 CONTINUE
101 CONTINUE 4 RETURN
HNG8.85
HNG9.1
HNG8.87
HNG9.2
HNG9.3
HNG8.90
HNG8.91
HNG8.92
HNG8.93
30
END
HNG8.9t
CARD NR. SEVERITY DETAILS
DIAGNOSIS OF PROBLEM
20
27
ZZ
_zz_
ARRAY NAME OPERAND NOT SUBSCRIPTED, FIRST ELEMENT WILL BE USED.
ARRAY NAME OPERAND NOT SUBSCRIPTED, FIRST ELEMENT U1LL BE USED.
-------
SUBROUTINE UVSTAR 7600-TfcOO OPT = 2 FTN 1.6+139/010 09 WAR IT 01.15.31* BUY PAfiE
00
1 SUBROUTINE UVSTAR( UV ,H ,N ) $ DIMENSION UV(1),H(1)
COMnON/CONTROL/HTZC 300 ),HTU( 300).HTV( 300).NNP( 20),NNCK 20)rNEr
lME,iMO,MAXDIM,nAXLAY,P10DEL,ROT
A,ZZ( 3),HLO( 1 ),HL( 3),ALPHA( 3),BETA( 3),LEV( 3 ) , A5RHOC 3 >, A5RH( 3 > , A3,
C DIMENSION ZLAYK(NAXDirn»7»( K-l ) ) (WHERE K - ACTUAL NO. LAYERS).
COWON // Z( 300), U( 300), V( 300), Zn< 300),ZST( 300)
A HGU< 300 ).HGV( 300)
K ,ZLAYK( 1 )
10 Ul ,XK( 300), YK( 300)
l^COH
5 Doioon=i,riE$ui=NiMn(n)$ip=NNp(n)$iF( IP .LT.iniGOTOioo
D090I=in.IP$ 11=1 $ IF(N.NE.l) GO TO 10
15 I3=I>1INO( I + 1,IPP ) * I2=I-NE S IF(I2.GT.O) GO TO 7
12=1 $ 11=13 * GO TO 20
7 H=I $ GO TO 20
10 I2 = P1AXO( I-l,Idm> $ I3=I+NE $ IFU3.LE.NO) GO TO 15
13=1 * 11=12 * GO TO 20
20 15 I1=P1AXO( 13-1. IPP + 1 )
20 DO 80 L=1,MAXLAY * JK=LEV(L)
80 ZST( JK + I ) = (UV( 11 + JK )+UV( I2+JK >+UV( I3+JK)+UV( I1+JK))».25
90 CONTINUE
100 CONTINUE *RETURN
25 END
HNG8.95
CONTRP1.1
HNG6.1
A6HNG8.1
HNG8.96
BLANKA.l
BLANKA.2
BLftNKA.3
BLANKA.1
BLANKA.5
HNG8.2
BARSTAH.8
HNG6.32
BARSTAR.10
HNG6.33
HNG6.31
HNG6.35
HNG6 .36
HNG8.98
HNGfc.38
HNG8.99
HNG8.101
BARSTAR.21
BARSTAR.22
BARSTAR.23
-------
SUBROUTINE BAR
7600-»7600 OPT=2
FTN 1.&+t39/0')0
09 PIAR 77 01.15.3U BKY PAGE
co
1
5
10
15
20
25
CARD NR.
17
17
18
19
20
21
SUBROUTINE BAR
COMMON /CONTROL/HTZC 300),HTU( 300)rHTV( 300),NNP< 20>rNNFI< 20),NE,
HNG8.105
CONTRM.l
1ME, NO, MAXDIM,Ii(t90) HNGS.IOA
C DIMENSION ZLAYK(MAXDIM*7»(K-1 ) > ( UIHERE K - ACTUAL NO. LAYERS).
COMMON // Z( 300), U(300),V< 300), ZW 300). ZST< 300)
A .HGU( 300KHGV< 300)
K ,ZLAYK(1)
U ,XK( 300), YK( 300)
l.COR
DO 30 Pl=l,niE 4 IM=NNr)(n) $ 1P=NNP(M) * IF(IP.LT.in) 60 TO 30
lPP=m*NE 4 InM=IPP-NE+l
DO 25 I = ln.IP 4 NP1=P1AXO( 1-1. imfl) $ NP=niNO( I + l.IPP)
IW=1-NESIF(P1P1.LE.O )inri=I $ nP=I + NE * IF(MP.GT.NO)MP=I
DO 20 L=l,mAXLAY * JK=LEV(L> $ JL=JK+I
ZZ = L $ IF(HTZ( I ).GE.ZZ) GO TO 1M *Zin(JL)=0. * GO TO 20
11 IF(HTZCNPl).LT.ZZ) NPI=I
IF(HTZCNP).LT.ZZ) NP = I
IF(HTZ
-------
SUBROUTINE BARUV
T600-*7600 OPT=2
FTN 1.4+139/010
09 MAR 77 01.15.3U BUY PAGE
oo
I SUBROUTINE BARUV(B,H,N> » DIMENSION 6(1), H(l)
COPinON /CONTROL /HTZ( 300),HTU< 300),HTV( 300)rNNP< 20 >,NNn( 20 > rNE .
lPlE,NO,MAXDIPl,nAXLAV(P10DEL,ROT
A,ZZ< 3),HLO( 1>,HL(3).ALPHA(3),BETA<3),LEV(3),A5RHO( 3 >, A5RHC 3 ), A3,
5 2,im.ip,mn,rrip.Nn,NP,iinM,ipp.JLrJK.lu(')SO)
C DIMENSION ZLAYK(nAXDIPt»7»(K-n) (UHERE K - ACTUAL NO. LAYERS).
COnMON // Z< 300J.UC 300),V< 300>,Zm 300>,ZST< 300)
A .HGU( 300).HGV( 300)
K ,ZLAYK(1>
10 U ,XK( 300), VK( 300)
l.COH
DO 30 M=1,P1E * H1=NNM(n) * 1P=NNP(C1) * IF(IP.LT.IPl) GO TO 30
DO 25 I = 1M,IP * Nn=P!AXO( I-1,IM) S NP=MINO< 1 + 1, IP )
Pin=I-NE$IF(PlM.LE.O)ni'l=I * flP = I+NE * IF(nP.GT.NO)MP = I
15 DO 20 L = 1,P1AXLAY $ JK=LEV(L) 4 JL=JK + I
IF(H( JD.GT.O. ) GO TO 15 * ZI-KJL>=0. * GO TO 20
15 IF(N.NE.O) GO TO 16
IF(H(M + JK).LE.O. ) m=I $ IF(H(MP+JK).LE.O. ) P1P=I » GO TO 17
16 IF(H(NM+JK KLE.O. ) NM=I * IF( H( NP + JK ) .LE .0 . ) NP=I
20 17 CONTINUE
ZN< JL ) = B( JL )»ALPHA( 1 )+( B< NR+JK )+B( NP+JK )+B( NPI+JK )+B( MP+JK ) )»
1BETA( 1 >
20 CONTINUE
25 CONTINUE
25 30 CONTINUE SRETURN
END
BARftOD.8
CONTRPI.1
HNG6.1
A6HNG6.1
HNG8.116
BLANKA.l
BLANKA.2
BLANKA.3
BLANKA.l
BLANKA.5
HNG8.2
HNG8.118
HNG8.119
BARSTAR.11
HNG8.120
HNG8.121
BARP10D .9
BARriOD.lO
BARMOD.ll
6flRP)OD.12
HNG8.121
HNG6.125
HNG8.126
BARSTAR.18
BAHSTAR.19
BARSTAR.50
-------
SUBROUTINE INPOUT 7600-7600 OPT=2
FTN 1.6+139/010
09 MAR 77 01.15.3U BKY PftGE
oo
1 SUBROUTINE INPOUTf A .LAYER, 1TY , IT, FACOUT )
COMMON /CONTROL/HTZ( 300).HTU( 300).HTV( 300).NNP( 20).NNN( 20).NE
lME,NO,MftXDIM,MAXL AY, MODEL, ROT
A,ZZ( 3),HLO( U,HL(3),ALPHA( 3),BETA( 3 > ,LEV( 3 ) , A5RHOC 3),A5RH( 3)
5 2,IW(500)
DIMENSION AC 1 )
DATAt 1REC=0),(M21=7777 7777B ), ( M12=7777B ), < LEN=500 >
l.U = 0>
C -LAYER=READ,+LAYER=URITE
10 PRINT 1000, IT, LAYER, 1TY,NE, ME, MODEL, IREC
1000 FORFiATf 10H INPOUT 101 10)
ID1 = SHIFTC IT, 36 ).OP. .SHIFT( IABSC LAYER), 12)
ID2=SHIFT( ITY,18).OR.SHIFT(NE,36 ). OR. SHI FTC ME, 21 ). OR. MODEL
KP=f1ftXO( l^JHINOf 3. 10-ITY ) )
15 DO 600 K = 1,KP * KK=( K-l )»MAXDIM
IF( LftYER.GE.O ) GO TO 100 $ FACIN=10 . /FACOUT» . 1
NW=501 S DO 300 1 = 1. NO * 1FC NU) .LT . LEN ) GO TO 150
IFCNW.EO.LEN.AND.J.NE.l ) 60 TO 150
C INPUT
20 100 BUFFER INC 1 . 1 )( JU< 1 ) , IU( LEN ) ) * IREC=IREC+1 $ NU=2
IF( UNITC 1 ) 1110,110,120
110 lERR^l $ GO TO 130
120 IERR=2
130 PRINT 1130, IUC 1 > , IU< 2 > , ID 1 , ID2 $ CALL IOERRC IERR, IREC ) * 60
25 1130 FORFIATC 6< 1X020) )
110 CONTINUE
IFC IUC 1 ).NE.ID1.0R.1U( 2KNE.ID2) GO TO 100
150 J=MOD( 1,5 ) + l 4 GO TO < 250, 210, 220, 230, 210 ) J
210 Nu=rao+i $ iuu=iu(NU) $ N=SHIFT< iuw.-36 > $ GO TO 300
30 220 N=SH1FT( 1WU,-12> $ GO TO 260
230 N=IUU.AND.Mi2 $ NU=NU + 1 * IUW=IUICNU>
N=SHIFT( N^J.OR.SHIFTC IUUI . 12 ) . AND .M12 * GO T0_ 270
210 N=SHIFT( 1WU,-21 ) $ GO TO 260
25 0 N = IUIUI
35 260 N=N.flND.M21
270 IFCN.GT.1000 OOOOB ) N = N.OR ..NOT. M21
300 ACKK + I )=FLOAT(N)*FACIN $ GO TO 600
C OUTPUT
100 NU=2 $ IU(1) = ID1 * IU(2) = ID2
10 DO 500 1=1, NO S N=A(KK+I )»FACOUT $ N=N.AND.M21
J=MOD< 1,5 > + l * GO TO (150.110.120.130.110)^
110 NW=NU+1 S IU(NU)=SHIFT(N,36 ) $ GO TO 170
120 IUKNU)=IWl>.OR.SHIFT( N.-12) $ NW=NUI+1 * N=N.AND.PU2
15 IU(NW) = SHIFTCN,18) * GO TO 170
110 1W(NU)=IU(NU).OR.SHIFT(N,21 ) $ 60 TO 170
150 IW(NU)=IU(NU).OR.N
170 IF(I.EQ.NO) GO TO 175 * IF( NU .LT.LEN ) GO TO 500
IFCNU.EQ.LEN.AND.J.NE.l ) 60 TO 500
50 175 BUFFER OUTU . 1 >( IU< 1 >, IU< NU>). $ IREC = IREC+1 $ NU=2
IFCUNITf 1 »500,180,190
180 IERR=3 * 60 TO 195
190 IERR=1
195 CALL IOERR( IERB, IREC ) » 60 TO 175
55 $00 CONTINUE
600 CONTINUE
RETURN
END
PR3.6
, CONTRM.l
HNG6.1
,A3,A6HNG8.1
UTIL.1
UTIL.5
UTIL.6
UTIL.7
UTIL.8
UTIL.9
UTIL.10
UTIL.ll
UTIL.12
UTILB.l
UTILA.2
PR3.7
UTIL.15
HNG6.39
UTIL.17
UTIL.18
UTIL.19
UTIL.20
UTIL.21
TO 100HNG6.10
HNG6.11
UT1LA.3
UTILA.1
UTIL.21
UTIL.25
UTIL.26
UTIL.27
UTIL.28
UTIL.29
UTIL.30
UTIL.31
UTIL.32
PR3.8
UTIL.31
UTIL.35
PR3.9
UTIL.37
UTIL.38
UTIL.39
UTIL.10
UTIL.ll
UTIL.12
UT1L.13
UTIL.ll
HNG6.12
UTIL.16
UTIL.17
UTIL.18
UTIL.19
UTIL.50
UTIL.51
UTIL.52
UTIL.53
-------
SUBROUTINE IOERR
7600-»T600
FTN M.6+139/010
09 MAR 77 01.15.314 BKY PAGE
00
00
1 SUBROUTINE 10ERR( L , IREC ) 4 D1P1ENSION LAB( 1 )
DATftf LAB=10HEOF ERR IN.10HPAR ERR IN.10HEOT ERROUT. 10HPAH ERROUT)
PRINT 1000, LAB(L>, IREC 4 IFCL-2) 1,3,2
5 1000 FORMAT! All , I 10 )
2 BACKSPACE 1 4 ENDFILE 1 4 BACKSPACE 1 4 IFCL.E0.3) GO TO 1
ENDFILE 1 4 BACKSPACE 1
3 IREC=IREC-1 4 K=K + 1 $ IF(K.LT.ZO) RETURN
1 STOP
10 END
UTIL.55
UTIL.56
UTILA.5
UTILA.6
UTILA.7
UTILA.8
UTILA.9
UTILA.10
UTIL.62
UTIL.63
SUBROUTINE XMIT 7600-7600 OPT = /! FTN 1.6+139/010
1 SUBROUTINE XPIIT( I A IB ,N, JA, JB ) 4DIC1ENSION 1 A( 1 ) , !B< 1 )4 I = J = 1
DO 5 NN=1.N 4 IB(J)=IA(1) 4 J=J+JB 4 1=I+JA
5 CONTINUE 4 RETURN
END
09 MAR 77 01.15.3U BKY PAGE 1
HNG5.28
HNG5 .29
HNG5.30
UTIL.66
SUBROUTINE CXPUT 7600-+T600 OPT=2 FTN 1.6+139/010
1 SUBROUTINE CXMIT( A,B ,C ,D,N, J A, JB , JC ) 4 DIMENSION A< 1 ),B( 1 ) ,C( 1 )
I = J = K = 1 4 DO 5 NN=1.N 4 IF( C( K ) .GT.D ) B( J ) = A( I)
K=K+JC 4 J=J+JB 4 I=I+JA
5 CONTINUE 4 RETURN 4 END
09 MAR 77 01.15.314 BKY PAGE 1
BC8.13
BC8.11
BC8.15
BC8.16
-------
SUBROUTINE PRTHAX 7600-7600 OPT=2
FTN 1.6*139/010
09 MAR 77 01.1S.3U BKV PAGE
00
1 SUBROUTINE PRTI"IAX< S.D, I S.NAPIEl ,LE,NC . IT )
DIMENSION S( l),D( 1 )
COnnON /CONTROL /HTZ( 300),HTU( 300),HTV( 300>,NNP(20 J.NNM 20),NE,
5 inE,NOrfnAXOtM,MAXLAY,nODEL.ROT
A,1Z( 3),KLO( n,HL( 3 ) , ALPHAC 3 ) , BETAf 3 > ,LEV< 3 > . A5RHO< 3 ), A5RHC 3 ), A3,
l,IK,IUU,JK,I.A,B,NOP,nF.ML,MnF,NN.NU,IU<188)
DATA(LU = 3> .< P1AXNAr>E = 8)
NU=(NC-1 )/NCPUI»l
10 NOP=(flE-l >/15*l * nF=-11
CALL XniT(NAI>F-l )»NE i nL=(ML-l)»NE S 1F(D( 1 KNE.tHNONE) BO TO 40
DO 50 J=1,NE * JK = J+PinF $ NN = flL + J
20 50 URITEC LU.2000 HS. Jr< S( K ),K = JKP' N.NE)
2000 FORMAT^ Al, IH.3X15F8.0 )
GO TO 100
60 DO 90 J = 1.NE S NN=10 $ DO 80 K = PinF , CIL f NE * JK = J + K * NN = NN+1
BMsJK-NE 4 lF60 TO 70 J IU(NN)=1H S GO TO 80
70 A = SQRT( X»X + V»V ) » B = AP10D( ROT-ATAN21 V X )»57 . 2957795 13, 360 .)
ENCODE* 9, 3000, 1UU )A,B $ IU< NN )=IUU .OR . 770000000000B
3000 FORriAT< F5 .O.F1 .0 >
30 80 CONTINUE
90 URITE(LU,tOOO)IS, J,< IU( K ),K = 1 1 ,NN )
MOOO FORPIflTI Alj It, 3X15A8 )
100 CONTINUE $ RETURN
END
CARD NR. SEVERITV DETAILS DIAGNOSIS OF PROBLEM
12 I CONTROL VARIABLE IN COnmON OR EOi! I VALENCEO^ OPTIH1ZATION HAY
16 I CONTROL VARIABLE IN COnnON OR EQU 1 VALENCED, OPTIMIZATION PIAY
HNG8.127
UTIL.68
UTIL.69
CONTRH.l
HNG6.1
A6HNG8.1
UTIL.71
UTIL.72
UTIL.73
UTIL.76
1)HNG6.H3
UTIL.77
UTIL.78
HNG8. 128
HNG5.35
UTIL .81
UTIL.82
UTIL. 83
UTIL. 81
UTIL. 85
UTIL. 86
UTIL. 87
HNG6.1H
HNG81 .8
HNG81.9
HNG81 . 10
HNG81 .11
UTIL. 92
UTIL. 93
UTIL. 91
HNG6.16
UTIL. 96
UTIL. 97
UTIL. 98
BE INHIBITED.
BE INHIBITED.
-------
LOAD PIAP.
LINK - BKY iOOO/TOOO 8.1
09 C1AR 77 01.1$.3$
PAGE
FL REQUIRED TO LOAD 36500
FL REOU1RED TO RUN 31100
INITIAL TRANSFER TO LAYERS
6150
BLOCK ASSIGNMENTS.
BLOCK
/CONTROL/
LAYERS
SETHUV
UVCftL
UVSTAH
BAR
BftflUV
INPOUT
IOERR
XMIT
CXP1IT
PRTP1AX
/STP.END/
/FCL.C./
/OLINLMT/
/QJOBSFL/
/OASAFLG/
/.QCOfflP./
/Q8.IO./
FTN4LIB
ATAN2
BflCKSP=
BUFIN=
BUFIO=
BUFOUT=
CIO =
COMIC=
DECODE=
ENCODE=
ENDF1L=
ERM$*
FECMSK=
FLT!N=
FLTOUT=
FdTAP=
FORSYS=
/QOVERL7/
FORUTL=
GETFIT=
GOTOER= .
INCOM=
INPC =
KOOER=
ADDRESS
100
2777
11065
11167
11350
11M53
11561
11673
12271
12354
12374
12123
13022
13023
13044
13045
13047
13050
13051
13315
13315
13415
13477
13560
13724
14001
14014
14075
14212
14356
14446
16134
16175
16353
16670
17262
20243
20244
20262
20325
20341
20623
21071
LENGTH
2677
6066
102
161
103
106
112
401
60
20
27
377
1
21
1
2
1
1
244
0
100
62
61
144
55
13
61
115
144
70
1466
41
156
315
372
761
1
16
43
14
262
246
467
FILE
LGO
LGO
LGO
LGO
LGO
LGO
LGO
LGO
LGO
LGO
LGO
FTN4LIB
FTN4L1B
FTN4LIB
FTN4LIB
FTN4LIB
FTN4LIB
FTN4LIB
FTN4LIB
FTN4LIB
FTN4L1B
FTN4LIB
FTN4LIB
FTN4LIB
FTN4LIB
FTN4LIB
FTN4LIB
FTN4LIB
FTN4L1B
FTN4LIB
FTN4LIB
FTN4L1B
FTN4LIB
FTN4LIB
KRAKER=
21560
435
FTN4L1B
-------
LOAD MAP.
LINK - BKY 6000/7000 8.1
09 ClftH 7T 01.15.35
PAGE
BLOCK ADDRESS
OUTCOM=
OUTC =
Rt)L =
RCO=
RDU=
SINCOS=
SORT
SVSAID =
SVS=
SVS=A1D
SVS=1ST
UNIT
UZERO. .
UITH =
UTL=
UTU =
//
NAME TEST BOX MODEL
GRID
LAYS
TIKE
FACT
TIDE
TIDE
TIDE
SORC
COP1P
INPOUT
INPOUT
TNPOUT
5 533333.3333
1 77777
600 600
15.011068613
1111 33165
5511 32565
2555 31665
222211
99
0
0
600
22215
22121
22656
22715
22710
23017
23113
23211
23212
23216
23255
23337
23121
23126
23156
23505
23573
0.
0
LENGTH FILE
201 FTN1LIB
235 FTN1LIB
37 FTN1LIB
23 FTN1L1B
107 FTN1LIB
71 FTN1L1B
16 FTN1LIB
1 FTN1LIB
31 FTN1L1B
7 FTN1L1B
62 FTN1LIB
65 FTN1LIB
2 FTN1LIB
30 FTN1LIB
27 FTN1L1B
66 FTN1L1B
5216
3.2E-6 12. .003
.99 1.022
2100 100 0 1200
.913035628.9811012 30.
117. 121. 359. 052. 1.0 2.0 3.0 1.0
117. 121. 359. 052. 1.0 2.0 3.0 1.0
117.
100
-1
-1
1
121. 359. 052. 1.0 2.0 3.0 1.0
1000000 100. 2.
1 5 5 99
2 5 5 99
2 5 5 99
2.
0
3
6
INPOUT
1200
9?
-------
WATER HEIGHT DEVIATION (CP1) LAYER
1200
1
2
3
1
5
1
2
3
1
5
INPOUT
INPOUT
1
1
1
1
1
0.
1
7i\82
8,112
•9,125
9:116
8;111
2 3
1. 1 .
1. 0.
1. 0.
1. 0.
1. 0.
CURRENT SPEED
2 3
10 1 181 13-. 181
10;150 11,158
9,130 8,132
8,111 7:118
8,109 6;11H
1800 1
2100 1
1
1 .
-0.
-1 .
-0.
0.
5
0.
-0.
-0.
-0.
-0.
AND DIRECTION ( CM/SEC, TRUE ) LAYER
1
9; 181
8,166
1,117
2;117
2, 106
2
2
5
2, 96
1,215
0;226
0,121
5 5 99 12
5 5 99 15
1 TII-IE= 1200
1
Z
3
1
5
1
2.
2.
2.
2.
1.
WATER HEIGHT
2 3
2. 2.
2. 1.
1. 1.
1. 1.
1. 1.
DEVIATION (CCl) LAYER
1
2.
1.
0.
1.
1.
5
0.
1.
1.
1.
1.
1 TIME= 2tOO
1
2
3
1
1
iO:i8i
10,160
10;138
10;125
CURRENT SPEED
2 3
11:183 11:182
11,161 10,167
10,112 10,115
10_jl26 9j_130
AND DIRECTION ( CM/SEC, TRUE ) LAYER
1
8: 181
8,157
8,112
6:128
5
5j 92
1; 97
3:101
1 TinE= 2100
10,120 10,120 8,125 5,120 3;101
-------
APPENDIX C
NON-LINEAR MODIFICATIONS
CONVCT AND CONVCU SERIES
The CONVCT and CONVCU modifications add non-linear terms to Phase II u
and v velocity computations. The procedure expands the logic in the BARUV
routine to calculate the non-linear terms when the routine is entered through
the special entry CONUV,- which sets the flag ICON. CONUV is called once for
each component, and the results are added to the previously computed velocity
component in array ZMi. To save the required work arrays in the computational
block in the main program, the UVCAL call is changed so the Vi values are
computed into the ZMi array rather than into the Vi array, and an XMIT is add-
ed to transfer ZMi into Vi at the end of the calculation. Control common is
expanded with the insertion of the constant A7 and variable CY before the 500
scratch positions. The computation of A7 has been added to the initialization
block. CY is set in CONUV to the maximum absolute value of the non-linear
component and printed in the output section as a debugging aid that has been
retained as a useful feature in the model.
Boundary conditions have been presented in (1). Boundary condition 1 is
applied by the default value at the land-sea boundary normal to flow and in-
teresected by tangential flow line. Boundary condition 3 is applied at open
normal boundaries by index substitution. Boundary condition 7 is applied
prior to the computation by scaling the terms for which the BARUV default
boundary condition 6 has applied index substitution at tangential open or
closed boundaries or normal boundaries not intersected by the tangential flow
93
-------
line. An indexing technique is used to accomplish the computation for both
the u and v non-linear term cases. The computation is symmetric or quasi-
symmetric in space depending on the influence of boundary conditions.
The names of these correction sets are in alphabetical sequence. Fur-
ther correction sets related to these computations would be named CONVCV,
etc. This naming technique will occur automatically if all related correc-
tion sets were preceded by *IDENT CONVCT due to renaming conventions of
UPDATE(6) in a duplicate name situation.
CONPRT, CONPRU AND CONPRV SERIES
The CONPRT, CONPRU and CONPRV modifications provide an output capability
for the non-linear term computations. There are two separate procedures:
printed output with other scheduled printed output is provided when CONPRT
is used alone, and output to the intermediate data set as type 10 data records
without printed output when all three modifications are used. Use of the
CONPRT feature alone modifies the format for print of all single arrays so
that normal printed display of layer surface deviations is in exponential for-
mat. An additional array of CZ (300) is inserted immediately following the
CY variable in control common. CZi must be modified to be at least MAXDIM*
(number of layers) when MAXDIM is changed and these modifications are used.
The initialization block is affected by the change in control common and by a
call to XMIT to initialize array CZi to zero. The u and v computational block
is modified by a call to INPOUT to save the array CZi on the intermediate data
set (or call to PRTMAX for CONPRT alone to print the array CZi) . The output
block is similarly modified by a call to INPOUT (or PRTMAX) . Routine BARUV
is expanded to define the CZi arrays as the computation of the non-linear
terms progresses. Naming conventions are similar to those employed for CONVCT
-Series corrections. These modifications may never be used without the CONVCT
94
-------
series.
LISTING FOR NON-LINEAR CORRECTION SETS
The following listing of the CONVCT and CONPRT series correction sets
is in card image form as originally entered into the library. The spacing
between correction sets is artificial and the sets are presented in a logi-
cal order rather than the order in which they were entered on the library.
These correction sets are applied to the basic Phase II program shown in
Figure 8. Since the correction sets are actually a part of the library the
user is not normally concerned with these decks. They are presented in this
report to document the complete model since it would be prohibitive to show
all possible FORTRAN compiler listings for the different combinations of
selectable options.
95
-------
-IDE NT CONVCT
'•'• I N S E h> T H N G ft , 1
« , A 7
"INSFwT HMGfl.17
A 7 = A 4 * . 5
MNSEPT HARMOD.l
CALL CONUV(U,HGU,-1)
J>OELFTE HARMOD.2
CALL RAPUV (V.HGVtO ) $ CALL UVC AL ( Z^ »-1 . ) f CALL CONUV(V,HGV»0)
00 260 L = 1,N"AXLAY $ JL=LEV(L)+1 % CALL X MI T ( ZM ( JL ) » V ( JL ) « NO , 1 , 1 )
260 CONTINUE
*^EFORE HMGfi.118
DIMENSION VU(2) S ICON = 0 * GO TO S * ENTRY CONUV <*, ICON =1 $ CY=0.
5 COMTINIIF
*DELFTF HAPMQD.12.HNGR.125
17 IF(ICON.EQ.1) GO TO 18
GO TO 20
18 VU(1)=VU (?) = !.
IF ( (MM.EO. I .OP.MP.EQ. I .AND.N.EO.O) .OR. (NM.EQ. I .OR.NP.EO. I .
».-!)> VU(2)=?.
IF ( (MOD(M+1,ME) .LE.2.AND.N.EQ.-1 ) .OR, {^00( 1+ 1»ME) .LE.?.AtML).N,EO.O)
*) VU (1 )=0.
VU(1-N) = <6(NM+JK)-8(NP + JK) )*VU(1-N)
CX = -A7* (fl (JL) «VU ( 1 )+ZST ( JL) *VU (2) ) S ZM ( JL ) -T.* (JL) +CX
*IDENT CONVCU
^DELETE CONVCT.13.CONVCT.16
IF(N.EQ.O) GO TO 181 S IF(NM.Eti.I.OP.NP.EQ.I) VU(?)=?.
IF(MOD(M+l.ME).LE.2) VU(1)=0. $ GO TO 18?
181 IF(MM.EQ.I.OR.MP.EQ.I) VU(2)=2. $ IF (MOD(I + 1»NE) ,LE . 2) VU(l)=n,
182 CONTINUE
*IDENT CONPRT
*OELETE CONVCT.1
B »A7,CY»CZ(800)
»INSERT BC8.7
IF (MOD (IT,IPMOD) .EQ.O.AND. IT.GE.IPO) CALL PKTMAX (CZ ( (JL-1) *N.O*1') «
*4HNONE»1H f42HCONVECTIV£ TERM IN U COMPUTATION FOR LAVERf1»42, IT)
-------
*INSFRT HNG8.54
CALL PPTMAX (CZ ( (L-l ) *MO+1 > ,4HNONE. 1H , 4?HCOMVFCT I VE TERM j N v COHP
*UTATION FOR LA YF R , 1 . 4? , I T )
*DELETE HNG8.1?!
J=(L-1)*NO+I S IF(H( JL) ,GT.O.)GOTO 15 $ ZM ( JL ) =CZ ( J ) =0 . $ GO TO ?0
CZ(J)=CX $ CY=AMAX1 (CY ,ABS (CX) )
799 IF (CY.LT.1000. ) GO TO P.OO !. PRINT 100ft, CY !», GO TO flOl
* INSERT HNG.?46
801 CONTINUE
*DELETE UTIL.R6
?000 FORMAT (Al, I4,3X15E8.2)
CONPRU
E CONPRT.l
B »A7,CY,CZ(300)
^INSERT HNGB.PO
CALL XMIT(-0.,CZ(JL+1) ,MAXDIM,0»1)
-^INSERT RCR.7
IF (MOD ( IT.MODOUT) .EQ.O.AND. IT.GE ,ITO)
* CALL INPOUT (CZ (K-NO+1 ) . JL»8,IT,CFAC)
»DELFTE COMPRT.2»CONPRT.3
^INSERT h!MG8.4Q
CALL INPOUT (CZ ( JL) t L t P • I T » CFAC )
*OELETE CONPPT.4fCONRRT.5
^DELETE CONPRT.fr
IF (CY.CST. 10. ) RRINT 100ft«CY
*OtLFTF COMPRT.7
^DELETE CONPRT.R
IF (H ( JL) .GT.O. ) GO TO IS $ ZM ( JL ) =CZ ( JL ) =0 . % GO TO ?n
^DELETE CONPPT.9
CZ(JL)=CX $ CY = AMAX1 (CY, ABS(CX) )
^DELETE CONPRT.10
^RESTORE UTIL.R6
*IDENT CONPRV
^DELETE CONPRU. 4
*CALL INPOUT (CZ (K -NO* 1 ) , JL» 1 0 » I T» CF AC )
^DELETE CONPRU. 5
CALL INPOUT(CZ (JL) tLtlOt IT, CFAC)
-------
APPENDIX D
LAYER DISAPPEARANCE MODIFICATIONS
FLAT, FLATA AND FLATB SERIES
The FLAT, FLATA and FLATB modifications allow layer disappearance cases
to be handled reasonably in the Phase II computation. The procedures main-
tain a closure status flag array based on a minimum allowed thickness, reset
lesser thicknesses to the minimum thickness, and reset the old velocity com-
ponents on closed sides to zero. For any closed side, the adjacent c;i
heights are recomputed without smoothing using the modified thickness and old
velocity fields. After the recomputation of £i height, new velocity compon-
ents are computed in the usual fashion.
The initialization section equivalences the flag field IZi to the ZSTi
array, which is available at this phase of the computation for this special
use. The LAYS card and interpretation section are modified to include a mini-
mum thickness value, HMIN, in columns 71-75 and to calculate the constant
HMIN4. SETHUV is modified to include an array index offset parameter so that
array Zi is used in the computation for initialization only. Thereafter, the
tidal and volume source specifications and thickness calculations are all
modified to use ZMi (rather than Zi).
The thickness computation block is most heavily affected by the FLAT
modifications. The bulk of the layer disappearance procedures follows the
call to CXMIT that is used to set the Hvi boundary condition. The flag field,
IZi, is initialized to display code zero (33B). For all sea points the average
98
-------
thickness must be greater than HMIN else the flat position is set to indi-
cate cell closure (30 B). If the individual thickness of a cell side is less
than HMIN, the flag field is incremented by an integer coded to show the par-
ticular side closed in the flag position except on full closure. The top side
of a cell is indicated by an increment of 2, the right side by an increment of
8, the bottom by 4 and the left by 1. These integers are independent so the
sides that are closed may be determined from the sum (TABLE D-l.) For each
closed side the normal velocity component is set to zero and the thickness if
less than HMIN is reset to HMIN. The boundaries of the flag field are ad-
justed based on the sequence of operation. Each cell is viewed only in terms
of its north and west sides so that the computation does not progressively
contaminate itself. The CXMIT routine is called to put the usual boundary
condition back on the Ui, Vi, HGUi and HGVi arrays. This precaution must be
taken since some of the recomputed £i may be on the boundary. The flag field
is used to control the recalculation of £i. Note that recomputation cannot
use the smoothed £i since the ZM array is overwritten during the calculation
of the required £i.
The non-availability of the smoothed £i has beneficial side effects: cer-
tain boundary configurations occur with the layer disappearance procedure which
would require specialized computation of the smoothed component, but the limit-
ed testing of such specialized smoothing computations did not produce signif-
icantly different results so no smoothed component is included. Instead the
old £i is used from array Zi for the recomputation. The array ZMi is trans-
mitted into array Zi after this computation by routine XMIT. The naming con-
ventions for these modifications are similar to those described for CONVCT
series, except that subsequent corrections must be named PLATA rather than
FLAT to obtain the correct sequence ident.
99
-------
TABLE D-l. TIDAL FLAT CODE CONVENTIONS
FLAG NO.
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16+
17
18++
ASCII
CHARACTER
0
1
2
3
4
5
6
7
8
9
+
-
*
/
(
)
X
L
OCTAL
33
34
35
36
37
40
41
42
43
44
45
46
47
50
51
52
30
14
CLOSURE STATUS
OPEN
W
N
NW
S
SW
SN
SNW
E
EW
EN
ENW
ES
ESW
ESN
CLOSED (FLAT CELL)
CLOSED (FLAT CELL)
CLOSED (LAND CELL)
+ Flag No 16 exists only temporarily before being replaced by Flag
No 17. it represents a partially closed cell being fully closed
by subsequent full closure of adjacent cells.
++ Flag No 18 exists only for printed display purposes.
100
-------
FLATIO AND FLATIP SERIES
The FLATIO series of modifications provide an output capability to dis-
play the results of the tidal flat closure during the layer disappearnace pro-
cedures. The general technique has been to save the flag array on the inter-
mediate data set, LU1, whenever a change of condition occurs in the particular
layer. An attempt was made to provide multi-layer generality as in the layer
disappearance procedure proper, but no testing of this function has been per-
formed. A subroutine DISPLE with several entry points was designed to read,
write and process the flag array as type 9 data records on LUl. Since the
velocity fields and height fields are modified by the layer disappearance
procedure, this procedure also implements the means of saving these modified
arrays on the intermediate data set as type 3 data records at the usual tape
save interval. At the conclusion of the run the intermediate data set file
is rewound and the character flag fields read in and displayed. A call to the
DSTART entry of DISPLE has been added to the initialization block to initial-
ize the display package. The thickness computation block contains a call to
routine DISPLE to process the flag field produced by the FLAT, FLATA series
modification. A call to INPOUT is added to save the modified arrays at data
save steps. The output block is expanded to incorporate a call to entry
DISPRT of routine DISPLE, which rewinds LUl and prints the flag fields saved
during the course of the run.
Routine INPOUT is modified to handle the special type 9 data records,
which unlike other arrays are packed as 6 bit display code characters rather
than 24 bit integer quantities. The record also contains in the data portion
just following the normal header two full word variables, the layer and the
time. Special control is incorporated to avoid the 24 bit compression and
101
-------
scaling usually provided by INPOUT.
Routine DISPLE (Block N;FLATIO.4,66)
Routine DISPLE provides the means of saving the layer disappearance flag
field as character arrays on the intermediate data set. The routine requires
two parameters, JL and IT, which represent the layer and the time, respective-
ly. Common/DLAY/ is introduced to save the time of last change, the layer and
the last condition of the character flag array. A comment in the program in-
dicates the required dimensionality for the members of common/DLAY. On the
main entry the IZi array is converted to a display code character array and
then compared to the old character array for the same layer. If they are the
same, the routine is complete and returns to the main program; otherwise,
the routine calls INPOUT for type 9 data record setting NO=1. When INPOUT is
complete NO is reset to its normal value. The time level and the character
array are moved into the space reserved for that level's old character array
and return is made to LAYERS.
Entry DISPIN provides a restart capability, which has been tested and
found to function properly. Its use requires special coding to recover the
record desired, because the character arrays are saved at irregular times.
Entry DISPRT is called at the end of the main program after the time loop
is complete. It provides the capability to read back in the saved type 9
records and restructures them into printable form using a DECODE* statement to
unpack the characters. One print line at a time is encoded and printed. The
process cycles through all arrays and the program ends with the end of file
error message printed by the IOERR routine called from INPOUT.
The actual form of the display of cell closure depends on the number of
layers and the number of rows and columns in the defining mesh. If the
102
-------
required array size is sufficiently small, up to the number of layers plus one
array displays may appear within the print line width. In any event at least
one array will appear within the width of a print line, and this single array
display is described.
The flag fields are printed in the form in which they are saved, i.e.,
as a series of two dimensional character arrays reflecting the changing char-
acter of the tidal flat region in change steps. The array includes the orig-
inal land areas as the double character "L" (flag 18 of TABLE D-l), the
closed tidal flat cells as double character "X" (flag 17) and the partially
closed cells as single character flags indicating closure status (flags 1-16).
Printed headings for each character array show the layer number for the array
and the associated computation time in seconds during the course of the model
run.
The DSTART entry computes constants required by the other calls to the
routine so must be called first. Comments in the program define the con-
stants.
LISTING FOR LAYER DISAPPEARANCE CORRECTION SETS
The following listing is the FLAT and FLATIO series correction sets in
card image form as originally entered into the library. The spacing between
correction sets is artificial and the sets are presented in a logical order
rather than the order in which they were entered on the library. The de-
scriptions rely on references to these listings and their relationship to the
listing of the basic Phase II program presented in Figure 8. These correc-
tion sets reside on the library and the means of using these sets are pre-
sented in TABLE 11.
103
-------
»IOENT FLAT
*INSERT HNGS.S
•t(IZtZST) $ DIMENSION 17(1)
*DELETE HNG8.14
18 DECODE(7*,1018,IC)HL,ALPHA,RHO,HMIN $ GO TO 1
MNSEPT HNG8.23
38 HMIN4=HMIN*4.
»DELETE HNG8.25
CALL SETHUV(-3*MAXDIM)
*DELFTE HNG8.32'
ZM(II)=ZM(II)-ZZ*A4 $ IF (JL.EQ.l)GO TO 50 % JL=JL-1 S GO TO 41
^DELETE HNG81.6
K=HTZ(I) $ IF(K.GE.JL.ANO.ZZ(JL).NE.O.) ZM(I+L)=ZZ(JL)
"DELETE MNG.134
140 CALL SFTHUV(O)
^INSERT BC8.6
DO 142 JL=1-MAXLAY $ L=LEV(JL)*1 * CALL XMIT(33B,IZ(L),NO»0»1)
14? CONTINUE
DO 151 M=1,MF $ IM = NNM(M) t IP=NNP(M) J, IF(IP.LT.IM) GO TO 151
IMM= (M-l)*NE+1
DO 150 I=IM,IP S IF(HT?(I) .LE.O.)GO TO 150
MM=I-NE $ IF(MM.LE.O)MM=I S NM = MAX0{I- 1 f IMM)
-------
150 CONTINUE
151 CONTINUE
DO 152 JL=1»MAXLAV $ K=LEV(JL)*NO
CALL CXMIT(HGU(K-NE),HGU(K),HGU(K),0.,NE,-1,-1,-1)
CALL CXMIT(HGV(K-1 ),HGV(K),HGV(K),0.,ME,-NE»-NE,-NE)
CALL CXMIT( U(K-NE), U(K),HGU(K),0.,NE,-1,-1,-1)
CALL CXMIT(V (K-l ),V (K),HGV(K),0.»ME,-NE,-NE,-NE)
152 CONTINUE
DO 161 M=1,ME * IM=NNM(M) $ IP=NNP(M) $ IF(IP.LT.IM) GO TO 161
IMM=(M-1)*NE*1 $ DO 160 I=IM,IP 1 IF(HTZ(I).LE.0.)GO TO 160
MM = I-NE $ IF (MM.LE.O)MM=I $ NM=MAX0(I-1.IMM) $ JL = HTZ(I)S ZZ=ZX=0,
153 L = LEV(JL) $ II=L*I S IF(IZ(I I)-33B) 155,157,154
154 I1=L+MM $ I2=L+NM $ IF(ZZ.EQ.O.) ZZ=ZX
ZX = ZZ = ZZ-A4*(HGU(I I)*U(I I)-HGU(I 1)*U(11) +HGV(I 2)*V(I 2)-HGV(II)*
155 Z(II)=Z(II)+ZX $ ZX = Z(II)-ZM(H) $ ZM(I I)sZ(11) S GO TO 159
157 ZM(I I)=ZM( II)+ZX
159 IF (JL.EO.l ) GO TO 160 5B JL = JL-1 $ GO TO 153
M 160 CONTINUE
g 161 CONTINUE
DO 171 JL=1,MAXLAY $ L = LEV(JL)+1 "B CALL XM I T ( ZM ( L ) , Z ( L ) , NO , 1 , 1 )
171 CONTINUE
*OELFTE HNG.249
SUBROUTINE SETHUV(KZ)
^DELETE HNGR.66
»IDENT PLATA
*DELETE FLAT.19
144 IF(I2.EQ.II) GO TO 146 * I F( IZ(I?) .NE.30B) I Z( I ?) = IZ(I 2) +4
^DELETE FLAT.21
146 IF(HTZ(MM) .LT.ZZ) GO TO 149 $ IF ( I Z(I I) .EQ.30B) GO TO 147
"DELETE FLAT.24
147 IF(Il.EQ.II) GO TO 149 ¥ IF(IZ(I 1) .NE.30B) IZ( 11)=IZ(I 1) +8
-------
MDENT FLATIO
«• INSERT HNG8.P3
CALL DSTART (L.ITS)
FLAT.4R
CALL DISPLA(JL.IT)
HNG.P47
REWIND 1 "K CALL OISPRT(1,IT)
"INSERT RAPSTAR.50
SUBROUTINE DISPLA(JL»IT)
"CALL CONTRL
*CALL PLANK
*» IW(500)
COMMON/DLAY/IFT(3)»IFL(3)»IFC(15M
C DIMENSION IFT(MAXLAY)«IFL(MAXLAY),IFC(K0«(NO.OF LAYERS+l))
DIMENSION IZ (1) »JZ( 1)
EQUIVALENCE (Z^(1) • JZ(1) ) , (ZST(1).IZ(1))
DATA (I8=1H )
L = LEV(JL) S K=l "5 KK = 1SO
M DO 100 I=1»NO«150 * IF(1*150,GT.NO) KK=NO-I + 1
§ IMsl+L 1 IP=IM*KK-1
ENCODE(KK,1000»IFC(K))(IZ(II),I1=IM,IP) * K=K+15
100 CONTINUE $ L=JL*KO
1000 FORMAT(150R1)
DO 140 K=1,KQ * IF(IFC(K),NE.IFC(K+L)) GO TO 150
140 CONTINUE $ RETURN
150 CALL INPOUTUFC(l) »JLt9,IT»XX)
CALL XMIT(IFC(1)»IFC(L+1),KO,1f1)
IFT(JL)=IT % IFL(JL)=JL
RETURN
c
ENTRY DISPIN $ IFP=0 % L=-JL $ K=JL*KO+1 * GO TO 300
C
ENTRY DISPRT
190 IFP=1
200 L=-IFP
300 KKrIT S CALL INPOUT(IFC(1) »L ,9»KK,XX) $ IF(IFP.EQ.O) GO TO 310
-------
IF fIFC (K) ,KOt1 ,1) * IFT(L)=KK $ IFL(L)=L
IF(IFP.EQ.O) RETURN $ IFP=IFP+1 $ GO TO POO
350 ?Z-JL=L * IT=KK % L=0 * K=l $ DO 3M KM=l,MZ * IF (KM .E Q. 1) GOTO 355
K=KM-I $ L=LF.V (K) +MAXDIM * ZZ = IFL $ K=K^KO*I
355 KK = 1BO * DO 360 I=1»NO,150 S> I F ( I + ISO . GT , NO ) KK = NO-I + 1
IM=I+L ^> IP=IM + KK-1
DECODE (KK, 1000 . IFC (K) ) ( J7 ( I I ) , I I = IN',I P )
-------
»INSERT UTIL.12
M=FACOUT
^BEFORE UTIL.15
IF(ITY.GT.7) GO TO 310
^BEFORE UTIL.34
310 BUFFER INdtl) (I W ( 1 ) , IW (LEN ) ) $ I REC= I hfEC + 1
IF(UNIT(l)) 350,320,330
320 IERP=1 $ GO TO 340
330 IEPR=2
340 PRINT H30»IW(1) ,IW (2) , ID1,ID2 $ CALL I OERR ( IERR , I REC ) <6 GO TO 310
350 ID2=SHIFT(IW(?),-48) * IF(ID2.NE.9) GO To 310
CALL XMIT ( IW (5) , A( 1) ,M.l ,1 ) <6 LAYER=IW(3) $IT=I*<4) * uO TO 700
*INSERT UTIL.35
IF ( ITY.GT.7) GO TO 510
^DELETE UTIL.51
500 CONTINUE * GO TO 600
510 CALL XMIT(A(1) , IW(5) »M,] ,1) $ NW=M
520 BUFFER OUT (1,1) ( I Wd) , IW(NW) ) -fi
M IF(UNITd)) 700,530,540
g 530 IERR=3 t GO TO 550
540 IEPR=4
550 CALL IOERR(IFRR,IREC) S GO TO 520
^INSERT uTiL.5?
700 CONTINUE
»IOENT FLATIP
^DELETE FLATI0.2
^INSERT FLAT.50
CALL DISPLF(JL,IT)
IF(^OD(IT,MDDOUT) .EQ . 0.AND.IT.GE.I TO)CALL INPOUT(Z(L) ,JL.3«IT»CFAC)
^DELETE FLATI0.4
SUHROUTINE DISPLE(JL.IT)
^DELETE FLATI0.21
150 N0=l "R CALL INPOUT ( I FC ( 1 ) , JL , 9 , I T , X X.) i, NO = NE»M£
^DELETE FLATI0.68
^INSERT UTILA.3
-------
IF(ITY.EQ.9) GO TO 160
'INSERT UTIL.24
160 ID2=SHIFT(IW<2),-48) $ IF(I02.NE.9) GO TO 100
CALL XMIT(IW<5)tA<1)»M,i,1) $ LAYER=Irt(3) 5 IT = IW(4) 5, GO TO 600
'DELETE FLATI0.69,FLATI0.75
'DELETE FLATI0.76
'BEFORE UTIL.37
IF ( ITY.E0.9) GO TO 405
'INSERT UTIL.37
405 CALL XMIT(A(1) , IW(5) ,M,j,1) $NW = M + 4 * IW(3)=LAYE& $ IW(4)=IT
GO TO 47B
'RESTORE UTIL.51
'DELETE FLATI0.77»FLATI0.83
'DELETE FLATI0.84
-------
APPENDIX E
.MONTE CARLO MODIFICATION
MONTE
The Monte Carlo modification uses a pseudo-random number generator to ef-
fect diffusion in an advected field of lagrangian marker points. HN Phase II
generated velocities are interpolated to the prior position of each marker
point, then modified by a random increment to determine the new position. The
version on the HN library now allows up to 1000 marker points, which are dis-
played through time to show the dispersion of each source array. In the pres-
ent code a single marker point is introduced at each time step simulating a
continuous source. The procedure operates as a part of Phase II and the
printed output feature is an integral part of the modification rather than a
separate option as in the convective and tidal flat extensions. The marker
point position data are not saved on the intermediate data set although a data
saving feature would fit within the present framework of data record types.
The two forms of printed output generated are described as the modification
to the Phase II blocks and additional routines are presented.
The initialization block contains the 1000 position arrays,,XS and YS,
that are inserted in control common immediately before the 500 word subroutine
scratch area and are used to keep track of the present position of each mark-
er point. The number of points followed depends on the number of SORC cards
and the number of time steps. If the number of SORC cards times the number of
time steps exceeds 1000 the XS and YS dimensions must be expanded. The SORC
card used with the Monte Carlo extension was modified to include 2 additional
110
-------
floating point values per card in order to specify the exact source point lo-
cation. This change required modification of array SORC and the DECODE pro-
cessing for the SORC card. The £ computation block contains the initializa-
tion of the position array, but the bulk of the Monte Carlo procedure follows
the thickness computation and is considered part of the thickness computation-
al block. When the layer disappearance procedure is included the Monte Carlo
procedure follows the tidal flat procedure. Each marker point considered for
a particular step is advanced by modification of the position arrays in accor-
dance with the interpolated and randomly incremented velocity and the time
step. The function WAYTUV does a simple bi-linear interpolation between the
nearest four like velocity components. Information is retained during the
position computation to provide a print line of statistics for each source
at each time step giving the x-y components of center of mass, distance to cen-
ter of mass from the source, horizontal width of the mass and standard devia-
tion of the marker points about the center of mass.
In addition to the statistical information, a second printout is produced
for each source at each print step showing the point positions as a two dimen-
sional character array in which empty rows are deleted. Each printed array
presents the essential portion of a variable dimension square window view of
the point distribution for each source card. Each element of the 50 row by
100 column array represents one five-thousandth of the entire area of the
square window. The dimensions of the window are determined by the extent of
the point distribution. The characters printed each represent a count of the
marker points falling within each sub-window array position at the print time.
Total row counts appear to the right of the printed output and total column
counts appear in an offset array format with each count right-justified
111
-------
immediately below its respective column. Empty rows are not printed so a
row identification of absolute floating N index and relative position within
the window is printed on the left. Similar absolute M index and window
relative information are provided above the columns. The absolute and rela-
tive information refer to the upper left hand corner of the print position.
The output block prints the final position of the random number generator for
use in restarting the model.
Function WAYTUV (Block O:MONTE.64,76)
The function WAYTUV returns the value of the interpolated velocity com-
ponent based on the input parameter list. The first parameter is the left-
most index in the x direction, the second parameter is the uppermost index in
the y direction, the third parameter specifies an offset which determines
whether a u or v computation will be made. The routine computes the bracket-
ing indices around the x (or y) position, then computes linearly interpolated
values in the x direction and finally linearly interpolates the resulting
values in the y direction.
LISTING FOR MOTE CARLO CORRECTION SETS
The following pages present the MONTE correction sets in card image
form as originally entered into the library as modifications to Phase II.
These correction sets reside on the library and the means of using these sets
are presented in TABLE 11.
112
-------
*IDENT MONTE
MNSERT HNGfl.l
* f XSUQOO) .YS(IOOO)
»DELETE HNG8.5
1 .TIDSPDU) , I FLOW (B»30) .FLOW { 4 » 30 ) » ISORC < 8 . 1 0 ) . SORC (4 , 1 0 ) .
*INSERT HNG.15
* « (NPOC = 61 )
»INSEPT HNGR.ll
POC=NPOC-1
"DELETE HNG.44
DECOOE(76,10?O.IC)(ISORC(I,ISNO),I=1,R),(SOPCU.ISNO).I=1.4)
*INSERT HNG.1?0
JL=0 $ DO Rl K=1,ISNO * NPTS=SORC( 1 «K)
IF (IT.NE.ISORC (7.K) ) GO TO 81
XSO=SORC(3»K) * YSO=SORC (4»K)
DO PO 11 = 1. NPTS 1> I = II+JL 5 XS(I)=XSO H> YS(I)=YSO
BO CONTINUE
81 JL = JL + N'PTS
*INSFPT FLAT. SI
IF (ISNO.FQ.O) GO TO 1^9 % JL=0
DO 198 K = 1,ISNO * NPTS=( IT-ISOPC (7.K) ) /TTINC + 1
IF ( IT.LT.ISOMC (7»K) .OR. IT.CiT. ISORC (H,K) ) ',0 TO 198
PC = 1./NPTS * XSfi=SORC(3.K) S YS()=SOKC ( 4 , K )
XMIN = YMIN = MAXDIM i, XM A X = Y M AX = X SUM= Y SlnVl = 0 .
DO IRQ 11 = 1. NPTS i I = IUJL
UH = PANF (FTN) -.S <* V « = H AMF ( F TN ) - . S
USx = WAYTL)V(XS(T)-.5.YS(I),0)+iJP
VSX=WAYTUV(XS(I)»YS(I)-.5«MAXDlM)+V'->
IF (UP. Nf- .USX.) XS ( I } =USX«A4 + XS ( I ) % ASU'V'= XS ( I ) » XSU^
IF (VP.Nf .VSX) YS ( I ) =-Vc,X*A4 + YS ( I ) * YSi )'' = YS ( I ) + Y SUM
X M A X = A M A X 1 ( X S ( I ) . X ry| A X ) * X M I N = A M I !M ] ( X S ( I ) , X M 1 <\' )
CONTINUE $ MP=XMAX + 1. t" NP = YMAX+1. f i^^-XMIN f,
O. ^ DO i«? II = I.NPTS ^. I = II*JL
SS = XS(I)-AX5 SS*S IG?X f SS = Yb(I)-AYS >S I G?Y = bS*bb + 5 I G?Y
182 CONTINUE % SIG?X=PC*SU-i?X % S I G?Y = PC*S I G ? Y
AX.S=AXS-XSn * AYS = AYS-YSO
-------
RVEC=SQRT(AXS*AXS+AYS»AYS)
RANG = AMOO( ROT- AT AN? ( A YS,AXS)»57. 29577951 3 »3,4)
IF(MOD(IT,IPMOD).NE.O) GO TO 198
PC=50./W $ NY=XMIN*1000. $ MX=
-------
•IDENT MONTE
•INSERT HNGB.I
* «XS (1000) fYS( inOO)
•DELETE HNG8.5
1 «TIDSPD(4) .IFLO* .OR. IT.GT. ISdRC(H,K) ) -.,0 TO 198
PC=1./NPTS * XSO=SORC(3»K) $ YSO=SOKC(A,K)
XMIM = YtviIN = MAXDIM S XM A X = Y M AX = X SUM= YS Uf" = 0 .
DD 180 II=1«NPTS * I=TI+JL
VSX=WAYTUV(XS(I)»YS(I)-.i:S»MAXOlM)+Vi->
IF (UP. Nf .USX) XS ( I ) =USX«A4*XS ( I ) % X SU>V'= XS ( 1 ) + X SU^
IF ( VP.NF . VSX) YS ( I ) =-VSX*A4 + YS ( I ) * YSi !'•• = Yb ( I ) + YSUM
XMAX=AMAX1 (XS ( I ) .X'^AX) * XMIN = AMIIM] (XS ( I ) iXN'lN)
YMAX=A^AX1(YS(I)«YMAX) t YMIN=AMIN1(YS(I),YMIN)
1HO CONTINUE $ MP=XMAX + 1« t. NP = YMAX+1. f i^N^XMIN f, iMvi = YMlN
M = ni = AMAXl ( XMAX-X^IN, YM/i X-Y^'IN) * IF(M.Fij.Q) M = »V = 1 .
AXS = XSL)M*PC ^ XMAX=X(^AX-AXS % AYS = YSUM*;-i-C ' Y^AX=YMAX-o,
SIG?X = SIG?Y=0. "> DO 18? 11 = 1, N^TS ^. I=TI*JL
SS = XS(I)-AXS •*• c:.IG?X=SS->SS*SIG?X * SS = Yb(I)-AYS ?S I G? Y=
18?. COMTTNUF f S I G?X=PC*S I G?X % S I G?Y=PC»S I G? Y
AX.S=AXS-XSO f AYS = AYS-YSO
-------
RVEC=SQRT( AXS»AXS + AYS*AYS)
RANG=AMOD(ROT-ATAN2(AYS,AXS>*57.295779513,360.)
PRINT llflO»NPTS,NM,NP,MM,MP,XMIN,XMAX,YMIN»YMAX,SIG2X,SIG2Y
* ,PVEC,PANG,AXS.AYS.UR,VR
1180 FORMATU5,8H POINTS.4 I 3»4H, X=?F7.3»4H, Y=2F7.3»5H, S?=2F7.3»
*4H, S=F7.3,1H,F5.1,4H, A=?F7.3,4H, P=?F6.4)
IF(MOD
-------
* INSERT BARSTAR.SO
FUNCTION WAYTUV (X, Y,NOFF)
*CALL CONTRL
*CALL BLANK
X1=M=X $ X1=X-X1 $ X2=1.-X1
Y1=N=Y $ Y1=Y-Y1
NP=MINO (N* ! »NE) $ N=MAXO(N,1)
$ M = MAX;0(M»1)
s MP=
-------
APPENDIX F
PHASE IV PROGRAM SERIES
THERMAL ADVECTION PROGRAM
Phase IV consists of the thermal advection routines which were built
from a version of the Pederson-Prahm advection .by method of moments and a
thermal treatment heat budget package originated by T. Laevastu. A major
difference between this and other phases of the HN system is the means of stor-
ing array elements by sequential column positions within sequential rows.
Phases I and II and the intermediate data set use the opposite convention of
sequential row positions within sequential columns. This change of conven-
tion requires modification to the I/O routine indexing both on input to
include conversion to the internal storage convention and on output to pro-
vide conversion back to the convention defined by the intermediate data set.
Program ADVECT (Blocks A-E:A2.1,A1.27)
The program ADVECT consists of blocks A-E as identified in Figure 6.
This routine acts as the control program not only for the advection by method
of moments but also to control the thermal treatment of the advected excess
heat quantity. Logical units are specified for the thermal advection routines
as follows:
Logical Unit Description
INPUT Punched card input file
OUTPUT Printed output file
1 (TAPE1) Intermediate data set file
2 (TAPE2) Thermal advection data set file
3(TAPES) = OUTPUT
116
-------
The blocks A-E described below are arbitrary modules of the program chosen
for convenience in description.
Advect Initialization (Block A:A2.1,36)
The program card specifies the logical units. Common /A/ includes gen-
eral control and some specific control parameter for the various subroutines;
blank common holds the major array variables associated with the advection
procedure and the heat/temperature/temperature anomaly array, T; common/1/
provides space required by array BDVECT for the advection by method of mo-
ments; common/ABC/ provides the major arrays associated with the thermal
process except the temperature array; and common/THERM/ provides central
values to the thermal process. Dimensionality and equivalence classes are
established for local routine use. Various constants are defined including
some repetitive definitions due to temporary modifications. The modification
TA800, for example, provides (when present) the capability to handle a 20 by
40 array. Without this modification only 300 position arrays of fewer than 30
rows are allowed. The master control cards are read in character mode and im-
mediately printed for verification. Control then transfers to processing the
specific card identified from the table of allowable card names, NAME. The
position of the card names in this array is used as an index for routing con-
trol as in the other phases. Some default values are provided where no value
is specified. Too many wind cards will not cause the routine to abort, but an
informative diagnostic is provided.
If ME is greater than or equal to IDIM an information diagnostic is print-
ed before early termination. Too many SORC cards cause immediate termination
without diagnostic. Unidentified cards are ignored to avoid having to provide
detailed recognition of cards not used in this phase: FACT, TIDE, FLOW, LAYS.
117
-------
All accepted cards return to process the next card except the COMP card, which
causes the execution of the rest of the routine. Additional constants are de-
fined and the first and second floating values from each SORC card are multi-
plied by the time step to allow the deviant usage of volume and concentration
(heat) per second rather than per step. Arrays are initialized by calls to
routine XMIT, and a restart procedure is invoked when ITS is not equal to zero.
The restart modification should be present only when the program is being start-
ed from a previous thermal advection data set. This modification causes type 20
data records to be read from LU2.
Thermal Treatment (Block B;A3.4,9)
From the beginning of this program block to the end of the output block the
processing occurs for each time step specified on the TIME card. The first time
through the loop the thermal treatment block is partially bypassed to avoid ad-
vection or thermal treatment in the first step. For all subsequent steps the
following procedures occur: the HN Phase II velocity fields are read from the
intermediate data set using INPOUT. If this is a print step, the velocity field
is printed by PRTMAX, the thermal treatment occurs in routine THERMAL, and the
output cell temperatures calculated in THERMAL are converted to a temperature
anomaly by interpolating the width factors produced in the last advection step.
For all steps including the first step the excess heat field is computed from
the temperature anomaly field. Then, for the first step only, a bypass occurs
to avoid the duplex advection procedure. Zero values in most arrays represent
ambient values rather than true zero values. This feature restricts the field
for computations and printed output to essential values required for computation
of non-ambient cases.
118
-------
Duplex Advection (Block C:A2.47,52)
Block C processing commences with three calls to routine XMIT. COMMON/1/
is used as temporary storage to swap the Q and V arrays, minimizing the changes
required in routine BDVECT. Routine BDVECT does the minimal advection compu-
tation for the excess heat component (zeroth moment only.) On return XMIT is
called again to reverse the swap. BDVECT is called again to advect the V ar-
ray containing variable V applying the full method of moments using the same
velocity field as in the first call for heat advection. This duplex advec-
tion block is bypassed for the first step. In the restart cases this prevents
duplication of the last step of the prior run. In the normal start the bypass
avoids advection of no quantity as the sources are applied subsequent to the
advection step.
Sources; Temperature Anomaly; Wind (Block D:A2.53,78)
At and after the time the earliest source is specified to start, all
source specifications are reviewed to determine whether they apply to the cur-
rent step. When a source is relevant it is applied. The start time of the ear-
liest source should correspond to the start time of the run on the TIME card.
If either center of mass indicator computed by BDVECT exceeds .5 in absolute
value at any point in the field, a warning indicating that the time step for the
process is too long is printed and the routine is terminated. The center of mass
is computed relative to the center of the cell in grid step units. The advected
temperature anomaly is computed from the advected excess heat and the advected
volume. Routine XMIT zeros the wind field, then the wind specifications are ap-
plied such that any wind scheduled to start at or prior to the current step ac-
tually begins to have effect in the subsequent step and carries its effect through
to the indicated end time. A wind intended to begin prior to the end of the
presently active specified wind will not take effect until the active wind stops.
119
-------
Wind cards must be in order of start times and their start times and end times
must be precise multiples of the run time step. Multiple wind subfields sched-
uled to start at the same time will take effect at the same time with the most
recent specifications determining in overlap cases.
Output Processing (Block E:A2.79,Al.27)
Data save steps are defined by the initial output time ITO and the inter-
val MODOUT from the TIME card. Routine INPOUT is called to place data record
type 20 on the thermal advection data set file, LU2, at each data save step.
Print steps are defined by the initial printed output time IPO and the modular
interval IPMOD from the time card. Routine PRTMAX is called to print several
arrays. If the current time exceeds the end time for the run, TE, the job
terminates; otherwise, the program loops for the next time step.
Routine PRTMAX (Block F:Al.28,61)
The major differences between PRTMAX and the routine of the same name in
Phase II are: for most arrays in the thermal advection, zero values are car-
ried in lieu of ambient values; rows consisting entirely of zeroes are not
printed; the level parameter is missing; and the method of indexing. Since
the internal storage conventions differ, the indexing to achieve printed out-
put is also different.
Routine THERMAL (Block G:Tl.9,TH.324)
All required constants and variables are either local or in»common. De-
fault values are provided for all constants, though for a particular case these
should be redefined to reflect the specific conditions desired. In the current
version, the ambient water temperature is assumed constant as is the effect of
solar heating. (The imposition of a diurnal cycle is such an individual
phenomenon it is left to be determined on a case by case basis.)
First the heated fluid temperature anomaly is interpolated to the cell
120
-------
temperature based on the distribution of the fluid within the cell and the am-
bient water temperature and the saturated and moist water vapor pressures are
defined for affected cells based on water temperature and relative humidity.
Effective back radiation (the net of oceanic and atmospheric longwave
radiation), convective transfer of sensible heat (conduction from the ocean
surface to the air molecules), and transfer of energy by evaporation (formation
of atmospheric latent heat) or condensation (recovery of heat from atmosphere)
are the heat budget components computed. The heat budget components and the
volume associated with the heated fluid determine the temperature decrement
for the step. If the cell temperature is adjusted to less than or equal to the
ambient temperature, the volume, temperature and associated center of mass
and width factors are set to zero. Several arrays are then printed by PRTMAX.
All of the arrays are printed in restricted format due to the convention of
maintaining zero rather than ambient values for unaffected cells or cells not
required in a computation. Return to the main program is made with the array
T=Q holding the value of the temperature in the cell.
Routine BDVECT (Block H:A2.97,A3.26)
A single parameter is required by routine BDVECT to distinguish between
calls for advection of the excess heat and for advection of the associated
fluid volume. The parameter is used specifically to bypass first and second
order moment computations for the heat advection, and to tailor the initiali-
zation processes for the difference in requirements. COMMON/1/ is used to hold
seven intermediate accumulation computations by row. Each accumulation array
holds space for present, prior and succeeding row position accumulations plus
space for adjacent outside column position accumulation. The computations are
done row by row with the accumulation arrays cycled to the next row setup after
transfer of gross accumulation for prior row into the advection array. During
121
-------
the cycling process, the accumulator positions for the new succeeding row are
set to 0. Quantities advected into the zeroth row or column, the ME+1 col-
umn or the NE+1 row are lost. The direction of the flow field restricts com-
putation to the four or fewer accumulation cells: the source cell, a diagonal
cell and two adjacent cells forming a square of side 2*DL. The advection sub-
routine code has only been modified to the extent required to handle the dual
problem. Much could be done to improve the speed if heavy use is anticipated.
Routine XMIT (Block I:A1.62,A1.65)
See Appendix B for description of routine XMIT.
Routines INPOUT, IOERR (Block J:A1.81,A .144)
Routine INPOUT for Phase IV is similar to the routine of the same name
used for Phases I and II from deck UTIL. Ultimately the same routine could be
used for all four phases with simple modifications applied for the special re-
quirements of each phase. The three major changes made to the Phase IV ver-
sion are: 1) when the velocity components are read from the intermediate data
set as type 2 or type 3 data records, the Ui values overlay the £i values
which are not required (to avoid future incompatibility this function was set
to occur when ITY<0; 2) when the data are transferred from and to internal ar-
rays the indexing is changed so internal array storage has the column index
varying faster; and 3) six arrays are processed in a single call to INPOUT us-
ing LU2 when type 20 is specified.
Routine IOERR has been modified to identify the logical unit in the error
message through an additional formal parameter in the call sequence.
LISTING FOR PHASE IV SAMPLE RUN
Following is a reproduction of a Phase IV run. The first page reflects
the UPDATE processing of the library according to the library directives set
122
-------
input on cards with listing option L=124 as shown in the SCOPE control card
example of TABLE 8. The FORTRAN compiler listing of Phase IV shows the UP-
DATE module identification and card number to the right. The run results are
presented immediately following the listing of routine THERMAL. The inpul;
master control cards are printed directly as input followed by a series of
arrays printed with empty lines deleted. The INPOUT lines appear as the type
20 restart fields are saved on LU2. The format for array print cannot print
negative numbers except as a row of asterisks. This phenomenon may occur in
the heat of evaporation field, indicating a net gain due to condensation at
the surface, and in the center of mass component fields, indicating negative
direction to center of mass.
123
-------
CURDS ENCOUNTERED IN INPUT
BKV 76/66 UPDATE 1.2/110-1H 09 MAR 77 01.15.35
PAGE 1
• • •*»
• COPIPILE AD.TH
MODIFICATIONS / DIRECTIVES
YANKSSS
YANK$$$
YftfJKSS*
YANK$S$
YANK$$$
YflNKSSS
YANKSS*
YANK$$»
YANK$$$
YANK***
YANK$$$
YP,NK$S$
YANK*$*
AD
AD
AD
AD
AD
AD
AD
AO
AD
TH
TH
TH
TH
»YANK
»YANK
• YANK
»YANK
»YANK
»YftNK
»YANK
»YAMK
• YBNK
»YANK
»YANK
»YANK
»YANK
• CALL
•CALL
• CALL
• CALL
•CALL
• CALL
•CALL
•CALL
• CALL
• CALL
• CALL
• CALL
• CALL
VARUINE
CONPRV
HNG7
VARCOHL
VARUIIND
VALDEZ
CONVCT.CONVCU
CONPRT,CONPRU
FLAT. PLATA
FLATIO.FLATIP
MONTE
HNJ
BKYPRT
A
B
ABC
THERPI
A
A
B
C
A
A
B
ABC
THERM
AVARU2
ACONV2
AHNG71
AVARC1
AVARU1
AVALDEI
ACONV1
ACONV1
AFLAT1
AFLAT1
AP10NT1
AHNIJ1
AHNIJ1
Al
Al
A2
A2
Al
Al
Al
Al
Al
T2
T2
T2
TZ
1
1
1
1
1
1
1
2
1
2
1
1
2
3
1
2
3
30
66
67
68
82
1
2
3
1
-------
PROGRAM MJVECT 7600-7600 OPT=2 FTN 1.6+139/010 09 MAR 77 01.15.36* BKY PAGE
1
5
C
c
c
PROGRAM ADVECTC INPUT, OUTPUT, TAPE1 ,TAPE2,TAPE3=OUTPUT )
CONVERTED TO CDC 7600 BY A .0 . STROUD.CSI - MARCH 1976
ROUTINE BASED ON J. HARDING PROGRAM FEDBACK
MODIFIED BY R. BAUERCOMPASS SYSTEMS, INC. FEB. 75
PROGRAMMED FOR CDC1601 UNDER CS1ROS
A2.1
A1.2
AD. 3
AD.1
AD. 5
COMMON /A/ ME.NE,ROT,DL,MAXDIM.MODEL,IDIM,NO,MINX,MINY.I'1AXX,MAXY,DTA.2
»,FRACTX,FRACTY IDIMH.NXPl A.3
10
15
20
COMMON CT< 300),FXT( 300 ),FYT( 300 ),RXTl 300),RYT( 300)
».T(300),VX(300).VY(300)
COMMON/1/1UK 960)
COMMON /ABC/TA( 300 ) ,EUI( 300 ),EA( 300 ), XK( 300>,YK< 300)
» ,0.6(300), OH( 300), QE(300)
COMMON/THERM/ IT, I TS , IPMOD ,TL ,CC .UU .K .C , EC rEK.K I. EK I. Q2T. DO. R
»,IPO
REAL K,K1
DIMENSION NAME(6) 1CH 9 ) r IC(1 ) , ISORCl 8r8 ) , SORCt 1 . B > .0.11) . V< 1 >
EO.UIVALENCE( ic 1( 2 ) , fc< 1 ) )
* (0,T),(V,CT)
DATft(NAME=1HNAME,1HGRID,1HTlME,1H$ORC,1HUlNO,1HCOMP),< IWND=0)
1,(MAXDIM = 300).( IDIM = 30),(MAXNAME=6 ) ( MAXSORC = 8 ),{ MAXWIND=1 )
2,( M INSORC = 99 999999 ) (MIMU I ND=9999 99 $9 ).(MINX = 1 ) (MINY = 1 ) ( 1SNO=0
COM. 2
COM. 3
Al .5
COM. 5
COM. 6
COM. 7
T3.28
THERM. 6
A2.1
A2.5
A1.8
A2.6
A2.7
A2.8
I A2.9
IC=5H( 1H1 »,(F = 1000. ),(TUO=10. ), ( TE = 1 .E-10 ) ,(CB= .956>7rRO=l .025 )A2.10
I— 25 »',(RPD=.01715329252) A2.ll
10 C UNITS C6S SYSTEM UIITH CALORIES( GRAM-CALORIES ) A2.12
Ln '
30
35
C
1
2
1000
1001
3
1
5
ROCB=HO»CB
PRINT 1C
REftD CONTROL CARDS
READ 1000, Id
PRINT1001,IC1 * DO 3
FORdftTC AMi7A10^A6 )
FORMftTC lXAt,7A10,A6 )
1F( ICH 1 ).EO.NAME( I )
CONTINUE * GO TO 1
GO TO (1,5,6,7,8,10)
DECODE(76,1001,IC)NE
IFCHE.LT. LDIM ) GO TO
1 = 1
) GO
,1
iLi
,MAXNAME
TO 1
DL,ROTANG,C3 * ROT=810 .-ROTANG
PRINT 900J $ GO TO 999
A2,
A2,
AD
Al.
AD.
Al,
Al .
AD.
AD
A2
A2.
A3,
.13
.11
.17
.10
.19
.11
.12
.22
.23
.15
.16
. 1
9005 FORPIAT(92H 1D1M MUST.GT.ME, COMMON/I/IU( 7»IDIM1 ),BDVECT COMMON/1/CA3. 2
10 »ACC(IDIM1 ),...,C(MAXO( 1,NO-7»IOIM1 ))) A3.3
10,01, FORFlftT( 213,7E10.3) AD_i_27_
6 DECOOEl76,1005,1C IMODOUT,1TO,ITS,ITE,ITINC,IPO,IPMOD AD.28
MODOUT=MAXO(ITINC,MOOOUT,!) $ IPMOD=MAXO(IPMOD,ITINC,1) A1.13
15 1005 FORMAT(16,7110) AD.31
7 ISNO=ISNO + 1 S IF(ISNO.GT.MAXSORC ) GO TO 999 AD.32
DECODE) 76,1007, 1C )( ISORC( 1 f ISN01,1 = 1,8 ),( 50RC( l.ISNO), 1 = 1,1 )
1007 FORMAT(1I3,212,2110,1E10.3) A2.18
MINSORC=M1NO(MINSORC,ISORC(7,ISNO)) S GO TO 1 Al.11
___J_1 F( 1UINO.LT.MAXUINO) GO TO 9 » PRINT 1008 $ GO TO 1 AJLJJL
1008 FORMAT!21HO»»»UNABLE TO PROCESS»»»/) A2.20
9 IUNO=IUNO+1 A2.21
D E CO DE( 76.1007.1C )( I|ilIND( 1, IUINO ), 1 = 1. 8 1, ( M1ND( I r 1UINO ), 1 = 1 ,H ) 6_2il2_
GO TO 1 ' A2.23
55 10 DECOOE(76,1005,1C) MODEL,IFIL,LAYER * NEM=NE-1 $ TEM=1.-TE A2.21
UlLftVER . Et.O) LAYER=1 A.L..15
^
DT=ITINC S MAXX=I"IE » MAXY=NE * NO=NE»ME * NXP1=ME+1 A1.16
-------
PROGRAM ADVECT 7600-7600 OPT=2 FTN 1.6+139/010 09 MAR T7 01.15.36* BKV PASE
to
cr>
60
65
70
75
80
65
90
95
100
10$
110
FRACTX=FRACTY=DT/DL 6 IDIni=IOIB»3
OZT=DT»DL»DL/< 21.*3600.»ROCB>
DO 19 J=1,ISNO * SORC( 1,J )=SORC< 1,J )»DT » SORC( 2, J > = SORC< 2, J )»DT
19 CONTINUE
TL=FRACTX 0 EKI=1./EK * KI=1./K
C3=C3»DT»10000. * MAXU1ND=IUNO » 1WNO=0 » ITU=ITS
CALL XniT .EQ .0 >
»CALL PRTMAX(VY,VXrlH .8HVELOC ITV , 8 . IT >
CALL THERMAL
DO 30 I = 1,NO $ P=RXT( I )»RYT< I > * IF(P.GT.TE) T< I > = ( T( 1 )-TUO >/P
C INTERPOLATE TO TEMP ANOMALY OF HEATED FLUID FROM CELL TEMPERATURE
30 CONTINUE
50 DO 60 1 = 1, NO $ IF(V(I).6T.TE) 0( I )=T(I )»V( I )»ROCB
C COMPUTE EXCESS HEAT IN FLUID BASED ON ANOMALOUS TEMPERATURE
60 CONTINUE 4 IF( IT .EQ . ITS ) GO TO 70
CALLXMITC 0,IU.NO, 1,1 )$CALLXMIT( V,0,NO,1,1 )*CALLXMIT< IU,V,NO,1,1 >
CALL BOVECT( 1 )
CALLXMIT< V, I U, NO, 1,1 HCALLXMITt Q,V,NO,1,1 )*CALLXMIT( IW,0,NO, 1 , 1 >
CALL BOVECTC 7)
70 IF( IT.LT.MINSOROGO TO 80 * DO 79 J = 1.ISNO
A3.1
A2.25
A2.26
A2.27
A2.28
A2.29
A2.30
A2.31
A2.32
AD. 39
A3. 5
A3. 6
A2.35
A2.36
A2.12
A3. 7
A3. 8
A2.16
A3. 9
A3. 10
A3. 11
A2.17
A2.19
A2.50
A2.52
A2.53
IF( ISORC(7,J ).GT.IT.OR.ISORC(8,J).LT.1T)GOT079 A2.51
LM=ISORC(5,J ) * LP = ISORC(6,J )$IF( LAYER .LT.LM .OR .LAYER .GT .LP )GOT079 AD .12
NM=ISORC( 1PJ ) * NP=ISORC( 2.J ) AD.H3
DO 78 N = NM,NP $ IP = (N-1)»ME * IM=IP + ISORC( 3, J ) $ IP = IP + 1SOHC( 1, J )
DO 77 l=ln,IP
C INSERT LOGIC FROM ROUTINE MERJE
AD. 11
AD .15
AD. 16
IF( ABSCFXT( I )).LE..5 .AND . ABS( FYT< I)).LE..5 ) GO TO 71 A1.18
PRINT 7000,IT,I,FXT( I),FYT( I) $ GO TO 999 AD.t8
7000 FORMATS POSSIBLE AOVECTION ERROR, TIME STEP TOO LONS* r2I8 ,2E10 .3)AO .t9
7t SU«=V( I ) + SORC(2,J )
0< I > = SUW>=0( I) + SORC( 1,J> $ IF(SUM.GE.TE) GO TO 75
CT(I) = SUM $ RXT(I)=RYT(I)=FXT(I)=FYT = SUM
FXT< I ) = FXT( I)»CONCSM $ FYT( I )=FYT( I )»CONCSM
RXT< n=l.-2.»ABS(FXT< I)) * RYT( J )=1 .-2 .»ABS( FYT( I ) )
77 CONTINUE
78 CONTINUE
79 CONTINUE
1F( IT.GE.IPO.AND.MOD( IT, IPMOD ) .EQ .0 )
» CALL PRTMAX(Q,tHNONE,lH ,20HADVECTED EXCESS HEAT, 20, IT)
80 DO 89 1 = 1, NO * IF(V(I).6T.TE)60 TO 88 * T(I)=0. « 60 TO 69
88 T( I ) = fl( I )/( VI I )»ROCB)
89 CONTINUE
90 IF< IT.NE.ITU) GO TO 110 6 IF( IUNO .EQ .0 ) GO TO 91
IF( IT.EO. IUIIND $ ITU=0
100 IUNO=IUNO+1 * IF< IUNO.GT.MAXU1ND) GO TO 110
IF( IT.GE.IUIND(6,IWNO» GO TO 100
ITW=IWIND(7,1UNO) $ IFflTU.GT.lT) GO TO 110
120 YW=(-ROTANG-90.-WIND< l.IUNO))»RPD
A2.56
A2.58
_AD.51
A2.59
A1.21
A1.22
AD. 56
AD. 57
AD. 58
A3. 12
A3. 13
A2.60
A2.62
A2.63
A2.61
A2.65
A2.66
A2.67
A2.68
A2.69
A2.70
TK=UIND(2,1UNO> $ TK=TK*TK*C3 « XU=COS(YU )*TK $ YU=SIN(YU )*TK A2.71
-------
PROGRAM ADVECT 7600->7600 OPT=2 FTN 1.6+139/010 09 P1AR 77 01.15.34* BKY PAGE
US Ntt=IU!ND( 1,IUNO> * NP = IUIND(2,1MNO>
DO 130 N = NClrNP $ IP = (N-1)»P1E $ ln=IP+IUIND< 3. IUNO >
1P = IP + IUIIND(1,1UNO) * DO 130 I=in,lP 6 XK< I >=XU * YKC I )=VU
130 CONTINUE
1F( TUNO.GE.PIAXUIND) GO TO 135
120 1F(IT.6E.IWIND(T,IUNO+1» GO TO 100
135 ITW=IWIND( 8,IUINO)
110 CONTINUE
A2.72
A2.73
A2.71
A2.75
A2.76
A2.77
A2.78
A2.79
1F(IT.GE.lTO.ftND.nODCIT.MODOUT).EQ.0)CALL 1NPOUT(V,LAYER.20,1T.F ) A2.80
1P(IT.LT.IPO.OR.mODtIT,lPnOD l.NE.O) GO TO 800 A2.81
125 CALL PRTP1AX(T.1HNONE.1H , 28HADVF.CTED TEMPERATURE ANOHALY ,28, IT ) A2._82
CALL PRTHAXCCT 1HNONE.1H ,11HCONCENTRATIONS.11,IT) AD.65
800 CONTINUE AD.66
999 CONTINUE » END A1.27
CARD NR. SEVERITY DETAILS DIAGNOSIS OF PROBLEM
67ICONTROL VARIABLE IN COMMON OR EQU1VALENCED, OPTInlZATION MAY BE INHIBITED.
-------
SUBROUTINE PRTC1AX 7600-7600 OPT = 2 FTN 1.6+439/010 09 MAR 77 01.15.36t BKV PAGE
1 SUBROUTINE PRTMAX! 5 0 , I S ,NAflEl ,NC , IT )
COMMON /A/ ME,NE,ROT,DL,flAXDIM,MODEL,IDIM(NO,fHNX,MlNY,rHAXX,l'lAXY
•,FRftCTX,FRACTY101M4,NXPl
^ »,NEH.TEMfTE,TWO
COMMON/I /IK . IWUI,JK,I,A,B,NOP,nF,ML,MMF,NN,NW,IU< 948)
DIMENSION S< I >,01 1 )
OATA
NW = !NC-l 1/NCPU+l
10 NOP=(flE-l )/15 + l $ MF=-14
CAJ-L XMIT! NAME1.1U NU^J^l ) SCALL XMIT11H , I W( NW+ 1 ) JIAXNAME-NU^ 0
DO 100 IK=1,NOP * r\F=flF + 15 $ ML = MINO( HF + 14 ,ME )
1F(NE.LT.25 .AND.mODC 1K,2).EQ.O.ANO.IS.EQ.1H ) GO TO 30
20 WRITE! LU, 1000 > < IU< I ) , I = 1 , MAXNAP1E ) . LEV . IT
15 C DEBUG
1000 FORNftT! 1H 20X6A10, 1 3, 7X5HTIME=1 10 )
30 URITE( LU . 1001 K I . I = 01F.nL )
1001 FOfl(1AT( //6X15I8 )
IF(D( 1 1.NE.4HNONE) GO TO 60
20 DO 51 J^l^ME * KJ = (J-»)»nE $ KK = KJ + ML S KJrKJ+MF
C DEBUG
DO 49 I = KJ,KK S IF( S( I > .NE .0 . > GO TO 50
49 CONTINUE S GO TO 51
50 WRITE! LU, 2000) I S,J,( S( JK ),JK=KJ,KK)
25 51 CONTINUE
C DEBUG
2000 FORMAT! Al , 11 , 3X 15E8 . 3)
GO TO 100
60 00 90 J=J_,NE * NN=10 $ 1J = (J-1)»ME * KK = 1J+ML S KJ=IJ+MF
30 DO 80 JK=KJ,KK * NN=NN+1
P1M=JX-1 * IF(Mn.LE.IJ) MMiJK * NM=JK-nF. S IF
-------
SUBROUTINE BDVECT T60(HT600 OPT=2 FTN t.6*t39/OtO 09 OAR 77 01.15.36* BKY PftSE
1 SUBROUTINE BDVECTt 10)
C MODIFIED TO SINGLE DIMENSION BV B. BAUER COMPASS SYSTEMS, INC.
COMMON/ A/ ME, NE. ROT, DL,MAXDIM, MODEL, IDIN, NO, KlINX,niNY,MAXX, WAXY
»,FRACTX,FRftCTY.IDlMt,NXPl
5 » NEM.TEfl.TE.TUO
COPIMON CT(300),FXT(300),FYT( 300),RXT( 300).RYT( 300)
• T(300),VX(300) VY( 300)
COitMON/l/CACC( 90) FXACC( 90).FYACC( 90).RXACC( 90>,RYACC( 90).
» F2XACC( 90),F2YACC( 90),C(1)
10 CALL XniTC .0,CACC,10.»IDIMt,0,l>
C MAIN CALCULATIONAL LOOPS FOLLOW. ENDING AT
C STATEMENTS LABELED 500 AND 800
DO 800 K = 1,NE S J = 2 * IM=( K-l )»ME * DO 500 1 = 1, ME $ IK = I + Ifl
IF(CT( IKKLT.TE) GO TO 500
15 IKM=MAXO( IK-1.1M+1)
SIG!»FRACTX
IKP1=C1AXO( IK-NErI )
SIGPiAY=(( .5+FYT( IK ))»VY( IK >+( .5-FYT( IK ) )»VY( IKM))»FRACTY
SIGinAY=-SIGnAY
20 C PLUS V DIRECTION IS NEGATIVE OF Y INDEX DIRECTION
SIGNX = SIGNY = 1.
IF( SIGP1AX.NE.O. ) SIGNX = SIGNC l.,SIGMAX)
IF( SIGMAY.NE.O. ) SIGNY = SIGN( 1. .SIGMftY)
11 = I + CSIGNX + 0.1) $ IF( SIGMAX.EO.O. ) 11 = I
25 Jl = J + (SIGNY + 0.1) * IF( SIGMAY.EQ.O. ) Jl = J
IF(Il.LT.l) I1=NXP1
PIJX = ( SIGNX*2.»(FXT( I K J + SIGMAX )+RXT( I K )-l . 0 )/( 2.»RXT( IK ) )
PIJY = (SIGNY»2.»(FYT( I K J + SIGdAY )+RYT( I K )-l . 0 )/( 2.»RYT( IK ) )
C NORflAL CONDITIONS
30 C IF PX,PY<0 OR PX,PY>1 MODIFICATIONS ARE MADE
NCOUNT=0 $ IF(PIJX.GE.TE.AND.PIJX.LE.TEM) GO TO 150
NCOUNT=1 S IF(PIJX.LT.TE) PIJX=0.
IFCPIJX.GT.TEM) PIJX=1.
150 IF(PIJY.GE.TE.ftND.PIJY.LE.TEM) GO TO 160
35 NCOUNT=1 S 1F( PI JY .LT.TE ) PIJY=0.
IF( PIJY .GT.TEC1) PIJY = 1.
160 C1IJT1=CT(I K )»PIJX»( 1.0-PIJY )
C2IJT1 = CT(I K )»( 1.0-PUX )»PIJY
C3IJT1 = CT( I K)»P1JX»PIJY
tO C1IJT1=CT(I K)»( 1.0-PIJX)»( l.-PIJY)
ij=< j-i )*JD:PI+I
iu=( j-i )*iDin+n » ui=( Ji-i )»IOIM+I » lui*< ji-i)»ioin+Ji
CACCCIl JJsCACCUl J)+C1IJT1
CACC( I Jl ) = CACC( I Jl 1+C2IJT1
15 CACCU1 Jl)=Cft'CC(Il JD+C3IJT1
CACC( I J )=CACC( I J )+C1IJTl
1F( IQ.EO.l ) GO TO 500
F1XT1 = (PIJX » RXT( I K) - 1.0) » 0.5 » SIGNX
50 F1YT1 = (1.0 - RYT< I K) * (1.0 - PIJY)) » 0.5 » SIGNY
F2XT1 = (1.0 - RXT( I K) » (1.0 - PIJX)) * 0.5 » SIGNX
F2YT1 = (PIJY » RYT( I K) - 1.0) » 0.5 » SIGNY
FXACCdl J) = FXACC(I1 J )+FlXTl»C 11 JT1
FYACCdl J) = FYACC(I1 J )+F 1YT1*C 1 1 JT1
55 FXACCd J1) = FXACC(I Jl )+F2XTl*C2IJTl
FYACCd J1)=FYACC(I Jl) + F2YT1»C2IJT1
FXACCdl Jl) = FXACCdl Jl) + F1XT1»C3IJT1
A2.97
MARAD.101
,DTA.2
A. 3
con.i
con. 2
COM. 3
A3. ^6
A3. 17
A2.98
AD. 115
A3.lt
A3. 15
A2.100
A3. 16
A2.102
A2.103
A2.101
A3. IT
A3. 18
AD. 121
A1.70
A1.71
AD. 121
AD. 125
A2.105
AD. 127
AD. 128
AD. 129
AD. 130
A2.106
A2.107
A2.108
A2.109
A2.110
A2. Ill
AD. 137
AD. 138
AD. 139
AD.ltO
6fij.lt!
AD. 142
AD. 113
AD.ltt
AD. It5
AD.lt6
A2.112
A2.113
AD.lt7
AD. It8
AD. It9
AD. 150
AD. 151
AD.152
AD. 153
AD.lSt
AD. 155
-------
SUBROUTINE BDVECT 7600-7600 OPT = 2 FTN 4.6+139/010 09 HAR T7 01.15.36* BUY PAGE
CO
O
60
65
TO
75
80
8$
90
95
100
105
110
FYACCdl Jl) = FYACCU1 J 1 > + F2YT1*C3IJT1
FXACCd J) = FXACC( I J) + F2XT1»C11JT1
FYACCd J) = FYACCd J) + F1YT1»C1IJT1
F2XACC(I1 J) = F2XACCdl J) + F1XT1»F1XT1»C1I JT1
FZYACCdl J) = FZYACCdl J) + F1YT1*F1YT1»C11 JT1
F2XACC( I Jl) = FZXACCd Jl) + F2XT1»F2XT1»C2I JT1
F2YACC( 1 Jl) = F2YACC(I Jl) + F2YT1»F2YT1»C21 JT1
F2XACC(11 Jl) = FZXACCdl Jl) + F 1 XT1 »F 1 XT1»C3I JT1
F2YACCdl Jl) = F2YACCdl Jl) + FZYT 1»F2YT1»C3I JT1
FZXACCd J) = F2XACCd J )+ F2XT1»F2XT1»C1IJT1
F2YACCCI J) = FZYACCd J )+F 1 YTl»F 1 YT 1 »C1 1 JT1
R1XT1 = (P1JX » RXT( 1 K)> * (PIJX » RXT( 1 K)>
R1YT1 = ((1.0 - P1JY) » RYTC I K)) • (( 1.0 - PIJY) » RYT( 1 K >>
R2XT1 = (d.O - PIJX) » RXTd K)) » (( 1.0 - PIJX) » RXT( I K))
R2YT1 = (PIJY » RYT( I K)) • (PIJY » RYT( I K»
RXACCdl J) = RXACCdl J )+ R1XT1» C1IJT1
RYACC(I1 J) = RYACCdl J ) + R1YT1» C1IJT1
RXACCd Jl)=RXACCd Jl) + R2XT1* C2IJT1
HYACC( I J1) = RYACC(I Jl)+ R2YT1» C2IJT1
RXACCdl J1)=RXACC(I1 J 1 ) + RlXTl»C3I JT1
RYACCdl J1) = RYACC(I1 J 1 )»R2YT1»C31 JT1
RXACCd J) - RXACCd J )+ R2XT1»C1UT1
KYACCd J) = RYACC(I J1+R1YT1 *C4IJT1
1F( NCOUNT.EQ.O ) GO TO 500
SIGP1X = SIGMAX+FXT( I K)
SIGMY = SIGP1AY + FYT( I K )
IF(PIJX.GT. 0.0001 ) GO TO 190
FXACCdl J) = FXACCdl J) - C1IJT1» F1XT1
FZXACCdl J)=F2XACCdl J)-CUJT1» F1XT1»F1XT1
FXACCd Jl) = FXACCd Jl) - C2I JT1»( F2XT1-S1GMX )
FZXACCd Jl)=F2XACCd J 1 )-C2I JT1»( FZXTl'F2XTl-SIGnX»SIGBX )
FXACC(I1 J1) = FXACC(I1 JU-C3IJT1* F1XT1
FZXACCdl Jl)=F2XACCdl J1)-C3IJT1» F1XT1»F1XT1
FXACCI I J )= FXACCl I J) - C1IJT1»< FZXT1 - S1GP1X)
F2XACC( I J) i FZXACCd J )-C1I JTl»( FZXT1«FZXT1-SIGMX»SISHX )
GO TO 192
190 IF(PIJX.LT. 0.9999) GO TO 192
FXACCdl J) = FXACCdl J )- CM JT1»< F1XT1-SIGMX +1.»SIGNX)
FZXACCdl J) = F2XACC(I1 J )-C 11 JTl»( FlXTl»F IXT1-
1 (SIGP1X - SIGNX) » (SIGP1X - 5IGNX))
FXACC( I Jl) = FXACC( I Jl) - CZ1JT1» FZXT1
FZXACCd Jl)-FZXACCd J1)-CZ1JT1» FZXT1»F2XT1
FXACC(I1 Jl) = FXACCdl Jl) - C3I JTl«( F 1XT1 - SIGMX +1.»SIGNX)
F2XACCdl Jl) =F2XACCdl J 1 )-C3I JT1»( F 1XT1»F 1XT1-
1 (5IGP1X - SIGNX) » (SIGMX - S1GNX»
FXACCd J) = FXACCd J) - C1IJT1» F2XT1
FZXACCC I J)=F2XACC
-------
SUBROUTINE BDVECT 7600-»7600 OPT=2 FTN 1.6+139/010 09 CIAR 77 01.15.36* BKY PAGE
115 191 IFIPIJY.LT. 0.9999) GO TO 500
FYACCU1 J) = FYACCCI1 J) - C1IJT1» F1YT1
F2YACCU1 J)=F2YACC<11 J )-C 1 1 JT1»( F1YT1*F 1YT1 )
FYACCU Jl) = FYACCC1 Jl) - C2I JT1*( F2YT1 - SIGHY +1.»SIGNY)
F2YACCCI J1)=F2YACC(I J 1 )-C21 JT1»< F2YT1»F2YT1-
120 1 (SIGPW - SIGNY) » ( SIGMY - SIGNY))
FYACCCI1 J1)=FYACC(I1 J 1 >-C31 JTl»r F2YT1-S1GRY +1.»SIGNY)
FZYftCCdl Jl )=F2YACC( 11 J 1 )-C3I JT l*( F2YT1»F2YT1-
1 (SIGfflY - SIGNY) » (SIGNY - SIGNY))
FYACCd J) = FYACCCI J) - CtIJTl* F1YT1
125 F2YftCC(l J)=F2YACC(I J)-CH1JT1» F1YT1»F1YT1
500 CONTINUE
IF(K.LE.l) GO TO 610 S IL1=( K-2 >»nE S 1J1=0
510 00 600 1=1. WE * 1L=1L1+I $ IJ=IJ1+I
C VALUES AT NEXT TIME STEP ARE CALCULATED
130 CT( I L) = CACCl 1 J)
IFC IQ.EO.l ) GO TO 600
CTREPL = CACCCI J )
IFC CTREPL.LT.l .E-10) GO TO 599 $ CTREPL=1 . /CTREPL
135 FX=FXT(I L)=FXACC(I J )»CTREPL
FY=FYT( I L)=FYACC(I J )»CTREPL
BXnAX = ABS(2.»FX-l.) $ RXM1N = ABS( 2 .»FX + 1 . )
IF( RXM1N.GT.RXMAX ) RXP1IN=RXMAX
RYP1AX = ABS( 2.»FY-1. ) S RYP1IN=ABS( 2.*FY + 1 . )
1MO IF(RYPlIN.GT.RYmAX) RYMIN = RYP1AX
RXT( IL ) = SQRT( RXACC( IJ )*CTREPL+12 .»( -FX*FX+F2XACC( IJ )»CTHEPL ))
RYT( 1L)=SQRT(RYACC( IJ )»CTREPL+12.»( -FY»FY+F2YACC( IJ )»CTREPL) )
IF(RXT(I D.GT.RXMIN) RXT( I L)=RXP1IN
IF(RYT(I D.GT.RYP1IN) RYT( I L )=RYP1IN » GO TO 600
115 599 FXT( IL )=FYT( IL )=RXT( IL)=RYT( IL)=0.
600 CONTINUE * IFC ILl-NO+C!E+nE ) 610.605.800
605 IL1=(NE-1 )*nE * IJ1 = 1DII1 $ GO TO 510
610 DO 650 1=1,10 * IL=( 1-1 )»IDint+IDIM+l
CALL XhlT( CACCC IL).CACC( IL-IDIM ) IDIP1»2 . 1 . 1 )
150 650 CALL XMITC 0 ., CACCC IL+IDIFI ),HE+1, 0, 1 )
800 CONTINUE $ RETURN $ END
AD. 213
ftD.21t
AD. 215
AD. 216
AD. 217
AD. 218
AD. 219
AD. 220
AD. 221
AD. 222
AD. 223
AD.22t
A3. 19
A3. 20
AD. 227
AD. 228
A2.11M
A2.115
AD. 229
A1.72
A1.73
A1.71
A1.75
AD. 237
A1.76
AD.2tO
A2.116
A2.117
AD. 215
AD. 216
AD. 217
A3. 21
A3. 22
A3. 23
A3. 21
A3. 25
A3. 26
-------
SUBROUTINE INPOUT 7600-«7600 OPT = 2
FTN 1.6+139/010
09 MAR 77 01.15.36* BKY PAGE
U)
NJ
1 SUBROUTINE INPOUT( A ,L AYER , ITY , I T , F ACOUT ) A1.81
COMI10N/A/ dE.NE,ROT.DLJ'1AXDIM.MODEL.IDIM,NOJMINX,f'lINY.nAXX.r)AXY,DTA.2
»,FRACTX,FRACTY,lDIfl4,NXPl
• NEM,TEtn,TE,TWO
5 COnMOD/1/ IUK 960)
DIMENSION ft( 1 ),IREC( 2)
DATAMREC = 0,0),( 024 = 7777 7 777B >,< ni 2=77778 >,( LEN=5 00 >,( J = 0 >
PRINT 1000, IT, LAYER, 1 T Y , NE ,ME .MODEL , I REC
10 1000 FORMAT! 10H INPOUT 10110)
FAC=FfiCOUT
1D1 = SHIFT( IT, 36 ). OR. SHIFTS 1ABS( LAYER ), 12 )
ID2=SH1FT( IABS( ITY ) ,48 ) . OR . SHIFT( NE . 36 ) .OR . SHIFT! CIE, 21 ) .OR .MODEL
KP=3 $ 1F( 1TV.E0.20) KP=6 * LU=1 S IF(KP.E0.6) LU=2
15 DO 600 K=1,KP 4 KK=( K-l >»MAXDIM+1
IFCK.EO. I .AND. ITY .E0.20 ) FACOUT=1 . E-9
C VOLUME RANGE ALLOWED- 1.E6 TO 8.E13
IF(LAYER.GT.O) GO TO 80 $ FACIN=1 , /FACOUT
NU=n21 S IFfK.NE.l . AND . ITY .LT . 0 ) KK=KK-MAXD1M t GO TO 90
20 80 NUI=2 a IU(1)=ID1 s IW(2) = 1D2
90 DO 500 1 = 1 NO * NX = (I-1)/NE $ NY=I-NX»NE $ KK 1 1 =( NY-1 )»nE+NX+KK
IF(LAYER.GT.O) GO TO 100 S IF( NU .LT.LEN ) GO TO 150
1F(NU.EO.LEN.AND.J .NE.l ) GO TO 150
C INPUT
25 100 BUFFER IN( LU , 1 K 1W< 1 ) IW( LEN ) ) $ !REC< LU ) = IREC( LU )»1 * NU=2
1FHINIT< LU ) 1110. 110,120
110 IERH=1 * GO TO 130
120 IERR=2
130 PRINT 1000, IU< 1 >.IUI( 2), 101.102 * CALL IOERR( 1ERR, IREC ,LU HGOT0100
30 110 IF( IU( 1 KNE.IDl.OR. IU( 21.NE.ID2) GO TO 100
150 J=MOD< 1,5 >+l S GO TO ( 250, 210, 220, 230, 240 ), J
210 NU=NU+1 « IWU=IU1(NU) * N=SHIFT( IUU.-36 ) $ GO TO 300
220 N=SHIFT< 1WU,-12) * GO TO 260
230 N=IUUI.ANO.M12 $ NU=NU + 1 $ IUU=IU(NU)
35 N=SHIFT( N, 12) .OR.SHIFTC IUU. 12 ) . AND .H12 * GO TO 270
210 N-SHIFTl IUW,-21 ) * GO TO 260
250 N=IWW
260 N=N.AND.M24
270 1F(N.GT.1000 OOOOB ) N=N.OR ..NOT. M21
10 300 A( KKII ) = FLOAT(N)»FAC1N » GO TO 500
C OUTPUT
100 CONTINUE
N=A( KKII )*FACOUT $ N=N.AND.H21
J = F10D( 1.5 ) + l $ GO TO ( 150,110,120,130,110)rJ
45 410 NU=NW+1 $ IW(NU)=SHIFT(N,36 ) $ GO TO 170
420 IU(NU) = IU(NU).OR.SHIFT
-------
SUBROUTINE INPOUT 7400-7600 OPT=2 FTN 1.6+139/010 09 MAR 77 01.15.36* BKV PAGE
500 CONTINUE A1.132
&00 FACOUT=FAC » RETURN » END A3.38
SUBROUTINE IOERR 7600-»7600 OPT=2 FTN 1. 6+139/010 09 HAH 77 01.15.36* BKY PAGE
1 SUBROUTINE IOERR( L . IREC,LU ) S DIHENSION LAB( t). IREC.( 1 ) A3. 39
DATA(lftB=10HEOF ERR IN,(OHPAR ERR INr10HEOT ERROUT,10HPAR ERROUT) A1.137
»,(K=0) A3.10
PRINT 1000,LAB(L).IRECtLUKLU * IF(L-2) H.3.2 A3.11
_| 1000 FOBHIBTC All.2110) A3.12
2 BACKSPACE LU 4 ENDFILE LU $ BACKSPACE LU A3.13
1F(L.E0.3> GO TO 1 « ENDFILE LU 4 BACKSPACE LU A3.11
3 1BEC(LU)=IREC(LU)-1 ? K=K-H $ IF(K.LT.ZO) RETURN A3.15
1 STOP A1.113
10 END A1.111
OJ
OJ
-------
SUBROUTINE THERMAL 7600-7600 OPT=2 FTN t .4+139/010 09 MAR T7 01.15.36* BKV PAGE
1 SUBROUTINE THERMAL
COMMON /A/ MErNE.ROTfDL.MAXDIn.MODEL.IDin,NO,nINX.mINY.nAXX.PIAXY.
»,FRACTX,FRACTY IDIM1.NXP1
».NEM.TEM,TE,TU6
5 COWr.ON CT( 300) FXT( 300),FYT( 300)rRXT( 300),RYT( 300)
».T( 300),VX( 300 i VY( 300)
COMMON/ ABC /TA( 360 >.EU( 300 ) ,EA( 300 ) XK( 300>,YK( 300)
» .OB( 300), OH( 300). QE< 300)
COMI«10N/THERM/IT,ITS,IPMOD,TL,CC,UU,K,C,EC,EK,KI,EKI,a2T,DD,R
10 »,IPO
REAL K,KI
COMP10N/1/IW(960>
DIMENSION TWm.Vm S EOUIVALENCE ( T,TW >,( CT, V )
DATA .(EK=.5 >,< C=.l >. IP=N»ME * IM=IP-ME+1
25 DO 100 I=IM,IP S MM=MAXOC I-l,in> * MP=MINO( 1+1 , IP )
NM=I-P1E S IF(N.EQ.l) NM=I * NP=I+ME * IF(N.EQ.NE) NP=I
IFCCTf I ) + CT(NP)+CT(MP ) + CT(MM)+CTCNM).LE.TE*5. > GO TO 99
C INTERPOLATE HEATED FLUID TEMP ANOMALY TO CELL TEMPERATURE
TUK I ):=T< 1 )»RXT( I )*RYT( I HTWO * TA( I >=TUO-2.
30 C AIR TEMP IS RESET AT EACH STEP TO AMBIENT TW -2C DEGREES
C TO SIMULATE WIND MIXING EFFECT
AK = 373. U/( TW( I ) + 273.16 )
AK1 = ( l.-AK )»7. 90298 S AK2=ALOG 10( AK )»5 .02808
AK3=1,-10.»»( ( 1 .-1 ./AK)*11.311) 6 AK1=10 .+»( -3 .19119»AK-1 . >-l .
35 EU( I )=10.»»( AKl + AK2+AK3»1.3816E-7+AK1»8.1328E-3)»10I3.216
IF(EA( I ).EO.O. ) EA< I)=.98»EW( I )
GO TO 100
99 EW( I ) = EA( I )=TA( I ) = 0.
100 CONTINUE
10 101 CONTINUE
CALL XMIT(EW(ME-1 ) , EW( ME ),NE,ME,ME )
TT=IT/3600. * TTT=EXP(-EK»TT) S TT=EXP( -K»TT )
T3=( l.-TTT)»EKI $ T2=( 1 .-TT )»K I
DO 301 N=1,NEM * IP=N»ME-1 S 1M=IP-ME*2
15 DO 300 I = IM,IP $ IF(TW( I ).LE.TUO+TE) 60 TO 300
J = I » IF(N.EO.l) J = J+ME * IF(I.EQ.IPI) J = J + 1
JJ = J-1 * IF( XK( J).LT.O. ) JJ = J + 1
JI=J+ME * IF(XK( J).LT.O. ) JI=J-ME
DIF =§( EW( JJ >-EW( J ) )»XK< J )+( EU( JI )-EU( J ) )»YK( J ) )»TL
50 EWEft=( EU( J )-EA( J ) )»TTT+( EC+DIF )»T3
EA( I )=EU( JJ-EWEA
GO TO 300
DIF=«TWt JJ >-TUI< J))»XK( J ) + (TW( JI )-TU< J ) >»YK< J ) )»TL
55 TWTA = (TW( J )-TA( J ) )»TT+( C+DIF )»T2
TA( I)=TW( J)-TWTA
T1.9
OTA. 2
A. 3
COM.l
COM. 2
COM. 3
COM. 5
COM. 6
COM. 7
T3.28
THERM. 6
T2.5
T3.1
T2.7
T2.8
T3.2
T3.3
T3.1
T3.5
TH.13T
TH.138
TH.139
TH.ltO
T2.10
T2.ll
T2.12
T3.6
T3.T
T3.8
T3.9
T3.10
T3.ll
Tl .It
T1.17
T1.18
T2.15
T3.12
T3.13
T2.U
T2.1T
T3.11
T2.18
T2.19
T2.20
T3.15
T2.22
T2.23
T2.2H
T2.25
T2.24
T2.27
T2.29
T3.16
T3.17
T3.18
300 CONTINUE T2.36
-------
SUBROUTINE THERMAL 7600-«7600 OPT=2 FTN 1.6+439/010 09 PIAR 77 01.15.36* BKY PAGE
60
C
C
C
C
65 C
C
70
C
75
C
C
80
301 CONTINUE
CALL XMIT(EA(NO-1 ),EA( NO ) .NE,-ME,-ME )
CALL XM1T(EA(NO-ME),EA( NO), ME ,-!,-! )
CALL XMIT(TA(NO-1 ),TA( NO ) ,NE , -PIE, -ME )
CALL XNIT( TA< NO-ME),TA< NO >f ME, -!,-!>
CALC
COMPUTES THE EFFECTIVE BACK RADIATION
COMPUTES THE CONVECTIVE TRANSFER OF SENSIBLE HEAT
COMPUTES THE TRANSFER OF LATENT HEAT BY EVAPORATION
DO 101 N=1,NE S 1P=N*ME $ IM=IP-ME+1
DO 400 I=IM,1P $ L=I
IF(V( I KGT.TE) GO TO 350
TA( n=BB( I) = QHCI)=OE( 1 )=TU( I ) = 0 . $ GO TO 400
UNITS OBfQHtQE»XK( I >+YK< I )»YK( I))».077+.26
QH(1) = 39.*SS*(TU(I )-TA( I))
OE(L)=53.9*SS»(EU( I )».98-EA< I))
TUTA=( QB( L ) + QH( L )+OE( L )-OSR )»Q2T/CT( 1 )
TUI( I )=TU( I )-TUTft
OUTPUT T(I) IS THE CELL TEMPERATUREFOR AFFECTED CELLS
ADJUSTED FOR SOLAR NET RADIATION (QSR) AND OTHER EFFECTS.
400 CONTINUE
T2.37
T2.38
T2.39
T2.40
T2.41
TH.278
TH.279
TH.280
TH.281
T2.42
T2.43
T3.19
T3.20
T3.21
T3.22
T2.48
T3.23
T2.50
T2.51
T3.24
T3.25
T3.26
T2.53
401 CONTINUE
IFC 1T.LT.IPO.OR.MODC IT, IPMOD ) .NE .0 ) GO TO 500 T3.27
CALL PRTMAX(EA,4HNONErlH .32HUIATER VAPOR PRESSURE OF AIR (MB). 32. T1.47
—
11.18
85 CALL PRTMAX(EU,4HNONE,1H ,4 1HEU-SATURATED UIATER V. P. AT SEA TEPIP T1.49
_ _ •(I-IB).41.IT) _ Tl .50
CALL PRTMAX(TA.4HNONE,1H ,15HAIR TEHPERATURE , 15 , IT ) T1.51
CALL PRTMAX(TU,4HNONE,1H , 20HTU-UATER TEMPERATURE, 20, IT ) T1.52
_ 1F( L.EO.l 1 GO TO 500 _ T2.56
90 CALL PRTMAX(QB,4HNONE,1H .47HQB-EFFECTIVE BACK RADIATION ADJ. FOR T1.54
»CLOUDINESS,47,IT) Tl .55
CALL PRTP1AX(OH.4HNONE.1H . 39HQH-CONVECTIVE TRANSFER BY SENSIBLE HET1.56
^TT,39,IT) T1.57
CALL PRTMAX(QE,4HNONE,1H ,4 1HQE-TRANSFER OF LATENT HEAT BY EVAPORAT2.57
2| _ »T10N.41 .IT) _ ; _ T1.59
500 CONTINUE T1.60
RETURN TH.323
END _ TH.321
SEVERITY DETAILS DIAGNOSIS OF PROBLEM
T~~THERE IS NO PATH TO THIS STATEMENT.
-------
LOAD MAP.
LINK - BKY 6000/7000 8.1
09 MAR 77 01.15.11
PAGE
CJ
CTl
FL REQUIRED TO LOAD 13300
FL REQUIRED TO RUN 35200
INITIAL TRANSFER TO ADVECT
13020
BLOCK ASSIGNMENTS.
BLOCK
/A/
l\l
/ABC/
/THERP1/
ADVECT
PHTMAX
XP11T
BDVECT
INPOUT
IOEBR
THERMAL
/STP.END/
/FCL.C./
/QLINLnT/
/OJOSSFL/
/QASAFL6/
/.QCOKP./
/Q8.IO./
FTN1LIB
ALOG
ATAN2
BACKSP=
BUF1N=
BUFIO=
BUFOUT=
CIO =
COC1IO =
DECODE=
ENCODED
ENDFIL=
ERFISS
EXP
FECC1SK =
FLTIN=
FLTOUT=
FPlTAP =
FORSVS=
/OOVERL?/
FORUTL=
GETFITs
GOTOER=
INCOPl:
INPC =
ADDRESS
100
»l5
2025
6565
6605
11131
15035
15055
16062
16515
16602
17121
17125
17116
17117
17151
17152
17153
17717
17717
20016
20116
20200
20261
20125
20502
20515
20576
20713
21057
21117
22635
22736
22777
23155
23172
21061
25015
25016
25061
25127
25113
25125,
LENGTH
25
1700
1510
20
5621
101
20
1005
133
65
622
1
21
1
2
1
1
211
0
77
100
62
61
111
55
13
61
115
111
70
1166
101
11
156
315
372
761
1
16
13
11
262
216
FILE
LGO
LGO
LGO
LGO
LGO
LGO
LGO
FTN1LIB
FTN1LIB
FTN1LIB
FTN1LIB
FTN1LIB
FTN1LIB
FTN1LIB
FTN1LIB
FTN1LIB
FTN1LIB
FTN1LIB
FTN1LIB
FTN1L1B
FTN1LIB
FTN1LIB
FTN1LIB
FTN1L1B
FTN1LIB
FTN1LTB
FTN1L1B
FTN1LIB
FTN1LIB
FTN1LIB
FTN1L1B
KOOER
25673
167 FTN1LIB
-------
LOAD HAP.
LINK - BUY 6000/TOOO 8.1
09 MAR TT 01.15 .11
PAGE
°L"i-» HJJH
BLOCK
KRAKER=
OUTCOd=
RDL=
R00=
SINCOS=
SORT
SYS=
SYS=AID
$VS=1ST
UNIT
UZEHO..
UTH =
UTL=
UTU=
XTQY5
jni-HLHi j.
ADDRESS
26362
27017
27223
27160
27517
27512
27651
27715
30013
30011
30050
30057
30111
30226
30230
30260
30307
30375
LENGTH
135
201
235
37
23
107
71
16
1
31
7
62
65
2
30
27
66
7
FILE
FTN1LIB
FTN1L1B
FTN1L1B
FTN1L1B
FTN1LIB
FTN^LIB
FTN1L1B
FTN1L1B
FTN1LIB _.__
FTN1LIB
FTN1LIB
FTN1LIB _ _ . . ...
FTN1L1B
FTN1L1B
FTN1LIB
FTN1L1B
FTM1LIB
FTN1L1B
30101
-------
00
00
NAIOE TEST BOX MODEL
GRID
LAYS
TIME
SORC
COUP
3
INPOUT
5 533333.3333 0. 3.2E-6 12. .003
1 TT777 .99 1.022
600 1200 1200 2100 600 1200 1200
3 3 3 311 1200 T77T7TT 1.76E8 1.30ET
99
ADVECTED EXCESS HEAT
12315
.0 .0 .286E+12.0 .0
1200 1 20 5 $ 99 0 0
ADVECTED TEMPERATURE ANOMALY
1 TinE= 1200
1 TIME= 1200
12315
3
3
INPOUT
INPOUT
INPOUT
I
2
3
1
5
2
3
1
5
2
3
t
$
3
1
.0 .0 .113E+02.0 .0
CONCENTRATIONS
12315
.0 .0 .258E + H.O .0
1800 -1 -2 5 5 99 0 6
1800 1 20 5 5 99 15 6
2100 -1-2 5 5 99 15 12
VELOCITY
12315
10,181 11;183 11;182 8;181
10,160 Ilf161 10:167 8:157 5; 92
10,138 10,112 10,115 8,112 1, 97
10,125 10,126 9,130 6,128 3,101
10,120 10:120 8:125 5:120 3, 101
WATER VAPOR PRESSURE OF AIR (MB)
«>
12315
.0 .0 .120E+02.120E+02.120E+02
.0 .120E+02.212E+02.121E+02.121E+02
.0 .120E+02.120E+02.119E+02.119E+02
« .120E+02.120E+02.119E+02.119E+02
EU-SATURATED WATER V. P. AT SEA TEI-1P ( PIB )
12315
.0 .0 .123E+02.123E+02.123E+02
.0 .123E+02.238E+02.132E+02.132E+02
.0 .123E+02.129E+02.123E+02.123E+02
.0 .0 .123E+02.123E+02.123E+02
AIR TEMPERATURE
1 2 3 1 $
.0 .0 .800E+OJ.800E+01.0
.0 .0 .800E+01 .800E+01.0
TU-UATER TEMPERATURE
1 TIHEE 1200
1 TIME= 2100
1 T1ME= 2100
1 TIME= 2100
1 TIME= 2100
1 TIME= 2100
-------
3
t
3
1
3
1
3
H
3
4
INPOUT
M
U)
3
4
3
4
.0
.0
.0
.0
.0
.0
.0
.0
.0
.0
.0
.0
.0
.0
.0
.0
1
.0
.0
1
.0
.0
1
.0
.0
1
.0
.0
2400
1
.0
.0
1
.0
.0
.203E+02.110E+02.0
.107E+02.934E+01 .0
QB-EFFECTIVE BftCK RADIATION ADJ. FOR CLOUDINESS
2315
.IMHE-t-OS. 157E + 03.0
.15TE+03.158E+03.0
QH-CONVECTIVE TRANSFER BY SENSIBLE HEAT
2315
.125E+03.31ME+02.0
.286E+02.213E+02.0
OE-TRANSFER OF LATENT HEAT BY EVAPORATION
2 3 1 5
.89TE+01 .220E+01.0
ftDVECTED EXCESS HEAT
2345
.6TTE+12.692E+11 .0
.856E + H .13TE + H .0
1 20 5 5 99 18 12
ADVECTED TEHPERATURE ANOdftLY
2315
.113E+02.108E+02.0
.109E+02.906E+01 .0
CONCENTRATIONS
2345
.613E+11.652E+10.0
.805E+10.154E+10.0
1 TIME= 2400
1 TIHE= 2400
1 TinE= 2400
1 T1HE= 2400
1 TinE= 2400
1 TIrtE= 2400
-------
TECHNICAL REPORT DATA .
(Please read Instructions on the reverse before completing)
i. REPORT NO.
EPA-600/3-78-074
'lENT'S
4. TITLE AND SUBTITLE
User Guide for the Enhanced Hydrodynamical-Numerical
Model
i. REPORT DATE
July 1978
6. PERFORMING ORGANIZATION CODE
7. AUTHOR(S)
A.D. Stroud and R.A. Bauer
8. PER
9. PERFORMING ORGANIZATION NAME AND ADDRESS
Compass Systems, Inc.
San Diego, California 92109
10. PROGRAM ELEMEI1
11. CONTRACT/GRANT NO.
68-03-2225
12. SPONSORING AGENCY NAME AND ADDRESS
Environmental Research Laboratory-Corvallis
Office of Research and Development
U.S. Environmental Protection Agency
Corvallis, Oregon
13 TYPE OF REPORT AND PERIOD COVERED
'Final 1975-1977
14. SPONSORING AGENCY CODE
EPA/600/02
15. SUPPLEMENTARY NOTES
This is the companion document to EPA-600/3-78-073.
16. ABSTRACT
This guide provides the documentation required for used of the Enhanced Hydrodynamical-
Numerical Model on operational problems. The enhanced model is a multilayer Hansen type
model extended to handle near-shore processes by including:
-Non-linear term extension to facilitate small-mesh studies of near-shore,
including river inflow dynamics;
-Layer disappearance extension to enable appropriate procedures in tidal
flat and marshy regions, as well as some down/upwelling cases;
-Thermal advection enhancement for treatment of thermal pollution cases
by method of moments coupled with heat budget procedures;
-Monte Carlo diffusion enhancement to deal with dispersion via statistical
methods and comparison to the method of moments experiments.
The guide includes a description of the model system, source code maintenance procedures
notes on implementation, functional descriptions of routines and option sets, listings
and example runs for a test data set.
17.
KEY WORDS AND DOCUMENT ANALYSIS
DESCRIPTORS
b.IDENTIFIERS/OPEN ENDEDTERMS
COSATI Field/Group
numerical methods marshy regions
thermal pollution operations research
monte carlo method
mathematical methods
near-shore
river inflow
tidal flat
08/C
12/A,B
Release to Public
19. SECURITY CLASS (This Report)"
Unclassified
20. SECURITY CLASS (This page)
Unclassified
21. NO. OF PAGES
152
22. PRICE
EPA Form 2220-1 (Rev. 4-77)
140
t, U. S. GOVERNf.ENT PRINTING OPFICE: 1978-797-714/233 REGION ,<,
------- |