**Today's Topics:**

- Fortran Code Generation by Computer Algebra Systems
- Temporary Address Change for Hans Mittlemann
- Models of the Human Head or Ear
- Transcendental Functions on the 80387
- Inverse Square Root of a Matrix
- International Conference on Applications of Supercomputers
- Floating-Point Indoctrination by Kahan

From: Heikki Apiola <APIOLA%FINFUN.BITNET@forsythe.stanford.edu>

Date: Sun, 3 Apr 88 20:37 O

I would like to know about existing computer algebra based (especially

MACSYMA) Fortran code generation packages and their availability.

There have been articles about GENTRAN by Paul S. Wang in J. of Symbolic

Computing.

Please reply to:

Heikki Apiola

e-mail: apiola@finfun.bitnet

Scientific Computing Centre

State Computer Centre

PO-BOX 40

SF-02101 Espoo

Finland

------------------------------

From: Hans D. Mittelmann <aihdm@asuacad.bitnet, na.mittelmann>

Date: 8 Apr 88 19:35 MST

My address from 1 May 1988 until 31 July 1988 will be:

Institut fuer Angewandte Mathematik

Universitaet Erlangen-Nuernberg

Martensstrasse 3

8520 Erlangen

West Germany

HJW@DERRZE0.BITNET or NA.MITTELMANN@NA-NET.STANFORD.EDU

--Hans Mittelmann

------------------------------

From: Jim Buchanan <ACAD8005%RYERSON.BITNET@forsythe.stanford.edu>

Date: Tue, 12 Apr 88 16:40:24 EST

I am looking for information on mathematical models of the human head

in general, or of the ear in particular. Does anyone have any

information in this area? Thanks in advance.

--Jim Buchanan

------------------------------

From: Warren E Ferguson <smu!ferguson@uunet.UU.NET>

Date: Tue, 12 Apr 88 15:31:27 CDT

I am thinking of having students in one of my classes create a partial

emulator of the INTEL 80387 chip. After looking at the spec sheets for the

80387, it appears that INTEL is using CORDIC's to compute their special

functions. Does anyone know of literature describing which CORDIC's INTEL is

using for these functions? Did members of the IEEE Binary Standard ever

suggest how the Standard could be used to compute the basic transcendental

functions?

[Editor's Note: See the announcement of Kahan's lectures at Sun

later in this Digest.]

Warren Ferguson, Southern Methodist University

na.ferguson

ferguson%smu.uucp@uunet

------------------------------

From: Cleve Moler <na.moler@na-net.stanford.edu>

Date: Tue, 12 Apr 88 22:14:26 PDT

I was asked recently if I had any software to compute the inverse

square root of a symmetric, positive definite matrix. All I had to

suggest was EISPACK and the obvious diagonalization algorithm.

Has anybody got anything faster? The setting involves lots of

matrices of small order, say less than 10-by-10.

--Cleve

------------------------------

From: [Three levels of forwarding]

Date: Fri, 08 Apr 88 14:05:16 -0800

Announcement of ASE 89

The International Conference on Applications of Supercomputers in

Engineering will be held 5-7 September 1989 (yes 1989) at Southampton

University, United Kingdom. This International Conference is convened

to increase the awareness of the potential of the new supercomputers

amongst scientists and engineers.

Papers are invited on the following topics, and the abstracts should reach

the Organizers by 30 December 1988. Topics include:

* computer architecture and current supercomputers

* use and development of transputers

* case studies in

hydrodynamics

groundwater

environmental sciences, etc.

* new parallel algorithms for finite differences

and finite and boundary element methods

* optimization problems

* structural and dynamic analysis including crash and impact

simulation applications in automotive engineering

* use of supercomputers in CAD and CIM

* benchmarks and comparison between different processors

* fluid dynamics and atmospheric science applications

* expert system applications

* future trends.

For further information, please contact:

Liz Newman

Conference Secretary

Computational Mechanics Institute

Ashurst Lodge

