Subject: Draft Report of Numerics Working Group

From: [email-masked] (Ronald F Boisvert)
Date: Fri, 12 Jun 1998 17:54:30 -0400

Attached is a draft of the report of the Java Grande Numerics Working
Group.  This draft was extensively revised by Marc Snir based on an
original that I put together after the last workshop.

Please submit any comments that you have to this mailing list ASAP.
We plan to revise this next week and circulate it more widely, still
as a draft.

(NOTE: Numerics Working Group members have already received a 
       copy of this.  I am simply increasing its circulation
       to the entire Forum.)

Ron Boisvert

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

                                Draft

                  REPORT OF THE NUMERICS WORKING GROUP

                             May 26, 1998

If Java is to become the environment of choice for high-performance scientific
applications, then it must provide, for large scale floating-point
computations, performance comparable to what is achieved in currently used
programming languages (C or Fortran).  In addition, it must have language
features and core libraries that enable the convenient expression of
mathematical algorithms. The goal of this working group is to assess the
suitability of Java for numerical computation, and to work towards community
consensus on actions which can be taken to overcome deficiencies of the
language and its run-time environment.  In this report, we present preliminary
findings of the working group.

I. CRITICAL JAVA LANGUAGE AND JAVA VIRTUAL MACHINE ISSUES

We begin by outlining critical issues that impede Java's effectiveness in
applications that are dominated by the use of floating-point arithmetic.
Unless these issues are satisfactorily resolved it is unlikely that Java will
be accepted by the numerical computation community.  This can impact the entire
Java enterprise by slowing the dissemination of high quality components for
solving commonly occurring mathematical and statistical problems. For each
issue we present solutions recommended by the working group. In selecting such
solutions, the working group has been careful to balance the needs of the
numerical community with those of Java's wider audience. The proposed solutions
require additions to the current Java and JVM design.  We have tried to
minimize the changes required in Java, relying on compiler technology, whenever
feasible. This minimizes the changes that affect all Java platforms, and enable
implementors to optimize for high numerical performance only in those
environments where such an effort is warranted.

ISSUE 1 : Complex arithmetic

REQUIREMENT:  The Complex field is an essential tool in the analysis and
solution of mathematical problems in all areas of science and engineering.
Thus, it is essential that the use of complex numbers be as convenient and
efficient as the use of floats and doubles.

The obvious solution is to develop a straightforward complex class with
methods for each arithmetic operation and use such objects as needed.  There
are several reasons why this approach fails.

  (a) The object overhead of complex methods makes them unacceptably
      inefficient.

  (b) The semantics of complex objects are different than those of floats
      and doubles.  For example, the = and == operators manipulate references
      rather than values.  Such differences lead to many errors.

  (c) Use of method calls for elementary arithmetic operations leads to
      inscrutable code which is very tedious to write and debug.
      Users would simply stay away.

The ideal solution is to add new base complex types to the language on par with
float and double.  This, of course, requires a significant change in the
language and in JVM to satisfy the needs of a relatively small community. In
particular, it requires the addition of a significant number of new opcodes.

An alternative solution acceptable to the working group requires that the
following actions occur:

  (1) A complex arithmetic package be developed and included as a core Java
      package, perhaps as a subpackage in Java.math.  Such a package will
      support assignment by value, and standard arithmetic operations and
      relations on complex values.

  (2) Use of these classes is made as efficient as float or double. This
      may require an extension to Java and JVM in support of lightweight
      classes; see following sections. This also requires the cooperation of
      compiler writers to use the opportunity provided to generate efficient
      code.

  (3) Operator overloading can be used to bind natural notation for
      arithmetic, logical and assignment operators to the methods of the
      complex classes; see following sections.

This alternate suite of changes requires fewer changes to Java and JVM, but
will require more compilation efforts for an efficient implementation.
Lightweight classes and operator overloading are general mechanisms which can
satisfy the needs of many groups for alternate arithmetic systems, such as
interval and multiple precision.

In addition to complex arithmetic, the Java.math library should be extended to
support complex transcendental functions.

The current proposal assumes that complex numbers are pairs of doubles. It is
deemed acceptable (at least initially) not to support complex numbers with
float components.

