# Venezuela I:

1967 – 1974

That is how my first Venezuelan period started in August 1967. I got a job at the Instituto de Cálculo and my wife got one at the Mathematics Department, with very generous salaries in strong currency. There was only a small group, in what was the germ of the Computer Science Department. A Licenciatura was started officially in 1968 using as a basis the recently published ACM Curriculum (I believe this was one of the first undergraduate programs in Computer Sciences to be created in the world; earlier, Manuel Sadosky had created a 3 years program in Scientific Computing in Buenos Aires).

By this time I was feeling some pressure to diversify and let deferred correction rest for a while. I thus retook the thread of my early work on least squares and pseudoinverses and started to study in earnest the problem of regularization of ill-posed least squares problems. This resulted in some publications:

“Stabilizing Linear Least Squares Problems”, Proc. IFIP 68, Supplement, p. 127 (1968)

Also in the area of nonlinear least squares I wrote another early convergence proof for the Gauss-Newton method:

After this detour I retook my program for applying deferred corrections to as many different kinds of problems as I could and started to look seriously into partial differential equations. This lead into two directions. In

I discussed the application of deferred corrections to nonlinear elliptic problems. In order to produce a computational algorithm I had to address the problem of generating automatically high order finite difference approximations to complicated multi-dimensional differential operators. This lead to a second profitable avenue of research that produced a number of papers on fast Vandermonde solvers and several generalizations of Vandermonde systems as they apply to numerical differentiation, interpolation and quadratures in one and several variables. This resulted in a number of publications with various collaborators:

These papers presented an O(n**2) method for the solution of linear systems of equations with Vandermonde and confluent Vandermonde matrices, which were known to be quite ill-conditioned. The method was not only faster, but also seemed to be impervious to ill-conditioning for larger values of the dimension than Gaussian elimination. Curiosly enough, the paper with Bjorck was entirely written via snail mail. We only met in person years later: another of Gene Golub’s magic touches.

Then I worked on several extensions to multi-dimensional problems:

The material in this paper could be useful to produce high order discontinuous elements in any dimension, but I have never pursue it.

These techniques were essential for the effective implementation of iterated deferred corrections that we were developing. Some aspects of this subject are now quite popular in Electrical Engineering applications and a number of generalizations have been developed. Still, the method we devised with Ake Borck seems to be unbeatable in terms of stability, accuracy and efficiency.

Nick Higham has provided recently a solid basis for what we had demonstrated in practice through the years.

Carl de Boor has looked into some of the multidimensional generalizations in his study of divided differences. He has also embraced our fast tensor product algorithms as the basis for fast tensor product spline manipulations and many other such applications, having written a paper with the same title, in which he explains more clearly our algorithm, without resorting to multilinear algebra notation and manipulations. Eric Grosse has extended those algorithms to the tensor least squares problem and I use them all in my current research.

Years later, Laura Perez (Universidad de Rio Cuarto, Cordoba, Argentina), during a one month visit to California, implemented the tensor product solution of systems of equations and least squares problems that I have used extensively to this day. Laura is now writing her thesis under my direction on control theory applied to hybrid vehicles using PASVA4 (see below).

A number of important developments occurred now concurrently. First I connected to the work of Hugo Scolnik on separable nonlinear least squares while visiting the Centre de Recherche Mathematique of the Universite de Montreal. In a contribution with Irwing Guttman and Hugo Scolnik, we generalized and clarified that work:

Then, as a consequence of a visit to the Fundación Bariloche in the Lake region of Southern Argentina, we extended further these ideas to what is now known as the Variable Projection method for separable nonlinear least squares problems:

This work was then extended and analyzed by a number of US and Swedish scientists and constitutes now an independent chapter in the area. It was shown that the elimination of linear variables not only reduced the dimensionality of the problem, but also caused the minima to be better defined. Thus, the same minimization algorithms would converge faster for the reduced problem than when applied to the unreduced one.