Ashurt

Southampton, S04 2AA

England

Telephone: (042129) 3223

Telex: 47388

Attn: COMPMECH

Fax: (042129) 2853

------------------------------

From: David Hough <dgh@Sun.COM>

Date: Fri, 8 Apr 88 17:47:25 PDT

Sun Microsystems presents a public series of lectures:

Computer System Support for Scientific and Engineering Computation

Prof. W. Kahan

University of California

The environment for computing

Comparisons of floating-point arithmetics in computers

Programmers' models of floating-point arithmetic

Comparison of computers' floating-point register architectures

Real elementary transcendental functions

Complex arithmetic

Exception handling

Numerical lapses in standard programming languages

A detailed syllabus follows this announcement.

FORMAT

The lecture series will start at 8:00 AM on Tuesday May

3, 1988, in the Training Room at Sun Microsystems building

7, 2700 Coast Avenue, Mountain View, CA. From US 101, take

the Amphitheater Parkway exit in Mountain View and turn left

on Garcia. Turn right on Marine Way, then right on Coast

Avenue to its end.

The first lecture will provide an overview of the sub-

ject matter to be discussed in the course, and the last lec-

ture will summarize the material presented previously and

survey advanced topics and directions for further research.

Approximately 24 lectures will be presented Tuesdays and

Thursdays through approximately July 29. Each lecture will

start at 8:00 AM and last for two hours.

Lectures are open to the public subject to space avail-

able, and persons with experience in some of the lecture

topics are encouraged to attend and contribute to the dis-

cussion. The amount of space available will be determined by

the number of pre-registrants, so persons desiring to insure

a place should register in advance with the Sun Training and

Development Department using the form below.

Persons with training or experience comparable to that

represented by a B. S. degree in Computer Science should be

able to understand most of the conclusions. The mathematical

background required for a deep understanding of the lecture

topics varies considerably, especially for some of the

proofs. Chapter 4 of Knuth's Seminumerical Algorithms exem-

plifies the presentation level - accessible to undergradu-

ates, but containing interesting mathematics at all levels.

Note that the scope of the course is much broader than

Knuth's book.

Written lecture notes will be distributed to regis-

trants about a week after each lecture. These notes may

subsequently be commercially published to supersede the 1973

Implementation of Algorithms lecture notes distributed by

NTIS.

The first and last lectures will be videotaped; regis-

trants will receive a copy of the first for their personal

use.

This lecture series does not provide academic credit.

Registrants interested in obtaining formal recognition of

participation may present an acceptable individual or group

project. Outstanding project reports will be included in

the published course lecture notes. An example of such a

project is to reimplement some of the 4.3 BSD libm functions

using fast table-driven techniques. Prof. Kahan will be

available after lectures for consultation on projects.

Questions pertaining to the technical scope and depth

of the course may be addressed to David Hough,

dhough@sun.com, (415) 691-2268. Questions pertaining to

registration procedures should be directed to Janet Rath,

jrath@sun.com, Sun Training and Development, (415) 494-6058

or 494-8008 x 6058.

VENDOR LITERATURE

Vendors of semiconductor devices, calculators, comput-

ers, and software products oriented toward scientific and

engineering computation are encouraged to contribute techni-

cal literature describing their products for distribution to

attendees. Please contact David Hough after 22 April to

determine quantities required.

REGISTRATION - general

The registration fee guarantees a place and provides

lecture notes and a videotape of the initial lecture of 3

May. Registration by 22 April 1988 is recommended, but late

registrants will be accommodated subject to available space

through at least the second class meeting, Thursday 5 May.

The registration fee is $1000, subject to reduction if more

participants register than projected.

Send a check for the registration fee, payable to Sun

Microsystems, with the Registration Form below.

Full-time students at accredited institutions may

register without fee; they will receive no videotape. Such

students should submit current evidence of their status in

lieu of a check.

REGISTRATION FORM

Please send the following information to:

Sun Microsystems Training and Development