ISSUE 2 : Lightweight classes

REQUIREMENT:  Implementation of alternative arithmetic systems, such as
complex, interval, and multiple precision requires the support of new objects
with value semantics. Compilers should be able to inline methods that operate
on such objects and avoid the overheads of additional dereferencing.  In
particular, lightweight classes are critical for the implementation of complex
arithmetic as described in issue 1.

A lightweight class is final.  It holds a value and supports deep assignments
and deep comparisons (who work on the object value, not its reference).

  o  a.assign(b), assigns to object a the value of object b.
  o  a.equals(b), tests that objects a and b have the same value.

Lightweight objects will usually support additional unary and binary
operators.

Note that there is no requirement for lightweight classes to be immutable (i.e.
instance variables need not be final).  Immutability would lead to unnecessary
copying, preventing the updating of array of objects in place, a key need of
the numerical community.

There are two alternatives to the design of lightweight objects.

  (1) Lightweight objects are new types of objects in the Java language.  They
      are explicitly declared as such.  A lightweight  object is always
      accessed by value, as if it had a primitive Java type. : e.g., the value
      of a complex variable is a complex number, not an object reference.
      Instances of lightweight objects can be assigned (using the assign
      operator) and can be compared (using the equals operator).  On the other
      hand, if v is an instance of a lightweight object, then v = null or v ==
      null are illegal expressions.  Lightweight objects are passed by value in
      method invocations.  Methods on such objects cannot be synchronized. It
      is expected that compilers will inline invocations of standard methods
      (assign, equals, ...) on such objects.

  (2) Lightweight objects are regular Java objects, and are accessed by
      reference.  Since they are final, and since the Java (back end) compiler
      has full knowledge of the semantics of the methods applied on these
      objects, it is expected that compilers will inline invocations of assign,
      equals and other predefined methods.

The first scheme is likely to lead to the best performance: one always saves
the storage required for a Java object: a complex will always require two words
of storage, no more. No additional referencing is needed to access an object
value.  Garbage collection for such objects is simplified.

On the other hand, this scheme seems to require significant changes in the Java
language and the JVM.  For example, the JVM instruction set does not support
method invocations that return non-scalar values. It does not support arrays
with entries that are not of a primitive or a reference type. Modifications
will also be required in the Java program verifier.

The second scheme is more dependent on compiler optimization techniques for
performance: a compiler will generally be able to inline invocations of
predefined final methods; however a lean storage layout, that holds data but no
object descriptor, can be generated only if the compiler can determine that the
object is not accessed by reference; changes in the garbage collector might
also be required in order to support a lean layout.

The second scheme has the added advantage that it supports both deep and
shallow assignments or comparisons on lightweight objects.  Thus, it provides
functionality equivalent to this achieved in C or Fortran by the use of
pointers. In this scheme, lightweight objects are always passed as reference
arguments; this provides more flexibility (e.g., allowing a method to return
multiple values).

Finally, the second scheme does not require changes in Java or JVM
specifications.

There are several alternative design points that should be evaluated.

   o  If the use of lightweight objects is restricted to predefined
      classes, such as Complex, then the inlining could be done by the
      front-end compiler.  However, such an approach does not extend to
      user-defined lightweight classes. It makes harder optimizations by
      back-end (dynamic or static) compilers.  Such a choice also impacts
      debuggers and other tools.

   o  An alternative approach to lightweight objects is to treat them
      as regular Java objects that must obey certain restrictions;
      e.g., no (reference) assignment; lightweight objects cannot be
      components of regular Java arrays or regular Java structures
      (they can be components of the special rectangular Java arrays,
      defined below). This approach would still require Java language
      extensions (lightweight classes have to be declared as such, and
      the restrictions need to be spelled out).  However, no changes
      are needed in JVM (beyond carrying in the class file an
      attribute that marks the class as lightweight).  Back-end
      compilers might be able to better optimize code using such
      lightweight objects, because of the added constraints.

ISSUE 3 : Operator overloading

REQUIREMENT: Usable implementation of complex arithmetic, as well as other
alternative arithmetics such as interval and multiprecision, requires that
code be as readable as those based only on float and double.

Operator overloading is the obvious solution to this problem.  Without it,
codes implementing complex arithmetic would be extremely difficult to develop,
understand and maintain.  Such a code will look very different than similar
code using real arithmetic, thus burdening library developers.  E.g., a simple
statement such as

a = b+c*d

will be replaced by

a.assign(sum(b, product(c,d))

Without operator overloading, a large portion of the scientific computing
community would choose to avoid Java as being too unfriendly.

Only a limited facility for operator overloading is necessary to fulfill this
requirement. If the first scheme is used, so that an object supports either
shallow or deep assignments, but not both, the assignment operator = will be
overloaded to signify deep assignment:  a = b is interpreted as a.assign(b)).
Similarly, a == b is syntactic sugar for a.equals(b), if a and b are
lightweight objects.  On the other hand, if the second scheme is used, then new
operators are needed for deep assignments and comparisons.  Thus a <- b is
syntactic sugar for a.assign(b) and a === b is syntactic sugar for
a.equals(b).  (It may be desirable to introduce new operators even if the first
scheme is used, so as to reduce confusion.  It may also be desirable to allow
'<-' to be used for assignment of primititive types, for consistency.)

Thus, if the second scheme is used, then

Complex c = new Complex(0.0, 0.0);
Complex d = new Complex(1.0, 1.0);
Complex e;
e = c;
c <- d;
System.out.println(e.real + "," e.imag);

prints:

1.0, 1.0

The arithmetic, assignment and logical operators must be extendable by
overloading, using their natural notation.  One may have a predefined naming
scheme for methods that overload existing operators: 'sum' for '+', 'product'
for '*', etc.  It is not necessary to admit the introduction of new operators
(beyond assign and equals). Thus, no new syntax is required for operator
overloading,  except that predefined operators apply to lightweight objects. An
expression of the form 'a+b' is merely syntactic sugar for sum(a,b) (and is
illegal if the 'sum' method is not defined on a and b).

There are several alternatives that should be evaluated.

  o  Operator overloading may be restricted to predefined lightweight
     objects (such as complex) or extended to user-defined lightweight
     objects, or extended to arbitrary classes. It is reasonable to
     couple operator overloading with lightweight objects: the
     language conveys the right intuition by using operators for
     'cheap' operations and method invocations for expensive
     operations.

  o  Binary operators, such as 'sum' can be restricted to the case were
     both operands are of the same type, or extended to operands of
     distinct types.  Type promotion is very natural in many cases
     (e.g., real*complex), and should be supported.

ISSUE 4 : Optimal use of hardware

REQUIREMENT: The high efficiency necessary for large-scale numerical
applications requires aggressive exploitation of the unique facilities of
local floating-point hardware.

The current insistence of bitwise reproducibility of results on all JVMs makes
it impossible to satisfy this requirement.  Efficient processing of Java
programs requires that compilers and JVMs provide the _option_ to

  (a) use IEEE extended arithmetic hardware anywhere in the computation

  (b) use the associative law to rearrange to order of computation

  (c) use possibly unsafe identities to eliminate computations

The use of IEEE extended arithmetic in intermediate computations can improve
the accuracy and reliability of numerical results.  Processor with hardware
support for IEEE extended should NOT be required to round intermediate
results.  This slows the computation and makes it less accurate.  There should
not be a requirement that such a store be forced on each assignment statement
in a user's program. Rounding should be done by the compiler when necessary.
For example, on a machine with extended precision registers, rounding should
occur only when registers must be spilled to memory.  On a machine with fused
multiply-add, a multiplication followed by an addition should always be
replaceable by a fused multiply-add.

The associative law can be used by optimizing compilers to reorder arithmetic
operations in order to make more efficient use of hardware.  Such optimizations
are crucial to improving the performance of numerical codes, and users should
have access to this technology. A typical example for such an optimization is
the use of associativity to execute a reduction in parallel.

On the other hand, bitwise reproducibility is important for code testing and is
needed in many environments. Users may want to ensure strict reproducibility by
enforcing the default Java model.  Programmers may want to disable unsafe
optimizations for selected codes in order to have better control on the
execution (e.g., guarantee better precision).

We do not feel that large changes to the Java specification are needed to
satisfy the needs of the numerical community on this issue.  In particular,
the current JVM specification for bitwise reproducibility can remain the
default behavior.  Instead, what is needed is

  (a) JVM flags that allow the _user_ to select efficiency over
      reproducibility at runtime.

  (b) a class and method modifier, StrictNumerics, which specifies that
      the given class or method must adhere to the more restrictive
      Java arithmetic specification, regardless of flags which may be
      set by the user.

The number and semantics of runtime flags should be left up to the JVM
developer.  The important features of this proposal are that

  o Strict Java semantics apply, by default.

  o JVMs may provide environment flags to overrule strict Java
    semantics.  These would be similar to compiler optimization options.
    Thereby, the user can decide to select efficiency over
    reproducibility at runtime.

  o Developers can shield critical segments of code where these
    relaxations should never occur, using the StrictNumerics attribute.

The existing proposals for LooseNumerics and IdealizedNumerics is unnecessarily
complicated -- we do not feel that code developers will want this level of fine
grain control.  The implementation mechanism would be the same as in these
proposals: the StrictNumerics attribute is carried in the class file, and
observed by back-end compilers.  Also, while optimizations may change the
numerical outcome of a computation, or even cause a Not a Number value to be
returned, rather than a regular value, the optimizations should still preserve
the "precise exception" model of Java: null pointer of index out of bound
exceptions should occur in the optimized code in the same state as they would
have occured in the unoptimized code.

ISSUE 5 : Rectangular arrays

REQUIREMENT: Operations on multidimensional arrays of base types must be easily
optimized.  In addition, the memory layout of such arrays must be known to the
algorithm developer in order to process array data in the most efficient way.

The performance of Java code can suffer from a deterioration of up to 25%
because of the lack of true rectangular arrays.  For native Java arrays code
generated for column traversal is less efficient because of pointer chasing.
Compiler elimination of run time tests for null pointers and out of bound
indices is harder if arrays can be jagged, or can change shape at run time.
More significantly, disambiguation is hard: even if two 2D arrays are not
identical they may still share a row. This forces compilers to generate
superfluous stores because of potential aliasing. Finally, a clearly defined
memory layout with guaranteed locality of data will allow developers to devise
algorithms which can be processed more efficiently.

We propose that standard Java classes be developed which implement
multidimensional rectangular arrays, and that these be included as subpackage
in Java.math.  These classes would store multidimensional arrays internally so
as to provide access that is as efficient as if the arrays were stored in a
canonical order (e.g., row-major).  The classes would support 1D, 2D, 3D, and
possibly 4D...7D arrays with Int, Long, Float, Double and Complex entries (a
different class is needed for each dimensionality and each element type --
since Java does not support templates).  The classes provide the following
methods.

  (a) Get and set to access and update an array entry.

  (b) Operations that correspond to Fortran 90 array intrinsics.
      In particular:

    (b.1) Operations to access the number of dimensions and
          the extends of an array.

    (b.2)  operations to reshape and transpose an array.

    (b.3) elemental conversion functions (e.g., the
          equivalent of REAL and AIMAG, that convert complex arrays
       into double arrays.

    (b.4) elemental transcendental functions

    (b.5) elemental boolean functions

    (b.6) array reduction functions (sum, minval, etc.)

    (b.7) array construction functions (merge, pack, spread, unpack)

    (b.8) array reshape function

    (b.9) array manipulation functions (shift, transpose)

    (b.10) array location functions (maxloc, minloc)

    (b.11) array scatter-gather and array scan operations (Fortran 95)

    (b.12) matrix multiply

Not all Fortran 90 and Fortran 95 operations are needed, upfront.  One can
likely do without elemental transcendental functions.

  (c) operations that correspond to array expressions (sum, scaling, etc.)

  (d) operations that create copies of or references to array sections.
      These operations allow one to copy subarrays (defined by
      subscript triplets or by vector subscripts) or to create
      references to such subarrays, thus supporting in place update of
      subarrays.  (As in Fortran 90, references to subarrays may be
      restricted to subarrays described by subscript triplets, so as
      to have succinct subarray descriptors.) A possible mechanism is to
      support the definition of index sets (or array shapes) and the
      extraction of a subarray defined by such an index set.

  (e) Operations to cast Java arrays into rectangular arrays, and
      vice-versa.

The array classes can be implemented with no changes in Java or JVM. However,
It is essential that the get and set methods be implemented as efficiently as
array indexing operations are in Fortran or in C.  We expect that inlining will
be used for this purpose, and that garbage collectors will recognize
rectangular arrays.  Multidimensional arrays are extremely common in numerical
computing, and hence we expect that efficient multidimensional arrays classes
will be heavily used.

Note that an array of complex entries need not be implemented as an array of
references (if Complex objects are regular Java objects), or as an array of
lightweight objects.  Rather, such an array can be implemented as an array of
doubles (with twice as many entries as the complex array).  The naive
implementation of the get method will access two double values and return a
(new) complex object -- better implementations will inline this code.
Additional methods will be provided to convert a complex array into a double
array (with twice as many entries), and vice-versa.

The inclusion of standard array classes in Java.math does not require any
change to the Java language. However, the use of explicit method invocation to
effect all array operations will significantly decrease the readability of Java
code, and incur the wrath of users. The introduction of a simple notation for
multidimensional arrays which maps to the standard array classes would make the
use of such arrays much more natural.  A multi-index notation, like a[i,j] to
refer to such array elements would be ideal.  This would allow statements like

   a.setElement(i,j,b.getElement(i,j)+s*c.getElement(k,l));

to be more naturally expressed as

   a[i,j] = b[i,j] + s*c[k,l];

Alternatively, one could reuse the bracket notation of Java, namely

  a[i][j] = b[i][j] + s*c[k][l].

The front-end compiler disambiguate the expression according to the type of a.
This requires changes in the Java language or (with the second alternative)
fancier operator overloading mechanisms.

Some alternatives that need be discussed:

  o Operator overloading may be applied to array arithmetic; e.g. A =
    B+C. This is nice, but not strictly necessary.

  o It would be nice to facilitate indexing operations by explicitly
    supporting triplet notation. This either implies new syntax, or
    fancy overloading of the indexing

We did not impose a strict requirement that rectangular arrays be stored in
contiguous memory in, say, row major order.  This for two reasons: (i) this
requirement would not have any semantic effect, since one cannot access a 2D or
3D array as if it was one-dimensional (we do not propose the equivalent of
Fortran 90 assumed-size arrays).  The requirement has only performance
implications; e.g., in place reshaping of a 2D array into a 1D array is
expected to be very fast, as no data copying is required. In any case,
contiguity is a significant requirement only within page boundaries: contiguous
pages are not necessarily contiguous in real memory.  Also (ii) a strict
requirement that arrays be stored contiguously would require changes in JVM.
Therefore the weaker requirement that access be as efficient as if the arrays
are stored in canonical order.

As for the storage order, one can follow two approaches.

  o There is a unique storage order, e.g., row major.

  o Arrays can be stored in distinct orders.  For example, storage order could
    be specified when the array is instantiated.  Possible choices would be (i)
    row major (C order), for better performance when native C methods are
    invoked; (ii) column major (Fortran order), for better performance when
    native Fortran methods are invoked; (iii) block major, for block oriented,
    recursive algorithms. A default, row-major layout would be used when users
    do not specify layout.

II. DEVELOPMENT OF CORE CLASSES AND INTERFACES FOR NUMERICAL COMPUTING

The numerics working group has agreed to begin the development of a variety
 of core numerical classes and interfaces to support the development of
substantial Java applications in the sciences and engineering.  The main
purpose of this work is to standardize the interfaces to common mathematical
operations.  A reference implementation will be developed in each case.  The
purpose of the implementation will be to clearly document the class and its
methods.  Although we expect these to be reasonably efficient, we expect that
highly tuned implementations or those relying on native methods will be
developed by others.  Also, the simple methods, such as get or set, will not
provide reasonable performance unless they are inlined, because the method
invocation overhead will be amortized over very few machine instructions.
Unless otherwise specified, we will initially only define classes based on
doubles, since computations with Java floats are less useful in numerical
computing.

The classes identified for first consideration are the following.  We expect to
have the first three fully developed this year, with the others to follow soon
after.

  (a) Complex

      This implements a complex data type for Java as described above.  It
      includes methods for complex arithmetic, assignment, as well as the
      elementary functions.

      Contacts: John Brophy, Visual Numerics
                Marc Snir, IBM

  (b) Multidimensional arrays

      This implements one, two and three-dimensional arrays for Java as
      described above.

      Contacts: Marc Snir, IBM
                Roldan Pozo, NIST

  (c) Linear algebra

      This implements matrices (in the linear algebraic sense) and operations
      on matrices such as the computation of norms, standard decompositions,
      the solution of linear systems, and eigenvalue problems.  A strawman
      proposal has already been developed here and will be released for comment
      soon.

      Contacts: Cleve Moler, The MathWorks
                Roldan Pozo, NIST
                Ron Boisvert, NIST

  (d) Basic Linear Algebra Subroutines (BLAS)

      These implement elementary operations on vectors and matrices of use to
      developers of linear algebra software (rather than to average users).
      This work will be done in conjunction with the BLAS Technical Forum.

      Contacts: Roldan Pozo, NIST
                Keith Seymour, University of Tennessee
                Steve Hague, NAG

  (e) Higher Mathematical Functions

      This includes functions such as the hyperbolics, erf, gamma, Bessel
      functions, etc.

      Contacts: Ron Boisvert, NIST
                John Brophy, Visual Numerics

  (f) Fourier Transforms

      This includes not only a general complex transform, but specialized
      real, sine and cosine transforms.

      Contact: Lennart Johnsson, University of Houston

  (g) Interval Arithmetic

      This implements an interval real data type for Java.  It includes methods
      for interval arithmetic, assignment, as well as elementary functions.

      Contact: Dmitri Chiriaev, Sun

  (h) Multiprecision Arithmetic

      This implements a multiprecision real data type for Java.  It includes
      methods for arithmetic, assignment, as well as elementary functions.

      Contact: Sid Chatterjee, University of North Carolina

The working group will review these proposals and open them up for public
comment.  It will also set standards for testing and documentation for numeric
classes.  It will work with Sun and others to have such classes widely
distributed.

III. OTHER ISSUES

The following problems were discussed by the forum, but no formal
position was taken.

  (1)  Alternative definition of the Java.math library of transcendental
       functions.  The current operational definition is imprecise and
       suboptimal (the functions are defined in terms of bitwise compatibility
       with a particular implementation).  Alternative definitions are (i)
       precise rounding -- result is as if computed in infinite precision
       arithmetic, next rounded; (ii) within fixed bound of precise result; or
       (iii) improved operation definition.  The first definition is very
       desirable if it can be achieved with acceptable performance overhead.
       The second weakens bitwise reproducibility.  Note that current Java
       implementations are not in strict adherence to this aspect of the Java
       standard: most JVMs use their native C math library.

  (2) Improved native interfaces between Java and Fortran.

  (3) Extensions to support multiple NaN values.  This seems to be
      already in the making.

IV. MEMBERS OF THE NUMERICS WORKING GROUP

The following individuals contributed to the development of this document
at the Java Grande Forum meeting on May 9-10 in Palo Alto, California.

  Ronald Boisvert, NIST, Co-chair
  John Brophy, Visual Numerics
  Sid Chatterjee, University of North Carolina
  Dmitri Chiriaev, Sun
  Geoffrey Fox, University of Syracuse
  Steve Hague, NAG
  Siamak Hassanzadeh, Sun
  Lennart Johnsson, University of Houston
  Roldan Pozo, NIST, Co-chair
  Keith Seymour, University of Tennessee
  Nik Shaylor, Sun
  Marc Snir, IBM
  George Thiruvathukal

The following additional individuals also contributed comments which helped
in the development of this document.

  Susan Flynn-Hummel, IBM
  Roger Golliver, Intel
  William Kahan, University of California, Berkeley
  Cleve Moler, The MathWorks
  Jose Moreira, IBM
  Sam Midkiff, IBM