It seems that this work was inspiration for algorithms for the detection of signals in antenna arrays and it also applies to the problem of best knot selection for the least squares fit of data by B-splines. A computer implementation that I wrote is embodied in the code VARPRO. VARPRO and its successors, written by J. Bolstad and R. LeVeque, and a Bell Labs. version by David Gay and Linda Kaufman, are still in wide use. Also, B. Rust was an early user at the National Bureau of Standards, where they wrapped VARPRO with an easy to use interface and detailed documentation.

In collaboration with Linda Kaufman we extended this method to problems with constraints:

Some applications to physics problems (Mossbauer spectroscopy) were carried out with my student Consuelo Maulino. Thirty years down the line we wrote a review paper and a couple of books related to this important topic, but more about that at the proper time.

A couple of unconnected contributions also occurred at this time. In my only collaboration with Seymour Parter, and as a consequence of a visit to Madison, we solved a problem in numerical bifurcation for ordinary boundary value problems of special type:

With my student Godela Scherer and I believe as a consequence of a very well done homework assignment we wrote:

that introduced the idea of combining multiple approaches (poly-algorithm) to solve a problem efficiently. In this case we produced a software implementation that competed well with the state of the art code available in EISPACK. Unfortunately we did not pursue it to the point of having a robust, publicly available code, so it eventually got lost and did not have the impact that it deserved. Anecdotically, this paper was refereed by Jim Wilkinson, who not only made some keen comments that helped improve the presentation, but also tested the code.

The work on deferred corrections for ordinary differential equations with boundary conditions took a large impulse when the seminal book on numerical solution of such problems by Herbert Keller came to light. This contained a general theory for the numerical solution of systems of two-point boundary value problems and allowed us to extend our methods and theory to that realm. The work of extension and implementation as computer software was done in collaboration with my student Marianela Lentini. The first contributions worked on uniform meshes:

This work had a disproportionate impact in the community, which until now had not seemed to pay much attention to the work on deferred corrections, despite the fact that the numerical results were quite spectacular. But now we had hit a nerve by producing a piece of fairly professional code that could be used by others and that solved a large class of problems of interest in many disciplines, and last but not least, which had no competitors.

Very rapidly we found some of the limitations of our approach, specially for singular perturbation problems. Thus adaptive meshes came to light quite naturally. I looked around in the literature for some pointers and found very little but ad hoc techniques that had been used with mild success in some special problems, noticeably, the work of a prominent Argentine Meteorologist, Eugenia Kalnay. Keller’s theory covered nonuniform meshes, so we were in solid ground there, but there was no hint on how to choose such meshes in a useful way.

Then I hit on some work from Carl de Boor on best meshes for spline interpolation that had originated with H. G. Burchard. Granville Sewell had written his Doctoral Dissertation at Purdue under R. E. Lynch, using some of those ideas in an early version of adaptive meshes for elliptic problems. This, combined with the serendipitous fact that Granville had moved to Caracas to work at Simon Bolivar University, led to a collaboration that produced one of the seminal contributions to the adaptive selection of meshes for boundary value problems:

It is in this paper where the term and concept of equidistributing meshes, that pervades the modern applications in ordinary and partial differential equations, was coined. From this work it was an easy step to produce a code with adaptive mesh capabilities, which, when combined with our variable order deferred correction techniques constituted a robust, versatile and efficient tool, in the same league with established initial value problem solvers (in modern terminology this was an h and p-adaptive method).

Versions of this code are still in routine use throughout the world as part of the IMSL, Harwell, NAG, na-net and other libraries of scientific subroutines. We have knowledge of successful applications in disciplines as different as mechanical engineering and biology.

These ideas were also the precursors of the h-p finite element techniques that were later introduced by Babuska and his school for ordinary and partial differential equations as well. It helped also to develop a competitor code, COLSYS, authored by J. Christiansen, R.D. Russell and Uri Ascher in Canada. We interacted a lot with Bob Russell during mutual visits to Vancouver and Caracas and shared with him and his students all our successful tricks, specially those related to mesh refinement, many of which were later incorporated in the collocation solver COLSYS. We even attempted to make some fair comparisons between the resulting codes that lead somewhat later to: