<TEXT N="BAD7210">
<BODY>
<DIV1>
<P><PB REF="00000001.tif" SEQ="00000001" RES="600dpi" FMT="TIFF6.0" FTR="TPG" CNF="894" N="1">
THE USE OF COMPUTERS IN CIVIL
ENGINEERING EDUCATION
Harold J. Welch
Department of Civil Engineering
The University of Michigan
June 1, 1962
This material is distributed by The Ford Foundation Project
on the Use of Computers in Engineering Education at The
University of Michigan. This report appears in the library
edition of the Final Report of the Project and is also issued
as a separate booklet. Similar "Curriculum Reports" for
other engineering disciplines are also available on request.
</P>
<P><PB REF="00000002.tif" SEQ="00000002" RES="600dpi" FMT="TIFF6.0" FTR="TOC" CNF="892" N="2">
Table of Contents
Page
I. Introduction C3
II. The Civil Engineering Department of The University of
Michigan: Faculty and Curricula c4
III. Introduction of Civil Engineering Students to the Use
of the Computer Cll
IV. Specific Uses of Computers in Various Fields of Study C16
V. Conclusions C20
VI. Example Problems Prepared by Project Participants C21
76 Use of the Computer as a Teaching Machine for Plane
Analytic Geometry C23
77 Optimum Regulation of a Reservoir for Power Production C35
78 Economical Size of Canal-Feeding Hydro-Power Plant C41
79 Most Probable Number Method for Determining Bacterial
Populations C47
80 Iteration of Sum of Principle Stresses at Interior
Points of Plane Rings C56
81 Elastic Buckling Load for Columns of Non-Uniform
Cross Section C67
82 A MAD Program for Truss Analysis C82
83 Stiffness Factors for a Flat Slab Column C97
84 Vibration of Beams on Spring Supports C104
VII. References C113
-C2 -
</P>
<P><PB REF="00000003.tif" SEQ="00000003" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="857" N="3">
USE OF COMPUTERS IN CIVIL ENGINEERING EDUCATION
ABSTRACT
The Civil Engineering Department of The University of Michigan has taken an active part
in computing work at the University with a total of seven faculty members serving as participants in the Project on the Use of Computers in Engineering Education during the past two-year
period. In addition, twenty-six civil engineering professors from other universities have
participated in Project activities for periods varying from one week to a full semester.
Most undergraduate engineering students at the University are introduced to computer work
through a required introductory course taught by Computing Center and Mathematics Department
personnel. The Civil Engineering Department, however, has elected to introduce its students
to computers via a departmental course, the second of the three required surveying bourses.
Details of the organization of this course are included as part of this report.
Digital and analog computer work has been assigned in more than 25 departmental courses
during the past year, giving students an opportunity to gain practice in the application of
computers in the solution of their engineering problems. The Civil Engineering curriculum, the
areas where computer activity has been strongest, and a sampling of opinion of faculty members
as to the effectiveness of computer usage in engineering instruction are included.
A selected set of example problems prepared by faculy Project participants is also
included. These may be considered as a supplement to the seventy-five example engineering
problems, including several related to Civil Engineering subject areas, which have been published previously by the Project.
I. INTRODUCTION
Civil Engineering, originally named to distinguish it from military engineering, bhas
always covered a wide field of engineering practice. Civil engineers plan, design and supervise the construction of roads, railroads, harbors, buildings, tunnels, waterways, bridges,
dams, airfields, canals, water supply and sewerage systems, and the many other facilities
necessary for public works and industrial development. They plan the conservation, utilization,
and control of water resources. They operate in the field of surveying and mapping.
In recent years, digital and analog computers have become important tools in many of
these areas. With this in mind, the Civil Engineering Department of The University of Michigan
has eagerly cooperated with the Ford Foundation Project on the Use of Computers in Engineering
Education.
Since the Project was initiated, seven members of the Civil Engineering Department at
Michigan and twenty-three visiting Civil Engineering faculty members have participated in the
Project's program. In addition, other members of the faculty have become familiar with computer
-C3 -
</P>
<P><PB REF="00000004.tif" SEQ="00000004" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="886" N="4">
philosophy and techniques as a result of weekly luncheon discussions and special lectures
sponsored by the Project.
It has been found that problem solving by computer methods adds to the depth, scope, and
logical understanding of many problems. It is believed that civil engineers of the future
will profit from an introduction to computers as part of their undergraduate training.
II. THE CIVIL ENGINEERING DEPARTMENT OF THE UNIVERSITY OF MICHIGAN: FACULTY AND CURRICULA
The courses which comprise civil engineering may be divided into seven general areas:
Construction, Geodetic, Hydraulic, Sanitary, Soil-Mechanics and Foundation, Structural, and
Transportation Engineering. These areas are referred to as technical options and represent
fields of specialization for the civil engineer. In addition, there are a few basic courses
dealing with legal problems, professional conduct, ethics, the scientific method, etc., which
are common to all fields of civil engineering.
Each of the general areas is described below. Members of the faculty primarily concerned
with each field and specific courses are listed. Courses which are required for all civil
engineering students are indicated with plus signs; courses in which the presentation is related
to or supplemented by computers are indicated with asterisks. A more detailed description is
given for those courses involving computer applications. Course numbers are given with each
course to provide a convenient means of referring to them later. The number of semester credit
hours is given in parentheses after each course.
A. Courses Common to All Fields of Civil Engineering
Lawrence C. Maugh, Professor of Civil Engineering,and Acting Chairman of the Department
of Civil Engineering
+ CE400 Contracts, Specifications, Professional Conduct, and EngineeringLegal Relationships (2).
CE500 Fundamentals of Experimental Research (2).
CE501 Legal Aspects of Engineering (3).
CE990 Civil Engineering Research (arranged).
B. Construction Engineering
The methods and techniques of modern construction; fundamental principles of construction
applicable to all types of engineering structures; business and legal principles of contracting
as applied in the field of construction.
Glen L. Alt, Professor of Civil Engineering (Advisor)
Frank E. Legg, Jr., Associate Professor of Construction Materials
+ CE350 Concrete Mixtures (1).
CE531 Cost Analysis and Estimating (2).
CE532 Construction Methods and Equipment (3).
CE533 Estimating Practice (1).
+ Course required for all Civil Engineering students.
</P>
<P><PB REF="00000005.tif" SEQ="00000005" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="890" N="5">
Use of Computers in Civil Engineering Education
C. Geodetic Engineering
Theory and practical applications of surveying and surveying techniques; theoretical and
applied geodesy, figure of the earth, local and extended control surveys; precise measurements
and the adjustment of observations; design and execution of municipal surveys, topographic
mapping projects, boundary surveys and land subdivision; industrial applications of surveying
techniques; research in fields of instrumentation, computation, mapping and photogrammetry,
and the problems of land surveying and route location.
Ralph Moore Berry, Professor of Geodetic Engineering (Advisor)
Harold J. McFarlan, Associate Professor of Geodetic Engineering
Harold J. Welch, Assistant Professor of Geodetic Engineering
+* CE260 Basic Surveying (3).
Use, care, and adjustment of basic surveying instruments; horizontal angle
measurement, leveling, taping; circular curves, grades and vertical parabolic
curves, profiles, earthwork; computation, significant figures, errors, desk
calculators, U. S. Public Land System.
+* CE261 Surveying Computations (3).
Principles of horizontal control, unified system of geometric computations
based on rectangular coordinates and the "Intersection Solution"; logical
synthesis of computations for complicated geometrical engineering applications; introduction, applications, and use of electronic computers; principles
of surveying astronomy, use of ephemeris and star catalog; elements of
photogrammetry.
+* CE362 Advanced Surveying Measurements (4).
Precise observational methods for triangulation and control levels; use of
theodolite and geodetic level, base-line measurements; electronic distance
measurements; field assessment of errors and adjustments; astronomical
observations for precise azimuth determinations; application of precise
control to layout of engineering projects.
CE505 Boundary Surveys (3).
CE560 Photogrammetry (2).
CE561 Geodesy (3).
* CE562 Geodetic Field Methods (2).
Reconnaissance for geodetic triangulation; special observing methods for
first-order horizontal and vertical control; Laplace stations and deflection
of the vertical; actual observations, reduction and adjustment of results in
an actual field situation. Term paper or project report required of each student.
* CE563 Adjustment of Geodetic Measurements (2).
Theory of least squares, applications to the adjustment of geodetic observations; arrangement for solution of complex adjustments on electronic computers;
actual solution of selected problems.
* CE564 Special Problems in Advanced Surveying.
Special advanced work can be provided for those who have received credit in
Advanced Surveying Measurements (CE362).
* CE565 Municipal Surveying (2).
Control surveys, methods and adjustments for use in municipal mapping and
administration, surveys for streets, utilities, property lines, tax maps,
subdivision control and development.
* CE960 Geodesy and Surveying Research.
Assigned work in geodesy, or other special field in surveying of interest
to the student and approved by the professor of geodesy and surveying.
* Course involves computer applications.
-C5 -
</P>
<P><PB REF="00000006.tif" SEQ="00000006" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="889" N="6">
Hydraulic Engineering
T'Sr application of the fundamental principles of hydraulics and hydrology to the
deve.ir pment of water power, flood control, drainage, and improvement of rivers and harbors,
and other hydraulic structures. Laboratory facilities and instruction are offered for students who wish to engage in research work in hydrology and hydraulics that will lead to
advanced degrees.
Ernest F. Brater, Professor of Hydraulic Engineering (Advisor)
Victor L. Streeter, Professor of Hydraulics
- CE420 Hydrology (3).
+-* CE421 Hydraulics (2).
Hydrostatic stability; orifices and weirs; Venturi meters; cavitation; pump
characteristics; flow in pipes and fittings; unsteady-uniform flow; steadynonuniform flow. Lecture, laboratory, and computation.
* CE522 Hydraulic Transients (2).
Introduction to water power development, including selection of type of
turbine; storage and pondage; surge in pipe lines; water hammer analysis;
digital programming of unsteady flow situations.
CE523 Flow in Open Channels (3).
* CE524 Advanced Hydraulics (3).
Two-dimensional potential flow; the flow net; percolation and hydrostatic
uplift; side-channel spillways; boundary-layer; hydraulic similitude;
hydraulic models; stilling pools.
C E526 Hydraulic Engineering Design (3).
Design of hydraulic structures such as diversion dams, head gates,
control works, silt traps, siphon spillways, side-channel spillways,
earth canals, and other structures involving accelerated flow, backwater,
hydraulic jump, sedimentation, and erosion.
CE623 Applied Hydromechanics (2).
CE825 Seminar in Hydraulic Engineering.
CE920 Hydrological Research.
* CE921 Hydraulic Engineering Research.
Assigned work in hydraulic research; a wide range of matter and method
permissible.
E. Sanitary Engineering
The planning, construction, and operation of water works, sewerage and drainage systems,
water-purification plants, and works for the treatment and disposal of city sewage and industrial
wastes; improvement and regulation of natural waters for purposes of sanitation; air sanitation;
principles and standards of ventilation.
Jack A. Borchardt, Professor of Civil Engineering (Advisor)
Lloyd L. Kempe, Professor of Civil Engineering and of Chemical Engineering
Eugene A. Glysson, Associate Professor of Civil Engineering
+* CE480 Water Supply and Treatment (3).
Sources of water supply, quality and quantity requirements, design fundamentals of works for development, collection, purification, and distribution
of water.
+* CE481 Sewerage and Sewage Treatment(2).
Requirements of residential and municipal sewerage systems, procedures for
the design and construction of sewerage and sewage treatment works.
-C6 -
</P>
<P><PB REF="00000007.tif" SEQ="00000007" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="895" N="7">
Use of Computers in Civil Engineering Education
CE580 Microbiology (3).
CE581 Sanitary Chemistry (2-3).
CE582 Sanitary Engineering Design (3).
CE583 Water Purification and Treatment (3).
CE584 Waste Water Treatment and Disposal (3).
CE585 Municipal and Industrial Sanitation (3).
CE586 Industrial Waste Treatment (2).
* CE587 Industrial Bacteriology (3).
Lectures and demonstrations to illustrate the application of microbiological
principles and techniques to industrial processes.
CE680 Microbiology II (2).
CE681 Advanced Sanitary Chemistry (3).
* CE682 Advanced Sanitary Engineering Design (3).
Functional design of sanitary engineering structures and typical plant layouts;
drafting room and field studies; preparation of design reports.
CE685 Special Problems in Sanitary Engineering I.
CE686 Special Problems in Sanitary Engineering II.
CE780 Public Water Supply (3).
CE880 Sanitary Engineering Seminar (1).
CE980 Sanitary Engineering Research.
F. Soil Mechanics and Foundation Engineering
Soil mechanics; evaluation of soil properties and environmental conditions in foundations
of earth-supported structures; mass stability in excavations and subsurface construction. Soil
engineering; use of soil characteristics and properties and soil classification in design and
construction of highways, railways, airports, and other surface facilities.
William S. Housel, Professor of Civil Engineering (Advisor)
Robert 0. Goetz, Instructor in Civil Engineering
Ulrich W. Stoll, Instructor in Civil Engineering
+ CE445 Engineering Properties of Soil (3).
CE543 Soils in Highway Engineering (2).
CE544 Airport Design and Construction (3).
CE545 Foundation and Underground Construction (3).
CE546 Soil Mechanics Laboratory (1).
CE946 Soil Mechanics Research.
G. Structural Engineering
The theory, design, and construction of structures, such as bridges, buildings, dams,
retaining walls, and reservoirs, involving the use of steel, reinforced concrete, and lumber;
the testing and utilization of soils in foundations and subsurface construction.
Lawrence C. Maugh, Professor of Civil Engineering, and Acting Chairman of the
Department of Civil Engineering (Advisor)
Glen V. Berg, Professor of Civil Engineering
Bruce G. Johnston, Professor of Structural Engineering
Leo M. Legatski, Professor of Civil Engineering
Robert B. Harris, Associate Professor of Civil Engineering
Wadi S. Rumman, Assistant Professor of Civil Engineering
-C7 -
</P>
<P><PB REF="00000008.tif" SEQ="00000008" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="895" N="8">
+* CE312 Theory of Structures (3).
Calculations of reactions, shears, and bending moments in simple, restrained,
and continuous beams due to fixed and moving loads; simple trusses with fixed
and moving loads; determinate frames; columns; tension members; girders;
introduction to design.
+* CE313 Elementary Design of Structures (3).
Design and details of simple beams, girders, columns, and trusses. Computations, drawing work, and laboratory experiments.
+* CE415 Reinforced Concrete (3).
Properties of materials; stress analysis and design of reinforced concrete
structures; introduction to prestressed concrete and ultimate strength analysis. Lectures, problems, and laboratory.
CE310 Structural Drafting (2).
CE511 Timber Construction (1).
* CE512 Advanced Theory of Structures (3).
Stresses in subdivided panel trusses; principle of virtual displacements and
virtual work; energy theorems; graphical methods; analysis of statically
indeterminate trusses and frames.
* CE513 Design of Structures (3).
Design of reinforced concrete and steel structures. Computations and drawing.
* CE514 Rigid Frame Structures (3).
Analysis of rigid frames by methods of successive approximations and slope
deflections; special problems in the design of continuous frames.
CE515 Prestressed Reinforced Concrete (2).
* CE516 Advanced Design of Structures (3).
Functional design of buildings; selection and analysis of structural elements;
reinforced concrete flat slab design; digital computer applications. Lectures
and computation laboratory.
CE517 Bridge Engineering and Design (3).
* CE611 Structural Dynamics (3).
Structural vibrations. Transient and steady state response to dynamic forces.
Response beyond the elastic range. Response to nuclear explosions. Earthquake forces. Structural response to earthquake. The response spectrum.
Seismic building codes and their relation to structural dynamics.
CE612 Structural Members (3).
* CE613 Structural Plate Analysis (2).
Stress analysis of flat plates loaded either in their plane or in bending.
Numerical analysis. Applications to special problems in flat slab construction.
* CE614 Advanced Problems in Statically Indeterminate Structures (3).
Continuous truss bents; hinged and fixed arches; rings; frames with curved
members; flexible members including suspension bridges; frames with semirigid connections.
CE615 Analysis and Design of Folded Plates, Domes and Shells (3).
CE616 Plastic Analysis and Design of Frames (2).
CE617 Mechanical Methods of Stress Analysis (1).
* CE618 Computer Analysis of Structures (3).
The analysis of beams, frames, trusses, and arches by high speed digital
computers. The general method of influence coefficients; matrix methods,
algorithm development, flow charts, programming. Students will solve a
sequence of problems on the high speed computer in the University
Computing Center.
CE810 Structural Engineering Seminar (1).
CE910 Structural Engineering Research.
-C8 -
</P>
<P><PB REF="00000009.tif" SEQ="00000009" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="893" N="9">
Use of Computers in Civil Engineering Education
H. Transportation Engineering
Transportation Engineering is further subdivided into Highway, Railroad and Traffic
Engineering. Highway Engineering involves the location, design, construction and maintenance
of various types of roads and streets, including materials, surveys, plans, specifications,
economics, and financing. Railroad Engineering involves the design, construction, and
operation of railroad properties, including metropolitan terminals, statistical analysis of
operating data, freight and passenger traffic, economics, financing, administration, and
regulation. Traffic Engineering involves methods of increasing the efficiency and safety of
traffic movement; street and off-street parking; traffic surveys, geometrical design of urban
and rural highways, traffic control devices, and other means of regulating and controlling the
use of highways.
John C. Kohl, Professor of Civil Engineering and Director of the Transportation
Institute (Advisor for Traffic Engineering)
Ward K. Parr, Associate Professor of Highway Engineering
Bruce D. Greenshields, Lecturer in Transportation Engineering
Donald N. Cortright, Associate Professor of Civil Engineering (Advisor for
Highway Engineering)
Clinton L. Heimbach, Lecturer in Civil Engineering (Advisor for Railroad Engineering)
+ CE370 Transportation Engineering (3).
CE550 Highway Materials (3).
CE551 Physical Properties of Concrete Masonry (2).
CE552 Bituminous Materials and Pavements (2).
CE570 Highway Traffic (2).
CE571 Traffic Engineering (3).
CE572 Highway Economics (2).
CE573 Highway Design (3).
CE574 Railroad Engineering (3).
* CE575 Terminal Design (3).
Design of railroad, highway, waterways, and airport terminals, joint terminals,
layout of the various types of yards and traffic facilities.
CE576 Economics of Railroad Construction and Operation (2).
CE652 Advanced Bituminous Materials and Flexible Pavement Design.
* CE670 Transportation Planning (3).
Analysis of supply and demand for transportation services, transport relationships to land use and other elements of regional and urban planning, and
planning techniques applied to transportation problems.
CE671 Advanced Highway Engineering (2).
CE672 Transportation (3).
CE673 Highway Transport (2).
CE674 Industrial Transport Management (4).
CE970 Highway Engineering Research.
CE971 Transportation Engineering Research.
-C9 -
</P>
<P><PB REF="00000010.tif" SEQ="00000010" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="871" N="10">
A list of courses taken by Civil Engineering students is given in the table below.
The table is arranged to indicate the semester in which the course is normally taken. Courses
in which the presentation is related to or supplemented by computers are indicated by an
asterisk (*).
LIST OF COURSES TAKEN BY CIVIL ENGINEERING STUDENTS
Total Hours
Semester Course Name and Number Hours for Semester
Analytic Geometry and Calculus I (Math 233) 4
General and Inorganic Chemistry (Chem 104) 4
1 Engineering Drawing (Engr. Graphics 101) 3 15
Freshman Composition I (English 111) 3
Introductory Speech (English 121) 1
Analytic Geometry and Calculus II (Math 234) 4
Descriptive Geometry (Engr. Graphics 102) 3
2 Freshman Composition II (English 112) 2 15
English Elective 2
General and Inorganic Chemistry (Chem 106) 4
Analytic Geometry and Calculus III (Math 371) 4
Mechanics, Sound, and Heat (Physics 145)5
3 Statics (Engr. Mech. 208) 3 16
Basic Surveying (CE260) 3
Concrete Mixtures (CE350) 1
Analytic Geometry and Calculus IV (Math 372) 4
Electricity and Light (Physics 146) 5
4 Mechanics of Material (Engr. Mech. 210) 4 17
Laboratory in Mechanics of Material (Engr. Mech. 212) 1
* Surveying Computations (CE261) 3
Summer * Advanced Surveying Me<- -urements (CE362) 4
Geology for Engineert.- Gieology 218) 4
Fluid Mechanics (Eng. i.-ich. 324) 3
* Theory of Structures::13) 3
5 Transportation Engirt - (CE370) 3 17
Hydrology (CE420) 3
(*) Advanced Mathematics:ive 3
Elective in Humanities...: Social Science 2
Dynamics (Engr. Mech. 343) 3
* Reinforced Concrete (CE415) 3
6 * Water Supply and Treatment (CE480) 3 17
* Hydraulics (CE421) 2
Engineering Properties of Soil (CE445) 3
Elective in Humanities or Social Science 3
Specifications and Contracts (CE400) 2
* Elementary Design of Structures (CE313) 3
7 Electrical Apparatus and Circuits (EE315)4 16
* Sewerage and Sewage Treatment (CE481) 2
English Elective 2
(*) Technical Elective 3
Thermodynamics and Heat Transfer (Mech. Engr. 333) 4
(*) Technical Electives 7
8 Modern Economic Society (Econ. 401) 3 17
Elective in Humanities or Social Science 3
* Course Related to or supplemented by computers.
(*) Student may elect course related to or supplemented by computers.
-C1C
</P>
<P><PB REF="00000011.tif" SEQ="00000011" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="885" N="11">
Use of Computers in Civil Engineering Education
III. INTRODUCTION OF CIVIL ENGINEERING STUDENTS TO THE USE OF THE COMPUTER
Surveying Computations (CE261) has been selected as the medium for introducing Civil
Engineering students to the use of the digital computer. This course is the second in a
series of three required courses in surveying and is normally taken in the second semester of
the sophomore year.
Civil Engineering students are not required to take the one-hour course in introductory
computing techniques (Math 373) which is required of students in most of the other fields of
engineering. Surveying Computations has, therefore, been designed as an effective substitute.
A. Preparation of Input for the Computer
The first half of the semester is devoted to applications of analytic geometry for plane,
rectangular coordinate computations. Emphasis is placed on five basic problem types. These
types, described in detail in Example Problem No. 76, have the following names: Forward
Computation, Inverse Computation, Line-Line Intersection, Line-Circle Intersection, and
Circle-Circle Intersection.
The students are shown the recommended computing methods, supported by a mathematical
derivation, for each problem type. After studying the assignment, listening to the lecture,
and assimilating the methods taught, each student is assigned a set of fifteen individual problems of each type. The student is instructed initially to perform the necessary calculations
for only the first five problems in each set and to prepare both the given information and the
solution on punched IBM cards.
These cards are then submitted as data to a computer program prepared specifically for
the course,called the "teaching machine program." This program, described in more detail in
Example Problem No. 76, has access to all of the given data and computes the correct results
for each of the students' problems. It compares all of the items on each of the students'
cards with the known, correct values. If any single item in any one of the five solutions is
found to be in error, the whole problem is considered incorrect. For each incorrect problem
the student is required to correct all of the errors, and also solve an additional problem in
the assignment. The original five problems (with errors corrected) plus the additional problems are then submitted in a second approach to the computer. Thus, for the second approach,
the student may be required to enter from six to ten problems. The same rules apply again, for
a third and final approach to the computer. The worst possible penalty would be that on the
final approach to the computer the student would be required to present solutions to all fifteen
of the assigned problems.
If all of the problems submitted in the final approach to the computer are not correct,
the student is penalized in his grade. On the other hand, if a correct solution is obtained
-Cll
</P>
<P><PB REF="00000012.tif" SEQ="00000012" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="863" N="12">
to all of the problems submitted on any given approach to the computer, the student is given
full credit regardless of whether success is achieved on the first, second or third trial.
It may be seen that the students quickly realize that it is to their advantage to obtain
the correct results for submission on the first pass through the computer. In addition to the
savings in their calculating time (on the desk calculator), by solving the first five problems
of each type correctly and finishing the assignment, it has been noticed that a friendly spirit
of rivalry generally develops between the students in which those who do not complete the assignment on the first pass receive frendly, but reslightly barbed comments on their proficiency.
The "teaching machine" theory seems to be working quite well for these problems. The
students quickly overcome their difficulties and by the time they have finished the five problem
types, they show a marked proficiency in the use of trigonometric tables and the desk calculator. In order to meet the requirements for tolerance in the teaching machine program, the students learn by experience the necessity for using significant figures properly in their
computations. For example, they are shown, by actually working problems both ways, that in an
intersection of a line and a circle, it is necessary to interpolate fractions of seconds in
order to obtain enough accuracy when the intersection angle is very small and the radius of the
circle is very large. These are things which cannot be easily taught in classroom lectures but
which can be illustrated very pointedly to the student during his laboratory work.
This three-semester-hour course is handled with a one-hour lecture per week and two
three-hour laboratory periods per week. The first hour of each period is devoted to a lecture,
followed by two hours of actual computation. During this time, the students solve problems in
the laboratory room using desk calculators. They are encouraged to ask the instructor for help
at any time on problems involving the theory of the computations, the mechanics of the desk
calculators, the manipulation of trigonometric tables, or an explanation of lack of accurate
checking of their completed problems.
The "teaching machine" idea serves also to introduce the students to the digital computer.
They learn to use the keypunch to punch their own cards before they are actually assigned to
write a program for the digital computer. They are strongly impressed with the capabilities of
the computer when they see how it is able to solve their problems and check them at a rate of
approximately two of their problems per second (on the IBM 709). They realize quite painfully
that these problems cost them sometimes as much as an hour or more of desk calculator work.
Because of this, when they are later permitted to begin programming for the computer, they do
not have the reaction so often expressed by students who are given elementary problems, "I could
have worked the problem faster on a desk calculator."
-C12 -
</P>
<P><PB REF="00000013.tif" SEQ="00000013" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="879" N="13">
Use of Computers in Civil Engineering Education
B. Geometric Analysis Using the Interpretive Language, AGILE
The writer has developed a computer program over the last several years which has recently
been named "AGILE" (Analytic Geometry Interpretive Language for Engineers). A number of small
programs to solve problems such as the five basic types described in Example Problem No. 76, have
been combined with a reading, storing, and scheduling master program. In this way, any geometry
problem (which may consist of a number of basic problems) can be solved with one pass through
the computer.
After learning to compute the basic problems already described, the students are taught
to recognize them when they occur in such applications as the design of a subdivision or a
highway interchange. The students are also taught to describe these problems in a logical manner,
using the AGILE language. A complete description of the AGILE language is being prepared. If
sufficient interest develops, a manual describing the language and its use will be made available.
C. Computer Programming
Beginning with the first Wednesday night of the semester, a series of three weekly lectures
of two hours each, on computer programming, is presented for the entire College of Engineering.
These lectures were presented in 1961-1962 by Mr. Brice Carnahan, Assistant Director of the
Project on Computers in Engineering Education. The lectures introduce the students to the
philosophy of computers, and specifically, to the MAD language. The students are not required
to attend these lectures, but are strongly urged to do so, being told that no substitute for
them will be presented in classroom lectures.
Beginning immediately after the series of night lectures, the students are assigned their
first problem which requires them to write a computer program for the course. This initial
program is for a forward computation (the first problem type listed in Example Problem No. 76 ),
and is handed to the students completely checked out and ready to run on the computer. They
are required simply to prepare cards for communicating with the computer. They are also asked
to punch data cards for a problem which is given to them. The form of the data cards is slightly
different than the form in which the data are presented. In addition to the punched cards, the
students are required to submit a flow diagram for the problem and a short discussion of the
problem, its solution, and the method which a user must follow in preparing data for a run on
the computer.
By this means the students are taken over a second, short step towards computer programming
by being introduced as painlessly as possible to the computer itself. The submission of a program which has already been written for them may sound like a wasted step but it has been found
that this is a step which, when combined with other complications, can become a very frustrating
experience to the neophyte.
-C13 -
</P>
<P><PB REF="00000014.tif" SEQ="00000014" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="887" N="14">
This first program is explained thoroughly to the students so that they understand what
is done by the computer as a result of each instruction in the program. Two very elementary
loops are included so that the student becomes familiar very quickly and painlessly with the
concept of iterative techniques on the computer. The students are permitted to submit this
program to the computer only once. It is felt that since everything has been handed to them
they should be graded on their ability to follow instructions.
The second computer program is for the inverse computation. In this case, they are given
a clear concept of the theory and the method of solution of the problem but are required to
make up their own instructions in MAD for the solution of the problem. Again, they are required
to submit a flow diagram, and a written report which includes an explanation of the method
necessary to use the program successfully on the computer. In the solution of this problem the
students are required to write their own Newton-Raphson iteration for the square root solution
of the distance. They are also asked to solve for the distance by using either the sine or the
cosine of the angle which has previously been computed from the arctangent.
This gives them their first introduction to the use of library subroutines, the formal
part of MAD programming including the loop concept, and also gives them a chance to compare the
relative accuracy by which they can obtain the distance by two completely independent methods.
In the third program which is assigned in this course, the students are instructed to
write a computer program which will calculate the volume of a borrow pit. In this problem,
the students are asked to bring in the "before" and "after" elevation readings over a grid pattern which is basically rectangular. They are asked to work with these two matrices, subtracting
the "after" readings from the "before" readings to obtain a new matrix called the "excavation."
This subtraction is done by the use of a double loop. Following this, the computations are
performed by means of some single and double loops. This helps the student become familiar with
methods of handling large sets of subscripted variables.
The fourth and last program is one which may vary considerably from semester to semester.
It is designed to offer a challenge to the students. It is intentionally changed radically from
one semester to the next in order to discourage any ideas of propagating information from one
class to another. Examples of the problems which have been given previously are as follows:
(1) Intersection of a straight line with a spiral envelope.
(2) Approximate adjustment of a triangulation quadrilateral.
(3) Progressive adjustment of interlocking traverses by the Bowditch (compass) Method.
(4) Crandall method of adjustment of a traverse.
At the time of this writing it is strongly believed that the inclusion of computer programming in this course has not detracted from the course work in any way. In fact, it is felt
that the blending image of the desk calculator, the digital computer, and surveying computation
methods, has resulted in a course which has materially strengthened the students' understanding
in each of the three fields mentioned.
-C14 -
</P>
<P><PB REF="00000015.tif" SEQ="00000015" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="870" N="15">
Use of Computers in Civil Engineering Education
Some of the leading exponents of computer philosophy and education feel that geometry
is one of the best vehicles for teaching computer ph-il-sophy. The student's ability to follow
graphically and pictorially every step in his solutionf makes it possible for him to concentrate
on the programming necessary to solve the problem. By the same token, his detailed computer
solution of the problem brings him into a very intimate relationship with the problem itself.
By the time he has successfully programmed a solution, he is completely familiar with every
aspect of the problem.
It has been found that many of the students are completely apathetic to computers and
the philosophy of computation. These students manage to pass the course by literally hanging
on for dear life. Other students in the middle category passively accept the work and perform
the assignments satisfactorily but with little imagination. Finally, we have the group who
become quite interested in computers and their possibilities and who go beyond the call of duty,
polishing up their programs to achieve more than the assignment requires. Certainly these
students are in the minority, but it is encouraging to notice that their number seems to be
increasing. It is believed that this enthusiasm is growing because the work is passed down to
succeeding generations by those who have run into advanced problems which required computer
applications. Also, it is quite apparent to the students that many interviewers from industrial
firms who come to the campus are now asking about the amount of computer work students have been
exposed to. These students, especially those who are doing well in their other subjects (primarily mathematics), are encouraged to sign up for Mathematics 473, Introduction to Digital
Computers, a comprehensive, senior-graduate level com coter course.
It is not our intention to force the majority of ->r students into computer technology,
but it is felt that all students should have some kno;wle-ge of computers and that some few students should have a specialized knowledge of computer programming techniques. In Surveying
Computations, it is possible for the students to decide which of these paths they wish to
follow. Upon completion of this course, no student will feel that a computer is a mysterious
black box which possesses magical properties. On the contrary, he will realize that it is most
simple for him to solve elementary problems and that it is not impossible for him to solve even
the most difficult problems by using it.
</P>
<P><PB REF="00000016.tif" SEQ="00000016" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="886" N="16">
IV. SPECIFIC USES OF COMPUTERS IN VARIOUS FIELDS OF STUDY
A. Geodetic Engineering
The required courses in surveying are designed to give the student a realistic appreciation for the theory of measurements, application of mathematics to geometric figures, logical
analysis, and systematic computation procedures. These courses provide six of the first seven
semester hours of contact between the Civil Engineering Department and the students.
The first course, Basic Surveying (CE260), does not lend itself to student use of computers. However, the instructor assigns a unique traverse computation problem to each student
and compares his answers to those provided by the computer. Other similar uses are being explorec
The second required course, Surveying Computations (CE261), has been described in detail
in Part III of this report.
The terminal required course, Advanced Surveying Measurements (CE362), is taught at Camp
Davis, Jackson, Wyoming. This course is primarily one engaging in precise control surveys using
theodolites, electronic distance measuring devices and invar tapes. Although no computers are
available in the vicinity of the camp which can be used for the course purposes, programs have
been written which will make it possible to send data from camp to the University Computing
Center and receive the results back in time to use them for the purposes of the course. At the
present time these programs are primarily concerned with reducing observation data, such as
Tellurometer observations, to yield the actual measured values from the observations. These
reductions, when performed by hand, are quite laborious. Even more important to the user, however, such reductions are filled with potential errors, which, because they deal with measurements, cannot be checked. The use of the computer for this purpose virtually eliminates all
chance of such blunders in reductions.
Geodetic Field Methods (CE562) consists of unique problems which vary with the interests
of the students and the program of the department. Usually there are several applications for
computer solutions of problems which occur in this course.
Adjustment of Geodetic Measurements (CE563) depends very heavily on digital computer
applications. Problems are programmed and solved by means of the Gauss-Seidel iterative technique. Matrix inversion and other methods of solving large systems of linear equations are
also applied to the computer.
Special Problems in Advanced Surveying (CE564). Problems in this course again are quite
variable and may or may not require computer applications.
Municipal Surveying (CE565) consists in part of several conducted tours to various private
and municipal or governmental agencies. Oneprivate computing firm, the Bureau of Surveys of
the City of Detroit, and the Michigan State Highway Department were three locations visited
where computers were being used. Students observed the computers in operation and were told
by the engineers in charge exactly what applications were being made in the specific fields of
their needs. -C16 -
</P>
<P><PB REF="00000017.tif" SEQ="00000017" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="891" N="17">
Use of Computers in Civil Engineering Education
Geodesy and Surveying Research (CE960) is another example of a course which may have
varying interests or needs. Certainly many students in the future who take this course will
be carrying on extensive computer research in the field of geodetic surveying.
B. Hydraulic Engineering
During the last three years a majority of the courses offered in hydraulics and hydraulic
engineering have been completely re-examined, and now require use of the digital computer as a
part of the course work.
The required course, Hydraulics (CE421), has in the past offered the student an opportunity
to program one or two problems on a voluntary basis. Three or four students, out of about 40,
usually wrote programs. Beginning with the spring semester (1962), however, two programs are
required of all students: the first the determination of critical depth, normal depth, water
surface profiles, and location of a hydraulic jump in a trapezoidal channel; and the second, the
determination of the best formula for a v-notch weir, from student data, by the method of
least squares.
The course in Hydraulic Transients (CE522), has evolved from a hydropower course, primarily
because of the great improvement in means of solving problems in unsteady flow. Problems in
surge and water hammer have been programmed, and the instructor has set up the differential
surge tank problem for class inspection. New methods for handling water hammer have been
developed by graduate students as a result of this course and the availability of the digital
computer.
The course in Advanced Hydraulics (CE524), relies heavily on the computer, and has used
it for such problems as relaxation methods for computing flow under a dam and for applications
of the momentum theory. One significant problem, the design of a side channel spillway, has
been given very thorough study, and a new method has been developed to optimize the design.
n
Two unknown parameters are needed: a and n in V = axn, for velocity at a cross section x
distance along the spillway for minimum total excavation. A "gradient method" is used, to
determine how a and n should be changed to move in the direction of minimum excavation.
The Hydraulic Design course (CE526), is built entirely on the use of the computer.
At present it consists of several programs relating to a project involving construction of a
gravity dam, mass diagram study of the inflow hydrograph, flood routing study, reservoir operation study for flood control, power and irrigation supply, and finally an economic study which
determines the height of the dam to be constructed.
It is believed that the students in the hydraulic option are much better prepared to meet
and adjust to the future demands of their profession as a result of their ability to think through
the logical processes required in a program. It also gives them a definite competitive edge in
securing interesting job assignments with greater growth potential.
-C17 -
</P>
<P><PB REF="00000018.tif" SEQ="00000018" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="874" N="18">
C. Sanitary Engineering
The digital computer is being used for network analysis in Water Supply and Treatment
(CE480), as an optional problem which has been worked out successfully by several students.
In Sewerage and Sewage Treatment (CE481), the design of a storm sewer system is at present being worked on as an optional assignment.
In Industrial Bacteriology (CE587), a problem has been developed requiring the solution
on the computer of the MPN (Most Probable Number) equation for determining'bacterial populations.
In Advanced Sanitary Engineering Design (CE682), the programming of the sewage network
analysis is to be described and discussed but no actual programming is anticipated. However,
the students may elect to do so if they desire.
D. Soil Mechanics and Foundation Engineering
The traditional aim of the courses in soil mechanics and foundation engineering at The
University of Michigan has been to train the student to view problems not only in terms of
mathematical formulations but also in relation to environmental associations. Particular
emphasis has been placed on framing problems in terms of simplified basic relationships which
can be conveniently handled by slide rule or desk calculator. One observes that the conscientious student develops a "feel" for the subject in this approach which is more readily adapted to
variations in field conditions. It is desired to maintain this emphasis on the practical
application of fundamental principles and, at the same time, to take full advantage of the
spectacular developments in computational technique with present day electronic computers.
In this respect, the computer makes possible more exacting solutions of certain soil
mechanics problems and eliminates compromising assumptions and approximations which have
heretofore been necessary. Consequently, the soil mechanics teaching and research staff have
re-examined a number of the more important problems that require extensive computations and
have developed computer programs for the following:
1. Analysis of the plate loading tests to determine the bearing capacity of soil.
2. Determination of pressure distribution from large loaded areas.
3. Static equilibrium of sheet pile retaining walls under a variety of nonhomogeneous soil
conditions.
Additional problems under investigation include computation of settlement in granular
materials from volume-density relationships, a general solution for critical stresses under
surface loading and the reproduction of pavement profiles for a variable length of base line
reference. The solution of these more complex problems has not advanced to the stage where they
can be utilized in undergraduate instruction. They are, however, being developed and used by
graduate students working on these specialized problems.
The objective in presenting the Buse of computers to students of soil mechanics is to
demonstrate the value of these more exacting solutions and familiarize the student with their
application to soil mechanics problems. It is felt that the real value of computer techniques
</P>
<P><PB REF="00000019.tif" SEQ="00000019" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="873" N="19">
Use of Computers in Civil Engineering Education
to soil mechanics lies in the rigorous organization of the problem in a sequence of logical
steps leading from fundamental principles and realistic assumptions to a useful practical
solution. It is contemplated that the use of computers in soil mechanics will be included
in the advanced design course.
E. Structural Engineering
Computers are used in the structures sequence both at undergraduate and graduate levels.
In Elementary Structural Design (CE313), some experiments have been carried on in the use
of the computer in simple design problems. This has not been entirely successful, largely
because the students have too little experience in the use of the computer in earlier work.
There is not sufficient time to include programming in the course lectures.
In Theory of Structures (CE512), and Rigid Frame Structures (CE514), the computer has been
used to solve problems in truss and beam deflections, using conventional procedures. The use
of computers in these courses is still on an experimental basis.
In Structural Design (CE513), the computer is used regularly each semester and one computer problem is assigned. In a typical semester, the students are asked to prepare section
modulus and shear connector tables for composite beams in accordance with the AISC 1962
Specifications. The problem is relatively simple but involves quite a number of logical decisions that require a more complete knowledge of composite beam behavior than would the comparable
hand solution for only a few special cases. There is nothing artificial about the problem.
It is the type of problem that one might expect to encounter not infrequently in practice.
In Structural Dynamics (CE611), the computer is regularly used for one or two problems
each semester. Ordinarily the problems involve numerical solution of the differential equations
of motion and constructing response spectrum curves for some given dynamic load.
In Structural Plate Analysis (CE613), one problem is assigned each semester in finding
the deflections of a plate by finite difference methods. The problem is taken up in class at
two levels; first the solution is obtained by hand for a coarse mesh, and later a finer mesh
solution is prepared using the computer.
Computer Analysis of Structures (CE618), is a course serving the specific purpose of
providing the background needed for using the computer efficiently in structural analysis.
Matrix methods are of special value in computer analysis of structures, so the course is built
on matrix techniques. Most of the course is devoted to formulating the theory of structures
in the language of matrices. Assigned problems include small, near-trivial problems where all
of the matrix operations are carried out by hand, and larger, non-trivial problems for which
the students program the matrix operations for the machine. The smaller problems serve the
purpose of familiarizing the student with the details of the procedures and also retaining a
physical interpretation for all of the various matrix operations.
-C19 -
</P>
<P><PB REF="00000020.tif" SEQ="00000020" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="876" N="20">
F. Transportation Engineering
The staff use of computers in the area of traffic and transportation is somewhat limited
at the present time, but will undoubtedly increase in the future. In CE575, Transportation
Terminal Design (highway, rail, water, and air), and CE670, Transportation Planning, the problems, insofar as possible, are placed in a local Ann Arbor context. Both the City of Ann
Arbor and Washtenaw County have acquired voluminous data on local transportation, land use, and
population. Much of this data is available on punched cards so that it can be sorted in any
manner desired. Thus in design problems the computer is used to summarize and compile data on
specific geographic areas interconnected by transportation facilities.
There are some iterative type solutions involved in layout and design of automobile parking lots and in the delineation of retail trade areas for shopping centers which lend themselves
quite nicely to computer application. At the present, a study is being made of the extent to
which it is desirable to require the student to use the computer in solving these problems.
Transportation planning problems quite often are examined on a systems analysis basis. While
the computer has not been used to any great extent in this field, there is a potential for
widespread usage.
To summarize, the various areas in traffic and transportation courses and the extent to
which it is desirable for the student to become involved with the computer are currently being
examined. Use at the present time, however, is largely experimental.
V. CONCLUSIONS
As in the case of every important innovation, it will take several years to evaluate the
ultimate benefits resulting from the introduction of computers into Civil Engineering education.
At present, opinions of the faculty are varied, as evidenced by the following sampling of
statements.
"The Civil Engineering Department at The University of Michigan has successfully initiated the use of computers into the curriculum and the faculty is gradually understanding the
value of adapting their work to accept and include these new techniques."
"The real benefits of using the computer come not from teaching the student how to program,
but rather from forcing him to make all of the logical decisions that are necessary to arrive at
a suitable algorithm for solving an engineering problem in a more general formulation than he
usually encounters. It is often' said that a good way to learn a subject is to teach it. In a
sense, the student is being put in a position of teaching the computer how to solve the problem
at hand. In order to accomplish this he must necessarily understand the problem more completely
than would be necessary for conventional hand solution, since all of the steps involving judgement that are at his disposal in hand methods must be analyzed and programmed when the machine
is doing the work."
-C20 -
</P>
<P><PB REF="00000021.tif" SEQ="00000021" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="889" N="21">
Use of Computers in Civil Engineering Education
"While it is agreed that computers are here to stay and are a valuable tool in the
civil engineering profession, we must be very careful not to sacrifice the teaching of fundamental engineering principles in order to teach the use of computers."
"The method of introducing Civil Engineering students at Michigan to computer techniques
by integrating it with the teaching of surveying computations seems to be an excellent approach."
"In introducing the use of computers, it is necessary to maintain the emphasis on the
practical application of fundamental principles and, at the same time, to take full advantage of
the spectacular developments in computational technique with present-day electronic computers.
The real value of computer techniques lies in the rigorous organization of the problem in a
sequence of logical steps leading from fundamental principles and realistic assumptions to a
useful, practical solution."
"We need to make fuller use of computers in all fields of Civil Engineering. In this
way our graduates will be much better prepared to meet and adjust to the future demands of
their profession as a result of their ability to think through the logical processes required
to develop algorithms for the computer."
VI. EXAMPLE PROBLEMS PREPARED BY PROJECT PARTICIPANTS
Each participant prepared solutions for example problems which he felt would be useful
in courses at the undergraduate or graduate level. Nine of these problems have been selected
for inclusion in this report and are listed below.
Number*
76
77
78
79
LIST OF PROBLEMS
Title
Use of the Computer as a Teaching Machine
for Plane Analytic Geometry
Optimum Regulation of a Reservoir for
Power Production
Economical Size of Canal-Feeding HydroPower Plant
Most Probable Number Method for Determining
Bacterial Populations
Iteration of Sum of Principal Stresses at
Interior Points of Plane Rings
Elastic Buckling Load for Columns of NonUniform Cross Section
A MAD Program for Truss Analysis
Stiffness Factors for a Flat Slab Column
Vibration of Beams on Spring Supports
Author
H. J. Welch
J. W. Howe
J. W. Howe
L. L. Kempe
W. Weaver, Jr.
K. H. Chu
K. H. Chu
R. B. Harris
S. S. Kuo
Page
C23
C35
C41
C47
C56
C67
C82
C97
C113
80
81
82
83
84
* These problems may be considered as a supplement to problems presented in
tions of the Project (3, 4, 5, 6) and are numbered accordingly.
previous publica
-C21 -
</P>
<P><PB REF="00000022.tif" SEQ="00000022" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="873" N="22">
All of these problems were solved on the digital computer using either the MAD or FORTRAN
language. These languages are described in References 1 and 2.
It should be noted that there are several other sources of example computer problems
dealing with Civil Engineering subjects. See, for example, Example Problem Numbers 8, 9, 10,
21, 22, 23, and 27 in Reference 3; Example Problem Numbers 48, 49, and 50 in Reference 4; and
Example Problem Numbers 63 and 64 in Reference 6. A number of problems prepared by Project
participants are abstracted in the Second Annual Report of the Project (4). Some of these
(listed above) have been included in this report. Most of the remaining problems are available
on microfilm from the Project.
-C22 -
</P>
<P><PB REF="00000023.tif" SEQ="00000023" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="817" N="23">
Example Problem No. 76
USE OF THE COMPUTER AS A TEACHING MACHINE
FOR PLANE ANALYTIC GEOMETRY
by
Harold J. Welch
Department of Civil Engineering
The University of Michigan
Course: Surveying Computations Credit Hours: 2 Level: Sophomore
Statement of Problem
In learning surveying computations, students are required to solve a number of problems
in plane, analytic geometry using a desk calculator and tables of trigonometric functions.
The problems consist of several examples of each of five basic problem types illustrated below.
1. Forward Computation (X2,Y2)
Given: X1, Y1, D1, 01 (X1,Y1)
i D1
Find: X2' Y2
01
2. Inverse Computation
(X2,Y2)
Given: X1, Y1, X2, Y2 (Xl,)
D1D
Find: Di, Q 1
3. Line-Line Intersection
Given: X1, Y1, G91 X2, Y2, 02
D D2
Find: X3, Y3, D1, D2 (
3 1(X' Y2)
4. Line-Circle Intersection
Given: X1, Y1, 91' X2, Y2, D2 x, Y3)
Find: X3, Y3 82, D1 2 ( Y2)
3-9 3 2D
-C23 -
</P>
<P><PB REF="00000024.tif" SEQ="00000024" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="879" N="24">
The Computer as a Teaching Machine
5. Circle-Circle Intersection
D \
Given: X1, Y1, D1, X2, Y2, D2 (X1 Y ),
Find: X3, Y3, 01, G2
The only difference between two problems of the same type is that the numerical data
are different. A computer program has been written to check the students' solutions to problems assigned in class and point out any errors. The program has been written to "penalize"
students who have incorrect solutions by assigning them additional problems. On the other
hand, students with correct solutions are "rewarded" by not being assigned additional problems.
Because of this feature of penalties and rewards, the program which grades the student problems
has been named the "Teaching Machine Program."
Use of the Teaching Machine Program
Each student is given an identification number and a set of basic data for each problem
type. He then uses his identification number to modify the basic data. For example, if the
student's number were 1234, selected distances and point coordinates would be changed by
4253.897
1.234
4255.131
Azimuths (line directions) would be changed according to
290~ 47' 18"
12 34
290~ 59' 52"
In this manner, all the students work essentially the same problem, but with slightly different data and results.
When the student has solved a problem to his satisfaction, he punches both his modified
data and his computed results on IBM cards. On the first card, the station number and the
coordinates of the first point are punched. The length and azimuth of the line between the
first two points are punched on the second card. This sequence is repeated until all of the
problem data havebeen given. Hence, for the Forward and Inverse Computations, each problem
requires three cards. For each Line-Line Intersection problem, five cards are required.
Examples of input cards are shown on the next page.
-C24 -
</P>
<P><PB REF="00000025.tif" SEQ="00000025" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="872" N="25">
Example Problem No. 76
A}, n Bra
Col 1 3 M id E2
655 3730.751 3636.576
Station Coordinates
Number
Col 26 38 42 45
1415.702 332 36 13
Length Azimuth
On the first run of each type of problem, the student must work 5 problems (3 for the
Circle-Circle intersection). In front of the data cards, he places a card bearing his name,
class number, computing center number, and the type number of the problems which follow. For
example, a name card might be
Col 1 21 32 39
/J Q STUDENT 1234 X106Z 4
Student Identifica- Computing Problem
Name tion Number Center Type
Number Number
The data for each type of problem must be submitted separately, with an individual name card
for the set of data for each type.
At execution, the program first reads the current date to be used in output headings.
Following this is the class roll. Each card in this deck contains a student's name and class
number, the number of problems due and the number of the run to be made for each of the five
types of problems. Next, the five sets of basic data are read and then the first student's
data is considered.
The class roll, the sets of basic data, and each set of student data is followed by a
card having an equal sign punched in column 1. Whenever this equal sign is read, control is
transferred to the portion of the program which handles the type of data expected to follow.
In processing student data, the information on the student's name card is read and
printed in a heading. Then the student's class number is checked to see if it is a legal
number. The corresponding run number is checked to see if the student is attempting an extra
run. Finally, the student's name is checked to see if it corresponds with the given class
-C25 -
</P>
<P><PB REF="00000026.tif" SEQ="00000026" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="888" N="26">
The Computer as a Teaching Machine
number. If any of these tests fail, an appropriate comment is printed and the program searches
for the name card of the next student. If the information is all legal, control is transferred
to the section of the program which handles the type of problem given on the name card. Each
of these five sections modifies the basic data in the same manner as the student and computes
the correct solution to the problem. These results are then used by the internal function CHECK.
CHECK reads the student's data and compares it with the computed data, with a tolerance
of.003 feet on distances and coordinates and 3 seconds on azimuths. The student's data are
then printed, with asterisks flagging any incorrect values. Whenever an error occurs, a
boolean flag is set, and at the end of the problem, the value of the flag is returned to the
main program. If the flag has been set, the main program increments an error counter. Then
the data for the next problem (of the same type) are examined, the correct solutions are
computed, and CHECK is again called. This process continues until the required number of problems for this sequence have been considered or until an equal-sign card is read. Extra
problems are ignored and simply scanned for an equal-sign card.
At the end of the set of data, if the error counter is not equal to zero, or if too few
problems have been submitted, the number due is increased. On the next run, the student must
resubmit all the correct problems, the problems which were incorrect or missing (hopefully
corrected by this time), and an extra problem for each incorrect or missing problem. This new
due number is printed following the data listing. If all the problems are correct, the due
number is set to zero and a comment is printed.
The run number is now incremented and the program goes on to the next student's data.
After all the problem sets have been processed, the updated class roll is punched for
use in the next run and printed for the information of the instructor.
Flow Diagram
A general flow diagram of the Teaching Machine Program is shown below.
C 1 - Forward Computation
C 2 - Inverse Computation
C 3 - Line-Line Intersection
C 4 - Line-Circle Intersection
C ) - Circle-Circle Intersection
-C26 -
</P>
<P><PB REF="00000027.tif" SEQ="00000027" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="772" N="27">
Example Problem No. 76
MAD Program
The Teaching Machine Program (which is written in the MAD language) is listed below.
H J WELCH 217 W ENG X106W 002 050 050 9C8
H J WELCH 217 W ENG X106W 002 050 050 9C8
$COMPILE MADEXECUTEPRINT OBJECTPUNCH OBJECTDUMP 629C8001
EXECUTE SETEOF.(LAST)
EXECUTE SETERR.(SRCH)
READ FORMAT RDATEDATEoATE1DATE2
THROUGH ROLL. FOR N=11tN*G1000
READ)- FORMAT RD1,EQUALSNAME(N),NAME1(-N) NAME2(N —CN-(-N-T-DUEIITN
1),DUE2(N) DUE3(N),DUE4(N),DUE5(N),RUN1(N)RUN2(N),RUN3(N) RUN
24(N) RUN5(N )
ROLL WHENEVER EQUALS.E.$=$,TRANSFER TO OUT
PRINT CO MMENT$4T00 MANY STUDENTS$
EXECUTE SYSTEM.
OUT THROUGH LFOR- Ji,,J- G25
READ FORMAT RDEQUALS,FROM(J),TO(J),FRX(J),FRY(J),DTO(J),D(J)
1.M(J) S(J)
L WHENEVER EQUALS*Ee$=$W1RANSFER TO NEXT3
RTINTCOfMNENT -'$4TO0 MUCH DATA FOR PROBLEM 1$
EXECUTE SYSTEM*
NEXT3 THROUGH L2,FOR J=l,1,J1G425
READ FORMAT RD2,EQUALSFROM2(J),FRX2(J)FR(J) (J)T02(J)*TOX2(J
1ltTOY2(J)
L2 WHENEVER EQUALS*E.$=$,TRANSFER TO NEXT4
PR I NT COMMENT $4TOO MUCH DATA O'R -PROBLEM —2$'
EXECUTE SYSTEM.
NEXT4 HTHROI -L3 FOR ~J1-i J G25
READ FORMATRD3,EQUALSFROM3(J),FRX3(J),FRY3(J),FD3(J),FM3(J),
1FS3(J) INT3(J) TD3(J),TM3 (J),TS3(J) T03T(J,TOX3(J) TOY3(J)
L3 WHENEVER EQUALS.Ea$=$,TRANSFER TO NEXT1
'P-INT —CO-MMENT$4TO'-1o-MUH-' DATA 'F'R- PR-OBLEM 3$ —'
EXECUTE SYSTEM.
NEXT 1 THiROU-GFI- 4 ' FOR-J=-1I -i"-J.G,25
READ FORMAT RD4,EQUALSFROM4(J),FRX4(J),FRY4(J),FD4(J),FM4(J)
-C27 -
</P>
<P><PB REF="00000028.tif" SEQ="00000028" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="832" N="28">
The Computer as a Teaching Machine
MAD Program, Continued
IFS4(J),INT4(J),DTO4(J),TO4(J),TOX4(J),TOY4(J)
L4 WHENEVER EQUALS.E.$=$,TRANSFER TO NEXT2
PRINT COMMENT$4TOO MUCH DATA FOR PROBLEM 4$
EXECUTE SYSTEM.
NEXT2 THROUGH L5, FOR J=1~1,J.G.25
READ FORMAT RDSEQUALStFROM5(J),FRX5(J),FRY5(J),DFM5(J),INT5(
1J),DTO5(J),T05(J),TOX5(J),TOY5(J)
L5 WHENEVER EQUALSE.$=$,TRANSFER TO NEWNAM
PRINT COMMENT$4TOO MUCH DATA FOR PROBLEM 5$
EXECUTE SYSTEM.
NEWNAM ERR=O
TER=06B
READ fORMAT RDNAMEEQUALSSNAME..,SNAME(2) SCN,PROB
WHENEVER EQUALS.E.$=$.TRANSFER TO NEWNAM
PRINT FORMAT TITLEBSNAME*.-SNAME(2),SCNPROB,DATE,DATE1
1,DATE2
THROUGH L1i FOR K=1,,K.E.NN
L1 WHENEVER SCN.E.CN(K),TRANSFER TO TSTRUN(PROB)
PRINT COMMENT $OINCORRECT CLASS NUMBER$
TRANSFER TO SRCH
TSTRUN(1) WHENEVER RUN1(K).L.4,TRANSFER TO TSTNAM
PRINT COMMENT $ORUN LIMIT EXCEEDED ON PROBLEM 1$
TRANSFER TO SRCH
TSTRUN(2) WHENEVER RUN2(K).L.4oTRANSFER TO TSTNAM
PRINT COMMENT $ORUN LIMIT EXCEEDED ON PROBLEM 2$
TRANSFER TO SRCH
TSTRUN(3) WHENEVER RUN3(K).L*4#TRANSFER TO TSTNAM
PRINT COMMLNT $ORUN LIMIT EXCEEDED ON PROBLEM 3$
TRANSFER TO SRCH
TSTRUN(4) WHENEVER RUN4(K).L.4,TRANSFER TO TSTNAM
PRINT COMMENT$ORUN LIMIT EXCEEDED ON PROBLEM 4$
TRANSFER TO SRCH
TSTRUN(5) WHENEVER RUNS(K).L.4,TRANSFER TO TSTNAM
PRINT COMMENT$ORUN LIMIT EXCEEDED ON PROBLEM 5$
TRANSFER TO SRCH
TSTNAM WHENEVER NAME(K)*NE.SNAME.OR.NAME1(K).NE.SNAME(1).OR.NAME2(K
1)oNE.SNAME(2)
PRINT COMMENT $OINCORRECT NAME$
SRCH READ FORMAT SEARCHEQUALS
WHENEVER EQUALS.E.$=$*TRANSFER TO NEWNAM
TRANSFER TO SRCH
END OF CONDITIONAL
TRANSFER TO C(PROB)
C(1) PRINT FORMAT TITLECPROBDUE1(K),RUN1(K)
THROUGH BIGLP1,FOR I=l,l,I.G.DUE1(K)
STA=FROM(I)
STA(1)=TO(I)
X=FRX(I)
Y=FRY(I)
DIST=DTO(I)+SCN/1000.
D1=D(I)
M1=M(I)+SCN/100
S1=S(I)+SCN-100*(SCN/100)
SEC=3600*Dl+60*M1+S1
RA=,48481368E-5*SEC
X(1)=FRX(I)-SIN.(RA)*DIST
Y(1)=FRY(I)-COS(RA)*DIST
TER=CHECK.(THRU1)
WHENEVER TER
ERR=ERR+1
TER=OB
BIGLP1 END OF CONDITIONAL
WHENEVER ERR*E.O
DUE1(K)=O
PRINT FORMAT OKPROB
TRANSFER TO RK1
END OF CONDITIONAL
DUE1(K)=DUE1(K)+ERR
TEST1 WHENEVER DUEl(K).G.15,DUE1(K)=15
PRINT FORMAT RESMITDUE1(K)
RK1 RUN1(K)rRUN1(K)+l
TRANSFER TO SRCH
THRU1 DUE1(K)=2*DUE1(K)+l+ERR-I
TRANSFER TO TEST1
C(2) PRINT FORMAT TITLECPROBDUE2(K),RUN2(K)
THROUGH BIGLP2,FOR I=1,,1I. G.DUE2(K)
STA=FROM2(I)
-C28 -
</P>
<P><PB REF="00000029.tif" SEQ="00000029" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="826" N="29">
Example Problem No. 76
MAD Program, Continued
X=FRX2(I)+SCN/1000.
Y=FRY2(I)
DX=TOX2(I)-SCN/1000.-X
DY=TOY2(I)-FRY2(I)
DIST=SQRT.(DX*DX+DY*DY)
SEC=ATN1.(-D —DY)/0.48481368E-5+0.5
STA(1)=T02(I)
X(1)=X+DX
Y(1)=TOY2(I)
TER=CHECK (THRU2)
WHENEVER TER
ERR=ERR+1
TER=OB
BIGLP2 END OF CONDITIONAL
WHENEVER ERR.E.0
DUE2(K)=0
PRINT FORMAT OKPROB
TRANSFER TO RK2
END OF CONDITIONAL
DUE2(K)=DUE2(K)+ERR
TEST2 WHENEVER DUE2(K)*G.15*DUE2(K)=15
PRINT FORMAT RESMITDUE2(K)
RK2 RUN2(K)=RUN2(K)+1
TRANSFER TO SRCH
THRU2 DUE2(K)=DUE2(K)*2+1+ERR-I
TRANSFER TO TEST2
C(3) PRINT FORMAT TITLECPROBDUE3(K),RUN3(K)
THROUGH BIGLP3# FOR I=1,lI.G.*DUE3(K)
STA=FROM3(I)
STA(1)=INT3(I)
STA(2)=T03(I)
X=FRX3(I)
Y=FRY3(I)
X(2)=TOX3-(I)
Y(2)=TOY3(I)
D1=FD3( I)
D2=TD3(I)
M1=FM3(I)+SCN/100
M2=TM3(I)+SCN/100
S1=FS3(I)+SCN-100*(SCN/100)
S2=TS3(I)+SCN-100*(SCN/100)
SEC=3600*D1+60*Ml+S1
SEC(1)=3600*D2+60*M2+S2
RA=.48481368E-5*SEC
RB=.48481368E-5*SEC(1)
SA=-SIN.(RA)
SB=-SIN.(RB)
CA=-COS.(RA)
CB=-COS.(RB)
DX=X(2)-X
DY=Y(2)-Y
DA=(DX*CB-DY*SB)/SIN.(RA-RB)
X(1)=X+DA*SA
Y(1)=Y+DA*CA
DB=(DY-DA*CA)/CB
DIST=DA
DIST(1)=DB
TER=CHECK.(THRU3)
WHENEVER TER
ERR=ERR+1
TER=OB
BIGLP3 END OF CONDITIONAL
WHENEVER ERR.E.O
DUE3(K)=O
PRINT FORMAT OK,PROB
TRANSFER TO RK3
END OF CONDITIONAL
DUE3(K)=DUE3(K)+ERR
-C29 -
</P>
<P><PB REF="00000030.tif" SEQ="00000030" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="817" N="30">
The Computer as a Teaching Machine
MAD Program, Continued
TEST3 WHENEVER DUE3(K)J.G15#DUE3(K)=15
PRINT FORMAT RESMITDUE3(K)
RK3 RUN3(K)=RUN3(K)+1
TRANSFER TO SRCH
THRU3 DUE3(K)=2*DUE3(K)+1+ERR-I
TRANSFER TO TEST3
C(4) PRINT FORMAT TITLECPROBDUE4(K),RUN4(K)
THROUGH BIGLP4~ FOR I1,1,I.G.DUE4(K)
STA=FROM4(I)
STA(1)=INT4(I)
STA(2)=T04(I)
X=FRX4(I)
Y=FRY4(I)
X(2)=TOX4(I)
Y(2)=TOY4(I)
D1=FD4(I)
M1=FM4(I)+SCN/100
S1=FS4(I)+SCN-100*(SCN/100)
SEC=3600*D1+60*Ml+S1
WHENEVER DT04(I).L.O.
DIST(1)=DT04(I)-SCN/1Q000
OTHERWISE
DIST(1)=DT04(I)+SCN/1000.
END OF CONDITIONAL
DX=X(2)-X
DY=Y(2)-Y
RA=.48481368E-5*SEC
SA=-SIN.(RA)
CA=-COS.(RA)
Q=DX*SA+DY*CA
H=DX*CA-DY*SA
F=SQRT.(DIST(1)*DIST(1)-H*H)
WHENEVER (DX*DX+DY*DY-DIST(1)*DIST(1) ).L*.0OR.DIST(1).L.*O
DIST=Q+F
OTHERWISE
DIST=Q-F
END OF CONDITIONAL
X(1)=X+SA*DIST
Y(1)=Y+CA*DIST
RB=ATN1l(X(1)-X(2).Y(1)-Y(2))
SEC(1)=RB/*8448i368L-5+.5
TER=CHECK*(THRU4)
WHENEVER TER
ERR=ERR+1
TER=OB
BIGLP4 END OF CONDITIONAL
WHENEVER ERR.E.O
DUE4(K)=O
PRINT FORMAT OK,PROB
TRANSFER TO RK4
END OF CONDITIONAL
DUE4(K)=DUE4(K)+ERR
TEST4 WHENEVER DUE4(K).G.15#DUE4(K)=15
PRINT FORMAT RESMITPDUE4(K)
RK4 RUN4(K)=RUN4(K)+1
TRANSFER TO SRCH
THRU4 DUE4(K)=2*DUE4(K)+1+ERR-I
TRANSFER TO TEST4
C(5) PRINT FORMAT TITLECPROBDUE5(K),RUN5(K)
THROUGH BIGLP5, FOR I=1,1,I*G.DUE5(K)
STA=FROM5(I)
STA(1)=INT5(I)
STA(2)=TO5(I)
X=FRX5(I)
Y=FRY5(I)
X(2)=TOX5(I)
Y(2)=TOY5(I)
DIST=DFM5(I)+SCN/1000.
DIST(1)=DTO5( )+SCN/1000.
DX=X(2)-X
DY=Y(2)-Y
DC=SQRT.(DX*DX+DY*DY)
F=(DC*DC+DIST*DIST-DIST(1)*DIST(1))/(2.*DC)
H=SQRT.(DIST*DIST-F*F)
S3=DX/DC
C3=DY/DC
X(1)=X+F*S3-H*C3
Y(1 ) Y+F*C3+H*53
RA=ATN1.(X-X(1) Y-Y(1))
-C30 -
</P>
<P><PB REF="00000031.tif" SEQ="00000031" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="815" N="31">
Example Problem No. 76
MAD Program, Continued
RB=ATN1.(X(1)-X(2),Y(1)-Y(2))
SA=-SIN.(RA)
SB=-SIN.(RB)
CA=-COS.(RA)
CB=-COS.(RB)
SEC=RA/ 48481368E-5+.5
SEC(1)=RB/.48481368E-5+.5
TER=CHECK.(THRU5)
WHENEVER TER
ERR=ERR+1
TEK=OB
BIGLP5 END OF CONDITIONAL
wHENEVER ERR.E.O
DUE5(K}=O
PRINT FORMAT OK,PROB
TRANSFER TO RK5
END OF CONDITIONAL
DUE5(K)=DUE5(K)+ERR
TEST5 WHENEVER DUE5(K).G.15,DUE5(K)=15
PRINT FORMAT RESMITtDUES(K)
RK5 RUN5(K)=RUN5(K)+1
TRANSFER TO SRCH
THRU5 DUE5(K)=2*DUE5(K)+l+ERR-I
TRANSFER TO TEST5
LAST PRINT FORMAT PDATEDATE,DATE1,DATE2
THROUGH DONE. FOR I=1,1,I.E.N
PRINT FORMAT PR1,NAME(I),NAME1(I),NAME2(I),CN(I),DUE1(I).DUE2
1(I),DUE3(I),DUE4(I),DUE5(I),RUN1(I),RUN2(I),RUN3(I).RUN4(I).
2RUN5(I)
DONE PUNCH FORMAT PCHNAME(I)tNAME1(I),NAME2(I CNCN(I)DUE1(I),DUE2
1(I),DUE3(I),DUE4(I)~DUE5(I)tRUN1(I)tRUN2(I)gRUN3(I)~RUN4(I),
2RUN5(I)
EXECUTE SYSTEM.
INTERNAL FUNCTION (HERE)
STATEMENT LABEL HERE
BOOLEAN FLAG
ENTRY TO CHECK.
CK=O
FLAG=OB
NWCD READ FORMAT NEWCADEQUALS,ST(1),ST(3),ST(5)
WHENEVER EQUALS.E.$=$,TRANSFER TO HERE
WHENEVER ST(1).E.STA(CK)
ST(2)=$ $
OTHERWISE
ST(2)=$***$
FLAG=lB
END OF CONDITIONAL
WHENEVER.ABS.(FST(3)-X(CK)).G..003
ST(4)=$***$
FLAG=1B
OTHERWISE
ST(4)=$ $
END OF CONDITIONAL
WHENEVER.ABS.(FST(5)-Y(CK)).G..003
ST(6)=$***$
FLAG=1B
OTHERWISE
ST(6)=$ $
END OF CONDITIONAL
PRINT FORMAT OUT2,ST(1)...ST(6)
WHENEVER CK.El.1AND.PROBI.L.3TRANSFER TO FINI
WHENEVER CK.LE.1
READ FORMAT NEXCAD,EQUALS,ST(7),ST(9)...ST(11)
WHENEVER EQUALStE.$=$,TRANSFER TO HERE
WHENEVER.ABS.(FST(7)-DIST(CK)).G..003
ST(8)=$***$
FLAG=1B
OTHERWISE
ST(8)=$ $
END OF CONDITIONAL
BIGSEC=3600*ST(9)+60*ST(10)+ST(11)
-C31 -
</P>
<P><PB REF="00000032.tif" SEQ="00000032" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="787" N="32">
The Computer as a Teaching Machine
MAD Program, Continued
WHENEVER.ABS*(SEC(CK)-BIGSEC) G.3
ST(12)=$***$
FLAG=lB
OTHERWISF
ST(12)=$ $
END OF CONDITIONAL
PRINT FORMAT OUT3,ST(7).**ST(12)
CK=CK+1
TRANSFER TO NWCD
END OF CONDITIONAL
FINI FUNCTION RETURN FLAG
END OF FUNCTION
EQUIVALENCE (FSTST)
BOOLEAN TER
BOOLEAN CHECK.
INTEGER STDATE.DATE1tNEQUALSNAMENAME1,NAME2,CNDUE3~DUE4,
1DUE5,RUN3,RUN4,RUN54JFROM3,FD3,FM3,FS3,INT3,TD4,TM3*TS3,T03,
2FROM4tFD4,FM4tFS4,INT4,T04,FROM5,INT5 T05 ERRSNAMESCNPROB,
3KSTAIDl1D2oMl1M21SlS2,SECoCKSTtBIGSEC
INTEGER DUElDUE2,RUN1,RUN2oFROMTODMS$FROM2,T02TDATE2
DIMENSION NAME(100),NAME1(100),NAME2(100),CN(100),DUE3(1 ),D
1UE4(100)tDUE5(100)~RUN3(100),RUN4(100),RUN5(100)tFROM3(25),FR
2X3(25) FRY3(25),FD3(25),FM3(25),FS3(25),INT3(25),TD3(25)
3TM3(25),TS3(25)*TO3(25) TOX3(25) TOY3(25)DFROM4(25),FRX4(25)~
4FRY4(25),FD4(25),FS4(25)tFM4(25),INT4(25),DT04(25),T04(25),
5TOX4(25),TOY4(25),FROM5(25),FRX5(25),FRY5(25),DFM5(25),INT5(2
65),DT05(25), T05(25)tTOX5(25)~TOY5(25)S5NAME(2)~STA(2),X(2),
7Y(2),DIST(1),SEC(1)ST(12),FST(12)
8,DUE1(100),DUE2(100) RUN1(100)RUN2(100)
DIMENSION FROM(25),TO(25),FRX(25),FRY(25),DTO(25),D(25),M(25)
1S(25),FROM2(25),FRX2(25),FRY2(25)~T02(25)~TOX2(25)>TOY2(25)
VECTOR VALUES RDATE=$3C6*$
VECTOR VALUES TITLEC=$13H-INTERSECTIONI3,I6,26H PROBLEMS DUE
1 FOR RUN NO-I3//*$
VECTOR VALUES PRI=$1H S,20t3C6,S6 tI4,S4,10I6 *$
VECTOR VALUES RDNAME=$C1,3C6,I,5S14l1*$
VECTOR VALUES TITLEB=$1H1S5,3C6,S5*I4,S5,7HPROBLEMI3,S52C6*$
VECTOR VALUES SEARCH=$C1*$
VECTOR VALUES OK=$1HOS5,45HYOU HAVE DONE ALL THE PROBLEMS FOR
1 ASSIGNMENT I4,13H CORRECTLY./45HOSEE YOUR INSTRUCTOR FOR
2FURTHER ASSIGNMENTS.*$
VECTOR VALUES RESMIT=$1HO5555HYOU HAVE AT LEAST ONE ERROR (I
1NDICATED BY ASTERISKS). /34HOCORRECT THE DATA AND RESUBMIT T
2HEI4,20H PROBLEMS NOW DUE.*$
VECTOR VALUES PCH=$S1,3C69151I2,9I3*$
VECTOR VALUES RD=$(iA2I4~3F10.3,S55313*$
VECTOR VALUES RD1=$C1*3C6t151I2.9I3*$
VECTOR VALUES RD2=$C1,Z(15v2F10.3)*$
VECTOR VALUES RD3=$Cl13.2F,.3I5,213,2I5t213,15,2F9.3*$
VECTOR VALUES RD4=$Ct113t2F9.3,I5,213,A,FlO.3>I6t2F9.3*$
VECTOR VALUES RD5=$Cl1I3t3F9.316,FlO.316,2F9.3*$
VECTOR VALUES NEXCAD=$C1,S25tF9.3#16,2I3*$
VECTOR VALUES NEWCAD=$C1,I4,2F10.3*$
VECTOR VALUES OUT2=$1H S55I4,C3,F10.3C3,F10.3sC3*$
VECTOR VALUES OUT3=$1H S37,F10.3,C3,I5,213IC3 *$
VECTOR VALUES PDATE=$1HltS19,65H PLANE GEOMETRY TEACHING
1MACHINE MODEL 3-62 PROGRAM 9C8 3C6//520,12HSTUDENT NAME,
258912HCLASS NUMBERS3,58HDUE1 DUE2 DUE3 DUE4 DUE5 RUN1
3RUN2 RUN3 RUN4 RUN5//*$
END OF PROGRAM
S DATA
-C32 -
</P>
<P><PB REF="00000033.tif" SEQ="00000033" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="817" N="33">
Example Problem No. 76
Data
The following is an abbreviated set of data which might be submitted to the Teaching
Machine Program. The class roll has been shortened to four names. Basic data for only one
problem type are shown and only two sets of student's results are shown.
$DATA
MAY 17, 1962
D J PIKE
J A FOX
H W ANDRES
S H HAMMOND
4116 3
1949 3
1144 3
27U3 3
3 3 3 3 1 1 1 1 1
3 3 3 3 1 1 1 1 1
3 3 33 1 1 1 1 1
3 3 33 1 1 1 1 1
Class Roll
147
148
161
233
168
167
234
99
85
18
19
21
45
36
46
89
254
162
234
227
228
161
26
29
15
20
16
41
44
1
5000.000
4569.374
3228.603
3325.440
3911.262
3803.7~6
3234.698
4896.323
4511.209
4303.953
3784.575
5110.979
4482.951
3539.884
4542.456
5000.000
4430.280
4458.202
4388.835
4533.392
4505.241
4344. 77
4957.408
4552.191
4984.691
4539.076
4 32.670
4623.0752
4589. 095
4538o.92
253.897
15'7 * 16
93.000
101. 16O
108. 887
103.256
114.2 8
1875.690
2165.140
908.394
7804284
Z89.450
099.568
543.067
65U. 158
341 1- 59
345 50 01
266 56 35
63 44 43
00 28 34
342 49 0U
176 56 34
157 39 59
171 17 46
324 45 23
272 27 17
157 27 17
145 12 50
272 54 44
224 43 41
Data for
First
Problem
Type
/
H W ANDRES
147 5000.000
1144
5000.000
TOO3A 1
255. 41 341 29 43
158.160 346 01 45
89 5080.946 4758.145
148 4569.374 4430.280
254 4607.958
161 3228.603
4276.799
4458.202
}
Data for
Other Four
Problem Types
Would be
Inserted Here
First Student's Solution
Seoond Student's Solution
94.144 267 08 19
162 3322.630 4462.996
D J PIKE
147 5000.000
89 5079.782
148 4569.374
254 4606.930
161 3228.603
162 3325.638
4118
5000.000
4574.630
4430.280
4273.564
a458. u2
4462.*16
T0c3P 1
258. 15 341 59 17
161.134 346 31 1i
97.113 267 37 53
-C33 -
</P>
<P><PB REF="00000034.tif" SEQ="00000034" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="806" N="34">
The Computer as a Teaching Machine
Computer Output
The following is the computer output which corresponds to the data shown in the
previous section.
H W ANDkES 1144 PRO1bLEM 1 MAY 17, 1962
INTERSECTION 1 3 PROirLEMS DUE FOR RUN NO. 1
147 5000.0O 500o0.uO0
255.041 341 29 43
89 50g0.946 47)8.145
148 4569.374 4430.280
158.160 346 1 45 Th
254 4607.5U 8 42 6. 99 tc
161 3228.603 4458.202
94.144 267 8 19
162 3322.630 4462.996***
YOU HAVE AT LEAST ONE ERROR (INOICATED BY ASTERISKS).
CORRECT THE DATA AND RESUBMIT THE 4 PROBLEMS NOW DUE.
D J PIKE 4118 PROBLEM 1 _ __ MAY 17, 1962
INTERSECTION 1 3 PROBLEMS DUE FOR RUJN NO. 1
147 5000.000 0.0o —o -
258.015 341 59 17 -
89 5079.782 4574.630***
148 4569.374 4430.280
161.134 346 31 19
254 46C6.930 4273.584
161 3228.603 4458.202
97.118 267 37 53
162 3325.638 4462.216. YOU HAVE AF LEAST ONE ERROR (INDICATED BY ASTERISKS).
CORRECT [HE DATA AND RLSUHtMIT THE 4 PROBLEMS NOW DUE.
This summary is printed for use by the instructor:
INVERSE CUMPUTATION TEACHING MACHINE MODEL 3-62 PROGRAM 9C8 MAY 17, 196
STUDEN NAME CLASS NUMBER DUE1 DUE2 DUE3 DUE4 DUE5 RUN1 RUiN2
D J PIKE 4118 4 3 3 3 3 2 1
J A F-'X 1949.3 3 3 3 3 1 1
H W ANDRES 1144 4 3 3 3 3 2 l
S H HAMMOND 2703 3 3 3 3 3 1 1
Number due on next run Run numb
for each problem --- ype —. f-or each
lis is returned
) first student.
iis is returned
D second student.
2
RUN3 RUN4 RUN 5
1 1 1 __
1 1 1
1 1 _ 1
1 1 1
er of next run
problem type. -
This is the updated class roll which is used as input on the next run:
UJ PIKe 411 4 3 3 3 3 2 1 1 1 1
J A FOX 1949 3 3 3 3 3 1 1 1 1 1
H W ANDRES 1144 4 3 3 3 3 2 1 1 1 1
S H HAMMOND 270 3 3 3 3 3 1 1 1 1 1
Critique
A complete discussion of the effectiveness of the use of the Teaching Machine Program is
given in the section "Introduction of Civil Engineering Students to the Use of the Computer,"
in the first part of this report.
-C34 -
</P>
<P><PB REF="00000035.tif" SEQ="00000035" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="869" N="35">
Example Problem No. 77
OPTIMUM REGULATION OF A RESERVOIR FOR POWER PRODUCTION
by
J. W. Howe
Department of Mechanics and Hydraulics
State University of Iowa
Course: Water Power Engineering II Credit Hours: 1-3 Level: Graduate
Statement of the Problem
It is desired to determine the best way of operating the turbines at a reservoir in
order to produce the most power in the period studied. The following physical relations are
assumed:
Reservoir volume = (Y)3'5/4 (acre feet)
Spillway discharge = 120,000 h3/2 (acre feet/month)
Turbine discharge = 12,000 H1/2 (acre feet/month)
Tailwater discharge = 21,000 Z3/2 (acre feet/month)
Reservoir spillway level at Y = 63 feet
Reservoir inflows: 12, 16, 90, 210, 47, 44, 27, 19, 45, 64, 31, 30, 10, 50, 199,
137, 59, 231, 307, 239, 75, 50, 34, 34 (thousand acre feet/month)
Determine the average power output per year if the turbines are shut down for various
percentages of time when the upstream water level is below spillway level.
Solution
A cross-sectional diagram of the reservoir is shown below to illustrate the nomenclature.
i |zH
i Q TW
QTURB Z
-C35 -
</P>
<P><PB REF="00000036.tif" SEQ="00000036" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="878" N="36">
Optimum Regulation of a Reservoir for Power Production
Initial levels Yl and Z1 at beginning of month.
Final levels YAV and Z2 at end of month.
VOLYAV = YAV3 5
VOLO = 633'5
VOLYI = YI3'5
QSTOR = (VOLYAV-VOLY1)/4
YAV = (YMAX + YMIN)/2 = level at end of month
AVHD = (YI + YAV)/2
STORAGE EQUATION:
INFLOW - QSTOR = OUTFLOW = QSPLWY + QTURB = QTW
INFLOW = Q(I) = QTW + QSTOR = QTOT (1)
Initial values of YMIN and YMAX of 20 and 68 feet are assumed and YAV determined.
Using YAV in determining terms of Eq. (1), if 1000 Q(I) > QTOT, Y should be increased. This
is accomplished by letting YAV replace YMIN. In the opposite situation, YAV becomes YMAX.
In 15 iterations YAV is within a few thousandths of the value it must have to satisfy Eq. (1).
It is assumed that the monthly discharge through the turbines or over the spillway can be
secured by using the average of the initial and final water levels.
In the case of a month having some spillway discharge in the month but not during the
entire month, due to rising or falling water levels, the spillway discharge for a full month
is multiplied by:
(VOLYAV - VOLO)/(VOLYAV - VOLYI) for a rising stage and by
(VOLYI - VOLO)/(VOLY1 - VOLYAV) for a falling stage
QZ+Z2 0.746X0.80 43560 Zi+Z2
Power production/month = QTURB (AVHD - 2 )x 7468 80 X 456 =.82 QTURB (AVHD Z-
List of Symbols
Q(I) Monthly flow, thousands of acre feet/month
Yl Upstream depth at beginning of month, feet
Y2 Upstream depth at end of month, feet
Z1 Downstream depth at beginning of month, feet
Z2 Downstream depth at end of month, feet
YMIN Minimum upstream depth at end of month, feet
YMAX Maximum upstream depth at end of month, feet
YAV Average upstream depth at end of month, feet
AVHD Average upstream depth during month, feet
-C36 -
</P>
<P><PB REF="00000037.tif" SEQ="00000037" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="667" N="37">
Example Problem No. 77
Lst of Symbols, Continued
PURB Monthly flow through turbine, acre feet/month
3TOR Volume change in storage in one month, acre feet
3PLWY Monthly spillway discharge, acre feet/month
DW Monthly outflow, acre feet/month
(I) Monthly power output, KWH/month
\POWN Annual power output, KWH/year
JLT Fraction of time turbines are operated when upstream water level is below
spillway level, dimensionless
Low Diagram
/^ ^ ( -- -- - / -. EXECUTE ZE1o o
ODIN( Q... (.Q(24) ^ UT tMV T< z 2 SUKC)-S<. - L I, \,p L;,>u(I
zz, s
TOTPOW=-O U(C,)...~Cu(4)
\^^ = 4 ______3__:T:, CUNT. ( = 1 "' 1. rO: V6Y.Y- 3$T. VOLYAV-VOLY1(
IK CO\NTrI,,CGUNT),, 2. yNT 'AOOOVV2 D W R3. *
-/ YMAW ^;D Y0 + ^ YAYV ________~~v,i:: Y1 +V
2.
-C37 -
</P>
<P><PB REF="00000038.tif" SEQ="00000038" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="730" N="38">
Optimum Regulation of a Reservoir for Power Production
Flow Diagram, Continued
KE \ 1 K>) Fv t AWM2K-TR(V AWMIN( o P(l-O. 20-QTaR(W-: YZ =2
SUM(H)- SU~~(trl tUM().-CVU<+
MULT
(7)-^ <?f ) ---^~NWANPTO T WN VOW - sum(.IH ^
MAD Program and Data
J. We HOWE S164A 005 010 000 STOREG
J. W. HOWE S164A 005 010 000 STOREG
SCOMPILE MAD EXECUTEDUMP
R
ROPTIMUM REGULATION OF A POWER DEVELOPMENT AT A RESERVOIR
RASSUMED RELATIONSHIPS, VOLUMNE OF RESERVOIR~Y.Pe3S5/4.
RQTW=21000*ZP* 1 5 OQTURB=12000*SORT (H)
RQSPLWYz120000*(Y-63) P* 1 5 POWER*0*H*e820 KWH/MO
RQUANTITIES IN ACRE FEET/MO. SPILLWAY CREST AT Y~63FT.
R
DIMENSION Q(24). SUM(4),CUM(4),P(24)
INTEGER Q#COUNT#ItJ4K#SUMMCUM
PRINT COMMENT S15
READ FORMAT DATA Q(1) * Q(24)
VECTOR VALUES DATA-$(1216)*S
PRINT FORMAT ECHO*Q(1)t l Q(24)
VECTOR VALUES ECHOIS1H 818*S
THROUGH OMEGA. FOR MULT'l*-lt1 MULT.LT.e.2
Y1l63.
Z1l2.64
Z2'3,
TOTPOW O
EXECUTE ZERO.(SUM( 1 ) * SUM(4 ).CUM( 1...CUM(4)
THROUGH ALPHAJ FOR I"11pIe. G24
YMIN"20
YMAXt68
THROUGH BETA*FOR COUNT"1, 1 COUNToG 15
YAV'tYMIN+YMAX)/2
AVHD (Y1+YAV)/2.
QTURB=12000**SQRT (AVHD-( Z1+Z2 ) /2 )
VOLYAV=YAV~P~3*5
VOLO63. P*3*5
VOLYV1Yl1Pe3e5
QSTOR (VOLYAV-VOLY1) /4
WHENEVER YAV GE.634
WHENEVER AVHD.Ga63
QSPLWYu120000*(AVHD-63..P.1.5
OTHERWISE
QSPLWY 120000 * YAV-63 )/2 ) P 1. 5* ( VOLYAV-VOLO )tVOLYAV1VOLY1)
END OF CONDITIONAL
QTW-QSPLWY+QTURB
OTHERWISE
WHENEVER AVHD.G~63o
QSPLWYr120000* ((Y1-63 ) /2 ) P 15*(VOLY1-VOLO) / (VOL1-VOLYAV
1)
QTW-QSPLWY+QTURB
OTHERWISE
OSPLWY-O
QTUR BQTURB*MULT
QTW-QTURB
END OF CONDITIONAL
END OF CONDITIONAL
QTOT-QTW+QSTOR
WHENEVER 1000*Q(I).G.QTOT
YMIN-YAV
OTHERWISE
-C38 -
</P>
<P><PB REF="00000039.tif" SEQ="00000039" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="792" N="39">
Example Problem No. 77
MAD Program and Data, Continued
YMAX=YAV
END OF CONDITIONAL
BETA Z2 (QTW/21000. ) P.667
PRINT FORMAT CHECK YAVtQTURBQSPLWYtQTOT
VECTOR VALUES CHECK$1HOFl0*2*3F100*$
VECTOR VALUES HWMAX(1)=65.64.564.,63.
VECTOR VALUES HWMIN(1)=56*,58*,60,61*
THROUGH ADD, FOR KltlvK*G.4
WHENEVER YAVEGE*HWMAX(K)
SUM( K)SUM ( K )+
OR WHENEVER YAV.LE*HWMIN(K)
CUM(K)=CUM(K)+1
ADD END OF CONDITIONAL
P(I) lQTURB*(AVHD(Z1 +Z2)/2 **820
Y1=YAV
Zl=Z2
ALPHA TOTPOW=TOTPOW+P (I)
ANPOWN TOTPOW*12 /24
OMEGA PRINT FORMAT POWRtMULTtANPOWNtSUM( 1)...SUM 4) CUM(1) **CUM(4)
VECTOR VALUES POWR$1lHO06H MULTFS5.2,12H ANNUAL KWH,Flll/
117H MONTHS ABOVE 65,I3/19H MONTHS ABOVE 6-45,I3/
217H MONTHS ABOVE 64,I3/17H MONTHS ABOVE 63I13/
317H MONTHS BELOW 56,I3/17H MONTHS BELOW 58,I3/
417H MONTHS BELOW 6013/17H MONTHS BELOW 61,I3*$
END OF PROGRAM
$ DATA
012 016 090 210 047 044 027 019 045 064 031 030
010 050 199 137 059 231 307 239 075 050 034 034
Computer
12
45
59
YAV
59.91
56.68
56.76
61,81
60*05
58*08
55,18
51.53
49.12
47.90
4-407
39.49
30.28
27.21
45 20
Output
16
64
231
92015
89510
88262
90324
91600
90126
88189
85512
82963
81411
79182
75318
68527
61909
70035
90 210
31 30
307 239
0
0
0
0
0
47
10
75
11997
15981
90001
209970
47029
44001
26972
18985
45002
64014
31008
30002
9993
49998
199002
44 27 19
50 199 137
50 34 34
I %
-C39 -
</P>
<P><PB REF="00000040.tif" SEQ="00000040" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="837" N="40">
Optimum Regulation of a Reservoir for Power Production
Computer Output, Continued
49040 80462
47T*4 81526
46s9s 84506
6444 9:1118 14'
6C9S 93039 1 9;
46t42 92167 26'
60*61 91934
58646 90575
55.94 a643
NrLT* 1.00 AWIIIAL KWH 424
MONTHS ABOVE 65 0
MONTHS ABOVE 64*5 0
QMNTHS ABOVE 64 1
MOSTHS ABOVE 63 2
MOTHS BELOW 56 12
MOTHS BELOW 58 15
O!HTHS BELOW 60 18
MONTHS BELOW 61 20
0 3
0
0;
785
229
0
0
0
15296 ~ 0
137016
59005
231023
107008
239057
74969
49989
34027
34014
Discussion of Results
Only two years (24 months) of data were used which required 2.2 and 0.4 minutes, respectively, of execution and compiling time. Because of space limitations, the printed output
from the computer is shown for only one percentage of turbine operation (i.e., the case of
full-time turbine operation). The results indicate the desirability of decreasing the range
of turbine operation examined and increasing the number of months studied. Eight years of
record could be studied in five minutes of machine time with this modification.
The results revealed an unexpected advantage to an operation involving considerable
drafting of the reservoir even at substantial reductions of head since such operation
produced the greatest amount of energy conversion.
Critique
This problem has been used in a second course in Water Power Engineering but, because of
the complexity of numerical work,,full solutions have not been required previously. Since the
computer makes a trial and error solution for the upstream and downstream water levels feasible,
the problem can now be fully solved and will be assigned to students in the future.
-C40 -
</P>
<P><PB REF="00000041.tif" SEQ="00000041" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="873" N="41">
Example Problem No. 78
ECONOMICAL SIZE OF CANAL FEEDING HYDRO-POWER PLANT
by
J. W. Howe
Department of Mechanics and Hydraulics
State University of Iowa
Course: Water Power Engineering Credit Hours: 3 Level: Senior and Graduate
Statement of the Problem
A canal having 1:2 side slopes is to feed a hydro-power plant. It is to have best
hydraulic efficiency (R = D/2) when carrying the mean discharge of 1500 cubic feet per second,
but capacity for the maximum discharge of 2000 cubic feet per second with a 2-feet freeboard.
Excavation is limited to 60 per cent of the maximum depth with the excavated material forming
berms on either side to provide the upper portion of the cross-section. Transverse slope of
the ground is assumed to be zero but any longitudinal slope of the canal is possible by slight
divergence from a contour. The following data are given: Canal roughness coefficient,
n = 0.025; excavation cost, $2.00/cubic yard; annual charge for maintenance, interest, and
taxes, 10 per cent; cost of power, $0.006/KWH; operation schedule, 16 hours/day, 300 days/year.
Slopes varying from 0.00020 to 0.00046 are to be used. Determine the dimensions of the crosssection for which the annual cost is least.
Solution
A diagram of a typical canal cross-section is shown below to illustrate the important
dimensions.
A = B X Y + 2Y2
__^ ^^\. R R A/P BxY + 2Y2 Y
R = A/P = B + 4.4722 x Y =
1^.z~ | VY B = 0.4722 X Y
For a series of slopes varying from 0.00020 to 0.00046 by 0.00002 increments, the Manning
formula is solved by a trial and error process. Since the channel is to have best hydraulic
efficiency at mean discharge (R = Y/2):
By the Manning formula
= 15 A R2/3 S1/2
-
1.5 (BXY + 2Y2)5/3 /2
N (B + 4.4722 X Y)23
-C41 -
</P>
<P><PB REF="00000042.tif" SEQ="00000042" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="843" N="42">
Economical Size of Canal Feeding Hydro-Power Plant
Substituting B = 0.4722Y
Q 1.5 (0.4722 y2 + 2y2)5/3 1/2 2.472 2)/
_(o.472 S (2.42
N (0.4722Y + 4.4722Y)2/3 (4.944 y)23
y = Q N)3/8 B Y 0.4722 (Q N)3/8
1.375 S3 1 1.375 S3/16
Hence, there is a particular bottom width for each slope which satisfies the condition R = Y/2
for the mean discharge.
The required depth for maximum discharge at each slope and corresponding bottom width
are next found by trial and error starting with given minimum and maximum values of Y and
using their mean as the depth. After 11 iterations the mean should have the correct value.
QMAX 1.5 [B X YAV + 2 x (YAV)2] 5/3 1/2
N (B + 4.4722 YAV)2/3
L = (QMAX X N)/1.5 X S1/2
R = [B X YAV + 2(YAV)2] 5/3 / (B + 4.4722 YAV)2/3
Whenever L > R Whenever L< R
YMIN = YAV YMAX = YAV
Depth = 0.6 x (Final value of YMAX + 2)
Annual excavation cost/foot length = C X INT X (Vol. of Exc.
C X INT X [B X 0.6 DEPTH + 2(0.6 DEPTH)2]
Ann. Cost of Power loss/foot length = QMEAN X S X H X D X 1 X ~74
= MEAN X S X H X D X
= QMEAN X S X H X D X ~475
Total Annual Cost = Annual chargefbr Excavation plus Annual Cost of Power Loss
List of Symbols
QMEAN average discharge rate, ft3/sec
QMAX maximum discharge rate, ft3/sec
A cross-sectional area, ft2
R hydraulic radius, ft
N roughness coefficient n in Manning formula
S slope of bottom in direction of flow, dimensionless
B bottom width of channel, ft
YMIN lowest depth of water, ft
YAV mean depth of water, ft
-C42 -
</P>
<P><PB REF="00000043.tif" SEQ="00000043" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="633" N="43">
Example Problem No. 78
List of Symbols, Continued
YMAX highest depth of water, ft
DEPTH total depth of canal with 2-ft freeboard
C cost of excavation, $/yd3
INT annual fixed charges (interest, taxes, maintenance), dimensionless
P value of power, cents/KWH
H hours of operation/day, hrs
D days of operation/yr, days
ANCSXC annual cost of canal, $/yr
ANCSPR annual cost of power lost, $/yr
ANCOST total annual cost, $/yr
Flow Diagram
V(PH De P r H (L) - lvl\\X2 - C^2.N r )AC. ^'^.^ -"^ ^ DePT*(x)EPT) A)bCSiR G ) DN 0 A P7.)
PPN(SA- 2 7 475
0 - ~ "' H() I\ ft ~ _c Irrri).bDFTH().2 D~r~r)) rI i () I XDr
0(
c- (f) ( A h.5Csr
ANCO'.)" =T q5&
&(s) S -1 5(14)
E', rI) - ( DP(' 4) (
SLO PES A >ANC>A(J) W(..AL14) 1 RT
\LU i/~S ~ANC& P-)" -P!F-itC'\ 4)
(ANCS Pi().. ~ NcPRL(, )
T,, — /'''
-C43 -
</P>
<P><PB REF="00000044.tif" SEQ="00000044" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="767" N="44">
Economical Size of Canal Feeding Hydro-Power Plant
MAD Program and Data
Jo We HOWE S164A 001 005 000 CANAL
J. We HOWE S164A 001 005 000 CANAL
S COMPILE MAD, EXECUTE, DUMPt PRINT OBJECT
R
RCANAL IN EARTH WITH 1/2SIDE SLOPES. MANNING ROUGHNESS N, BEST
RHYDRAYLIC RADIUS (0.4722*Y FOR QMEAN),OMEAN AND QMAX IN CFS.
REXCAVATED DEPTH 0.6(YMAX+2). COST OF EXCAVATION C $/CUYDe
RFIXED CHARGE INT, POWER COST P CENTS/KWH, OPERATION H HOURS/
RDAY, D DAYS/YR~ ANY SLOPE AVAILABLE BY DIVERGENCE FROM
RCONTOUR. TRANSVERSE SLOPE ASSUMED ZERO*
R
DIMENSION S(14),B(14) ANCSXC(14),ANCSPR(14) ANCOST(14),
1DEPTH(14)
INTEGER IJ
PRINT COMMENT $1$
READ FORMAT DATAQMEANOMAXNCINTPHD*S(1).*.S(14)
VECTOR VALUES DATA=$2F1O.OFlO~.3,3F8*22F8.0/( 7F105)*S
PRINT FORMAT ECHOQMEANQMAX.NCINTPHDS( 1)**.4S( 14)
VECTOR VALUES ECHO = $1H 2F10s0PF1O.3s3F82t22F8.0/
1 (1H07OFIO5)*$
THROUGH SLOPESFOR Isl, l IG,14
B(I) =.4722*( MEAN*N )P 375/( 1375*S(! ) *P 187'5)
YMIN=11~
YMAXuS'18
THROUGH TRIALSFFOR J="11.tJ.G.11
YAV=tYM+YN+YMAX)/2
Ls(QMAX*N)/( 15*SQRT (S( I )
Ra(B I)*YAV+2*YAV*YAV) P ~ 1667/(B(I )+4 4722*YAV ),P *667
WHENEVER Le G.R
YMIN=YAV
OTHERWISE
YMAX=YAV
TRIALS END OF CONDITIONAL
DEPTH (I) -YMAX+2
ANCSXC(I)~C*INT*(B(I)*6+*DEPTH( I)+*72*DEPTH(I)*DEPTH( ))/27
ANCSPR(I)=QMEAN*S(I)*H*D*P/1475
SLOPES ANCOST(I) mANCSXC(I)+ANCSPR (I)
PRINT COMMENT $4$
PRINT FORMAT RESLTSS(1)*..S(14),B(1)*.*.B(14),DEPTH(1)*.4DEPT
1H(14) ANCSXC(1) *l ANCSXC(14) ANCSPR(1) **ANCSPR(14) ANCOST(l)
2 ~ ANCOST(14)
VECTOR VALUES RESLTS = S1H 7F10.5/1H07F10.5/(1HO7F10O2)*S
END OF PROGRAM
$ DATA
1500. 2000. *025 2.00 ~10 *60 16. 300.
00020 *00022 *00024.00026 *00028 @00030.00032
*00034 *00036 *00038.00040 *00042.00044.00046
-c44 -
</P>
<P><PB REF="00000045.tif" SEQ="00000045" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="801" N="45">
Example Problem No. 78
Computer Output
1500
000 20
*00034.00020.00034
6.60
5.98
17 74
16.25
2.20
1.84.59
1.00
2,79
2.84
2000.00022.00036.00022.00036
6.48
5.91
17.46
16.10
2.13
1.81
*64
1.05
2.77
2.86
*025.00024.00038.00024,00038
6.38
5.85
17.21
15.96
2.07
1,77
*70
1.11
2*77
2,89
2.00,10 *60 16 300,00026,00028 *00030,00032.00040,00042 *00044 00046
*00026
*00040
6.28
5.80
16.99
15.83
2.01
1074
*76
1.17
2.78
2,92. 00028.00042
6.20
5.74
16.78
15.70
1,96
1.72
*82
1.23
2.78
2.95.00030.00044
6.12
5.69
16.59
15,58
1.92
1.69.88
1.29
2.80
2.98.00032. 00046
6.04
5,65
16.42
15.47
1X88
1o66.94
1,35
2.82
3.01
}
}
)
Input Data
Slopes
Bottom
Widths
Depths
Annual Cost
of Canal
Annual Cost of
Power Lost
Total
Annual Cost
*
Discussion of Results
A plot is shown below of the total annual cost as a function of the canal slope. It is
seen that a minimum cost of $2.77 occurs for a slope of 0.00023.
E-4
Co
0
0
g
14I
$3.00
2.90
2.80
2.70
0.0002 0.0003
0.0004
CANAL SLOPE
The problem involves 0.3 minutes of compiling time and 0.7 minutes total time.
Critique
This is a problem regularly given to students in a Water Power course. Without the use
of the computer, the students can obtain a satisfactory but laborious solution using special
tables. Plotting of computed results is usually necessary for interpolation between points.
-C45 -
</P>
<P><PB REF="00000046.tif" SEQ="00000046" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="898" N="46">
Economical Size of Canal Feeding Hydro-Power Plant
The use of the computer for obtaining a trial and error solution replaces the tables, makes
it feasible to compute many more points and affords an easy way to secure solutions for
different original data with little additional work.
The problem, though short, illustrates the advantage of the computer approach to an otherwise laborious problem and provides a good exercise for students who are learning to program.
</P>
<P><PB REF="00000047.tif" SEQ="00000047" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="857" N="47">
Example Problem No. 79
MOST PROBABLE NUMBER METHOD FOR DETERMINING BACTERIAL POPULATIONS
by
L. L. Kempe
Sanitary Engineering
The University of Michigan
Course: Industrial Bacteriology Credit hours: 3 Level: Senior or
Graduate
Statement of Problem
Write a MAD program that will solve, using the Newton-Raphson method, the MPN (Most
Probable Number) equation. The input of the program will be N, ai, ni, Pi. The output
should be X, the most probable number of bacteria, and P, the probability associated with X.
Solution
The MPN equation is based on the extinction dilution method for finding the number of
viable bacteria in a solution. The development of the MPN equation is well documented ( )
and will only be reviewed here.
Assuming a negative exponential probability of growth, the probability of bacteria
growing in n tubes of culture medium innoculated with "a" milliliters of a solution containing X bacteria per mililiter is
nj -aX)p /-aX) n-p
p:(n-p)! (1 - ea ) (ea)n
If several different dilutions are used, the probability of having Pi tubes show growth
out of a total of n. tubes innoculated with a. ml. each is
n n ni.p
i=l Pi!n ei-p i) (1)
This gives an expression for the probability of ni out of pi tubes showing growth when the
solution contains X bacteria per ml. Now an expression will be found for the most probable X
given a set of ai, ni, Pi. To do this, it is noted that if P (the probability) is a maximum,
then In P will also be a maximum. Taking the natural logarithm of (1)
ln p = i l [ln (p i-Pi - (ni-Pi)aiX + Pi ln (l-e-aiX. (2)
Thus
ln p= ai (nia -ni) =0, (3)
dX i=l i a)X i
</P>
<P><PB REF="00000048.tif" SEQ="00000048" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="887" N="48">
Method for Determining Bacterial Populations
which is an implicit equation in X and gives the most probable number (MPN) of bacteria per
ml. of a solution used in an extinction dilution.
The approach to the problem solution is to use the Newton-Raphson method to find the X
that satisfies the MPN equation. Then, given this X, the probability associated with it is
calculated using equation (2). The probability will give a measure of the reliability of the
experiment. A low probability would say that the experimental results were not reliable.
Several extra features were incorporated in the solution. These include the ability
to do a single calculation or a table of calculations. An automatic estimation of the initial
value of X is also included. Any number of dilutions from 3 to 9 may be used for input data.
The output format is automatically adjusted to the number of dilutions.
This program is much more elegant than required, but illustrates the power of a reasonably
short computer program. Further information is included as remarks in the program.
List of Symbols
Problem Notation Program
N number of different dilutions used NBR
ai number of ml. used for innoculation A(I)
ni total number of tubesinnoculated N(I)
Pi number of tubes innoculated which show growth P(I)
X most probable number of bacteria X
P(%) probability associated with X PPER
Important additional variables:
G = MPN function
GP dG
dX
DELX = Ax
Flow Diagram
-c48 -
</P>
<P><PB REF="00000049.tif" SEQ="00000049" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="738" N="49">
Example Problem No. 79
Flow Diagram (continued)
P R I "T
X \PPEK= (e- oo; EX P i L L L (,'LI T ' - p, oe ) ) /
- pp re
MAD Program and Data
SCOMPILE MADePRINT OBJECTPUNCH OBJECT*EXECTEDUW4P ' PN 001
R
PROGRAM FOR DETERMINING BACTERIAL POPULATION
RIJING THE MPN EQUATION
R
R
RCODED BY W.J. SANDERS 4/3/62
R
R
R
R THIS PROGRAM USES A METHOD DEVELPED BY HALVORSOI AND
RZIEGLER (Jo BACT.,25*101-121(1933)) AND PROGRAMMED FOR THE
RIBM 650 BY NORMAN AND KEMPE (BIOTECHNOLOGY AND BIOENGINEERING
RII-2,l57-1-3(1960))e* IT INVOLVES THE COMPUTATION OF THE
RiOST PROBABLE NUMBER OF VIABLE BACTERIA IN A SAMPLE OF
RMATERIAL. INPUT DATA ARE THE RESULTS OF AN EXTINCTION
RDILUTION EXPERIMENT AND ARE ARRANGED AS FOLLOWS
R
R
RCARD TYPE I
R
R COL
R
R 1-4
R 5-7
R 8-10
R 11-14
R 15-22
R 23.24
R 41 42
R
R
RCARD TYPE II
R
VAR
CHECK
HOW
NBR
RAT
Al
P(i)
N(1)
All)
25.26 P(2) *~. 39,40 P(9)
43,44 N(2) *.. 57,58 N(9)
R
R
1-8
9-16
A(2) ** 65-72 A(9)
-C49 -
</P>
<P><PB REF="00000050.tif" SEQ="00000050" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="753" N="50">
Method for Determining Bacterial Populations
MAD Program and Data (Continued)
R
R.Rt~ST NOTAT ION FOLLOWS THAT OF HALVORSON AND ZIEGLER.
R NBR - NUMBER OF DILUTIONS
R gN(! NUMBER OF SAMPLES OF VOLUME A(I) ML. (.l Is.t1,*. NBR)
R P(I NUMBER OF FAILURES OUT OF N(I) TRIALS* FAILURE
R MEANS G6ROWTH IN TUBES
R
R
R
RTHERE ARE TWO METHODS OF IN PUT OF THE A (I) THEY DEPEND
RUPON THE SETTING OF THE VARIABLE HOW*
R
RHOW a REG
fR THIS OPTION IS USED WHEN A REGULAR DILUTION RATIO
R (E.G* DECIMAL) S USED I IS QE IT I LY NECESSARY TO PROVIDE
R THE FOLLOWING ADDITIONAL DATA
R
R Al Al ) A
R RAT ~ DILUTION RATIO
R
R THE FORtULA A lI) a Al*RAT*P.(Il-) IS THEN USED TO
R COMPUTE THE A(l), CARD TYPE II IS NOT NEEDED.
R
RHOW ~ VAR
R THIS OPTION IS USED WHEN THERE IS NO REGULAR DILUTION
R RATfO* VARIABLES Al AND RAT ARE IGNORED AND AN
R ADOITIONAL CARD OF TYPE II IS REQUIRE)D
R
R
R
RTHERE IS ALSO THE OPTION OF DOING A SINGLE CALCULATION OR THE
RPRINTING OF A TABLE OF CALCULATIONS* THIS DEPENDS UPON THE
RSETTING OF THE VARIABLE CHECK.
R
RCOECK B (BLANKS)
R A SINGLE CALCULATION IS MADE USING THE VALUES OF N(il
R AND P(1 ) ON THE CARD.
R
RCHECK 3 TABL
R A TABLE IS PRINTED4 ONLY P(1) IS READ. THIS IS TAKEN
R AS THE UPPER LIMIT FOR ALL OF THE N(I: AND P(CI. ALL
R FOLLOWING NUMBERS ON THE CARD ARE IGNORED.
R
R
R
RRESTRICTIONS ON VARIABLES
R NBR nUST BE AN INTEGER BETWEEN 3 AND 9
R N(1) AND P(I) ARE INTEGERS
R Alt RAT. AND THE All) MUST BE FLOATING POINT.
R (IF THERE IS DIFFICULTY IN OBTAINING A SOLUTION,
R THE All) SHOULD BE RESCALED TO A VALUE CLOSER TO.0O)
R
R
R
RSAMPLE INPUT CARDS ARE
R
R VAR 3 8 87 888
Rl. l.01
R
RTABLREG 3 *. 1 l 8
-C50 -
</P>
<P><PB REF="00000051.tif" SEQ="00000051" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="756" N="51">
Example Problem No. 79
MAD Program and Data (continued)
R
R
R
RTHE OUTPUT OF THE PROGRAM CONSISTS OF A LISTING OF THE INPUT
RVARIABLES PLUS THE VARIABLES X AND P. THE VARIABLE X IS THE
RMOST PROBABLE NUMBER OF BACTERIA PER ML* OF THE
RSOLUTIONF P IS THE PERCENTAGE PROBABILITY THAT THE
RCOMBINATION OF THE N(I) AND P(1) WOULD HAVE RESULTED IF X
RWERE THE TRUE NUMBER OF BACTERIA.
R
R
RIT SHOULD BE NOTED THAT ALL CALCULATIONS THAT YIELD A P
RLESS THAN *01 ARE SUPPRESSED WHEN A TABLE IS PRINTED.
R
R
R
DIMENSION A( 10)P(10 )PI(10) *N(10) *NI(10)
INTEGER CHECK, HOW, NBR.PI *INICOUNTtRET.KEY
BOOLEAN ZEROSeTABL
VECTOR VALUES BCD( 3 )S 36 S $4S8SS 60S5$72S S84S,$965.$108$
EXECUTE SETERRe START)
R
RREAD IN CARD TYPE I
R
START READ FORMAT IN. CHECK HOW,* NIBR.RAT.:Al*PI( 1})..**PI(9)NI( 1)
1. e*,N (9)
VECTOR VALUES INMSC4,C3t13,F4.3.F8*6t181I2*S
R
RCHECK NBR
R
WHENEVER 1BR.L 3 OR NBR.G4 9
PRINT CO4MENTS4NBR NOT IN RANGES
TRANSFER TO START
R
RREGULAR DILUTI ON RATIO
R
OR WHENEVER HOW*E*SREGS
THROUGH OMEGA FOR Il1 t I*G.NBR
OMEGA A(I)"Al*RATP ( I-1)
PRINT FORMAT REG6NBRtAlRAT
VECTOR VALUES REGaSx20H4NUMBER OF GROUPS IS.13t914 A(1) IS.
1F8.4*2H */18HODILUTION RATIO ISF6.4.2H **S
R
RVARIABLE DILUTIONS
R
OR WHENEVER HOW E.:,SVAR$
R
RREAD IN CARD TYPE Ii
R
READ FORMAT IN1l A( 1):.:..A(NBR)
VECTOR VALUES INl$=9F8.2*S
PRINT FORMAT VAR*A:( l)+,eA(NBR)
VECTOR VALUES VAR=S2H4A*S3,9E12 3*S
OTHERWISE
PRINT COMMENTS4INPUT ERRORS
EXECUTE SYSTEM.
END OF CONDITI ONAL
R
RPOSITION OUTPUT FORMAT
R
SLIDE(2)aBCD(NBR)
OUT(2 )BCD(NBR)
PRINT FORMAT SLIDE
VECTOR VALUES SLIDE=$1H /1H+tS9,S,1HX.S695HP 0/O*S
R
RSINGLE CALCULATION
R
WHENEVER CHECK~E.S S
TABL-OB
RET-O
-C51 -
</P>
<P><PB REF="00000052.tif" SEQ="00000052" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="711" N="52">
Method for Determining Bacterial Populations
MAD Program and Data (continued)
R
RFLOAT N(i- AND P(I
R
THROUGH PSI.FOR I-s1-.l: Ie6*BR
PSI p(-)PI(I)
PRINT FORMAT NOUT.*NI ( 1* )*NI(NBR)
VECTOR VALUES NOUTS$2H NS399t I17'5S5)*$
TRANSFER TO 6RIND
R
R TABLE
R
OR WHENEVER CHECK*,E.4STALS
TABL1 B
RET-1
P a PI1l)
EXECUTE SPRAY I(P,*N(l) 1 *N(NBR. tP:( 1 ).is,.*P( NBRi )
PRINT FORMAT NOUT.PI (1)
PRINT COMNENTSOS
R
RCALCULATION OF P( I
R
START( 1 ) THROUGH CHIFOR KEY:NBR'-1 KEY.L. 1
P('KEY) wP(KEY)-1.
R
RALL P(I) O I.E~!* FINISHED
R
WHENEVER ZEROS. (P*NBR)AMD.-P ( 2 ). *E 0*.AND P ( 1 )1.E.0 *
ITRANSPER TO START
WHENEVER P KEY) *GE 0.
R
RFIX P( I FOR PRINTING
R
THROUGH PHI IFOR I 'I * I4.&.NBR
PHt PI(I)-P(I)
TRANSFER TO GRIND
END OF CONDITIONAL
CHI P (KEY) P
OTHERWISE
PRINT CO#IENTS4INPUT ERRORS
EXECUTE SYSTEM.
END OF CONDITIONAL
R
RCALCULATION OF INITIAL GUESS FOR X
R
GRIND THROUGH UPSe FOR I=l1.1IE.GNBR
WHENEVER P( I)**0&.AND*P( I *L N( I)
X.f( At NBR/2+1 ) * SP ( 1l.1/t;(N ( I I.-*A ( I )
TRANSFER TO TAU
UPS END OF CONDITIONAL
TAU DELX=1.
R
RNEWTON-RAPHSON SOLUTION FOR X
R
THROUGH ALPHAt FOR COUNT'l lCOUNT*e.T&50,OR..ABS*DELX.L...005
GP-O0
THROUGH BETA. FOR I1.l. I.'.NBR
EXPON11*-EXP* (-A (1)*X)
sG& A 1P 1 I /EXPN-N( I ))
BETA GPGP-A ( I *A tA)*P(1 t ( l-EXPON)/ EXPONEXPON )
DELXG/GP
X-X-DELX
R
R(SOLUTIO UNSTABLE FOR X LESS THAN 0)
-C52 -
</P>
<P><PB REF="00000053.tif" SEQ="00000053" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="753" N="53">
Example Problem No. 79
MAD Program and Data (Continued)
R
ALPHA WHENEVER X*L*0*.X=.000000001
R
RITERATION LIMIT EXCEEDED WITHOUT REDUCING DELX TO LESS THAN
R 60005
R
WHENEVER COUNT wGe 50
PRINT FORMAT CMPTBLtPI (1)-...ePI (N R)
PRINT FORMAT XCED
VECTOR VALUES XCED=S1H,S113913HITER EXCEEDED*$
TRANSFER TO START(RET)
END OF COND.ITIONAL
R
RCtOPUTATION OF P
R
LOGEPtO
THROUGH GAMIMA FOR I l-.l ItG:NBR
GAMMA LOGEPaLOGEP+ELOG ( B I NO t 14(N (Z t P ( I I i -A (.) * (N I(:.)-P (I Xt. X
1+P( I )*ELOG. ( 1,-EXPe (-A( I *X)
PPER= OO *EXP (LOGEP)
WHENEVER PPERsL. O01.AND*TABLtTRANSFER TO START(RET)
PRINT FORMAT CMPTBLtPI ( 1):*PI(NBR)
VECTOR VALUES CMPTBL=$1H /2H+P,S3,9(I7tS5)*S
R
RPRINT ANSWERS
R
PRINT FORMAT OUT:*X#PPER
VECTOR VALUES OUTSIH *S4.S,E10*3tF7.2*$
TRANSFER TO STARTIRET)
END OF PROGRAM
$COIPILE MAD-PRINT OBJECT PUNCH OBJECT ZEROS001
R
RPROGRAM FOR CHECKING TO SEE IF P(2),.**..*P(NBR) ARE ZERO
R
R
R
EXTERNAL FUNCTION (PNBR)
INTEG4ER I NBR
ENTRY TO ZEROS*
THROUGH ALPHA FOR 1=21.-1I. eGNBR
ALPHA WHENEVER P( I,.G.0.tFUNCTION RETURN OB
FUNCTION RETURN lB
END OF FUNCTION
END OF PROGRAM
$COMPILE MADtPRINT OBJECTTPUNCH OBJECT
R
RPROGRAM FOR COMPUTING THE BINOMIAL COEFFICIENT
R
R
R
BINOMOOl
ALPHA
BETA
$DATA
VAR
1.
TABLREG
EXTERNAL FUNCTION (UP.DWN)
ENTRY TO BINOM.
FACT=1l
THROUGH ALPHAt FOR TOP=DWN+1 t,1,TOP.G4UP
FACT=FACT*TOP
THROUGH BETA. FOR BOT=2*l18BOT*G (UP-DWN)
FACT=FACT/BOT
FUNCTION RETURN FACT
END OF FUNCTION
END OF PROGRAM
3
1 *01
3 *1 1.
887
888
3
-C53 -
</P>
<P><PB REF="00000054.tif" SEQ="00000054" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="771" N="54">
Method for Determining Bacterial Populations
Computer Output
A
N
P
100E 01 100OE-01 1OOOE-02
8 8 8
8 8 7
X P 0/0
*208E 03 39.27
NUMBER OF GROUPS IS 3* Al ) IS
DILUTION RATIO IS.1000
N 3
P 3 3
P 3 3
P 3 3
P 3 2
P 3 2
P 3 2
P 3 2
P 3 1
P 31
P 3 1
P 3 1
P 3 0
P 3 0
P 2 3
P 2 3
P 2 2
P 2 2
P 2 1
P 2 1
P 2 1
P 2 0
P 2 0
P 2 0
P 3
P 1 2
1,0000.
x
X
2 *110E 03
1 *462E 02
0 *240E 02
3.292E 02
2.215E 02
1,149E 02
0 *933E 01
3 *159E 02
2 *115E 02
1 *749E 01
0.427E 01
2 *636E 01
1 *3-5E 01
1 *360E 01
0 *286E 01
1 *276E 01
0 *211E 01
2 *268E 01
1 *205E 01
0 147E 01
2 *199E 01
1 *143E 01
0 9-18E 00
0 *157E 01
1 *154E 01
P 0/0
44.44
42 77
36 59
*24
2.47
12. 51
32 82.03
*65
6.56
37 43
*16
3*10
*02
*22
*17
2 32
*02
*63
11 96
*02
1*12
31.93
*03
*03
-C54 -
</P>
<P><PB REF="00000055.tif" SEQ="00000055" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="870" N="55">
Example Problem No. 79
Computer Output (continued)
P 1 2 0 *114E 01 63
P 1 1 1 *112E 01 18
P 1 1 0 *7i6E 00 6.45
P 1 0 1 *72E 00 *62
P 1 0 0 *35E 00 39*20
P 0 2 0 *619E 00 *16
P 0 1 1 *11E 00 *05
p 0 1 0 *305E 00 3.36
Discussion of Results
The output shown above was designed to be illustrative rather than useful. In practice
more than three tubes would be used in each dilution group.
Compiling time for the program was 1.1 minutes. Computing the results shown above
took 0.5 minutes.
Critique
This problem is very appropriate for student solution. It gives the student an opportunity
to become familiar with numerical solution of complex equations as well as an opportunity to
become familiar with the use of the MPN equation. Solution of the MPN equation by hand is
extremely tedious, so that the tables prepared by the instructor's program will be useful
in later work.
References
1. Halvorson, H. 0. and Ziegler, N. R. J. Bact., 25, 101-121 (1933).
2. Norman, R. L. and Kempe, L. L., Biotechnology and Bioengineering, II,
No. 2, 157-163 (1960).
-C55 -
</P>
<P><PB REF="00000056.tif" SEQ="00000056" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="885" N="56">
Example Problem No. 80
ITERATION OF SUM OF PRINCIPAL STRESSES
AT INTERIOR POINTS OF PLANE RINGS
by
William Weaver, Jr.
Department of Civil Engineering
Stanford University
Course: Experimental Stress Analysis Credit Hours: 3 Level: Graduate
Statement of the Problem
Write a computer program which will evaluate the sum of principal stresses at
interior points of plane rings for which the boundary stresses are known from photoelastic analysis. The recommended procedure is to iterate the Laplace difference
equation applied to an x-y grid network of approximately square elements obtained by
the conformal transformation of the R-theta polar network using the relationships:
x = In (R)
y = theta
The program should be flexible enough to accommodate a whole ring or any segment of a ring for which the stresses on the radial boundaries are known to be equal
by symmetry. Divide the central angle of the polar network into 5~ increments. Make
the elements of the transformed rectangular grid as nearly square as possible. The
number of radial increments will depend upon the ratio of the outer and the inner
radii of the ring.
Limit the number of storage locations required by the grid to about 3,000.
Iterate the sum of principal stresses to within 3 decimal point accuracy. Limit the
number of iteration cycles to 50.
Solution
The Laplace difference equation for the sum of principal stresses at a lattice
point 0 and its four neighboring points 1, 2, 3, 4 (see Fig. 1) may be written:
S3-S0 SO-S1 S2-S0 S0-S4
Ax3 Axl AY2 AY4
+ = o (1)
Ax1 + Ax3 AY2 + AY4
2 2
-C56 -
</P>
<P><PB REF="00000057.tif" SEQ="00000057" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="843" N="57">
For the present problem let Ax1 = Ax3 = PAY2 = PAY4
Solving for SO in equation (1),
1 2 5\ P
S - (s + + + S4) (2)
0 2(p2+1) (2 3 2 (S2 4)
+1) +2+ / + s(2
in which, S = Sum of principal stresses at a point
p = Ratio of x-intervals to y-intervals.
Equation (2) may be used for the iterative evaluation of S at the lattice points
of the transformed grid.
In the photoelastic analysis of plane rings the principal stresses at the inner
and outer edges are determined experimentally at even intervals of the central angle.
A practical minimum interval is 5~, or r/36 radians.
By conformal transformation:
x = in (R); y = 9
Let Ay = AG = T /36. Define M as the number of angular increments, i.e., increments
in the y-direction. M is given by
M total angle of ring
AG
Thus, for a full ring and AG = 5~, M = 3600/5~ = 72. For a quarter ring
M = 90~/5~ = 18.
In order to obtain an exactly square grid, the number of subdivisions in the
radial direction must be:
xN-x0 ln(RN)-ln(RO) 36
- = - ^ - -A in(RN/RO)
N1 = /36 = /36 -
in which,
R = Inner radius of ring.
RN = Outer radius of ring.
x = Transformed inner radius.
xN = Transformed outer radius.
In general, N1 will not be an integer, but an approximately square grid may be
obtained by computing the actual number of intervals N used in the radial direction
as the next higher integer value of N1.
-C57 -
</P>
<P><PB REF="00000058.tif" SEQ="00000058" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="857" N="58">
Iteration of Sum of Principal Stresses
N = Xln(RN R)] Rounded up
The dimensionless ratio, p, may be computed from:
Ax N1 36 R
P = y - = n(RN/R) (4)
Radii of interior network circles are determined by the relationship:
(x ) (xo+jAx) ln(RN/R) j/N= R(RN/R)
Rj = e = e R(e )N = Ro(RN/RO) (5)
in which, Rj = Radius of jth network circle.
j = Radial index.
Assuming that a storage space of about 3,000 locations is available for the
M x (N+1) matrix of lattice points required for the solution,
M X (N+l) _ 3000. (6)
For a full ring the maximum value of M is 3600/50 = 72, and the value of N is
given by substitution of equation (3) into equation (6).
72{ [ l(N Rounded u + I1 3000 (7)
Equation (7) yields the maximum values of N and the corresponding ratio of radii.
(N)max = 41; (RN/R)max 35.0
The minimum non-trivial values for these parameters correspond to a network of
only two intervals in the radial direction.
(N)min = 2; (RN/RO)min 1 19
Figures 2 and 3 show the polar and transformed lattic networks for a quarter
ring having RN/R0 = 4, M = 18, and N = 16. This configuration, accompanied by a set
of hypothetical boundary stresses, constitutes the first of two trial problems submitted
to the computer. The second problem, for which no figure is given, consists of a
quarter ring having RN/R0 = 4/3, M = 18, and N = 4.
The first row of the transformed grid may be omitted because the stresses must
be the same as in the last row.
-C58 -
</P>
<P><PB REF="00000059.tif" SEQ="00000059" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="673" N="59">
Example Problem No. 80
Flow Diagram
M,)lINNER (i) P/REST I NDEX(
R R-)R) a
TITLE A INN^ER(,(M) vINE(I)
ReRNeUTER(1) \I>M E I>M / iT~u(r'.- *UTER(e )
N= Nf + f F t F,2 P PZ F I zZO) R T)= R e\ Ri -b
5 (Il,J ) = INNFR(I) " I TER'
z= If1= J = (FJ/IF.NNER(I) N 5WESWEPO 5EP=5WEEP+t
' \,,>l - eUTER(Z)) J COUNT '= 0
TEMP= FI(S(1,J-I)+S(i,J-+I))
+FZ(S(M, J) +S (,))
--- FE~R Tre) --- D
Special Output Form when N < 10:
INNER(I) = Boundary stress at inner edge.
OUTER(I) = Boundary stress at outer edge.
Fl = Weighting factor from first term in equation (2).
F2 = Weighting factor from second term in equation (2).
FJ = J expressed as a floating point number.
FN = N expressed as a floating point number.
-C59 -
</P>
<P><PB REF="00000060.tif" SEQ="00000060" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="839" N="60">
Iteration of Sum of Principal Stresses
MAD Program
W. WEAVER D046N 005 005 000
W. WEAVER D046N 005 005 000
$COMPILE MAD, EXECUTE
R
RTHIS PROGRAM IS INTENDED FOR USE IN EVALUATING THE
RSUM OF PRINCIPAL STRESSES- AT INTERIOR POINTS OF
RPLANE RINGS ANALYZED EXPERIMENTALLY BY THE
RMETHOD OF PHOTOELASTICITY* THE PROCEDURE CONSISTS OF
RTHE ITERATION OF THE LAPLACE DIFFERENCE EQUATION
RAPPLIED TO AN X-Y GRID NETWORK OF APPROXIMATELY
RSQUARE ELEMENTS OBTAINED BY THE CONFORMAL
RTRANSFORMATION OF THE R-THETA NETWORK USING
RTHE RELATIONSHIPS
R X=LN(R)
R Y=THETA
RTHE PROGRAM IS APPLICABLE FOR A WHOLE RING OR
RANY SEGMENT THEREOF,THE RADIAL BOUNDARIES OF
RWHICH LIE MIDWAY BETWEEN ADJACENT RADII OF SYMMETRY.
RTHE CENTRAL ANGLE IS DIVIDED INTO FIVE DEGREE
RINCREMENTStRESULTING IN M=72 INCREMENTS FOR A
RFULL RINGM=36 INCREMENTS FOR A HALF RINGM=18
RFOR A THIRD RINGETC. THE NUMBER OF RADIAL INCRERMENTS DEPENDS UPON THE RATIO OF THE OUTER AND THE
RINNER RADII OF THE RING. IN ORDER TO LIMIT THE NUMBER
ROF STORAGE LOCATIONS REQUIRED BY THE GRID TO ABOUT
R3,000, THE MAXIMUM RATIO OF RADII PERMITTED IS 35.
RTHE MINIMUM NON-TRIVIAL RATIO IS 1.19.
RTHE CALCULATIONS ARE CARRIED TO 3 DECIMAL PLACE ACCURACY
ROR 50 ITERATION CYCLES, WHICHEVER OCCURS FIRST.
R
PRINT FORMAT TITLE
VECTOR VALUES TITLE=$1H1,29HITERATION OF SUM OF PRINCIPAL/
11HO.42HSTRESSES AT INTERIOR POINTS OF PLANE RINGS*$
R
RREAD AND ECHO CHECK DATA
R
START READ FORMAT DATA1,M,INNER(1)...INNER(M)
READ FORMAT DATA2,RO.RN,OUTER(1)...OUTER(M)
VECTOR VALUES DATA1=$I2,S18S5F10.3/(7F10.3)*$
VECTOR VALUES DATA2=$2F10.45F10,3/(7F10.3)*$
THROUGH PRESET, FOR I=0,1,I.G.M
PRESET INDEX(I)=I
PRINT FORMAT DDATA1
PUNCH FORMAT DDATA1
VECTOR VALUES DDATA1-$HO5,4HDATA515t5HINDEXtS555HINNERt55
1,5HOUTER*$
THROUGH ALPHA, FOR I=1,lI.*G.M
ALPHA PRINT FORMAT DDATA2,INDEX(I),INNER(I),OUTER(I)
THROUGH ALPHA1i FOR I=1l1,t.G.M
ALPHA1 PUNCH FORMAT DDATA2,INDEX(I),INNER(I),OUTER(I)
VECTOR VALUES DDATA2=$1H,S27,I2,2F10.3*$
R
RCOMPUTE ITERATION FACTORS
R
Nl(36./3.14159)*ELOG.(RN/RO)
NN1l+1
Pz(36.*ELOG*(RN/RO))/(3.14159*N)
F11.o/(2.*(P*P+1.))
F2zP*P*F1
R
RCOMPUTE INTERMEDIATE RADII
R(USING FLOATING POINT NF=N AND JF=J)
R
NF=N
THROUGH BETA, FOR J=0,1,J.G N
JF=J
BETA R(J)xRO*EXP.((JF/NF)*ELOG.(RN/RO))
R
RSET STARTING VALUES FOR S USING LINEAR APPROXIMATION
RTO ESTIMATE INTERIOR VALUES BASED ON BOUNDARY VALUES,
R(USING FLOATING POINT NF=N AND JF=J)
R
DIMS(2)}N+1
THROUGH GAMMA, FOR Il.tltI.G.M
THROUGH GAMMA, FOR JOltJ*G.N
JFzJ
-c6o
</P>
<P><PB REF="00000061.tif" SEQ="00000061" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="801" N="61">
Example Problem No. 80
MAD Program, Continued
GAMMA S( IJ)=INNER( I )-(JF/NF)*( INNER( I )-OUTER( I))
R
RINITIALIZE COUNTERS FOR NO. OF CYCLES AND ACCURACY
R
SWEEP=O
R
RITERATE FIRST ROW*USING VALUES FROM FIRSTLASTAND
RSECOND ROWS.
R
ITER SWEEP=SWEEP+1
COUNT-O
THROUGH DELTA, FOR J=l1lJ.G.(N-1)
TEMP=Fl*(S(1,J-1 )+S(1 J+ )+F2*(S(MJ)+S(2,J))
WHENEVER.ABS (S(1,J)-TEMP )L.O. OO1COUNTCOUNT+l
DELTA S(1,J)=TEMP
R
RITERATE INTERMEDIATE ROWS
R
THROUGH EPSIL, FOR I=2,1.I.G*(M-1)
THROUGH EPSIL, FOR J=11,iJ.G.(N-1)
TEMP=Fl*(S( IJ-1)+S( IJ+1))+F2*(S( I-1J)+S(i+1,J))
WHENEVER.ABS.(S( I,J)-TEMP).L.0*001COUNT=COUNT+l
EPSIL S(IJ)=TEMP
R
RITERATE LAST ROWIUSING VALUES FROM FIRSTLASTAND NEXTRTO-LAST ROWS.
R
THROUGH ZETA. FOR J=ltlJ.G,(N-1)
TEMP=F1*(S(MJ-1)+S(MJ+1))+F2*(S(M-1,J)+S( 1 ))
WHENEVER.ABS.(S( MJ)-TEMP ) L.O0OO1COUNTzCOUNT+l
ZETA S(MtJ)=TEMP
R
RTEST FOR NO.OF CYCLES AND ACCURACY
R
WHENEVER SWEEP.GE.50,TRANSFER TO QUIT
WHENEVER COUNT.L*(M*(N-1)),TRANSFER TO ITER
QUIT WHENEVER N.L.10
TRANSFER TO ANS2
OTHERWISE
R
ROUTPUT FORM FOR N GREATER THAN OR EQUAL TO 10.
R
PRINT FORMAT FORM1,R(Q)..*R(N)
VECTOR VALUES FORM1=$1HO,11HVALUES OF R/(7F10.4)*$
PRINT FORMAT FORM1A.M,N
VECTOR VALUES FORM1A=-1HO,15HITERATED VALUESS10,2HMI,12,S10O
12HN=*I2//*$
PRINT FORMAT FORM1BS(10,)..*S(MtN)
VECTOR VALUES FORM1B=$(7F10.3)*$
PRINT FORMAT FORMlCSWEEPCOUNTM*(N-l)
VECTOR VALUES FORM1C=$1HO~14HNO, OF CYCLES=,SllI2/1HO023HN0.
1 OF ACCURATE VALUES=,I4/1HO,23HNO. OF ITERATED VALUES=,I4*$
END OF CONDITIONAL
TRANSFER TO START
R
RSPECIAL OUTPUT FORM FOR N LESS THAN 10.
R
ANS2 PRINT FORMAT FORM2,MNR(O)...R(N)
VECTOR VALUES FORM2=$IHO,15HITERATED VALUES,S102HM=,I2,S1l0
12HN=:I2/1HOS3'6HVALUES /1HOQS522HOFS7,2HR=,S3,lOF10.4*$
PRINT FORMAT FORM2AINDEX(O),.*INDEX(N)
VECTOR VALUES FORM2A=$1HOS4,5HTHETAS3,7HINDICES,
1S6,I1,9I10*$
THROUGH ETA. FOR I=1lltI.G.M
ETA PRINT FORMAT FORM2B,5*INDEX(I),INDEX(I),S(IO)...S(I*N)
VECTOR VALUES FORM2B=$55tI3~S7,I2,S3,1OF1O.3*$
PRINT FORMAT FORM2C9SWEEPCOUNTM*(N-1)
VECTOR VALUES FORM2C=$1HO,14HNO. OF CYCLES=S11,I2/1HO,
123HNO. OF ACCURATE VALUES=vI4/lHO,23HNO. OF ITERATED VALUES=
2I4*$
TRANSFER TO START
INTEGER M,NIJ9SWEEPCOUNTINDEX
DIMENSION INNER(72),OUTER(72),INDEX(72),R(42)
DIMENSION S(3000,DIMS)
VECTOR VALUES DIMS=2,1,42
END OF PROGRAM
-C61 -
</P>
<P><PB REF="00000062.tif" SEQ="00000062" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="827" N="62">
Iteration of Sum of Principal Stresses
Data
$DATA
18 -5000 -3222 -1800 -457 613
1711 2641 3416 3800 4000 3800 3416
2641 1711 613 -457 -1800 -3222
15000 60000 3000 7 -2697 -5002 -6896
-8005 -8985 -9625 -9900 -10000 -9900 -9625
-8985 -8005 -6896 -5002 -2697 7
18 -5000 -3222 -1800 -457 613
1711 2641 3416 3800 4000 3800 3416
2641 1711 613 -457 -1800 -3222
45000 60000 3000 7 -2697 -5002 -6896
-8005 -8985 -9625 -9900 -10000 -9900 -9625
-8985 -8005 -6896 -5002 -2697 7
Computer Output
DATA INDEX INNER OUTER
1 -5 000 3.000
2 -3.222 0.007
3 -1.800 -2.697
4 -0O457 -5.002
5 0.613 -6.896
6 1.711 -8.005
7 2.641 -8.985
8 3.416 -9.625
9 3.800 -9.900
10 4*000 -10.000
11 3.800 -9.900
12 3 416 -9.625
13 2.641 -8.985
14 1.711 -8.005
15 0.613 -6.896
16 -0.457 -5.002
17 -1.800 -2.697
18 -3.222 0.007
VALUES OF R
1.5000 1.6358 1.7838 1.9453 2.1213 2.3133 2.5227
2.7510 3.0000 3.2715 3.5676 3.8905 4.2426 4.6266
5.0454 5.5020 6.0000
ITERATED VALUES M=18 N=16
-5.000 -3.075 -2.290 -1.982 -1.920 -1.999 -2.160
-2.365 -2.584 -2.790 -2 953 -3.034 -2.977 -2.688
-1.988 -0.471 3.000 -3.222 -2.498 -2.049 -1.860
-1.852 -1.961 -2.142 -2.362 -2.595 -2.816 -2.998
-3 108 -3.098 -2.899 -2*404 -1.464 0.007 -1.800
-1.642 -1.545 -1 556 -1.666 -1.851 -2,085 -2.348
-2.620 -2.883 -3.119 -3*304 -3.412 -3.411 -3.271
-2.994 -2.697 -0.457 -0.725 -0.933 -1.153 -1.406
-1.691 -2.001 -2.325 -2.654 -2.979 -3.290 -3.578
-3.836 -4.063 -4.278 -4.543 -5.002 0.613 0.132
-0.307 -0.717 -1.114 -1.507 -1.901 -2.297 -2.693
-3.089 -3.484 -3 882 -4.290 -4.727 -5.234 -5 897
-6.896 1.711 0.945 0*288 -0.293 -0.824 -1.321
-1 799 -2.266 -2.730 -3.197 -3.674 -4.174 -4.713
-5.318 -6.031 -6.909 -8.005 2.641 1.649 0.807
0.079 -0.566 -1.154 -1.706 -2.238 -2,762 -3 291
-3.840 -4.425 -5.068 -5.799 -6.659 -7.699 -8.985
3.416 2.202 1.210 0.370 -0.363 -1.022 -1.632
-2.213 -2.785 -3.363 -3 966 -4.614 -5.332 -6.149
-7.103 -8.241 -9*625 3.800 2.529 1.461 0.554
-0.233 -0.936 -1.583 -2.196 -2.798 -3.407 -4.044
-4.730 -5.491 -6.357 -7 360 -8.533 -9.900 4.000
2.653 1.549 0.618 -0.188 -0*905 -1.564 -2.188
-2.800 -3.420 -4.068 -4 767 -5.544 -6.425 -7 443
-8.628 -10.000 3.800 2 530 1.463 0.557 -0.230
-0 931 -1.577 -2.190 -2.791 -3.401 -4. 038 -4.724
-5.487 -6.354 -7.358 -8.532 -9.900 3.416 2.204
1.214 0.376 -0.355 -1.013 -1.621 -2.202 -2.773
-3.351 -3.955 -4.604 -5.323 -6.142 -7.098 -8.239
-9.625 2.641 1.652 0.812 0.087 -0.556 -1.142
-1.692 -2.222 -2.746 -3 275 -3.825 -4.411 -5.056
-5.790 -6.653 -7.696 -8.985 1.711 0.948 0.294
-0.285 -0.812 -1.308 -1.783 -2.249 -2.712 -3.178
-062 -
</P>
<P><PB REF="00000063.tif" SEQ="00000063" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="790" N="63">
Example Problem No. 80
Computer Output, Continued
-3.657
0.613
-2.280
-5.227
-1.397
-3.564
-1.640
-2*609
-2 991
-1.957
-3.094
-4.158
0.135
-2.675
-5.893
-1.680
-3.824
-1.542
-2.872
-2.697
-2.137
-2.896
-4.700
-0.302
-3.071
-6.896
-1.987
-4.054
-1.551
-3.108
-3.222
-2.357
-2.402
-5.308
-0.708
-3 467
-0.457
-2.310
-4.272
-1.660
-3.294
-2.497
-2.590
-1 462
-6*024
-1.102
-3.866
-0.722
-2.639
-4.540
-1.843
-3.403
-2.047
-2.811
0.007
-6.905
-1.494
-4.276
-0.928
-2.964
-5.002
-2. 076
-3.404
-1.858
-2.993
-8.005
-1.886
-4 716
-1.146
-3.275
-1.800
-2.337
-3.267
-1.849
-3.103
NO. OF CYCLES=
50
NO. OF ACCURATE VALUES= 211
NO. OF ITERATED VALUES= 270
DATA
INDEX
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
INNER
-5.000
-3.222
-1.800
-0.457
0.613
1.711
2.641
3.416
3.800
4.000
3.800
3.416
2.641
1.711
0.613
-0.457
-1.800
-3.222
OUTER
3.000
0.007
-2.697
-5.002
-6.896
-8.005
-8.985
-9.625
-9.900
-10.000
-9.900
-9.625
-8.985
-8.005
-6,896
-5.002
-2.697
0.007
ITERATED VALUES
M=18
N= 4
VALUES
OF R 4.5000 4.8356
THETA INDICES 0 1
5 1 -5.000 -2.999
10 2 -3.222 -2.524
15 3 -1.800 -2 023
20 4 -0.457 -1 551
25 5 0.613 -1 132
30 6 1.711 -0.675
35 7 2.641 -0.248
40 8 3.416 0.114
45 9 3.800 0.326
50 10 4.000 0.414
55 11 3.800 0.326
60 12 3.416 0.113
65 13 2.641 -0.248
70 14 1.711 -0.674
75 15 0.613 -1.130
80 16 -0.457 -1 549
85 17 -1.800 -2 022
90 18 -3.222 -2 524
5.1962
2
-1.644
-1.845
-2.226
-2.607
-2.901
-3.039
-3.092
-3.086
-3.064
-3.050
-3.064
-3.086
-3.092
-3.039
-2.899
-2.606
-2.225
-1.844
5.5836 6.000
3
-0.018
-1.042
-2.430
-3.723
-4.775
-5.463
-5.974
-6.297
-6*447
-6.493
-6.447
-6.297
-5.974
-5.462
-4.774
-3.722
-2.429
-1.042
4
3.000
0.007
-2.697
-5.002
-6.896
-8.005
-8.985
-9.625
-9.900
-10.000
-9.900
-9.625
-8.985
-8.005
-6.896
-5.002
-2.697
0.007
11
NO. OF CYCLES=
NO. OF ACCURATE VALUES= 54
NO. OF ITERATED VALUES= 54
-C63 -
</P>
<P><PB REF="00000064.tif" SEQ="00000064" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="865" N="64">
Iteration of Sum of Principal Stresses
Discussion of Results
The first problem submitted to the computer involved 270 iterations in each
cycle. After 50 cycles, 211, (or about 78%)of the values were accurate to three
decimal places. Because N was large, values of R and S were printed sequentially
according to the general format labeled QUIT.
The second problem involved 54 iterations in each cycle, and all values were
accurate to three decimal places after 11 cycles. Because N was, in this case, less
than 10, values of R and S were printed in rows and columns corresponding to Figure 3,
according to the special format labeled ANS2.
Critique
This problem is very well suited for digital computation. The analysis of plane
rings is a recurring problem in photoelastic experiments, and the iteration of the sum
of principal stresses at interior points is an essential, though tedious, step in the
total analysis. The fineness of the lattice network is not significantly limited
by the storage space of a computer of reasonable size. A more important factor in
limiting fineness is the practical difficulty of measuring boundary stresses accurately
at the inner edges of plane rings for small increments of the central angle. Increments of 5~ are therefore used as a practical minimum in establishing the lattice
network.
x
yi 4
- Y4
) 1 Axl 0 Ax3 t )3
AY2
2
---------- ~~~~ -------- - - 2- - --
Figure 1.
Lattice Point 0 and its Four Neighboring Points 1,2,3,4
-C64 -
</P>
<P><PB REF="00000065.tif" SEQ="00000065" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="619" N="65">
Example Problem No. 80
o o
O ')
/c~~~~~ 0
/3
/0
7
S
Figure 2.
4
3
2
/o
0
+R~0
0 2 3, 5' G 7 8 9 / /3 2 13 /6
I /- /."6,0 7
Figure 2.
Quarter Ring: RN/R = 4
M = 18 N = 16
-C65 -
</P>
<P><PB REF="00000066.tif" SEQ="00000066" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="804" N="66">
Iteration of Sum of Principal Stresses
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16
N
0
0
1
2
3
4
5
6
7
I
i 8
9
10
11
12
13
14
15
16
17
M 18
0
00
50
10~
15~
200
25~
300
350
400
450
50o
550
60~
65~
70~
75~
800
850
90~
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16
Figure 3.
Transformed Quarter Ring
RN/R = 4 M = 18 N = 16
-C66 -
</P>
<P><PB REF="00000067.tif" SEQ="00000067" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="879" N="67">
Example Problem No. 81
ELASTIC BUCKLING LOAD FOR COLUMNS OF NON-UNIFORM CROSS SFCTION
by
K. H. Chu
Civil Engineering Department
Illinois Institute of Technology
Course: Advanced Structural Theory and Design Credit hours: 3 Level: Seniors
or First Year
Graduates
Statement of the Problem
The differential equation for elastic buckling of a column as shown below is
as follows:
P
IY
2
d Y M = P*Y P*_y
dX E*IA E*IA E*IO*I (1]
where X is the coordinate along the longitudinal axis of the column, Y is the
deflection of the column, M is the bending moment, P is the axial load, E is the
modulus of elasticity, IA is the moment of inertia of the cross-section about the
Z axis normal to the XY plane at the point under consideration, IO is the same
at a certain reference point and I = IA/IO.
The buckling load, PCR, is defined as the smallest value of P which will
satisfy the above equation with appropriate boundary conditions. The value of
PCR is to be determined by an iterative procedure using a digital computer.
Solution
Equation (1) is analogous to the following equation
=W (2)
dX
where W is the distributed load on a beam whose longitudinal axis coincides
with the X axis and M' is the bending moment due to W. W is positive when it
-C67 -
</P>
<P><PB REF="00000068.tif" SEQ="00000068" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="863" N="68">
Elastic Buckling Load for Columns of Non-Uniform Cross Section
is acting upward. M' is positive when it causes compression in the upper fiber
Qf a beam.
If YQ is a solution of Eq. (1), then using - as elastic weight W, the
bending moment obtained will be equal to
M' = YR - dX*dX + A*X + B (3)
where M' is renamed as YR and A and B are arbitrary constants determined by the
boundary conditions.
But, from Eq. (1) PCR
YQ -EI * YR (4a)
Therefore YQ _ PCR
YR -E*IO (4b)
or PC-R= - *E*IO (4c)
The procedure of obtaining the buckling load then is as follows: The
column length L is divided into N equal spacings* of length A. Assume an initial
set of deflectionsYO(J) for each point at X=J*A. Let I(J) be the I at the
same point. Then take - O(J) as the distributed elastic weight W at the point
X=J*^. Obtain the bending moment Y(J) at each of the equally spaced points. If
the ratio
YO(J) =C for 0 < J < N (5)
Y(J)
where C is a constant (within a tolerable limit) then
PCR = C*E*IO (6)
If not, then take the Y(J)'S as the new values for YO(J)'S and repeat the
procedure. Note that the real deflection Y(J) due to P=l is given by the
following equation:
Y(J) =P* Y = (J) (for P = 1) (7)
E*IO E*I0
* Spacings do not have to be equal, but computations will be simpler if equal
spacing is used.
-C68 -
</P>
<P><PB REF="00000069.tif" SEQ="00000069" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="869" N="69">
Example Problem No. 81
The problem is essentially an eigenvalue problem and it can be proved
mathematically that the procedure will always converge to the eigenvalue which
corresponds to the smallest value of P, i.e. PCR.
Newmark (reference 1) has shown that the numerical integration of
equations (1) and (2) can be facilitated by the following method. If W is a
distributed load, then we may pass a parabola through three points of W, say
W(1), W(2) and W(3), which are equally spaced at a distance A and convert the
distributed load to concentrated loads by the formulas as shown below.
CiCW I C2,W(3)
CW() CW(2) W(3)
CW() = [7*W(1) + 6*W(2) - W(5) * (8a)
CW(2) = [(1) + 10*W(2) +W(5)] * (8b)
CW(3) = [7*W(3) + 6*W(2) -W(l)1 (8c)
The shear and moment due to a system of concentrated loads at equally
spaced points can be easily obtained without further comment.
The solution to the problem consists of writing a program using the above procedure
which will determine the elastic buckling load of a column with simply supported
ends and with either uniform or non-uniform cross-section. For columns of nonuniform cross-sections, the cross-sections may be varying continuously or discontinuously
with steps occurring at any of the equally spaced points along the column. If
the discontinuity occurs between two of the neighboring equally spaced points,
the column will be treated as if the cross-section varies continuously between
these two points.
It should be noted that the solution technique, with slight modification,
can be applied to a cantilever column with one end free and the other end fixed.
For example, the cantilever column and the simply supported column as shown below
will each have the same buckling load.
-c69 -
</P>
<P><PB REF="00000070.tif" SEQ="00000070" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="843" N="70">
Elastic Buckling Load for Columns of Non-Uniform Cross Section
(b)
4 LI i LI ]
I L r-...
The following suggestions may be helpful in programming the solution. Let us
define the following quantities:
KY(J) Y(J) *E*I Y( (9)
2 N) 2 2 "' 2 N
L2*YO() L2*YO( L2*YO() N2*YO(()
E*IO*Y() Y(N) 2*Ky ( N) KY ()
T N2O(J) (11)
(CPCR) *KY(N)
Let "EPS" be a preassigned tolerance value. Then the condition for constancy of the ratio 0) (Eq. (5)) may be restated as follows:
Y(J)
IT - 11 < EPS and T >0 (for all N-l points)
(with 0 <c J<c N ) (12)
The formula for the critical load (Eq. (6)) may be rewritten into the
following form
PCR = (CPCR)* E*IO (13)
L2
If the number of points satisfying the above condition (Eq. (12)) is Q,
then if Q < N - 1, we shall take the KY's as the new values of YO's and repeat
the process. In order to prevent the machine from doing an endless number of iterations, a maximum number of cycles should be specified.
Moreover, in order to prevent the value of YO's from becoming very large after
h/ ^
many iterations, we may take the new YO's in terms of the maximum value of YO
* Note that N is the number of divisions and "integer division" is assumed
for evaluating the quantity N
2'
-C70 -
</P>
<P><PB REF="00000071.tif" SEQ="00000071" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="814" N="71">
Example Problem No. 81
or some value of YO which is close to the maximum, say the YO at a point near
the midpoint (J = ) of the column. In other words, we shall define a new
variable YO, such that
YO(J) = ^ N for initial values (14a)
YO( 7a
YO(J) = KY(J) for subsequent iterations (14b)
KY ()
In Eq. (14b) N
KY(j)
KY(J ) KY( ) = N (15)
YO(7 2
The KY's are obtained by taking
YOJ = W(J) (16)
EW(J) (17)
A
and
KY(J) = - where M' is the bending moment (18)
caused by EW(J)
From the definitions of YO and KY, we have
yo(N) TO(N) / To(N) _ o(N)
KY(O ) / (N KY( 2
KY( N K'(N)/ O() N KY (N9
Therefore, the formulas for CPCR and T may be rewritten as follows
N2*YO (N)
CPCR= (20)
KY(2)
T N2*YO(J) (21)
(CPCR)*KY(N)
-C71 -
</P>
<P><PB REF="00000072.tif" SEQ="00000072" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="778" N="72">
Elastic Buckling Load for Columns of Non-Uniform Cross Section
If the moment of inertia varies discontinuously at any of the equally
spaced points, then one way of solving the problem is as follows: Let the number
of points with a step change in moment of inertia be M. For M >0, we will have two
sets of index numbers for the division points as shown below.
M= 2
I i -
I I I
I i: I
I j
I
I.Dumm;n.ex.S 0, / - I 1
/ Dumry Ind4ex S o / 2{ 3
I I I I I I, '
2 Indx No. Set ( O0 2 3 4 6? 8 9 /0
- 1-119. 9- -. - I -
3 srpo(s)o I 3 I I 7 I I
4 Yo Yo(o) Y O() YY 3) YoW4) Yoro) Yo() o(7)J Y YOS) YO) Y.I.. t t ' I _
5 Index No. Set (Z) 2 3 4 5 6 7 8 9 /O // /2
6 STPP(S) -/ 8 1 I 83
7JI (o) 1(O ) 1(2) 1i)(4)) 11(6) 1(7) (I() 1(0) I
8 Y1 Y(o6) Y/(O/ YI(29 YtI) Y o()YI)Y ) Y/( Y/2)
q Ini-tirl Yalu/S (OW YOve) Y) YbO) Y(q) Y)Y)Y (1/0)|
0 VaO/us of YJ YO(o) Yo(l) YOZ() ) YYt ) YOf(4) YO(5) YO(6) Y)71YiK8) Yo{q) YO(IO)
after frst mo)
IVa/es of mv Y() YaO) Yoi O(2) "YO() YO(4) YO(5) W)( Vt7 t7) Y Y(9) YO(IO)
Iaf t e r — d, I v -__ _ _ _ _ _ _ _ _I_ I
The index number set (1) will have one number for each division point while
the index number set (2) will have two numbers at each step. The I's are stored
according to the index number set (2). (See line 7 of the figure.) At each step
in I, the index number in set (1) of the point with the step is stored in the
vector* "STPO". (The two ends are also considered as steps.) For example, the
index number 3 will be in the location for STPO(l). (See line 3 of the figure.)
The values of YO(J) are first stored in the vector Yl(J) with J corresponding to the index number set (1). (See line 9 of the figure.) Then the step numbers,
each corresponding to the index number in set (2) on the left side of a step,
are computed and stored in the vector "STPP". For example, the index number 8 will
be in the location STPP(2). (See line 6 of the figure.) Then starting from the
right end, the value of Yl(J) is put into the location of Yl(J+l) until the index
number of STPP(l) is reached. (See line 10 of the figure.) The process is then
repeated until STPP(2) is reached and so on until STPP(M). (See line 11 of the figure.)
* The "vector" is a linear array of numbers as defined in modern algebra
and used in computer language.
</P>
<P><PB REF="00000073.tif" SEQ="00000073" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="887" N="73">
Example Problem No. 81
Next, the concentrated elastic weights (EW) at the ends are computed. Then
the elastic weights at two neighboring points at each step and the elastic weights
between the steps are computed. Then the two concentrated elastic weights
at the two neighboring points of each step are combined into one value and the values
EW(J) are moved back to the locations with index number set (1). The procedure is
similar to the one for Yl(J) described above but in reverse.
With the elastic weights known, the procedure will be the same as that for
columns of continuously varying moment of inertia.
Table of Symbols
In the following, the definitions of the variables which have been introduced previously and will be used in the program are recapitulated. Some new
variables needed for the programming are also defined.
N number of divisions
M number of division points with step change in moment of inertia
EPS an assigned tolerance value
STPO(O)...STPO(M + 1) Index numbers of points with steps corresponding to
index number Set (1). The two end points are also
considered as steps.
I(O)...I(M + N) The ratio of the actual moment of inertia (IA) about
an axis normal to the plane of bending to the moment
of inertia at a reference point (IO). The I's carry
the index numbers of Set (1) for continuous variation
and they carry the index number of Set (2) for discontinuous variation at the division points.
YO(O)...YO(N) The ratio of the initial deflection at a point YO(J)
to the same at a reference point near the center of
the span YO(N). The YO's initially read in are
assumed. For subsequent iterations, YO(J) =
The YO's carry the index numbers of Set (1). K()
STPP(O)...STPP(M + 1) Index numbers of the steps each corresponding to
the index number in Set (21 on the left side of a
step.
Y1(O)...Yl(M + N) YO's carrying the index numbers of Set (2).
EW(O)...EW(M + N) Concentrated elastic weights due to the distributed elastic load - Y(J) with division length A as
I(J)
a unit. These elastic weights carry the index
numbers of Set (2).
EWNN refers to the set of elastic weights EW(O)...EW(M + N).
-C73 -
</P>
<P><PB REF="00000074.tif" SEQ="00000074" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="867" N="74">
Elastic Buckling Load for Columns of Non-Uniform Cross Section
EW(0)...EW(N) Similar to the weights in the set EWNN except that
the weights carry index numbers of Set (1). The two
weights at each stepping point of the Set EWNN have
been combined into one.
EW refers to the set of elastic weight EW(0)...EW(N).
RL, RR respectively left and right reactions due to EW
SW shear due to EW
A2*Yo(N)
KY(0)...KY(N) deflections with E*I as a unit. They are equal to the
bending moments (with division length A as a unit) due to the
concentrated elastic weights EW.
N2*YO(N)
CPCR =
KY(2)
2
T N2*YO(J)
(CPCR) *KY(N)
CYL = number of cycles of iterations
The data to be read in are as follows:
N, M, EPS
STPO(O)...STPO(M + 1)
I(O)...I(M + N
YO(O)...YO(N) (assumed)
The data read in should be printed out immediately as an echo check.
The output should be as follows:
CYL the final number of iteration when convergence is obtained.
CPCR the coefficient of the critical load.
The critical load is understood as given by Eq. 13
PCR = (CPCR)* E*I
L
Where E is the modulus of elasticity, IO is the moment of inertia at a reference
point and L is the length of the column as defined previously.
In order to check the results by manual computations, the following intermediate steps shall be printed out for each iteration:
Yl(0)...Yl(N + M)
EWNN tW(O)...EW(N + M)]
EW CEW(o...EW(N
KY(O)...KY(N)
-C74 -
</P>
<P><PB REF="00000075.tif" SEQ="00000075" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="351" N="75">
Example Problem No. 81
Flow Diagram
0= O,. G.(+l =5TP(S) TPO(5) -, 1, $.
S, I 152 \. _______, \ PRIN T — E______, rw~o - r 7*Y( 6 Y/ 7 VW
— > ~V J I G + * a') 24
~B J= N+S-L zI ) (O) ~-) -J
\\ EW(NN I(NN+) ~N+N 6*YI-(N)-I) Y.rNN-2
L. 571 P(S) TP(N+M) L3W(NN)= -G+ IZ4
-- -N- Zi-(u f) Z (NN-Z))
PP/T / L \ UWUS-PS)Tu, + I (U+I) + |(U) F
| 3r PP( —Ir I, m<r7oo v Y(2)]on / /-= — S r*pp(s) -
~EW(-=- (71; I(J)G;~ - I(J) EW(U) — =* WvU) + 6
/J53 pP(s) L
=0, I,$.G. U +(J- /OY/(.) yI(J+I) (?/
Ul OPP(S+/) 'T. &. =-STPPTS+I)4
LI(J-)^ ^^ )70 i Tg)-/
PRINT ----- / r^ \ KYW = —L P
(=, I,) 7S.G. -+ Eu - = U+z, I
EW(-N) JE G.(NAW-$+I)/ - EW
0co) 0c,) zO.
E w (0).. r 7 my!2,'o N + (oYOoNmI).YO.(N-'2)1, T G.N-
E-W( EW(N)- N) +r xON - N-Z).
Ew()-. ) 0 Yo) - +f., J. C7. F W(J.
L7 L8ysw(YI —
R RL= -99 T= 0, I, XC. N RL=RL- EW(1) 5W(o)= -TL = o, Z 5L8 I
L 9 KY JS o.- PRIN T C, ='/? )tlm(6XWO -T 0, X, It YO') + s wC Ti)J-~j.N- YO (NKZ)
pR/iNr
( <f)A4t c J CYL- START
CPcCR
-C75 -
</P>
<P><PB REF="00000076.tif" SEQ="00000076" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="779" N="76">
Elastic Buckling Load for Columns of Non-Uniform Cross Section
MAD Program
KUANG-HAN CHU D028N 001 010 000
KUANG-HAN CHU D028N 001 010 000
SCOMPILE MAD, EXECUTE
R ELASTIC BUCKLING LOAD FOR COLUMNS
R OF NON-UNIFORM CROSS SECTION
DIMENSION YO(50),I(100),EW(100) SW(50),KY(50),
1 STPO(50),STPP(100),Y1(100)
INTEGER N,JCYLQMNN SU U1 STPO,STPP
START READ FORMAT DATA1,NMEPS
VECTOR VALUES DATA1=$2I4,F7.5*$
READ FORMAT DATA2,STPO(0')**STPO(M+1)
VECTOR VALUES DATA2=$(18I4)*$
READ FORMAT DATA3,I(0).**I(N+M)
VECTOR VALUES DATA3=$(8F9.3)*$
READ FORMAT DATA4,YO(0) **YO(N)
VECTOR VALUES DATA4=$(8F9.3)*$
PRINT FORMAT ECHO1,N,M.EPS
VECTOR VALUES ECHO1=$S1H1S993HN =,I4,S3~3HM =I3,5S4~
1 4HEPSzF6.4*$
PRINT FORMAT ECH02,STPO(0)* *STPO(M+1)
VECTOR VALUES ECHO2=$1HO,59,8HSTEPS AT/(S10,1215)*$
PRINT FORMAT ECH03,I(0)***I(N+M)
VECTOR VALUES ECH03=$1HOQS9o10HRELATIVE I/(S10, 6F10.3)$*
PRINT FORMAT ECH04,YO(0)*o.YO(N)
VECTOR VALUES ECH04=$1H0S9,2HYO/(S10, 6F10.3)*$
CYLz1
Al WHENEVER M.E*0,TRANSFER TO A2
THROUGH L1,FOR J-=11,J.G*N
L1 Y1(J)=YO(J)
THROUGH L2#FOR S=01,[S.G.M+1
L2 STPP(S)=STPO(S)+S-1
THROUGH LS2,FOR S=1,lS.G.M
THROUGH LS2,FOR J=N+S-1,-1,J.L.STPP(S)
LS2 Y1(J+1)=Yl(J)
PRINT FORMAT CHECK1,Y1(0).**Yl(N+M)
VECTOR VALUES CHECK1=$1HO,59,2HY1/(S10. 6F10.3)*$
NN~N+M
EW(O)=-(7*Y1(O)/I(O)+6*Yl(1l)/I(1)-Yl(2)/I(2))/24
EW(NN)=-(7*Y1(NN)/I(NN)
1 +6*Y1(NN-1)/I(NN-1)-Yl(NN-2)/I(NN-2))/24
STPP(M+1)=NN
THROUGH L3~FOR S=1a1,S.G.M
UzSTPP(S)
EW<U+1)=-(7*YI(U+1)/I(U+l)
1 +6*Y1(U+2)/I(U+2)-Y1(U+3)/I(U+3))/24
L3 EW(U)=-(7*YI(U)/I(U)
1 +6*Y1(U-1)/I(U-1)-Y1(U-2)/I(U-2))/24
THROUGH LS3,FOR S=0,1,S.G.M
U=STPP (S)
U1lSTPP(S+1)
THROUGH LS3,FOR J=U+2t1,J.G.(U-1))
LS3 EW(J)z-(Y1(J-1)/I(J-1)+10*Y1(J)/I(J)+Y1(J+1)/I(J+1))/12
PRINT FORMAT CHECK2,EW(O).**EW(NN)
VECTOR VELUES CHECK2=$1HOS9,4HEWNN/(S10, 6F10.3)*$
THROUGH L4tFOR S=1lS.G.MM
U=STPP(S)-S+1
EW(U)=EW(U)+EW(U+1)
THROUGH L4,FOR J=U+2,19J.G.(NN-S+1)
L4 EW(J-1)=EW(J)
PRINT FORMAT CHECK3,EW(0).**EW(N)
VECTOR VALUES CHECK3=$1HOS9,2HEW/(S10, 6F10.3)*$
TRANSFER TO A3
A2 EW(O)m-(7*YO(0)/I(0)+6*YO(1)/I (1)-YO(2)/I(2))/24
EW(N)=-(7*YO(N)/I(N)+6*YO(N-1)/I(N-1)-YO(N-2)/I(N-2))/24
THROUGH L5 FOR J1,1,J.G~.(N-1)
L5 EW(J)m-(YO(J-1)/I(J-1)+lYO*O(J)/I(J)+YO(J+1)/I(J+l))/12
A3 RR=O
THROUGH L6~FOR J*0,1JeGeN
L6 RR~RR-EW(J)*J/N
RL'-RR
THROUGH L7,FOR Jx01[J.G.NN
-C76 -
</P>
<P><PB REF="00000077.tif" SEQ="00000077" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="772" N="77">
Example Problem No. 81
MAD Program, Continued
L7 RL=RL-EW(J)
SW(O)=RL
THROUGH L89FOR J=01,J.G.(N-1)
L8 SW(J+1)=SW(J)+EW(J)
KY(O)=O
THROUGH L9,FOR J=Oi,JJG. (N-1)
L9 KY(J+1)=KY(J)+SW(J+1)
PRINT FORMAT CHECK4,KY(0). *KY(N)
VECTOR VALUES CHECK4=$1HOS9 2HKY/(S10, 6F10.3)*S
Q=O
CPCR = N*N*YO(N/2)/KY(N/2)
THROUGH LlOFOR J=l11,J.G.(N-1)
T=N*N*YO(J)/(KY(J)*CPCR)
L10 WHENEVER (oABS (T-1) L.EPS)IAND (TGG*0 ),Q=Q+1
WHENEVER Q.E.(N-1)
TRANSFER TO A4
OTHERWISE
CYL=CYL+l
END OF CONDITIONAL
WHENEVER CYL.G.20
PRINT FORMAT NGD
VECTOR VALUES NGD=$1HQOS9.10HYO NO GOOD*$
TRANSFER TO START
OTHERWISE
THROUGH LllFOR J=1,1~J.G.N
L11 YO{J)IKY(J)/KYIN/2)
TRANSFER TO Al
END OF CONDITIONAL
A4 PRINT FORMAT CYNCYL
VECTOR VALUES CYN=$iHOS9.15HNO. OF CYCLES =tI4*$
PRINT FORMAT CRLtCPCR
VECTOR VALUES CRL=$lHO S9,15HCOEFF. OF PCR =tF7.2*$
TRANSFER TO START
END OF PROGRAM
Data
$DATA
10 2 *005
0 3 7 10
*2 2 2. 2 *2 1. 1. 1 1
1. *2 *2 *2 *2
0. *4 *6 *8 *9 1. *9.8
*6 *4 0.
10 0 *005
O 10
*25 *36.49 *64.81 1..81 *64
*49.36.25
0. *4 *6.8.9 1..9.8.6.4 0.
12 0 *005
0 12.2.2 *2.2 1. 1. 1. 1.
1 *.2 *2 *2,2
0..4.6.8.9.95 1..95.9.8 *6.4 0.
-C77 -
</P>
<P><PB REF="00000078.tif" SEQ="00000078" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="773" N="78">
Elastic Buckling Load for Columns of Non-Uniform Cross Section
Computer Output
N a 10 M a
STEPS AT
0 3
RELATIVE I
0.200
1.000
0.200
2
EPS=0.0050
7 10
0.200 0.200 0.200 1.000 1.000
1.000 1.000 0.200 0.200 0.200
YO
Y1
EWNI
EW
KY
Yl
m
00 0.000 0.600 06 0.,800 0.900
0.900 0.800 0.600 0.400 0.000
1.000
0.000
1*000
0.000
-0.375
-0.983
-0.375
0.400 0.600 0.800 0.800 0.900
0.900 0.800 0.800 0.600 0.400
-1.917 -3.000 -1.833 -0.417 -0.900
-0.900 -0.417 -1.833 -3.000 -1.917
-0.375 -1.917 -3.000 -2.250 -0.900
-0.900 -2.250 -3.000 -1.917 -0.375
0.000 8.558 15.200 18.842 20.233
20.233 18.842 15.200 8,558 0.000
-0.983
20.725
0.000
1.000
0.000
EWNN
-0.363
-0.996
-0.363
0.413 0.733 0.909 0.909 0.976
0.976 0.909 0.909 0.733 0.413
-2.026 -3.607 -2.157 -0,468 -0.973
-0.973 -0.468 -2.157 -3.607 -2.026
EW
KY
Yl
-0.363 -2.026 -3.607 -2.624 -0.973
-0.973 -2.624 -3.607 -2.026 -0.363
0.000 9.728 17.429 21.524 22.995
22.995 21.524 17.429 9.728 0.000
-0.996
23.493
0.000
1 000
0.000
EWNN
-0.363
-0.996
-0.363
0.414 0.742 0.916 0.916 0.979
0.979 0.916 0.916 0.742 0.414
-2.034 -3.646 -2*177 -0.470 -0.975
-0.975 -0.470 -2.177 -3.646 -2.034
EW
-0.363 -2.034 -3.646 -2.647 -0.975
-0.975 -2.647 -3.646 -2.034 -0.363
-0.996
KY
0.000
23 162
NO. OF CYCLES
9.801
21.689
a 3
17.568 21.689 23.162
17.568 9.801 0.000
23.661
COEFF. OF PCR a 4.23
</P>
<P><PB REF="00000079.tif" SEQ="00000079" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="785" N="79">
Example Problem No. 81
Computer Output, Continued
N 10 M 0
EPS=0 0050
STEPS AT
0 10
RELATIVE I
0.250
0.810
0.360 0.490 0.640 0.810
0.640 0,490 0.360 0.250
YO
KY
KY
KY
0.000 0.400 0.600 0.800 0,900
0.900 0.800 0.600 0.400 0,000
0.000 5.104 9,180 12.039 13,662
13.662 12.039 9.180 5.104 0.000
0,000 5.257 9.570 12*587 14e289
14.289 12.587 9.570 5.257 0,000
0.000 5.243 9.555 12.574 14,277
14.277 12,574 9.555 5.243 0,000
1.000
1.000
14.171
14.804
14.793
NO. OF CYCLES = 3
COEFF. OF PCR = 6.76
N = 12 M 0
EPS=0.0050
STEPS AT
0 12
RELATIVE I
0,200
1.000
0,200
YO
0.000
1.000
0.000
KY
0,000
31.267
0,000
KY
0.000
32.735
0.000
KY
0,000
32.778
0,000
NO. OF CYCLES
COEFF. OF PCR
0.200 0.200 0.200 1,000 1.000
1.000 1.000 0.200 0.200 0.200
0.400 0.600 0,800 0.900 0.950
0.950 0.900 0.800 0,600 0.400
11.183 20,450 26.717 29.325 30.771
30.771 29.325 26.717 20,450 11.183
11.604 21.445 28.056 30,756 32.236
32.236 30.756 28.056 21.445 11.604
11.610 21*470 28.095 30,799 32.280
32.280 30.799 28.095 21.470 11.610
3
= 4.39
-C79 -
</P>
<P><PB REF="00000080.tif" SEQ="00000080" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="841" N="80">
Elastic Buckling Load for Columns of Non-Uniform Cross Section
Discussion of Results
Three example problems are shown. The first example is a column shown below
which is the same as shown in the previous figure. With 10 division points and
values of EPS and YO as shown in the program, the procedure converges in 3 cycles
and gives the answer of (CPCR) = 4.23. The coefficient given by exact solution of
the differential equation is 4.22 (see reference 2).
IO.2 I I= 0.I
~I= O'.IA -I I= /O.Z
I i I
0 / 3 4 5 6 7 8 9 /
The second example is for a column as shown below. The coefficient determined
by the numerical method is 6.76, as compared with the value 6.72 which is interpolated
from a table given in reference 3.
o *^ ~ - -~ t ^ 6 4 ^^_,
0.2 0.'6 o. 4 64 0 0 0/ -4f 049 0
The third example for the column shown below is the same as the first, except
that it uses 12 divisions and any discontinuity of moment of inertia between two
neighboring division points is ignored and the column is treated as if it were a
column with moment of inertia varying continuously. The value of (CPCR) obtained
is 4.39 as compared with the true value 4.22.
L, = o. < { g sL - I r I 0.2
~.* **= ( ~ ~ ~ 0. 1X.0.0
0 2 3 4 5 6 7 8 9 / / 12
It should be noted that although the examples shown are for columns with
moment of inertia varying symmetrically with respect to the mid-point of the
column, the program is perfectly general and is applicable to columns with
moment of inertia varying unsymmetrically.
-c80 -
</P>
<P><PB REF="00000081.tif" SEQ="00000081" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="884" N="81">
Example Problem No. 81
Critique
The program is a natural one for the computer since the computer is most
suitable for any iterative procedure which requires us to do the same computation process over and over again.
References
1. N. M. Newmark,"Numerical Procedure for Computing Deflections Moments
and Buckling Loads, ASCE Transactions,1943.
2. S. Timoshenko, Theory of Elastic Stability.
3. S. Timoshenko & J. M. Lessels, Applied Elasticity.
</P>
<P><PB REF="00000082.tif" SEQ="00000082" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="885" N="82">
Example Problem No. 82
A MAD PROGRAM FOR TRUSS ANALYSIS
by
K. H. Chu
Civil Engineering Department
Illinois Institute of Technology
Course: Structural Theory Credit hours: 3 Level: Junior
Statement of the Problem
The problem is to write a program which will determine stresses in members
of trusses where the method of joints is applicable. The problem is of repetitive
type which is well adapted for computers. Also it is of fundamental importance
because if it is used as a subroutine, then a series of problems in truss analysis
can be solved easily; for example, influence lines for trusses, deflection of
truss by virtual work, stresses in indeterminate trusses, etc.
The method of joints is chosen instead of setting up a matrix and solving a
system of simultaneous equations because there will be too many zero elements in
the matrix. (The number of simultaneous equations the computer can solve by elimination is
usually limited, say to 30.) Therefore, there is a need to develop a program which will
not only conserve storage space but also will do the problem as we do it manually,
i.e., the machine would select the joint to be solved and proceed to a joint which
can be solved next.
Solution
Let us consider the truss shown on the next page. We will number the joints
and members in any fashion and set up an arbitrary but convenient coordinate
system.
-C82 -
</P>
<P><PB REF="00000083.tif" SEQ="00000083" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="865" N="83">
yI
The data to be read in are as follows:
N the total number of joints
M the total number of members
EPS1 a small number to guard against the case IDI| 0 (explained later)
JT(1)...JT(N) index number of the joints
JX 1...JX N X and Y coordinates of the joints
JY 1... JY NI
FXO(l)...FXO(N) the X and Y components of the force acting on each joint. Positive
FYO(l)...FYO(N) if the component is in the positive direction of the coordinate
axis, zero if there is none.
UNMB(1)...UNMB(N) the number of members meeting at each joint. Their stresses
are unknowns to be determined. For example, in the figure, the
value of UNMB(3) will be equal to 5.
MB(1...MB M
JN( 1...JN(M
JF 1..JF(M
the index number of the members and the corresponding index
number of the joints at the ends of each member. It does not
matter which end is referred to as the near joint (JN) and which
end is referred to as the far joint (JF). For example in the figure
the value of MB(15) will be 15 and the values of JN(15i and
JF(15) may be 5 and 7 respectively or 7 and 5 respectively.
-C83 -
</P>
<P><PB REF="00000084.tif" SEQ="00000084" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="878" N="84">
A MAD Program for Truss Analysis
So far as reactions are concerned, we need to read in the following data:
JRO The joint number where the reaction of unknown magnitude and direction is
applied. For example in the figure the value of JRO is 5.
R1 This value is equal to 1. The purpose of reading in this number is for
extending the program to the case of finding reactions along 3 given lines
of applications. The reactions will be identified as FR(Rl), FR(R2),
FR(R3) with R1, R2, R3 being equal to 1, 2, 3 respectively. For the sample
problem under consideration, it serves no useful purpose except to identify
the reaction with given line of application as FR(R1) or FR(1).
JR1 The joint number where the reaction of given line of application is applied.
For example, in the figure the value of JR1 is 6.
LXR1, LYR1, LR1 These are the numbers which specify the given line of application
in the following way. The line of action is assumed to be a
vector directed away from the joint and length of the vector
(LR1) may be taken as unity or any other value. The X and Y
components of this vector are denoted as LXR1 and LYR1
respectively.
In order that we understand the correlation of the various data which have
been read in, a concrete example is shown below. This is the same truss as shown
earlier except that the dimensions and loadings are specified and the origin of the
coordinate axes is at joint 1 and the X axis coincides with the lower chord. The
data for JT, JX, JY, FXO, FYO and UNMB are given in Table 1, and the data for MB, JN
and JF are given in Table 2. The data for the reactions are given in Table 3. In
the tables, there are many items which have not been defined. They will be explained
later., (s) ~(~) FR( FXRI
O/WXX
FYRO
4()30'
Having read in all the data, the first step is to compute the following
quantities:
LX(l)...LX(M) where L's are the length of the members and LX and LY are
LY(l)..LY(M) respectively the components of the vector of length L which is
L(1)..L(M) directed from the joint (JN) toward the joint (JF). The
computed values of LX, LY, and L for the members of the truss
shown are given in Table 2.
-c84 -
</P>
<P><PB REF="00000085.tif" SEQ="00000085" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="782" N="85">
Example Problem No. 82
Next we shall compute the reactions with the definition of terms as follows:
SFX, SFY, SMO
FR(R1)
FXR1, FYR1
are respectively the sum of the X components of the external
forces, the sum of the Y components of the external forces and
the sum of moments of the external forces about the joint (JRO).
is the force of the reaction of known line of action at the
joint (JR1). It will be positive if it is directed away from
the joint. This quantity is denoted as a subscripted variable
with the intention of writing a subroutine (or external function)
for computing 3 reactions with given lines of action. The
method of section is a possible application of this subroutine.
are respectively the X and Y components of FR(R1). The signs of
FXR1 and FYR1 are positive if they are in the positive directions
of the coordinate axes.
FXO(JR1), FYO(JR1)
FXRO, FYRO
FXO(JRO), FYO(JRO)
The external forces acting on the joint (JR1) must be
modified to include the reaction components just computed.
are values at the joint (JRO) corresponding to those as
defined for the joint (JR1).
The answers for the values of FR(R1), FXR1, FYR1, FXRO, FYRO for the example
shown are given in Table 3. The values of FXO, FYO at joints 5 and 6 after the
reactions have been computed are as shown by the circled values in
of Table 1.
lines 7 and 8
Table 1
lJne Dummy Index I Z 3 4 5 6 7 8 9
/ JT(I) /_ 2 3 4 5 6 7 8 9
2 IX(I) o. 30. 60. 90. /2 /20. 90. 60. 3 0.
3 JY() o. o.. o. o. 60. 40. 40. 40.
4 FXO() 0....0.. 0. 0000. 0.
5 FYO(I) -/0. -10. -/0.. 0. O. 0. 0o. O
6 UNMB1) 2 3 5 3 3 2 5 3 4
7 Af rSoAy/ng FXO 0. 0. 0. 0. '-5050.. 0 0. 0.
8 +hereactions FYO -/0. -/0. -/0. O4 0. ~ o. o. o.
9 After I UNMAB @. 5 3 3 5 3 (3
/0 solv/ng the FXO0, A750 0. 0 -5'C1. o0.. O. 0, 756
// | first Jointj FYr -/o. /0.00 -/o. - o/. t40. O. 0. 0.aO
Table 3
Dqtq Answers
JRo 5 FR(R/) 50o.
TR / 6 FXR 50.
ZXRI.00 FY I o.
LYRI 0.00 FXRO -50.
LRI. oo0 FYIqO 40.
Table 2
LinelDummY Index J I 2 3 4 5 6 7 8 /O /I 1/2 3 4 /5
/ M/B(J) 2_ 3 4. 5 6 7 /o // /Z./3 /4 /5
2 TJN() / 2 3 4 5 6 7 9 2 3 3 45
3 JF() 2 3 478 5 6 9 7 78
4 LX(J) 30. 30.. 30.3. 0. -30. -30. -30. -30. 0. -30. o. 30, 0. -30.
5 LLy(J) o. 0. 0. 60. -Zo. 0. 0. -40a +4040. +4-40. 40. 40. 40.
6 L(7) 30. 130. 30. 60. 06.o6 30. 30. 50. 40. 50. 40. 50. 40. 50.
7 j (J ) I /
8 Va1ves SeLX(J) 30. u 30.
for the joint'.
9 Es.-forYi LY(J), ______0 -. 40.
10 After LY)9 M I_____
I// so/ving SF(J) 9999 __999
_Z the first FX C) 0.0 _ 7__ 50 __1
/3 jointf F-YJ) -7.5 f/0O0
/4 __ F (J) -750o ___ _ /2.5_
-C85 -
</P>
<P><PB REF="00000086.tif" SEQ="00000086" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="880" N="86">
A MAD Program for Truss Analysis
Next we shall compute the stresses in the members. We shall use the following
symbols:
JT(I) is the index number of the joint to be solved as identified by
the dummy index I.
MS(l), MS(2) are the index numbers of the two members meeting at JT(I).
F(MS(l)), F(MS(2)) are the two forces in the members to be solved. They are
positive if they are directed away from the joint JT(I).
JTS(MS(1)), JTS(MS(2)) both store the value of JT(I) in order that in the
output, we know the forces in the members MS(1) and
MS(2) are the results of the solution of that
particular joint.
FX(MS(1)), FY(MS(l)) are the components of the forces acting on JT(I)
FX(MS(2)), FYMS 2 in the members MS() and MS(2). They are positive
if they are along the positive direction of the
coordinate axes.
IS(1), IS(2) are respectively the joint numbers of the other ends
of the members MS(l) and MS(2) which meet at the
joint JT(I).
FXO(IS(1)), FYO(IS(1)) The external forces FXO, FYO acting on the joints
FXO IS(2), FYO IS(2)) IS() and IS(2) shall be revised in order to include
the components FX(MS(1)), FY(MS(1)), etc. as known
forces. It should be noted that the signs of
FX(MS(1)), FY(MS(1)), etc. should be reversed since
they are originally considered as applied at the
other ends of the members (i.e. at the joint JT(I)).
The joint equations can be written very easily. In writing the joint
equations, it should be noted that the joint under consideration JT(I) is taken
as the near joint (JN). The forces are acting on the joint JT(I) and shall be
considered as positive when they are directed away from the joint JT(I). We
should also guard against the case where the denominator of the joint equation
is nearly equal to zero. This is really the case when the joint is unstable.
The remaining problem is how the machine can tell which joint is to be
solved first and which to be solved next. The following suggestions may be
helpful. Have the machine examine a list of UNMB. If there is a joint with
UNMB=2, the index number I of the joint JT(I) is picked up. First look into the
list for JN. If at the index number J, the value of JN(J) is equal to JT(I), then
the index number J is picked up and stored in the location MS(l). Notice that this
number J identifies the member MB(J) (See Table 1) which is renamed as MS(l),
since this is the first member to be solved for at the joint under consideration.
Then,in order to prevent the number in JN(J) and JF(J) from being picked up in future
examination of the list, we change the values of JN(J) and JF(J) into some unreasonable number, say 9999.
-C86 -
</P>
<P><PB REF="00000087.tif" SEQ="00000087" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="886" N="87">
Example Problem No. 82
Then look further down the list of JN. If we find that at another value of
J (say j'), JN(J') is equal to JT(I), then we store J' in MS(2) and the value
of JF(J') of the number MB(J') in IS(2). Also we change the values of JN(J')
and JF(J') into 9999.
If we examine the list of JN and find either none of the values or only
one of the values is equal to JT(I), then we must examine the list of JF to
find either both of the members MS(l) and MS(2) or only the member MS(2). The
procedure is similar to the one described above with the following exceptions.
Suppose that we find the value of JT(I) in the location of JF(J). We should
store the value of JN(J) in the location of IS and we should reverse the signs
of the values LX(MB(J)) and LY(MB(J)). The sign change is necessary because
the joint equation is written for the joint JT(I) which is considered as the
near joint (JN). Since we find the value of JT(I) in the list of far joints
(JF), the vector which represents the positive direction of the member MB(J)
must be reversed. (Note that the value of JT(I) is identical with the index
number I and the value of MB(J) is identical with the index number J.)
Having found the two members meeting at the joint JT(I) ready to be solved,
we end the search. We may set up a counter P to see whether P reaches the
value of 2. Also we would like to update the number in the list of UNMB for
the two members meeting at JT(I). The value of UNMiB(I) corresponding to the
joint solved JT(I) should be reduced by 2. Each of the two values of UNMB
corresponding to the far ends of the 2 members solved as identified by the
joint numbers corresponding to IS(1) and IS(2) shall be reduced by 1.
After going through the above procedure of searching and updating, we can
go through the following steps: (1) Applying the joint equations to find the
forces in the two members MS(l) and MS(2), (2) Obtaining the components of these
forces, and (3) Modifying the forces applied at the joints at the far end of the
member to include the components obtained as known forces.
As a concrete example of the above procedure, we examine the list of UNMB
in Table 1 and find that UNMB=2 at I=1. Thus JT(1) is the joint to be solved.
Next we examine the list for JN in Table 2 and find that JN(J)=JT(I)=I=l at
J=l. Thus we have MS(l)=l, noting that MS(l) really represents the member
MB(1). We set IS(l)=JF(1)=2 and change the values of JN(1) and JF(1) to 9999.
-C87 -
</P>
<P><PB REF="00000088.tif" SEQ="00000088" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="889" N="88">
A MAD Program for Truss Analysis
Looking further down the list of JN, we do not find any value equal to JT(I)=I=l.
When we look into the list of JF, we find JF(J)=JT(I)=I=1 at J=9. Thus
MS(2)=MB(9)=9. We set IS(2)=JN(9)=9 and change the values JN(9) and JF(9) to
9999 (See Table 2, lines 10 and 11). Since the number JT(l) is found in the list
of JF, we change the signs of LX(9) and LY(9) (See lines 8 and 9 in Table 2).
We will stop the search since we have already found the two members to be
solved. We then reduce the value of UNMB(JT(l)) by 2 and the values of
UNMB(IS(1)) and UNMB(IS(2)) each by 1 (See circled numbers in line 9 of Table 1).
Since the joint solved is JT(l), the values of JTS(l) and JTS(9) for the
members MB(1) and MB(9) respectively are equal to JT(l)=l. (See line 7 of Table 2).
The forces F(MS(1)) and F(MS(2)) are then determined by applying the joint
equations to the joint JT(1). (See line 14 of Table 2.) Their components with
respect to the joint JT(1) are then determined. (See lines 12 and 13 of Table 2.)
The forces FXO and FYO applied at joints IS(1) and IS(2) are then modified by
adding the force components of MS(1) and MS(2) with sign reversed to the external
forces acting at these joints. (See lines 10 and 11 in Table 1).
After the first joint has been solved, we re-examine the list of UNMB
(line 9 of Table 1) and repeat the whole procedure. We shall stop the process
when all stresses in all the M members are solved. Since after each application
of the joint equation, the total number of unknowns is reduced by 2, we may set
up a counter, say Q, which will be increased by 2 after each cycle. Since M
may be either even or odd, finally we may have (M-Q) equal to 0 or 1. If it is
equal to zero, the stresses in all the members have been determined. If it is
equal to 1, we shall find the remaining member using a procedure similar to the
one described previously. The force in the remaining member can easily be
determined. It should be noticed, however, that in order to maintain accuracy,
the force should be determined from the component with absolute value larger than
the other.
We would like to print out the data as originally read in. They are:
JT, JX, JY, FXO, FYO, UNMB, MB, JN, JF, JRO, R1, JR1, LXR1, LYR1 and LR1. The
answers to be printed out are: FR(R1), FXR1, FYR1, FXRO, FYRO, FXO(JRO), FYO(JRO),
FXO(JR1), FYO(JR1), LX, LY, L, JTS, FX, FY and F.
-c88 -
</P>
<P><PB REF="00000089.tif" SEQ="00000089" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="857" N="89">
Example Problem No. 82
With the detailed instructions as given above, the programming becomes
more or less straight forward. From statement labeled "Start" to the statement
just preceding the through loop L5, the data are read in and printed out in thespecific form as stated in the instructions.
The values of LX, LY and L are computed in the loop L5. The fifteen statements
following the statement SFX=O are essentially for determining the reactions. The
following quantities are computed: SFX, SFY, SMO, FR(R1), FXR1, FYR1, FXO(JR1),
FYO(JR1), FXRO, FYRO, FXO(JRO) and FYO(JRO). The formula for SMO may be written in
the following form:
SMO = {[JX(I)-JX(JRO] *FYO(I)- EJY(I)-JY(JRO)] *FXO(I) (1)
With computer language, in a THROUGH loop for I, the equation for SMO becomes
SMO = SMO + [JX(I)-JX(JRO)] *FYO(I)- [JY(I)-JY(JRO)] *FXO(I) (la)
The formula for FR(R1) is written in the following form:
FR(R1) = SMO*LR1 (2)
[JY(JR1)-JY(JRO)] *LXR1- LJX(JR1)-JX(JRO)J *LYR1
The searching for the members to be solved and updating various items by
the procedure as described in the instructions are performed from twenty-seven
statements after and including the statement labeled C3.
The joint equations are solved and the quantities F(MS(1)), FX(MS(1)),
FY(MS(1)), FXO(IS(1)), FYO(IS(1)), F(MS(2)), FX(MS(2)), FY(MS(2)), FXO(IS(2)),
and FYO(IS(2)) are computed from the statement (D =...) to the statement
(FYO(IS(2)) =...). The joint equations are written in the following form:
-FXO(JT(I)) * LY(MS(2) + FYO(JT(I)) * LX(MS(2))
L(MS(2)) L(MS(2)) (3)
F(MS(1))
D
LY(MStl)) LX(MStl/t
FXO(JT(I)) * LMS - FO(JT(I)) * LMS 1(
F(MS(2)) (4)
D
where
LX(MS(1)) * LY(MS(2)) LY(MS(1)) * LX(MS(2)) (5)
L(MS(1)) * L(MS(2)) L(MS(1)) * L(MS(2))
-c89 -
</P>
<P><PB REF="00000090.tif" SEQ="00000090" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="704" N="90">
A MAD Program for Truss Analysis
The counter Q registers the number of members already solved. If only one
member is left to be solved, (M-Q) = 1, the joint with UNMB=1 is searched out by the
loop LL4. The member to be solved is searched out by the loops LL5 and LL6.
The force in the member is found by the statement starting from the statement
labeled C6 to the statement END OF CONDITIONAL just before the statement
labeled C4.
The rest of the statements are for printing out the answers.
Flow Diagram
START
READ N, M,EP51
JTO) * ~ ~ 3JTN)
JX(O) ~~* JX(N)
JY('). Y(N)
FXO(1). FXO(N)
FYO()-. FYO (N)
UNMB(I).. UN(MB(N)
MB (0).. * M6 (M)
J.TN (1).... JN(M)
JF(l) ~ ~. rF (M)
Print
N,, A/,
EPS1
lRead T'RO 0Prin -. YR 0 L__ LX)= SYX(JF(7))- TJX(.TN())
(s _Rea d TR/,Pri ntL J/R o R LR LY(7)=YY(7F(J))- Y(XN(J)) - /
RI R.TR, JLXRI,LYR LRI R, LYRY LR JRI RILYRI, LR.M
L,(. =J LX(J)-x L*YX(J) t LY(J) LY(J),~~~~~~~Z
5SFX=u LNI SFX= S SFX+ f -XuI ()
( SFY= 0 L -,LN SFY=SFY+ FYo(ir) ASMo=o See E. (/a)in textf or SMO
5ee Eq.(Z) in tex- for FR(RI)
FXR I FA (R 1) * LXR I/LR /
FYKI = FR(RI) *LYRI/LRI
FXO ( JR) -= FXO(J-RI) -- FXRI
F YO (R I) = FYO (IR/) J FYR I I
GH
FXRO = - (FX- FXRI) -
FYR O = - (SFY - FYfI F) Pr'n t R I, FR(RI), FXF,/, Fi /, FX'Rc, FYRO
FXo (JRo) = FAX (J7Ro) FXRO JRO, FXO(FZRo) cJRo, FYO( ), /JRI, FXO(JRI), TRI, FY (JRI) E
FYO(JRO) = FYO (JRo) t FYR OQ,
-C90 -
</P>
<P><PB REF="00000091.tif" SEQ="00000091" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="705" N="91">
Example Problem No. 82
Flow Diagram, Continued
=-O!LL1 LL.
= QOI~. GN. 'NP48E P= LO =.
IT
U NMB (JT(I)) = See E'.(3) in fext for F(MS(i))
UN ) MB (Jrr))-Z ee2 Eq.(5) i I< EF E X(t))= FS(x )
7-UNW - tr SnFY(Ms I))= F(MSW()*Ly(LfM))
/NVM (I5(0->/ -Pr/intf L(MS(i))
UNMB (IS(2)) I UNSfrA8LE rAr) FO( )
IUT rsM(Z) ())T) 5 fxo(L5()) - FX(.MS())
JTS(MS( = ))=(I) / =FYO(I ()=
JTS(M())= Jr(l)FYO(IS())- FYS())
See Eq. (4) in textfor F(MS(Z))
* FX( 2 = F(45(2))* LX(MS(Z))
L(MS(Z))
FY(S(z2)) = F(Pf.52)) LY(MS z))
FXO(IS(z))= L(MS(z))
FXO (IS())- FX(M"(z))
rFO( Y(IS()= ) -
FYo (ISz)) -Fy(s5(z))
MAD Program
KUANG-HAN CHU D028N 001 003 000
KUANG-HAN CHU D028N 001 003 000
$COMPILE MAD, EXECUTE
R A MAD PROGRAM FOR TRUSS ANALYSIS
DIMENSION JT(50),FXO(50),FYO(50),UNMB(50)tJX(50) JY(50)
DIMENSION MB(100) JN(100) JF(100),LX(lOO) LY(100) L(100)
DIMENSION JTS(lOO) FX( 100)F Y( 100),F(lUO IS(2),MS(2) FR(3)
INTEGER JT,UNMB,JNJFMBJTS,IS,MSPQ,I,J,M,NJRO,JR1,R1
START READ FORMAT DATA1,N,M,EPS1
VECTOR VALUES DATA1i$2I4,F8.5*$
READ FORMAT DATA2,JT(1)*..JT(N)
VECTOR VALUES DATA2z$(184)*$
READ FORMAT DATA3,JX(1)...JX(N)
VECTOR VALUES DATA3=$(8F9.3)*$
READ FORMAT DATA4,JY(1)...JY(N)
VECTOR VALUES DATA4=$(8F9.3)*$
READ FORMAT DATA5,FXO(1)..*FXO(N)
VECTOR VALUES DATA5=$(8F9.2)*$
-C91 -
</P>
<P><PB REF="00000092.tif" SEQ="00000092" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="794" N="92">
A MAD Program for Truss Analysis
MAD Program, Continued
READ FORMAT DATA6,FYO( 1).FYO(N)
VECTOR VALUES DATA6S$(8F9-2)*$
READ FORMAT DATA7tUNMB(1)***UNMB(N)
VECTOR VALUES DATA7-$(1814)*$
READ FORMAT DATAllMB(l) 1.MB(M)
VECTOR VALUES DATA1l=$(1814)*$
READ FORMAT DATA12,JN(1),**JN(M)
VECTOR VALUES DATA12=$(1814)*S
READ FORMAT DATA13tJF(1),**JF(M)
VECTOR VALUES DATA13=$(1814)*$
PRINT FORMAT ECHOlN#MtEPSl
VECTOR VALUES ECHOls=$H1S7t16HNO. OF JOINTS, I2.54,
1 17HNO. OF MEMBERS =- I3vS57HEPS1 =,F8.5/*S
WHENEVER N.G.5
THROUGH L1,FOR I5=1 5*I.G*N
WHENEVER (I+4)*GE.N
TRANSFER TO Al
OTHERWISE
PRINT FORMAT ECH02 JT(I)**.*JT(1+4)
VECTOR VALUES ECH02=$lHOS7,2HJTtS2, 5I10/*$
PRINT FORMAT ECH03#JX(I)..JX(I+4)
VECTOR VALUES ECH03=$S892HJXS5 5F10.3/*$
PRINT FORMAT ECH04,JY(I)**.JY(I+4)
VECTOR VALUES ECH04=$S8,2HJYS5 5F10.3/*$
PRINT FORMAT ECH05,FXO(I).**FXO(I+4)
VECTOR VALUES ECH05=$S8,3HFXO,549 5F10*2/*$
PRINT FORMAT ECH06#FYO(I)***FYO(I+4)
VECTOR VALUES ECH06=$S8,3HFYOS4. 5F10.2/*$
PRINT FORMAT ECH07,UNMB(I)**.UNMB(1+4)
VECTOR VALUES ECH07=$S8,4HUNMB, 51l0/*$
L1 END OF CONDITIONAL
OTHERWISE
I1l
END OF CONDITIONAL
Al PRINT FORMAT EECH02. JT(I)***JT(N)
VECTOR VALUES EECH02=$S892HJTS2, 5I10/*$
PRINT FORMAT EECH03,JX(I)*..JX(N)
VECTOR VALUES EECH03=$S8.2HJX,55 5F103/*$
PRINT FORMAT EECH04,JY(I)**.JY(N)
VECTOR VALUES EECH04=$S8,2HJY55, 5F10*3/*$
PRINT FORMAT EECH05,FXO(I)**.FXO(N)
VECTOR VALUES EECH05O$S8,3HFXOS4, 5F10.2/*$
PRINT FORMAT EECH06,FYO(I)e**FYO(N)
VECTOR VALUES EECHO06-$S83HFYOS4, 5F10.2/*$
PRINT FORMAT EECH07,UNMB(I)**.UNMB(N)
VECTOR VALUES EECH07=$S8t4HUNMB, 5I10/*$
WHENEVER M.G.5
THROUGH L2,FOR JIl, 5,J.G.M
WHENEVER (J+4)GE.M
TRANSFER TO A2
OTHERWISE
PRINT FORMAT ECHOllMB(J)...MB(J+4)
VECTOR VALUES ECHOll$lHOS722HMBS2, 5I10/*$
PRINT FORMAT ECH012,JN(J).*.JN(J+4)
VECTOR VALUES ECH012=$S8,2HJN,52, 5I10/*$
PRINT FORMAT ECH013,JF(J)**.JF(J+4)
VECTOR VALUES ECH013=$S8,2HJFtS2. 5I10/*$
L2 END OF CONDITIONAL
OTHERWISE
J1l
END OF CONDITIONAL
A2 PRINT FORMAT EECOlIMB(J)***MB(M)
VECTOR VALUES EECOll=$lHO,57,2HMBS2o 5I10/*$
PRINT FORMAT EEC012,JN(J)*..JN(M)
VECTOR VALUES EEC012=$S8,2HJNS2, 5I10/*$
PRINT FORMAT EEC013#JF(J)**.JF(M)
VECTOR VALUES EEC013=$S8,2HJFS2, 5I10/*$
READ FORMAT ROoJRO
VECTOR VALUES RO$I4*$
READ FORMAT RARlJR1,LXR1~LYR1,LR1
VECTOR VALUES RA=$I4,S29I4,3F10*.5*
PRINT FORMAT RROJRO
VECTOR VALUES RRO=$S8,4HJRO,I3/*$
-C92 -
</P>
<P><PB REF="00000093.tif" SEQ="00000093" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="789" N="93">
Example Problem No. 82
MAD Program, Continued
PRINT FORMAT RR1,JR1,LXR1,LYR1,LR1
VECTOR VALUES RR1=$S8,4HJR1 I 3/S8 4HLXR1 S3 F10.3
1 S5,5HLYR1,F10.3#S5,5HLR1,F10.3/*$
THROUGH L3,FOR J=-lltJ.G*M
LX(J)=JX(JF(J))-JX(JN(J))
LY(J)=JY(JF(J))-JY(JN(J))
L3 L(J)=SQRT.(LX(J)*LX(J)+LY(J)*LY(J))
SFX=O
SFY=O
SMO=O
THROUGH LN1,FOR I=11PI.G..N
SFX=SFX+FXO(I)
SFY=SFY+FYO(I)
LN1 SMO=SMO+(JX(I)-JX(JRO))*FYO(I)-(Jy(I)-JY(JRO))*FXO(I)
FR(R1)=SMO*LR1/{(JY(JR1)-JY(JRO))*LXR1
1 -(JX(JR1)-JX(JRO))*LYR1)
FXR1=FR(R1)*LXR1/LR1
FYR1=FR(Rl)*LYR1/LR1
FXO(JR1)=FXO(JR1)+FXR1
FYO(JR1)=FYO(JR1)+FYR1
FXRO=-(SFX+FXR1)
FYRO=-(SFY+FYR1)
FXO(JRO)=FXO(JRO)+FXRO
FYO(JRO)=FYO(JRO)+FYRO
PRINT FORMAT ANS4R1i>FR(R1)tFXRltFYR1.FXROFYRO
VECTOR VALUES ANS4=$1HOS7t3HFR(9I291H),S1lF1O.2,S5S4HFXR1,
1 S1tF102,S5F4HFYR1S1,lF1Oe2/S8,4HFXROS3tF1lO2,55$4HFYRO~S1l
2F10.2/*$
PRINT FORMAT ANS5,JROFXO(JRO),JROFYO(JRO),
1 JR1,FXO(JR1) JR1,FYO(JR1)
VECTOR VALUES ANS5=$S8,4HFXO(,I2t1H),F10.2,S3,4HFYOC(I291H)t
1F10.2/S8,4HFXO(.I21H) F1O.2,S3~4HFYO(,12,1H) F10.2*$
Q=O
C3 THROUGH LL1,FORI=,1 1I.GGN
LL1 WHENEVER UNMB(I)*E.2,TRANSFER TO C1
Cl P=0
THROUGH LL2,FOR J=11J.*G.M
WHENEVER JN(J)*E.JT(I)
P=P+1
IS(P)=JF(J)
MS(P)=J
JN(J)=999
JF(J)=9999
WHENEVER P.E.2.TRANSFER TO C2
LL2 END OF CONDITIONAL
THROUGH LL3OFOR J=11t~J.G.M
WHENEVER JF(J).E.JT(I)
P=P+1
IS(P)=JN(J)
MS(P)=J
JN(J)=9999
JF(J)=9999
LX(J)=-LX(J)
LY(J)=-LY(J)
WHENEVER P.E.2tTRANSFER TO C2
LL3 END OF CONDITIONAL
C2 UNMB(JT(I))=UNMB(JT(I))-2
UNMB(IS() )=UNMB(IS(1))-1
UNMB(IS(2))=UNMB(IS(2))-1
JTS(MS(1) )JT(I)
JTS(MS(2))=JT(I)
D=LX(MS(1))*LY(MS(2))/(L(MS{(1)*L(MS(2)))
1 -LY(MS(1))*LX(MS(2))/(L(MS(l))*L(MS(2)))
WHENEVER.ABS.(D).L.EPS1
PRINT FORMAT UNSTBL
VECTOR VALUES UNSTBL=$S8,22HTHIS TRUSS IS UNSTABLE*S
TRANSFER TO START
END OF CONDITIONAL
-C93 -
</P>
<P><PB REF="00000094.tif" SEQ="00000094" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="790" N="94">
A MAD Program for Truss Analysis
MAD Program, Continued
F(MS(1))=(-FXO(JT(I)L*LY(MS(2))/L(MS(2))
1 +FYO(T(I))*LX(MS(2) )/L(MS(2)))/D
FX(MS(1))=F(MS(1))*LX(MS(1))/L(MS(1))
FY(MS(1))=F(MS(1))*LY(MS(1))/L(MS( 1))
FXO(IS(1 )=FXO(IS(1))-FX(MS(1))
FYO(IS(1) )FYO(IS1) )-FY(MS(1))
F(MS(2))-(FXO(JT(I))*LY(MS(l))/L(MS(l))
1 -FYO(JT(I))*LX(MS(1))/L(MS(1)))/D
FX(MS(2) )F(MS(2))*LX(MS(2))(2/L(MS))
FYMS(2)) =F(MS(2))*LY(MS(2) /L(MS(2))
FXO(IS(2) )FXO(IS(2))-FX(MS(2))
FYO(IS(2) )FYO(IS(2))-FY(MS(2))
Q=Q+2
WHENEVER (M-0).G61,TRANSFER TO C3
WHENEVER (M-Q).E.OTRANSFER TO C4
THROUGH LL4,FOR IltlI.G.N
LL4 WHENEVER UNMB(I)*E.1,TRANSFER TO C5
C5 THROUGH LL5,FOR J=1,1,J.G.M
LL5 WHENEVER JN(J).E.JT(I)JTRANSFER TO C6
THROUGH LL6,FOR J=1,1,J.G.M
WHENEVER JF(J).*EJT(I)
LX(J)=-LX(J)
LY{J)=-LY(J)
TRANSFER TO C6
LL6 END OF CONDITIONAL
C6 FX(J)=-FXO(I)
FY(J)=-FYO( I)
JTS(J)=JT(I)
WHENEVER *ABS(F Y(J))*GE.*ABS^(FX(J))
F(J)=FY(J)*L(J)/LY(J)
OTHERWISE
F(J)=FX(J)*L(J)/LX(J)
END OF CONDITIONAL
C4 WHENEVER M.G*5
THROUGH L4OFOR J=l1 5tJ*.GM
WHENEVER (J+4) GE.M
TRANSFER TO A3
OTHERWISE
PRINT FOPMAT ECHO11,MB(J).**MB(J+4)
PRINT FORMAT ANS1,LX(J)*..LX(J+4)
VECTOR VALUES ANS1=$S892HLX#S5, 5F10.3/*$
PRINT FORMAT ANS2,LY(J)**.LY(J+4)
VECTOR VALUES ANS2=$S8,2HLY#55t 5F10.3/*$
PRINT FORMAT ANS3.L(J)...L(J+4)
VECTOR VALUES ANS3=$S8,1HL,56, 5F10.3/*$
PRINT FORMAT ANS11,JTS(J)***JTS{J+4)
VECTOR VALUES ANS11=$S8,3HJTSS1s 5110/*$
PRINT FORMAT ANS129FX(J).**FX(J+4)
VECTOR VALUES ANS12=$S8,2HFXS.5, 5F10.2/*$
PRINT FORMAT ANS13,FY(J).**FY(J+4)
VECTOR VALUES ANS13=$S8,2HFYS55S 5F10.2/*$
PRINT FORMAT ANS14tF(J)**.F(J+4)
VECTOR VALUES ANS14=$S8,1HF9S6, 5F10.2/*$
L4 END OF CONDITIONAL
OTHERWISE
Jul
END OF CONDITIONAL
A3 PRINT FORMAT EECO11lMB(J)***MB(M)
PRINT FORMAT AANS1*LX(J)*.,LX(M)
VECTOR VALUES AANS1=$S8,2HLX~S5 5F10.3/*$
PRINT FORMAT AANS2,LY(J)*..LY(M)
VECTOR VALUES AANS2=$S8,2HLY,55, 5F10.3/*$
PRINT FORMAT AANS3tL(J)*..L(M)
VECTOR VALUES AANS3=$58,1HLS6, 5F10.3*$
PRINT FORMAT AANS119JTS(J).**JTS(M)
VECTOR VALUES AANS11=$S893HJTSS1* 5I10/*$
PRINT FORMAT AANS12,FX(J).**FX(M)
VECTOR VALUES AANS12=$S8,2HFX9S5t 5F10.2/*$
PRINT FORMAT AANS13,FY(J)...FY(M)
VECTOR VALUES AANS13S$S8,2HFYS59 5F10.2/*$
PRINT FORMAT AANS14$F(J).**F(M)
VECTOR VALUES AANS14$S8,1HFS6, 5F10.2*$
TRANSFER TO START
END OF PROGRAM
-C94 -
</P>
<P><PB REF="00000095.tif" SEQ="00000095" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="774" N="95">
Example Problem No. 82
Data
$DATA
9
1
2
1
1
2
5
1
15
2
0O
30.
0.
O.
40.
0O
O.
-10,
3
3
2
2
3
3
4001
4 5
30.
0.
0O
-10.
3 3
3 4 5
4 5
5 6
6 7 8 9
60. 90*
120. 120*
90. 60.
40. 40.
0. 0. 0. 60.
0. 0. 00. 0. 0.
5
3
3
4
-10. -10.
2 5 3 4
6 7 8 9
6 7 8 9
7 8 9 1
0.0
O. 0. 0O 0O
10
2
9
1.0
11
3
9
12
3
8
13
3
7
14
4
7
15
5
7
6 1.0
Computer Output
NO. OF JOINTS = 9
JT 1
JX 0*000
JY 0*000
FXO 0.00
FYO -10.00
UNMB 2
JT 6
JX 120.000
JY 60.000
FXO 0.00
FYO 0.00
UNMB 2
MB 1
JN 1
JF 2
MB 6
JN 6
JF 7
MB 11
JN 3
JF 9
JRO 5
JR1 6
LXR1 1.000
FR( 1) 50.00
FXRO -50.00
FXO( 5) -50.00
FX3( 6) 50.00
MB 1
LX 30.000
LY 0.000
L 30.000
JTS 1
FX -7.50
FY -0.00
F -7.50
MS 6
LX -30.000
LY -20.000
L 36.056
JTS 6
FX -50.00
FY -33.33
F 60.09
MB 11
LX -30.000
LY 40.000
L 50.000
JTS 3
FX 15.00
FY -20.00
F -25.00
NO. OF
2
30.000
0.000
0.00
-10.00
3
7
90.000
40.000
0.00
0.00
5
2
2
3
7
7
8
12
3
8
LYR1
FXR1
FYRO
FYO( 5)
FYO( 6)
2
30.000
0.000
30.000
2
-7.50
-0.00
-7.50
7
-30.000
0.000
30.000
7
-22.50
0.00
22.50
12
0.000
40.000
40.000
3
0.00
0.00
0.00
MEMBERS =
3
60.000
0.000
0.00
-10.00
5
8
60.000
40.000
0.00
0.00
3
3
3
4
8
8
9
13
3
7
0.000
50.00
40.00
40.00
0.00
3
-30.000
-0.000
30 000
4
45.00
0*00
-45.00
8
-30.000
0.000
30.000
8
-22.50
0.00
22,50
13
-30.000
-40.000
50.000
7
-22.50
-30.00
37.50
15 EPS1 = 0*00100
4 5
90.000 120.000
0.000 0.000
0.00 0.00
-10.00 0.00
3 3
9
30.000
40.000
0.00
0.00
4
4 5
4 5
5 6
9 10
9 2
1 9
14 15
4 5
7 7
LR1 1.000
FYR1 0,00
4 5
-30.000 -0.000
-0.000 -60.000
30.000 60.000
5 6
45.00 0e00
0.00 33.33
-45.00 -33.33
9 10
30.000 0.000
40.000 40.000
50.000 40,000
1 2
7.50 0.00
10.00 10.00
12.50 10.00
14 15
0.000 -30.000
40.000 40.000
40.000 50.000
4 5
0.00 5.00
10o00 -6.67
10.00 -8.33
-C95 -
</P>
<P><PB REF="00000096.tif" SEQ="00000096" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="889" N="96">
A MAD Program for Truss Analysis
Critique
In the writer's opinion, the program is well suited for instructive purposes. The problem can be separated into three parts. In the first part, the
students can be asked to write a program (or draw a flow diagram) for computing
the lengths and components of the length vector of the members. In the second
part, the students can be asked to write a program for computing the reactions
and modifying the external forces at the joints where the reactions are applied.
There are also possibilities of writing the program for computing length or
reactions or both as an external function.* The students may be asked to write
a program for computing three reactions along three given lines of action and
applied at three different joints of the truss.* This program, if it is written
as an external function, may be applicable to trusses which can be solved only
by method of sections.*
Then the students can be asked to write a program for computing two unknown
forces in members meeting at a joint using the method of joints and modifying
the external forces applying at the joint at the other ends of the members. The
process of searching for the joint to be solved and updating the various items
involved is a rather complicated procedure. Detailed instructions as given in
the text seem to be necessary. However, this complicated procedure is very
instructive since it dramatically illustrates the computer's ability to handle
problems which are essentially logical rather than computational. Often a human
being is scarcely aware of the logical detail which is involved where he solves
a problem manually.
* Due to lack of time, these problems have not been worked out by the writer.
-C96 -
</P>
<P><PB REF="00000097.tif" SEQ="00000097" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="876" N="97">
Example Problem No. 83
STIFFNESS FACTORS FOR A FLAT SLAB COLUMN
by
Robert B. Harris
Department of Civil EngnQring
The University of Michigan
Course: Advanced Design of Structures Credit hours: 3 Level: Graduate
Statement of Problem
Estimate the slab thickness, drop panel thickness, and the column diameter
for a typical panel in your flat slab structure. Use the ACI methods for this
approximation.
Using these results along with the story height, study the effect of column
capital size upon the stiffness coefficients of the column and the carry-over
factors.
Program your solution for the digital computer using the MAD language.
Solution
Stiffness coefficients for members with variable cross sections may be found
by expressing the relationship between unit end moments applied to a member and the
rotations of the end tangents. If the ends are defined by a and b respectively,
and any end moments Mab and Mba are applied simultaneously, the rotation of tjie end
tangents a and Gb must be,according to reference 1:
Ga = Mabal - Mbaa2
b -Mabbl + MbaSb2
Dal and Pbl are the angle changes due to a unit moment applied to end a and Pa2
and,b2 are the angle changes due to a unit moment applied at end b. By the reciprocal
theorem,angle Pbl is equal to Pa2'
Solving the above equations for Mab and Mba:
b2 + bl
Nab k= a -+ b
Mba -. a A b
-C97 -
</P>
<P><PB REF="00000098.tif" SEQ="00000098" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="784" N="98">
Stiffness Factors for a Flat Slab Column
in which
-A.= alb2 - (bl)2
EI
If I is the minimum moment of inertia then K = - and:
11
Pal
1 KA
C bl
c2 - K~-;
c b2
3 = KA_
The following sketch defines lengths used.
The following sketch shows the conjugate beams loaded with the appropriate
diagrams from which the values of al bl and b2
EI diagrams from which the values of p B P an< Bb? can 1be 2 oun<
B
-C98 -
</P>
<P><PB REF="00000099.tif" SEQ="00000099" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="807" N="99">
Example Problem No. 83
Assuming circular elements in the capital, the moment of inertia varies as
the fourth power of the diameter, therefore
I(H-Q) ^ D and
(H-T) C
(H-Q), D
) (DY)4
I(y) (Dy)4
In order to control the accuracy of the solution a variable was introduced called
DEL which represents the size of increment AY. This was assigned the value of 0.2
in the sample program. This was particularly necessary since rectilinear integration was used.
The following summation equations are the result of the application of the
conjugate beam method and the letters appearing in brackets below the various
terms are the letters used in the program for that quantity.
n = H-Q
iEI6A1 =
n = A
2
n = H-T 2
Yn
Y AY + D4 I -- bAY
n =/ l (Q (DY)
n = (H-Q A )
-" -
Al A
2EIfBi = (H - Yr
n^n
=
2
B
n = H-T
4 --
) Yn AY + D T
n = (H-Q4S )
2
(H - Yn) Y
4(DY
(DY)
I I I L
Bl
I I
G
G
I
K
n = H-Q
H2IB2 =,
b y
aY
n "-2
n = H-T
(H - Yn)2 AY + D4
n = (H-Q+y-)
2
(H - Y )
(DY)
I I L
B2
I L
F
E
-099 -
</P>
<P><PB REF="00000100.tif" SEQ="00000100" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="773" N="100">
Stiffness Factors for a Flat Slab Column
Letting X = H2EI, then
A1 = Al
DB1 Bl
PB2 = B2
PB2 = i
1
r
A =PAlB2 - Pl2
A = 5A1lB2 - 5B1
(Al) (B2)
X2
- B1l'
where U = (A1)(B2) - B12
U
X2
X Al
K U
(Al)X2
XKU
Coefficient at bottom, C1 =
X HEI = H
K EI
H
Hence,
Coefficient at top,
H3Al
1 U
H3B2
C3 U
For the carry-over factors the coefficient C2 is needed.
C2
Carry over - top to bottom
Carry over - bottom to top
HBl1
U
B1
B2
B1
A1
Flow Diagram
(Start) —* Print IRead - HC = 0 /Through
Title Data BACK
D,T,H For C
D,1.,HC.E.H/2
Calculate A = 0 /Through
) HC, Q E =E = 0 Y SUMA
G = 0 For Y
DEL = 0.2.l,01YG.H
-C100 -
</P>
<P><PB REF="00000101.tif" SEQ="00000101" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="786" N="101">
Example Problem No. 83
MAD Program, and Data
R. B. HARRIS S164E 003 020 000
R. B. HARRIS S164E 003 020 000
$COMPILE MAD, EXECUTE
R
R PROGRAM FOR STIFFNESS COEFFICIENTS FOR A COLUMN WITH 90 DE
1GREE CAPITALS
R
START PRINT FORMAT TITLE
R
R HEAD INGS
R
VECTOR VALUES TITLE=$1H1,S23 26HCOLUMN STIFFNESS STUDY
1 //S26,19H90 DEGREE CAPITAL ///S5,3HCAPS3,3HCOL~S4,3HCAP,
2S3,4HSLABS5,3HCOLS4t3HTOPS6,3HBOTS5,5HCARRY*S4t5HCARRY/
3S5 4HOIAM,S2 4HDIAM S2,5HDEPTHS2,5HTHICKS2 6HLENGTHS2,
45HCOEFF,S4,5HCOEFFS3,7HTOP-7BOT S2,7H0T-TOP*$
R
R INPUT
R
READ FORMAT DATA D, T. H
VECTOR VALUES DATA = $3F10.3*$
HC = 0
THROUGH BACK. FOR C a Dl HCHE.H/2
HC = C/2 - D/2
Q - T + HC
R
R INTEGRATING STRAIGHT LINE PORTION
R
A = 0.0
E - 0.0
G = 0.0
DEL = 0.2
THROUGH SUMA, FOR Y = *1~*2,Y*G.(H-Q)
A = A +(Y.P.2)*DEL
E = E +((H-Y).P.2)*DEL
SUMA G = G +((H-Y)*Y)*DEL
R
R INTEGRATING CURVED PORTION
R
DY = 0.0
B 0.0
F a 0.0
K = 0*0
DEL = 0.1
THROUGH SUMB~ FOR Y = (H-Q+.05).S1,Y*G.(H-T)
DY = DY + D + 2*(C-D)*(Q-H+Y)/HC
B = B +(D.P.4)*DEL*Y*P.2/DY*P.4
F = F +(D.P.4)*DEL*(H-Y).P.2/DY.Ps4
SUMB K = K +(D.P.4)*DEL*(H-Y)*Y/DY#P.4
R
R STIFFNESS FACTORS AND OUTPUT
R
Al = A+B
B1 = G+K
B2 = E+F
U = A1*B2-B1.P*2
CBOT = A1*H*P.3/U
CTOP = B2*H*P.3/U
COBT = B1/A1
COTB = B1/B2
BACK PRINT FORMAT COEFFt C~ Dt HC. THt, CTOPCBOT.COTB.COBT
VECTOR VALUES COEFF = $F8*1,4F7,1,4F9F4*$
TRANSFER TO START
END OF PROGRAM
SDATA
1*0 0,000 10.0
1.0 1,000 10.0
-C101 -
</P>
<P><PB REF="00000102.tif" SEQ="00000102" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="821" N="102">
Stiffness Factors for a Flat Slab Column
Computer Output
COLUMN STIFFNESS STUDY
______..__. _.___ 90. DEGREE CAPITAL
CAP COL CAP SLAB COL TOP OT CARRY CARRY
DIAM DIAM DEPTH THICK LENGTH COEFF COEFF TOP-BOT BOT-TOP
1.0. ~~1.0.U.. 160.0 4.0012 4.0012.5002.5002
2.0 1.0..5.0 10.0 4.6159 4.1473_.4985 _.5548
3.0 1.0 1.0.0 10.0 5.9587 4.4210.4878.6575
4.0 1.0 1.5.0 10.0 7.1344 4.6274.4765.7346
5.0 1.0 2.0.0 10.0 9.4565 4.9705.4538.8634
6.0 1.0 2.5.0 10.0 11.5429 5.2332.4359.9616
7.0 1.0 3.0.0 10.0 15. 1811 5.6761.4056 1.1278
8.0 1.0 3.5.0 10.0 19.7202 6.0218.3838 1.2568
9.0 1.0 4.0.0 10.0 28.0289 6.6153.3492 1.4796
10.0 1.0 4.5.0 10.0 36.1126 7.0907.3253 1.6568
11.0 1.0 5.0.0 10.0 54.0806 7.9271.2889 1.970 -CULUMN STIFFNS S STUDY
90 DEGREE CAPITAL
CAP COL CAP
DIAM1 DIAM I-EP r
1.0 1.0.0
2.0 1.0.5
3.0 1.0 1.0
4.0 1.0 1.5
5.0 1.0 2.0
6.0 1.0 2.5
7.0 1.0 3.0
8.0 1.0 3.5
9.O 1.0 4.0
10.0 1.0 4.5
11.0 1.0 5.0
SLAB COLL TOP BOT CARRY CARRY
THICK LENGTH COEFF CUEFF TOP-BUT bOT-TOP
1.0 10.0 6.C930 4.4461.4866.6669
1.0 10.0 7.1344 4.6274._ 4765.7346
1.0 lu.0 9.4565 4.9705.4538.8634
1.0 10.0 11.5429 5.2332.4359.9616
1.0 10.0 15.7811 5.6761.4056 1.1278
1.0 10.0 19.7202 6.0218.3838 1.2568
1.0 10.0 28.0289 6.6153.3492 1.4796
1.0 10.0 36.1126 7.0907.3253 1.6568
1.0 10.0 54.0806 7.9271.2889 1.9707
1.0 10.0 72.7273 8.6220.2643 2.2293
1.0 10.0 117.4918 9.8889.2276 2.7045
Discussion of Results
Capital diameter, column diameter, capital depth, slab thickness, and column
length must all be put in using the same units. In the sample these are all in
feet and are printed out with the results for checking. The program limits the
capital depth to one-half the column length and increments by one-half a length unit.
Coefficients obtained give an accuracy of four significant figures. This
accuracy may be improved by making the integration increments smaller in the SUMA
and SUMB iterations.
-C102 -
</P>
<P><PB REF="00000103.tif" SEQ="00000103" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="891" N="103">
Example Problem No. 83
Critique
This problem was selected as representative of the use of the computer in making
repetitive computations which would take considerable time to do longhand. In a sense,
the problem is synthetic. If, however, such information were available with the detail
obtained, it would materially help in preliminary design and checking of flat slab
construction.
The problem was designed to serve as an introduction to the use of a more complicated program prepared by the instructor, the use of which would allow the student
to make adjustments in columns, slabs, drop panels, and capitals, and exercise the
broadest judgement in design.
As an introduction, the student reaction was quite good and little difficulty
was experienced in programming. About one week and a half (3 laboratory periods)
were used in this activity, but actual results were not turned in for some time after
this since mechanical delays with the Computing Center prolonged the submission date.
This problem served the purpose for which it was intended and in addition firmed
up in the students' minds the effectiveness of the column capital in flat slab construction.
Reference
1. Maugh, L. C. Statically Indeterminate Structures, page 133.
-C103 -
</P>
<P><PB REF="00000104.tif" SEQ="00000104" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="872" N="104">
Example Problem No. 84
VIBRATION OF BEAMS ON SPRING SUPPORTS
by
Shan S. Kuo
Department of Civil Engineering
Yale University
Course: Structural Dynamics
Credit Hours: 2
Level: Senior, Graduate
Statement of the Problem
Write a computer program to determine the natural
shapes for a beam on spring supports as shown in Figure
5 4 3 2
k
frequencies and the modal
1. The beam is idealized
1
kII
L h h h h
X = 4h
<,~~~~
Figure 1. Beam on spring supports
as a multi-mass system. Each of the two end masses is taken as one half of an interior
mass, and the sum of all masses is equal to the total mass of beam m. In order to
study the effect of spring constant k, we take c R/ M H as a parameter where
R = 2k/m = "rigid body"frequency
6H =7 2 EI/m = fundamental frequency of simply supported beam
The program should be so written that any number of lumped masses can be handled.
Solution
We are concerned with the following matrix equation:
X = FMX
60 2
where
F = flexibility matrix due to the spring action and the beam action
M = mass matrix, a diagonal matrix
-C104 -
</P>
<P><PB REF="00000105.tif" SEQ="00000105" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="820" N="105">
X = displacement of lumped mass
O = frequency
It should be noted that the matrix FM is not symmetrical; therefore,
Jacobi's method cannot be applied to this case directly. The following is a general
procedure for solving this problem using a digital computer.
(a) Formulate the flexibility matrix F. It is considered here as the sum of three
action matrices, one from the beam action and the other two from the spring
action. In connection with the beam action matrix, we take
kh T 4 R 2
EI 2(n+l)3 RH
where n = number of interior masses, or number of total masses -2.
(b) Formulate the mass matrix M.
(c) Formulate the matrix M1/2 FM1/2
(d) Make use of the Jacobi's method to find the eigenvalues and the eigenvectors
of the matrix found in the step c.
(e) The natural frequencies are then equal to (eigenvalues)1/2
(f) The true eigenvectors are equal to M1/2 times the computed eigenvectors.
The steps (a) and (b) are incorporated in the subroutine "INPUT", whereas the
steps (c) through (f) are carried out by the subroutine "RITE'.
Flow Diagram
44/ 7/ \- Lp^ /^7/ \-5TA 42 97 CALL C<L
Pox
[NIN -F/,- - wr < J / )- (T'wwHe-))CO
C/ I )= O7 STp
N~y Npz# XJ7 N 21/ 141 AA(X)~& / K a<ct) o. \/
-C105 -
</P>
<P><PB REF="00000106.tif" SEQ="00000106" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="461" N="106">
Vibration of Beams on Spring Supports
Flow Diagram, Continued
hi; 2o/ \ /;l\ V CTOR /7 ps 2V E R02
/3 y r,:^ <:l r
(xTJl + + /&CT, J P- ) ^ L53\ ^j / - 3 X4,J NpZT
rL(x Cz,J) =JI1A P= /, 2 j/, NPZ T CrA, - -t) Tzp/, t f)/>
4- V C-C7T^ (I), -)
-c106 -
</P>
<P><PB REF="00000107.tif" SEQ="00000107" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="825" N="107">
Example Problem No. 84
Flow Diagram, Continued
FORTRAN Programs and Data
SHAN S. KUO (YALE UNIV)
SHAN S. KUO (YALE UNIV)
$COMPILE FORTRAN. EXECUTE
D003N 9
DOO3N 9
004 027 200
004 027 200
SPRING 6
SPRING 6
C
C
C
C
C
C
VIBRATION OF BEAMS ON SPRING SUPPORTS SHAN S. KUO 7/29/61
DIMENSION FLEX(5,5),ZMASS(5,5) TEMP(5,5),VECTOR(5t5),
1 XX(5),AA(5) WRWH(5)
COMMON NSPANBETA
NFIRST IS THE FIRST N AND NLAST IS THE LAST N
READ INPUT TAPE 7,12, NFIRSTNLASTNRWRWH
12 FORMAT(312)
READ INPUT TAPE 7,14, (WRWH(I1),I1=1,NRWRWRH)
14 FORMAT(6F6.3)
DO 471 N=NFIRSTNLAST
SPAN=N+1
DO 471 J=1iNRWRWH
BETA=48 *70455*WRWH(J)*WRWH( J ) / (SPAN*SPAN*SPAN)
WRITE OUTPUT TAPE 6,223tWRWH(J)
WRITE OUTPUT TAPE 5,223,WRWH(J)
223 FORMAT(67H1VIBRATION OF BEAM ON SPRING SUPPORTS. A PROBLEM OF SHOC
1K ISOLATION/64HORIGID BODY FREQUENCY / FUNDAMENTAL FREQ. OF HINGED
2-HINGED BEAM=F5.2)
CALL INPUT(FLEXZMASSVECTORTEMP XXAA)
CALL RITE (FLEX9ZMASSVECTORTEMP)
471 CONTINUE
-C107 -
</P>
<P><PB REF="00000108.tif" SEQ="00000108" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="790" N="108">
Vibration of Beams on Spring Supports
FORTRAN Programs and Data, Continued
SCOMPILE FORTRAN
C
C INPUT MATRICES SHAN S. KUO 7/29/61 8/5/61
C
SUBROUTINE INPUT (FLEX#ZMASSVECTOR,TEMPXXAA)
DIMENSION FLEX(595) ZMASS(5,5) VECTOR(5,5),XX(5) AA(5),TEMP(5.5)
COMMON N SPANBETA
EN=N
NM1=N-1
NP1=N+1
NP2=N+2
N2'N*2
C
C FORM THE INFLUENCE MATRIX* II AND JJ FOR MASS NUMBERING
C
DO 16 I=-1N
EI=I
IPl=I+1
XX(IP1)=EI
16 AA(I)=EI
XX(NP2 )SPAN
XX( 1 ) 0
DO 15 II=lN
DO 15 JJ=1iN
IP1=II+1
IF(XX(IP1)-AA(JJ) ) 17t17 18
17 BB=SPAN-AA(JJ)
FLEX(II JJ)= (BB/(6*SPAN))*(XX(IP1)*XX(IP1)*XX(IP1)
1 -(SPAN*SPAN-BB*B8)*XX(IP1))*(- 1)
GO TO 15
18 BBASPAN-AA(JJ)
FLEX(IIJJ)=-(BB/(6**SPAN))*(XX(IP1)*XX(IP1)*XX(IP1)
1 -(SPAN*SPAN-BB*BB)*XX( IP)
2 -(SPAN*(XX( IP1)-AA(JJ) )*(XX( IP )-AA(JJ) )*(XX( IP1)-AA(JJ) )/BB))
15 CONTINUE
C
C SHIFT FLEX MATRIX DUE TO BEAM ACTION. BORROW VECTOR
C
DO 201 I=1,NP2
DO 201 J=1iNP2
201 VECTOR(ItJ)=O0
DO 202 I=1tN
IP1=I+1
DO 202 J=lN
JP1=J+1
202 VECTOR(IP1,JP1)=(FLEX(I J))*BETA
WRITE OUTPUT TAPE 6, 777, BETA
WRITE OUTPUT TAPE 5. 777, BETA
777 FORMAT(7HOBETA =F6.3)
C
C FLEX IS NOW ONE OF THE 3 MATRICES TO FORM THE FINAL MATRIX FLEX
C
DO 203 I=1,NP2
EI=I
DO 203 J=1,NP2
203 FLEX( I J)=(SPAN+1.-EI)/SPAN
C
C BORROW TEMP, ONE OF THE 3 MATRICES TO FORM THE FINAL FLEX MATRIX
C
DO 204 I=1,NP2
EI=I
DO 204 J=1NP2
EJ=J
204 TEMP (IJ)=(EJ-1.)*(SPAN+2.*(1.-EI))/(SPAN*SPAN)
C
C HAPPY FINAL FLEX
C
DO 205 I=1,NP2
DO 205 J31,NP2
205 FLEX(I J)=FLEX( IJ)-TEMP (IJ)+VECTOR(I J)
-clo8 -
</P>
<P><PB REF="00000109.tif" SEQ="00000109" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="821" N="109">
Example Problem No. 84
FORTRAN Programs and Data, Continued
C
C FORM MASS MATRIX
C
DO 114 II=1iNP2
DO 114 JJ=1.NP2
IF(II-JJ) 112,111,112
111 ZMASS(IIJJ)=l./(EN+1.)
GO TO 114
112 ZMASS(IIJJ)sO.
114 CONTINUE
ZMASS( 11)=O05*ZMASS(1,1)
ZMASS(NP2,NP2) 0.5*ZMASS(NP2,NP2)
RETURN
$COMPILE FORTRAN
C
C WRITE OUTPUT SHAN S. KUO 7/29/61
C
SUBROUTINE RITE(FLEXZMASSVECTORTEMP)
DIMENSION FLEX(5,5),ZMASS(5,5),TEMP(5,5),VECTOR(5,5)
COMMON N,SPANBETA
K=N+2
DO 157 I=1,K
WRITE OUTPUT TAPE 5,158tI,(FLEX(IJ), J1lK)
157 WRITE OUTPUT TAPE 6,158.I,(FLEX(IJ), J=lK)
158 FORMAT(23H FLEXIBILITY MATRIX ROWI3/(1H 5F10.5))
DO 701 I=1,K
WRITE OUTPUT TAPE 5711,I*(ZMASS(I9J),J=lzK)
701 WRITE OUTPUT TAPE 6,711,I,(ZMASS(IJ)*J=1~K)
711 FORMAT(16H MASS MATRIX ROWI3/(1H 5F9.5))
C
C FORM TEMP MATRIX = SQ. RT. OF MASS MATRIX
C
DO 401 I=1iK
DO 401 J=1,K
IF(I-J) 412,411,412
C IN THE UNIV. OF MICHIGAN USE SORT NOT SQRTF
411 TEMP(I,J)=SQRT (ZMASS(I,I))
GO TO 401
412 TEMP(I,J)=0O
401 CONTINUE
C
C FORM H MATRIX = TEMP TIMES FLEX TIMES TEMP
C STORE FLEX TIMES TEMP INTO VECTOR ( TEMPORARY USE)
C
CALL MATMPY(FLEXKKTEMPKVECTOR)
CALL MATMPY(TEMPKK,VECTOR,KFLEX)
DO 402 I=1,K
WRITE OUTPUT TAPE 5,403,I,(FLEX(IJ),J=1,K)
402 WRITE OUTPUT TAPE 6,403,I,(FLEX(IJ),J=1,K)
403 FORMAT(13H H MATRIX ROWI3/(1H 5F10.5))
C
C EIGN IS A SUBROUTINE AVAILABLE IN THE UNIVERSITY OF MICHIGAN C.CO
C
S=EIGN(FLEXK,1,VECTOR(1,1),F)
IF (S-1.5) 13,13,14
13 WRITE OUTPUT TAPE 6,12,5
12 FORMAT(15HOBAD RETURN S =F10.6)
GO TO 777
14 IF(S-2.5) 15.15,16
15 WRITE OUTPUT TAPE 6,11,F
11 FORMAT(15HOSCALE FACTOR =F10.6)
C
C GET TRUE VECTORS
C
16 DO 404 I=1,K
404 TEMP(I,I)=1./(TEMP(II))
C
C USE MASS AS A TEMPORARY STORAGE
C
CALL MATMPY(TEMPKKVECTORKZMASS)
DO 405 I=1tK
WRITE OUTPUT TAPE 5,406*I,(ZMASS(IJ),J=1,K)
405 WRITE OUTPUT TAPE 6,406,I,(ZMASS(IJ),Ju1lK)
406 FORMAT(16H TRUE VECTOR ROWI3/(1H 5F10.5))
-C109 -
</P>
<P><PB REF="00000110.tif" SEQ="00000110" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="748" N="110">
Vibration of Beams on Spring Supports
FORTRAN Programs and Data, Continued
C
C
C
TRANSFORM THE EIGENVALUES TO FREQUENCIES
DO 150 I=1,K
150 FLEX(II)=SQRT (l./FLEX4I,Il)
DO 177 I=1,K
WRITE OUTPUT TAPE 5,178,IFLEX(I.I)
177 WRITE OUTPUT TAPE 6,178,I,FLEX(II)
178 FORMAT(28H FREQUENCY IN SQ. RT. OF K/MI6/(1H F30.6))
777 RETURN
C
C MATRIX MULTIPLICATION
C
SUBROUTINE MATMPY(ANoMBLC)
DIMENSION A(5,5) B(5,5) C(5,5)
DO 75 I=1,N
DO 75 J=1tL
C(I,J)=0.
DO 75 K=1,M
75 C(IJ)=C(IoJ)+A(I K)*B(K#J)
RETURN
$ DATA
030302
0.5 4,
Computer Output
VIBRE TIO t OF BERM i-OH SFRIHG SUPPORTS. FPROE:LEt OF SHOCK ISOCLATIOH
RIGID BODV FREQUEHC'.'
BETT = 0. 190
FLE>',IBILITYV MFPTRIX ROWd
1.00000 0.75000
FLEXIBILITY M1RTRIX FRO
0.75000 0. 76769
FLEXIBILITY MATRIX OWl
'.50000 0.6744
FLEXIBILITY MATRIX ROWil
2. 25000. 485 9:
FLE IBILIITY MATRIX ROW
Q. 00. 25C00'
MRSS MATRIX ROW 1
0.12500 0. 0.
M SS MATRIX ROW 2
0. 0.25000 0.
MASS MATRIX ROW 3
0!. 0. 0.
MiSS MATRIX ROW 4
O. O. 0.
MASS MRTRIX ROW 5
0. 0. 0.
H MATRIX ROW 1
0.12500 0. 1325
H MPTRI ROW 2
O. 13 258 0. 1 9192
H MPTRIX ROW 3
0.08839 0.160
H M'TRIX ROW 4
0.04419 0.12150
H MRTRIX ROW 5
0i.. 0.0441
TRUE VECTOR ROW 1
1.60933 -1.18746
TRUE VECTOR ROW 2
0.83966 1.1 3797
TRUE VECTOR ROW 3
-0. 00000 -0. 00000
TRUE VECTOR ROW 4
-0. 83966 -1.1 3797
TRUE VECTOR ROW 5
-1.60933 1.18746
FUDARMENTPL FREQ. OF HINHED-HINGED EEAEl= 0.50
1
0.67440 0.48598
3
0.75367 0.67440
4
0.67440 0.76769
5
0.50000 0. '5000
O. 0.
0. 0.
25000 0...,..-._,
0.25000. 0.
0.
0.25000
0.50000
0.75000
1.00000
0. 0.12500
0.08839
0. 1 6860
0. 18842
0.16860
0.08839
-1.71950
0.16622
0.99402
0.16622. 04419 0.
0. 2150 0.0441.I6860 O.08839
-L. 191 92 -... J-L1c 325 -,:.i
1325 0. 12500
0. 33 96 0. 58978
1.02312 -0. 9208
1.10045 1.34 199
1.02312 -0. 96208
-1.71950 0. 3396 0.58973
-C110 -
</P>
<P><PB REF="00000111.tif" SEQ="00000111" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="588" N="111">
Example Problem No. 84
Computer Output, Continued
FREQUENCY IN SQ. RT. OF K/M 1
2.29S4......_.._________ ----------. --- —-------------
FREQUENCY IN SQ. RT. OF K/M 2
-_ ~_ -.. -13- 8.55_16D-_______ --- —---------— _
FREQUENCY IN SQ. RT. OF K/M 3
5_ 91 5994
FREQUENCY IN SQ. RT. OF K/M 4
1. 294614 __ —___ _ --------------------------------
FREQUENCY IN SQ. RT. OF K/M 5
— ___ ___-_ _ _ ___ _24. 901JiS4 --- —------------— _______- — __________
_ _v IBRIT.I.Ot _OE.._BERMJ_ S._ PRI NG S_. PTS PBLEM SUP RTS_ P B LEE HO CKLSO _LLO N _
RIGIn BODn FPFRUFNCY / FlINDAMFNTrt FPFO- nF HTNffn-TMTNGFn P.FQM- 4. 0
BET..I2._76 __ --- —---- _.7 --- —------ --------------
FLEXIBILITY MATRIX ROW 1. —._ Il....D._,-_7QSt-__ — O5.. Q0o ---_ Qno --- _-______________ --- — --------------
FLEXIBILITY MATRIX ROW 2. 75noo 9. 7S71 n 11. 6, 1 4i. 7 47775 nrt. ^nnn
FLEXIBILITY MATRIX ROW 3..-__-____Q5Q-__O.. LL._6_J_4. - __ J73_45 _1~6J- ___._-__.LQ ------------------— ____ — ___
FLEXIBILITY MATRIX ROW 4. —_ _ __2..50O_..-.- 7.-_47775___1 --- L46____O675n — __ -- ---------— 5LO-
FLEXIBILITY MATRIX ROW 5
0c 0. 250-Inn n. snnn n- 70innn 1 nnnn
MASS M.TRIX ROW 1
_....... 0 s- _10.. -.....0 ---— a --------------------— ______________
MASS MATRIX ROW 2
-_-....... O....25Q....-__-___._.Q._._.-.._____,____ _______ -------
MASS MATRIX ROW 3
0. n..o25n00 n n,
MASS MATRIX ROW 4
-. -. _..............-____Q.__________ _____________ --- —----------— __ ------
MASS MATRIX ROW 5............. _......... _.. J2.. 25Q- ---— ______________ _ -----
H MATRIX ROW 1
0. 12500 0.13258 0.08839 0.0441 - 0.
H MATRIX ROW 2.13258.6 -- -...2.4...92_2. 536_. ------- _,_ --- —-__1_-_ t-9
H MATRIX ROW 3
_. 0.08839 2. 91 536.4_.!i1 Z -S6._ --- —------- ____________
H MRTRIX ROW 4
__. 04419 1 86944 2.9153.?-4q98 II 1 9.q
H MATRIX ROW 5
0. 0.0 449....... 0.889....... _.1-325.. -L-L2.l-= — -------------
TRUE VECTOR ROW 1
1.96434 0.37601 0.....0517_-. 7-.L7 r79.J4.5..__0_____ --- —--
TRUE VECTOR ROW 2
-0. 26588 _ -. 8- i 0. 07.,o 2 -? -7nin -n. 7A ---TRUE VECTOR ROW 3.2. 0.000000 -0. 000000. 1. 283. — _ -l -D6.36__ t.l.............. - -------
TRUE VECTOR ROW 4
0. 26588 -1. 38-00 1.00732..7.. _ 0 10 -.._..3-,i _..._...........
TRUE VECTOR ROW 5
o -.96434....-_O_.360_ n-_ n.....5.l_7EL 1_ - 1.- 3.7979, 44 r:,9~
FREQUENCY IN SQ. RT. OF K/'M 1
3, 041763.. _.-._...... _.....
FREQUENCY IN SQ. RT. OF K/M 2
1,305481............... —. --- —---
FREQUENCY IN SQ. RT. OF K-NM 3
0.. 345544........
FREQUENCY IN SQ. RT. OF K/.M 4
2. 530636
FREQUENCY IN SQ. RT. OF K M 5
3.407844
-Clll
</P>
<P><PB REF="00000112.tif" SEQ="00000112" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="752" N="112">
Vibration of Beams on Spring Supports
Discussion of Results
The natural frequencies are equal to the computed values times the square root
of k/m. The modal shapes are normalized in the subroutine BIGN according to the
following relation:
n+2. m. i 2= m j = 1, 2...,n+2
Figure shows the computed normal modes for two different alues.
Figure 2 shows the computed normal modes for two different &6R/aH values.
1 44 RH
i1 = 1.295 k
M2 =2.293
2 V
1 = 0.346 -
1 m.37 () = 1.305
2 m
is 1.377
t.ct~ ooc - '
3 = 5.916 k
3 ~~m
M3 = 2.531
3 Vm
/
\ I /
- 1. oS
I 9sll
4 = 13.855
a5 = 24.902 m
5
4) = 3.042 0
5=g 3.408 i
Figure 2. Normal modes and frequencies
__________ c R H = 0.5
----- - 4.0
-C112 -
</P>
<P><PB REF="00000113.tif" SEQ="00000113" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="892" N="113">
Example Problem No. 84
Critique
A problem of this type is of great importance in the design of underground
structures to isolate shock loads. Many textbooks on vibration do not mention
Jacobi's method. With the advent of high speed digital computers, this method
becomes very important. Its algorithm is simple and all the eigenvalues and eigenvectors are computed at the same time. Further extensions of this problem include
the computation of bending moment matrix, shear matrix, reaction coefficients, and
the moment-time relation for a given shock function.
VII. REFERENCES
1. Organick, E. I., A Computer Primer for the MAD Language, University of Houston, Computing
and Data Processing Center, 1961.
2. Organick, E. I., A Computer Primer for the FORTRAN Language, University of Houston,
Computing and Data Processing Center, 1961.
3. Project on the Use of Computers in Engineering Education, First Annual Report, The
University of Michigan College of Engineering, Ann Arbor, August 2b, 1960.
4. Project on the Use of Computers in Engineering Education, Second Annual Report, The
University of Michigan College of Engineering, Ann Arbor, December 15, 1961.
5. Westervelt, F. H., Use of Computers in Mechanical Engineering Education, Project on
the Use of Computers in Engineering Education, The University of Michigan College of
Engineering, Ann Arbor, June 1, 1962.
6. Wilson, R. C., Use of Computers in Industrial Engineering Education, Project on the Use
of Computers in Engineering Education, The University of Michigan College of Engineering,
Ann Arbor, January 27, 1962.
-C113 -
</P>
<P><PB REF="00000114.tif" SEQ="00000114" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="100" N="114">
</P>
<P><PB REF="00000115.tif" SEQ="00000115" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="100" N="115">
</P>
<P><PB REF="00000116.tif" SEQ="00000116" RES="600dpi" FMT="TIFF6.0" FTR="UNSPEC" CNF="200" N="116">
)r^
teo
N
4 _~s~sY
</P>
</DIV1>
</BODY>
</TEXT>