Attention: Janet Rath - Floating-Point Lectures

MS A6-05

2550 Garcia

Mountain View, CA 94043

Your name:

USENET or ARPANET address for electronic mail, if any:

Mailing address:

FULL-TIME STUDENTS ONLY: Your institution, department, next degree and projected date:

Telephone number during business hours:

Highest academic degree and major field:

Circle any of the following languages of which you have a working knowledge:

Ada APL BASIC C C++ Fortran Modula-2 Pascal PL/I

Name computers, if any, of which you have working knowledge of

assembly-language programming:

Name computers, if any, of which you have a working knowledge of the floating-point

architecture, such as number representations, rounding methods, etc.

What do you want to get out of this course?

SYLLABUS

The following details indicate the scope of the lec-

tures and point to references. While the focus will be on

general principles, particular areas of interest expressed

by attendees will also be discussed.

Lecture Plan

Lecture number Topic

1 Prospective overview from an elementary level (videotaped)

2-21 See details below

22-23 Student project reports

24 Retrospective overview from an advanced level (videotaped)

1. Environment for computing

Gresham's Law: "the bad drives out the good"; similarly the fast

drives out the slow even if the fast is wrong.

Computer Users: sophisticated but numerically naive.

Many layers of software from diverse sources.

Portable software in the public domain: LINPACK, EISPACK, ...

Commercial Libraries: IMSL, NAG, ... [Rice 1983]

Standards: for programming languages but little else.

Proprietary rights, copyright, patent, trade secrets.

L. J. Comrie's trick

Execute-only code, copy protection.

Liability and exculpation: law of torts

Reliability means conformity to reasonable expectations

Marketing claims vs. absolute and negligent liabilities

Obligation to correct errors; change vs. compatibility

2. Comparisons of floating-point arithmetics in computers

Conventional floating-point arithmetics: [Sterbenz 1975]

Radix, Word size, Range/Precision trade-off

Why is binary best? When is decimal best?

Commercially significant machines:

IBM 7090/7094, 360, 370, PC's

Univac 11xx

CDC 6600, Cyber 1xx, 205

Burroughs 6500

DEC PDP-11, VAX

Crays

Hewlett-Packard 3000

Calculators by Hewlett-Packard and Texas Instruments

IEEE Standard 754:

Chips: Intel 80x87, Motorola 6888X, ATT 32106

National 32081, Weitek 1164/5, AMD 29027,

Fairchild Clipper, Analog Devices 3212/22

BIT, Inmos T800, ...

Micros: IBM PC's, Sun-3 and Sun-4, Apple Macintosh II, ...

Mainframes: ELXSI 6400, HP Spectrum

Nonconventional representations of real numbers:

Logarithmic (Radix 1)

Variable exponent width - why are they all bad ideas? Demmel 1987

Matsui-Iri 1983, Hamada 1987

Symmetric level-index arithmetic, Clenshaw-Olver 1987

Floating-slash rational arithmetic: Kornerup-Matula 1987

Interval arithmetic: Moore, Herzberger and Alefeld

what it is and what good it is, briefly

Multi-precision schemes:

Conventional double and extended precisions

Numerical Turing language - Hull 1987

MP package for Fortran - Brent

MACSYMA Bigfloat - Fateman

Kulisch-Miranker paradigm, IBM ACRITH package 1983

Kahan's indefinitely extended IEEE 754/854

Implementation of algebraic operations +-*/ rem sqrt:

Fast carry propagation: fast prefix calculation, VLSI

Floating-point add/subtract, shifters, normalizers

Multiplication: iterative, Booth re-encoding,

Wallace tree of carry-save adders

Division: higher radix [Taylor]; preconditioned

iterative using fast multiply, can be grubby [Cray]

can be cleaned up and tested [Kahan 1987]

Why keep (division time)/(multiplication time) < 4

Interruptible remainder

Square root by hardware, by software, or without division

Rounding in theory: guard-round-sticky bits

rerounding without double rounding (C. Lee 1988)

carry-free roundings (chop, jam, ROM)

Rounding as practiced on commercially significant machines

DPROD, or multiply-add primitives

Conversions: of precisions, float-integer, BCD-binary

Exceptions: Traps/Faults, Flags, Defaults, Presubstitution,

Invalid operation, Division by zero, Over/Underflow,

as practiced on commercially significant machines

partial overflow on Cray

partial underflow on CDC Cyber 1xx

gradual underflow in IEEE 754/854

NaNs, infinities, and reserved operands

Range extension via over/underflow tallies [Kahan 1966]

Precision doubling, and further extensions [Pichat 1972]

Ultra-wide multiplication (Karatsuba, Cook, Schoenhage)

3. Models of floating-point arithmetic for programmers

Purpose of models: promote correctness proofs and portability

Axioms: non-categorical axiomatization is too loose

32 by Wijngaarden 1966

20 by Brown 1981

Axioms: over-axiomatization is too slow

Kulisch's ringoids, semi-morphism, optimality

Descriptions:

crude inequalities, Wilkinson 1959

parameterization, Sterbenz 1973

environmental inquiries in Fortran, ADA, ...

Prescriptions:

Computer manufacturer's manuals

IEEE standards 754/854

Tests for correct floating-point arithmetic:

Test batteries: the IEEE 754 test-vectors tape

Test patterns: Schryer 1980, NAG 1987

Inquiries and consistency tests: MACHAR Cody 1980, PARANOIA Kahan 1986

P-adic tests for correctly rounded sqrt: Kahan 1966

P-adic tests for correctly rounded division: Kahan 1987

Identities that persist despite roundoff

Weak monotonicity

Cancellation:

p-q is exact if 1/2 <= p/q <= 2

When does integer m = (m/n)*n rounded?

When does x = (x+y) - y rounded ?

When does sqrt(x*x) = |x| ?

Avoiding drift:

during conversion, Matula 1968

during iteration; Knuth-Reiser proof

Harry Diamond's theorem

Error analysis:

Conventional "forward" analysis (Hotelling)

"Backward" analysis (Turing, Givens, Bauer, Wilkinson)

Combined approaches: the web of relationships

Misconceptions corrected:

Rounding errors are NOT infinitesimal

Rounding errors are NOT random variables

Unstable computations are often free from cancellation

Cancellation is often a GOOD thing in equation solving,

extrapolation, divided differences

Simulation is easy, prediction difficult - removable singularities

Backward error analysis is an explanation, not a license:

accuracy not transitive under composition

Precision does not limit accuracy (Dekker, Pichat)

Why double+ precision is a good rule of thumb:

double zeros, normal equations, divided differences

What "unstable", "ill-conditioned", "ill-posed" mean; how to

help an ill-posed problem get well

Interval arithmetic vs. Running Error-Analysis

Examples:

When arbitrarily high precision might not help

Triangle with given side-lengths

Quadratic and Cubic equations: betrayed by books

Polynomial equations in general

Linear systems of equations and eigenproblems

Numerical quadrature and infinite series

Trajectory problems and boundary value problems

4. Comparison of computers' floating-point register archi-

tectures

One accumulator: IBM 7094, DEC PDP 10, GE/Honeywell 635

Stack-top: B5500, HP 3000

Small stack: HP calculators, Inmos T800, Intel 80x87

indefinitely long stack intended for 8087

Effect on procedures' argument passing and register saving

Extended width:

internally in HP calculators

Kulisch-Miranker's hidden Super-Accumulator

in registers: IBM 7094, DEC PDP-10, GE 635

IEEE 754: 8087, 68881, WE 32106

Orthogonal registers: IBM 370, CDC 6600, NS 32081, Cray,

DEC VAX, ELXSI 6400, Weitek, AMD, Analog Devices, etc.

Vector registers: Cray; simulated vectors in CDC Cyber 205

Fast vector arithmetic without vector registers.

Influences on programming languages' expression evaluation:

Precisions available

Mixed-precision expressions:

Fundamentalist's Fortran for IBM 370, DEC VAX

Widest available for C on PDP-11; 8087; 68881

"scan for widest" seems best [Corbett]

Conflict between compiler's pseudo-machine and actual registers

5. Real elementary transcendental functions

Argument reduction for log, exponential, and trig functions

How to simulate 2/pi to infinite precision, or not.

CORDIC algorithms = argument reduction carried to its limit

Polynomial, rational, and other "best" approximation

The Remez algorithm, with compensation for roundoff

Interpolation in BIG tables: IBM 370, DEC VAX VMS

J. C. P. Miller's trick to suppress roundoff: Gal

"Kernel function" approach; monotonicity achieved thereby

Special properties for special functions:

Weak monotonicity: EXP(X+Y**2) >= EXP(X)

Cardinal values: cos(0) = 1, |cos(x)| <= 1 .

Huge arguments: cos(3083975227) = -.00000000007486714572...

Removable singularities: does 0**0 = 1 ?

Unremovable singularities: ln(0)=-infinity; ln(-1) is NaN.

Accuracy claims and proofs: the Table-Maker's Dilemma

Testing: Cody 1980, Liu 1987

6. Complex Arithmetic

Mixed Real and Complex expressions

Multiply/Divide despite exceptions in concurrent environments

Different approaches in infinity

The sign of zero: functions like sqrt on slitted domains

Accurate elementary transcendental functions

Peter Tang's complex Remez algorithms

7. Exception handling

An Exception is not an Error unless handled badly

Why we must neither always ABORT nor always merely CONTINUE

Retrospective Diagnostics resolve the dilemma.

Let's not get enmeshed in the operating system's kernel.

Concurrent/parallel/pipelines machines restrict the options:

how can we pause and resume during debugging?

Infinity and NaN/Indefinite/Reserved operand allow CONTINUE option.

Presubstitution, and some of its applications.

Flags permit recovery from ill-presubstituted exceptions.

Circumventing over/underflow in extended products/quotients:

Kahan's 1966 implementations on 7094 and 360.

Barnett's 1987 implementations on VAX and Sun-3.

The problem of Scope:

solved without language support by Apple SANE 1986

with minimal language support:

APL's localization of system variables

proper language support achieves minimal run-time cost:

diminish incidence of forgotten SAVE/RESTOREs

avoid unnecessary SAVE/RESTORE in leaf subprograms.

Comparisons with:

"ON condition" in PL/I and Basic

"drop-through" in ADA.

"Error-Handlers" in DEC VMS Fortran.

Signals thrown and caught by Unix.

Augmented data-structures for Data-Flow.

8. Numerical lapses in standard programming languages:

Misguided attempts:

to define the semantics of arithmetic

Brown's model in ADA, Fortran 8x, ANSI C

Ambiguous environmental inquiries

to define occurrences like 1.0/0.0 and 0.0**0.0 as Errors.

to ignore integer overflow

Running roughshod over parentheses (C)

Moving code out of loops despite side-effects.

Other would-be optimizations like

misuse of register variables causing schizophrenia

misuse of "0-x" instead of "-x", or

"-(x-y)" instead of "y-x" , and

" x-y > 0" instead of "x > y".

Pascal's fixed array-dimensions, lack of double precision,

inadequate provision for precompiled libraries

C's awkward multi-subscripted arrays inside procedures

Fortran's unprotected argument lists

Inability to overload arithmetic operators to introduce:

complex arithmetic

matrix arithmetic

interval arithmetic

Awkward treatment of

"Eureka" problem (Break)

"Lost search-party" problem (compare Exception Handling)

Mixed-precision expressions (Corbett)

Re-running a program with higher precision

Mixed real-complex expressions

Inexact constants, inaccurate conversions

Output too wide for its field

What else should compiler writers know about floating point:

Farnum's 1988 paper

------------------------------

End of NA Digest

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

-------