Can Computer Algebra be Liberated from its Algebraic Yoke ?
An Algorithmic Approach to Functional Analysis from a Functional Approach to Computing
Rémi Barrère
Université de Franche-Comté, ENSMM
26, chemin de l’Epitaphe, F-25000 Besançon, France
Fax : 33(0)381 809 870
E-mail : rbarrere@ens2m.fr

Abstract

computations or to interpret purely numerical results: the
lack of symbolic representation prevents the user from
identifying patterns among models, hence from finding
relevant generic models.
The purpose of the paper is to show that symbolic
manipulations are still worthwhile when associated with
functional methods. Here, "functional" refers to the
cooperation of function theoretic methods and the functional
programming paradigm. Indeed, most function space
methods are "constructive" in the sense that they basically
generate symbolic approximations. Moreover, the resulting
algorithms turn out to be well adapted to the functional
paradigm. So functional programming can be regarded as
the natural computational counterpart of functional analysis.

So far, the scope of computer algebra has been needlessly
restricted to exact algebraic methods. Its possible extension
to approximate analytical methods is discussed. The
entangled roles of functional analysis and symbolic
programming, especially the functional and transformational paradigms, are put forward. In the future, algebraic
algorithms could constitute the core of extended symbolic
manipulation systems including primitives for symbolic
approximations.

1 Introduction
In recent decades, the extensive development of numerical
techniques entailed a lower use of traditional approximate
analytical methods, such as series expansions or functional
iteration. Admittedly, computer algebra is currently reviving
the analytical approach, and active researches in that area
bring about new ideas and let expect future progresses (see
for instance the Journal of Symbolic Computation or the
proceedings of the ISSAC conferences). The underlying
algorithms yield symbolic expressions that are generally
more compact, easier to manipulate and more meaningful
than numerical data.
However, computer algebra techniques are more or less
restricted to exact solutions, so they often fail in the case of
physical or engineering problems. Indeed, mathematical
modelling commonly leads to equations that have no known
closed-form solution. In that case the general trend consists
in resorting to numerical techniques. Unfortunately, these
have drawbacks such as the difficulty to link different

2 Functional and rule-based programming
The algorithms described below were drafted in
Mathematica, which consistently combines a computer
algebra system with a mutiparadigm programming language
built on top of a rule-based core. Its sophisticated patternmatching capabilities facilitate the writing of concise
symbolic programs, the syntax of which mimics that of
mathematics [31] (see also the critical reviews by
Fateman [14, 15] and Maeder’s presentations [19, 21]).
Actually, Mathematica results from a symbiosis of the
functional and rule-based (transformational) paradigms: the
user defines functions which are internally processed as
transformation rules by a rewriting mechanism. As a
functional language, Mathematica can be thought of as a
descendant of Lisp; as a transformational one, it can also be
regarded as a descendant of Snobol with a noticeable
influence of Prolog.
Born in the late fifties with Lisp in a context of artificial
intelligence, the functional paradigm really took off in the
eighties thanks to novel design techniques such as graph
rewriting [26]. With a functional language, a computation is
processed as a function evaluation. The theory of recursive
functions and the λ-calculus constitute the mathematical
foundation of the functional paradigm.

1

insight than intricate exact solutions or bulky numerical
results. Then, approximating can be a way of simplifying.
A few fundamental methods from functional analysis
(e.g., fixed points, perturbation expansions or variational
formulations) constitute the theoretical background for this
processes. The present investigation focusses on some of
their algorithmic counterparts, namely successive
approximations, Newton's method, perturbation series and
the Galerkin procedure.
A comment might be worth stressing here in order to
dismiss a common misinterpretation: function theoretic
methods are inherently symbolic (precisely because they
apply to function spaces), although they are usually thought
of as numerical. In particular, they require algebraic
facilities, hence rely on the algebraic core of computer
algebra.
These methods have been used long, even before the
emergence of computing, but they are tackled here as
generic tools for the programming of symbolic approximate
methods. Papers referring to this approach are still scattered
about the scientific literature, which testifies its structuring
power in scientific computing has not been fully recognized
yet.

Born in the mid sixties with Snobol, the
transformational paradigm remained more or less hidden
behind the scenes until it was acknowledged as a genuine
paradigm on the occasion of researches in programming
paradigms [27]. Rule-based languages emphasize the
transformation of symbolic expressions according to their
form, so a computation amounts to the restructuring of an
expression by term rewriting. The Post-Markov theory of
algorithms constitute the mathematical counterpart of this
paradigm.
These languages enable an abstract programming style
thanks to high level constructs such as patterns, rewrite rules
or higher-order functions that avoid classical control
structures. They proved to be a sound basis for computer
algebra systems (e.g., Reduce, Macsyma), so the primitives
of the language can be combined with the mathematical
operators, thus leading to a concise programming style
where computing and programming tend to blend into a
single activity [5].
By combining the functional and transformational
paradigms, Mathematica facilitates the recourse to
denotational semantics either as a tool for program
modelling or more generally for translating relational
algebraic formulations into transformational algorithmic
ones (Knuth-Bendix algorithm or its variants).
These considerations are far from superfluous. Computer
scientists have become aware of the role of "programming
languages as thought models" [29], which is of particular
interest in scientific computing. They now take into account
that "the language in which a programmer thinks a problem
will be solved will color and alter, in a basic and
fundamental way, the fashion in which that programmer will
develop an algorithm" [11]. In particular, thanks to higherorder functions, Mathematica facilitate the implementation
of functional methods, which lend themselves well to
symbolic programming.

4 Fixed points and successive approximations
Principle and programming of the method
Although it is often restricted to existence and uniqueness
demonstrations [22], the method of successive approximations also turns out to be valuable for computing
approximate (occasionally exact) symbolic solutions. The
method is valid in (complete) metric spaces, so it applies to
numerical as well as functional or even geometrical
problems. Its statement can be briefly recalled; as long as f
is a contraction mapping (f satisfies a Lipschitz condition
with k<1) for some distance:
(1)
d[f(x),f(y)] ≤ k d(x,y) ,

3 Role of symbolic approximations

we can solve the equation f(x)=x, i.e., we can find the
(unique) fixed point of f, by computing the limit:

Symbolic approximations combine the advantages of
symbolic computations and approximate methods. As
opposed to numerical techniques, they enable the
manipulation of literal expressions with symbolic
parameters, thus entailing a better understanding, hence a
better mastery of the underlying physical systems.
Actually, physical or engineering applications do not
necessarily require exact solutions nor accurate
approximations. Indeed, most (if not all) models are based
on assumptions, simplifications or approximations, so it
does not necessarily make sense to seek exact solutions nor
very precise approximations.
Moreover, exact solutions occasionally have such a
complexity that approximations may be more practical.
Sometimes, even rough representations may give more

on

x ∞ = lim f ( x )
n→∞

(2)

where fon denotes the nth iterate. Moreover, provided some
initial approximation xo is in the contraction domain of f,
any iterate
xn= fon(xo)
can be considered as an approximation to the fixed
point x∞. Such approximants are sometimes called
continued-functions or iterated functions. Above all,
functional languages have an operator ("Nest" in
Mathematica) that computes the nth iterate, so the program
comes down to using a primitive of the language. When f is
an operator in a function space, its algorithmic counterpart
is a higher-order function.

2

5 The Newton method

The example below shows the first approximation to the
solution to the equation u’+(1+ϕ2)u/2=ϕ, with u(0)=0 and
ϕ=a sin(wt), drawn from a model of an electric motor. The
second approximation enabled the identification of a
periodic steady-state behaviour, which was validated by
comparison with a numerical approximation.

Principle of the method
Although it is commonly thought ofas a numerical
technique, Newton’s method applies to Banach spaces, so it
basically solves functional equations in a symbolic way.
This was mentioned by Kantorovich [17] and to some extent
used by Bellman [7] who called it quasilinearization, for it
relies on solving a sequence of linear equations. Since then,
the nascent idea seems to have been overwhelmed by the
numerical tidal wave. It should not be likened to the
particular use of Newton’s method for series
computations [18].
Besides, optimization problems are connected with it
since finding an extremum to g amounts to seeking a zero to
f=grad g (hence the idea of approximate symbolic
optimization). The following example is worth mentioning
for it is often used as a test for numerical algorithms: the
(exact) extremum (1,1) of the Rosenbrock function
f(x,y) = (x-1)2+p (x2-y)2 turns out to be reached in two
steps, from the initial approximation (x,y)=(0,0).

x[t_]:= a Sin[w t]
operator[u_][t_]:= Integrate[
Exp[-(t-τ)/2]*(x[τ]-u[τ]*x[τ]^2/2),
{τ,0,t}
]
Nest[operator[#][t]&, 0 &[t], 1]
4aw
2a (2 w Cos[tw] - Sin[tw])^2
----------------------- – ------------------------------------------------------t/2
2
1 + 4w 2
E
( 1 + 4w )

Fixed point equations can be written in different ways
(e.g., f(x)=x <=> x=f-1(x) ), which is of practical importance
for finding a form where f is a contraction with the smallest
possible value of k. For instance, it can be necessary to turn
a differential equation into an integral equation.

Rosenbrock[p_,x_,y_]:=(x-1)^2+p(x^2-y)^2
equations={D[Rosenbrock[p,x,y],x]==0,

Convergence considerations

D[Rosenbrock[p,x,y],y]==0}

The major restriction regards the Lipschitz condition, since
the convergence is limited to the contraction domain of f.
When it is assured, the method leads to a linear
convergence, i.e., d(xn+1,x∞) ≤ k d(xn,x∞). Nevertheless, the
non-linear Aitken-Shanks transformation improves the rate
of convergence [9].

NewtonASolve[equations, {{x,y}, {0,0}, 2},
Simplify->True]
{{y -> 1, x -> 1}}

Here is the schematic statement of Newton’s method: as
long as f has an invertible Fréchet derivative grad f, and the
initial approximation x0 is sufficiently near the solution
(which can be more precisely stated for specific problems in
connection with convexity properties), a sequence of
approximate solutions to the equation f(x)=0 can be
computed by means of the recurrence formula:

x n + 2 x n – x n2 + 1
a n = a ( x n ) = ---------------------------------------------x n + 2 + x n – 2x n + 1
1
= x n + 2 + -------------------------------------------------------------1
1
------------------------------- – -----------------------xn + 2 – xn + 1 xn + 1 – xn

(3)

x n + 1 = x n – grad f(x n)

It (the second form) is generally used in numerical
contexts, although there is no impediment for its use in a
symbolic context. In the case of successive approximations,
its efficiency results from the fact that it extrapolates the
(asymptotic) geometric behaviour of the sequence. In that
case, the algorithm (called Steffesen’s method) leads to the
following expression, the programming of which is
straightforward in a functional style.

o(n + 2)

on

o(n + 1)

[ f ( xn ) ]

grad f(x n) x n + 1 – x n + f ( x n ) = 0
2

(5)

In other words, x∞ is the fixed point of x-[grad f (x)]-1[f(x)].
Since it is a functional iteration, the procedure can also be
implemented in a functional language by means of an
iteration operator, but it requires further programming for
computing the Fréchet derivative. In the case of differential
or integral equations, the recurrence formula is rather
written in the following form, which actually designates a
sequence of linear equations to be solved:

o2

f ( xn ) xn – f ( xn ) 2
s n = s ( x n ) = --------------------------------------------------o2
f ( x n ) + x n – 2f ( x n )

–1

(6)

The Prager model for the deflection of a beam [7] gives
an example of a non-linear Sturm problem :
u’’+(a u’)2+1=0,
u(0)=u(1)=0.
Below,
the
first
approximation is computed (with u0(x)=0 as initial
approximation) and is shown together with the exact
solution (in solid line).

(4)

f
( x0 ) f ( x0 ) – f
( x0 )
= ---------------------------------------------------------------------------------------o(n + 2)
on
o(n + 1)
f
( x0 ) + f ( x0 ) – f
( x0 )

3

The example concerns the common case of a quadratic
functional, which leads to the usual second degree form.

NewtonDSolve[{u''[x] + (a u'[x])^2 + 1 == 0,
u[0]==0, u[1]==0},
u[x], x, {0, 1}]

x
2

Discussion and future directions

x2
2

{u[x] -> --- - ----- }
2
The quadratic convergence, namely ||xn+1-x∞|| ≤ k ||xn-x∞|| ,
is an advantage of Newton's method; this generally
compensates its greater computational complexity.
However, Fateman [13] produced examples for which the
quadratically convergent algorithm is only marginally faster
and sometimes slower than a linear one. Nevertheless, in
most circumstances, quasilinearization is known to provide
with a wider interval of convergence than successive
approximations. A drawback may be the difficulty for
computing the inverse of the Fréchet derivative. This
problem especially arises with (non-linear) differential
equations, for the Fréchet derivative generally leads to
(linear) differential equations with variable coefficients.
Nevertheless in that case, approximations to grad f may
lead to satisfactory results [7], which means for instance, the
possibility of combining the Newton method with the
Galerkin one (section 7) for solving the associated linear
problems.

Figure 1: Prager’s model for the deflection of a beam.

Programming the Fréchet derivative
According to the circumstances, the Fréchet derivative (or
tangent linear map) can be expressed by means of the
ordinary derivative, the jacobian matrix or the gradient
operator. In this last case, it would be unreasonable to aim at
an algorithm computing any Fréchet derivative. But, a quite
general program can be designed for processing the
functional expressions commonly encountered when
solving differential or integral equations. It is based on a few
usual derivatives and the following formulas [22], where
FD[f] denotes the Fréchet (or functional) derivative of f with
respect to ϕ :
FD [ f (ϕ(x), ϕ'(x), …, ϕ

(n)

(x), x) ] [ h(x) ]

n+1

∑ ∂i f ( ϕ(x), ϕ'(x), …, ϕ

=

(n)

(x), x )h

(i – 1)

(x)

6 Perturbation expansions

i=1

FD
=

∫Ω f (ϕ(x), ϕ'(x), …, ϕ

(n)

(x), x)dx [ h(x) ]

∫Ω FD [ f (ϕ(x), ϕ'(x), …, ϕ

(n)

Principle of the method and functional formulation

(7)

Only regular perturbations are discussed here. Perturbation
series arise from equations of the form:
f(x, ε) =0
(8)
where f is an operator between Banach spaces, which means
the algorithm is valid whatever the precise nature of the
equation. It is assumed f(x,0)=0 can be solved and the
parameter ε is supposed to be small, so f(x, ε) =0 can be
viewed as a perturbed model of f(x,0)=0; likewise f(x,0)=0
can be viewed as an approximate model of f(x, ε)=0.
The perturbation procedure is related to the implicit
function theorem [16]: if f(x, ε) is a Ck mapping between
Banach spaces, and there is an x0 satisfying f(x0,0)=0 such
that the partial (Fréchet) derivative ∂xf(x0,0) is an invertible
linear map, then there is a unique solution to the implicit
equation f(x, ε)=0 in the neighborhood of ε=0, given by
x=ϕ(ε), where ϕ is a Ck mapping.
The procedure for regular perturbations [6] consists in
seeking a formal solution in the form of a power series in
epsilon:

(x), x) ] [ h(x) ] dx

Pattern-matching and transformation rules facilitate the
programming of such formulas. Here is a short extract of the
program.
FDerivative[
Integrate[expr_,{x_,a_,b_}],u_[x_]
][h_]:=
Integrate[
FDerivative[expr,u[x]][h],
{x,a,b}
]/;FreeQ[a,x]&&FreeQ[b,x]

A subsidiary application of the Fréchet derivative is
worth mentioning in this context. It regards the computation
of variational formulations associated with potential
operators in the case of conservative systems [25].

x=g(ε)=∑ εk xk
(9)
where the xk are the unknowns. The method of
indeterminate coefficients gives rise to an infinite sequence
of equations that can be solved recursively:
fk(xi,i≤k)=0
(10)
In practice, a truncated series is computed by solving a finite
system. Obtaining a functional iterative formulation is a bit

FDerivative[
Integrate[
1/2 u’[x]^2 - f[x] u[x],
{x,a,b}
], u[x]][v[x]]
Integrate[u’[x] v’[x] - f[x] v[x],{x,a,b}]

4

7 Weak solutions and the Galerkin procedure

tricky, but this can be achieved by way of the accumulation
operator "Fold" in Mathematica.
As an illustrative example, here is a perturbative
approximation to the Bernoulli differential equation
u’+k u+e u2=0, u(0)=1.

Principle of the method
The Galerkin method regards boundary-value problems,
which typically arise from continuum physics models; here,
we will limit ourselves to linear equations. Most of these
models can be placed into one of the three following classes:
elliptic (or Laplace) equations for static phenomena, spectral
(or Helmholtz) problems for stationary phenomena, time
evolution (parabolic or hyperbolic) equations in the case of
dynamic behaviour (diffusion or wave propagation). This
classifying turns out to be useful from the algorithmic point
of view. The method is closely related to weak (or
variational) formulations, which provide with both
theoretical results and approximation techniques [8, 30].
Schematically, starting with the initial differential
formulation with its boundary conditions :

theSeries=PerturbationADSolve[{
D[u[t],t]+k u[t]+e u[t]^2==0,
u[0]==1
},u[t],t,{e,1}]
e
{u[t] -> e-kt + --- (e-2kt - e-kt)}
k

Approximate resummation
Perturbation expansions tend to generate slowly convergent
series or even divergent ones. It has been a common practice
to improve the convergence by means of Padé approximants,
which in fact constitute a resummation procedure [9]. The
idea of Padé summation is to replace a truncated power
series by a rational function whose first terms in the series
expansion match the given series. Its advantage comes from
the fact that its computation involves only algebraic
manipulations. Although Padé approximants have been used
essentially from a numerical point of view, they can be
implemented in a symbolic way. Moreover, they are known
to work well, even beyond their proven range of
applicability [1].
In the previous example, the Padé approximant P0/Q1
restores the exact solution from the second approximation.

Au=f
(11)
(respectively : Au=λu, or Au=Dtu), the approximation is
obtained in two stages. First, a weak formulation leads to an
integral expression based on a bilinear operator, i.e., a
functional inner product :
∀v, <Au,v>=<f,v>, where : <u,v>= ∫Ω u v dΩ

(12)
(respectively : <Au,v>=λ<u,v>, or <Au,v>=Dt<u,v>). Then
an approximate solution un is searched in the form of a
linear combination of n (previously chosen) basis functions :
un=∑ ci wi
(13)
In other words, un is the projection of the exact solution onto
the finite dimensional subspace spanned by the system of
functions wi. The test functions v and the basis functions wi
are supposed to belong to a convenient Sobolev space. This
leads to a finite set of algebraic (respectively proper value,
or differential) equations.

Needs["Calculus‘Pade‘"]
Pade[u[t]/.theSeries,{e,0,0,1}]//Simplify
k
--------------------------------------kt
– e + ( e + k )E

∑ ci 〈 Awi, w j〉

Padé approximants are at the heart of a more general
resummation method suggested by Bergeron and
Plouffe [10] in the context of combinatorics. An extension
of their heuristic procedure to perturbation series was
suggested in [3].

= 〈 f , w j〉 , j ∈ { 1, n }

(14)

i

respectively :

∑ ci 〈 Awi, w j〉
i

∑ ci ( t ) 〈 Awi, w j〉
i

= λ ∑ c i 〈 w i, w j〉
i

= D t ∑ c i ( t ) 〈 w i, w j〉
i

Piecewise linear or piecewise polynomial approximations can be obtained by the finite element technique; this
is a common numerical approach. Nevertheless, other kinds
of basis functions can be chosen, according to the shape of
the domain, orthogonality considerations or possibly
experimental data. On the one hand, purely numerical
solutions are obtained, by solving large systems; no
symbolic parameter is available, but intricate boundaries are
possible. On the other hand, approximate symbolic solutions
are reachable, by solving small systems; analytical
expressions with symbolic parameters can be manipulated,
but intricate boundaries are difficult to handle.

Future work
The perturbation method often leads to non uniform
(singular) perturbations, and non-linear equations may
generate bifurcation phenomena. Appropriate computational
techniques have been settled to overcome those
difficulties [24, 28]. Their implementation as generic
operators could be investigated. More generally, a full
implementation of the perturbation method would facilitate
the symbolic treatment of homogeneization or boundary
layer problems.

5

Functional formulation

8 Discussion

The procedure is entirely determined by inner products and
a truncated system of functions, so it lends itself well to the
functional style. In fact, according to whether the problem is
elliptic, a spectral problem or an evolution equation, three
generic operators are required, the arguments of which are
the inner products and the basis [2]. Here is the first one (the
others being similar), where "Outer" (generalized tensor
product) generates the left-hand matrix, whereas "Map"
generates the right-hand vector. Less generic programs are

Symbolic approximations
The aforementioned operators for symbolic approximations
can be used like built-in ones. For instance, they can be
nested or combined with each other, so that intermediate
results can be "piped" from (or to) each other. Actually, they
extend the capabilities of the software environment, where
most computations tend to be entirely processed.
In practice, this kind of computation requires a
permanent control from the user, which is a characteristic of
most computer algebra manipulations. In particular, the
importance of a judicious choice of some parameters (e.g.,
initial approximations or basis functions) must be stressed,
so an intuitive knowledge of the solutions can be of great
service.
These methods are (or can be made) "adaptive" in the
sense that the approximations can be refined step by step; so
the user can choose any intermediate between rough cheap
approximations, and more precise but more expensive ones.
If need be, these can be implemented by means of a
stream [5, section 3.2.5].
As opposed to numerical computations, symbolic
manipulations pave the way for the discovery of relevant
patterns among the models and their solutions, which enable
in turn the emergence of novel mathematical abstractions.
However, these methods still rely on computer algebra
capabilities, the current limitations of which sometimes
hinder the whole computation. In particular, like algebraic
manipulations, analytical ones tend to generate intricate
expressions; in such cases, only low-order approximations
are practically reachable.

given in [8] or [23].
GalerkinSolve[a_,l_,basis_List]:=
LinearSolve[
Outer[a,basis, basis],
Map[l, basis]
].basis

However, this simplicity has a drawback: constructing
the inner products, which characterize each particular
problem, is up to the user and may be a more or less difficult
task. The worthwhile counterpart is the possibility of
choosing a purely symbolic, a purely numerical or a mixed
symbolic-numerical solution. Transferring this choice on the
definitions of the inner products makes the procedure
"GalerkinSolve" fully generic.
In the exemple, the finite element technique
(approximation with 3 hat functions without interpolation)
was qualitatively compared with the symbolic Galerkin
method (3 trigonometric functions) for the 1-D Laplace
equation: -u’’=exsin 5x with homogeneous boudary
conditions.
n  3

n  3

Symbolic programming
Functional languages enable the manipulation of higherorder functions which are particularly useful for
implementing function theoretic algorithms. Also, patternmatching facilitates the design of generic operators that
apply to classes of arguments (equations in this context).
This actually means that more or less general mathematical
methods can be translated into generic operators. Broadly
speaking, symbolic languages facilitate the expression of
scientific knowledge and the programming of mathematical
abstractions.
For historical reasons, programs written in symbolic
languages have obtained a reputation for lack of
performance. Part of this results from the high-level of
abstraction that is available in such languages and from the
background operations that they require. However, in the
eighties, research in functional languages resulted in new
design techniques leading to better performances [26].
Besides, although procedural programming is known to be
more efficient for numerical computations, functional

Figure 2: finite elements vs. symbolic Galerkin.

Possible developments
The Galerkin procedure constitutes the springboard for a
mixed symbolic-numerical treatment of continuum physics
equations. As an illustration of this expanding approach,
Dasgupta [12] resorted to the symbolic capabilities of
Mathematica to generate shape functions for concave
quadrilateral elements by means of Padé approximants.
In the case of nonlinear problems, the finite element
technique is quite usually associated with a numerical
version of Newton’s method. Similarly, the Galerkin
procedure could be associated with the symbolic version of
Newton’s method or with the perturbative approach.

6

References

programming often turns out to be more efficient for
symbolic ones [20].
Also, the spreading of software environments that
combine a computer algebra system with a programming
language arouses new activities in mathematical modelling
and scientific computing where the development time
reaches the same order of magnitude as the time of use; so it
can no longer be neglected, and benchmarks should now
take this into account. Furthermore, under the influence of
computer algebra and symbolic programming, scientists
become aware of the role of programming languages as a
medium of scientific knowledge, so in some circumstances,
the form of a program may take priority over its efficiency.
Finally, the referential transparency of (purely)
functional programs facilitates their mathematical
description. Actually, the programs and their descriptions
tend to become two complementary aspects of a single
formulation: a transformational version and a relational one.
Moreover, the emergence of very high-level languages that
implement the main three theories of computation, thus
enabling the direct expression of mathematical knowledge,
tends to extenuate the traditional distinction between
algorithms and programs [4].

[1] BAKER G.A. Essentials of Pade approximants,
Academic Press, 1975.
[2] BARRÈRE R., CARMASOL A. Approximate solutions of
boundary-value problems using computer algebra
software, in Keränen V., Mitic P. (Eds.) Mathematics
with vision, Computational Mechanics Publications,
1995, 45-52
[3] BARRÈRE R. A Heuristic Method for the Resummation
of Perturbation Series, in Keränen V., Mitic P.,
Hietamäki A. (Eds.) Innovation in Mathematics,
Computational Mechanics Publications, 1997, 61-68
[4] BARRÈRE R. The Structuring Power of Mathematica in
Mathematics and Mathematical Education, IMS’99
(International Mathematica Symposium), Hagenberg,
August 1999 (http://south.rotol.ramk.fi/~keranen/
IMS99/ims99papers)
[5] BARRÈRE R. Mathematica, Vuibert, 2002.
[6] BELLMAN R. Perturbation Techniques in Mathematics,
Physics and Engineering, Holt, Rinehart, Winston, 1964.
[7] BELLMAN R., KALABA R. Quasilinearization &
Nonlinear Boundary-Value Problems, American
Elsevier, 1965.
[8] BELTZER A. Variational and Finite Element Methods,
Springer-Verlag, 1990.
[9] BENDER C.M., ORSZAG S.A. Advanced Mathematical
Methods for Scientists and Enginers, McGraw-Hill,
1978.
[10] BERGERON F., PLOUFFE S. Computing the Generating
Function of a Series Given Its First Few Terms,
Experimental Mathematics, 1(4), 1992, 307-312
[11] BUDD T. Multipadigm Programming in Leda,
Addison-Wesley, 1995.
[12] DASGUPTA G., TREIL J. Maxillo-Facial Frame: Finite
Element Shapes, The Mathematica Journal 8(2), 2001,
235-246.
[13] FATEMAN R.J. Series Solutions of Algebraic and
Differential Equations: A Comparison of Linear and
Quadratic Algebraic Convergence, in ACM-SIGSAM
1989 Proceedings of the International Symposium on
Symbolic and Algebraic Computation, ACM Press, 1989
[14] FATEMAN R.J. A Review of Mathematica, J. Symbolic
Computation 13, 1992, 545-579.
[15] FATEMAN R.J. Symbolic Mathematics System
Evaluators, in Lakshman Y.N., Ed, Proceedings of the
1996 International Symposium on Symbolic and
Algebraic Computation, ACM Press, 1996, 86-94.
[16] KEENER J.P. Principles of applied mathematics,
Addison-Wesley, 1988.
[17] KANTOROVICH L.V., KRYLOV V.I. Approximate
Methods of Higher Analysis, Interscience, 1958.
[18] LIPSON J.D. Newton’s method: A great algebraic
algorithm, in Jenks R.D., Ed, Proceedings of the 1976
ACM Symposium on Symbolic and Algebraic
Computation, 1976

9 Conclusion
When equations have no known closed-form solution,
symbolic approximations may constitute an alternative to
numerical techniques, thanks to the symbiosis of the
expressive power of symbolic languages and the
computational power of function space methods.
On the one hand, these turn out to be a uniform
conceptual approach to the study of functional equations:
they provide not only potent theoretical tools (existence,
uniqueness and convergence theorems), but also efficient
symbolic approximation procedures. In that sense, these are
constructive methods.
On the other hand, the resulting algorithms lend
themselves well to the functional and transformational
programming paradigms. They are quite easily implemented
as generic operators in a computer algebra and symbolic
programming compound environment.
To sum up, not only functional analysis is the natural
theoretical background for symbolic approximation theory,
but also functional programming turns out to be its natural
computational counterpart. Together, they pave the way for
an algorithmic approach to functional analysis supported by
a functional approach to computing. In the future,the
integration of symbolic approximation algorithms into
computer algebra systems might contribute to extending the
restrictive connotation of computer algebra to the wider
notion of analytical manipulation.
A set of experimental Mathematica packages and notebooks,
including the code and a documentation is freely available
(http://macmaths.ens2m.fr/research/pages/symbolic.html).

7

[19] MAEDER R. Programming in Mathematica, AddisonWesley, 1991
[20] MAEDER R. The Mathematica Programmer II,
Academic Press, 1996
[21] MAEDER R. Computer Science with Mathematica,
Cambridge University Press, 2000
[22] MILNE R.D. Applied Functional Analysis, Pitman,
1980.
[23] MURA T., KOYA T. Variational Methods in Mechanics,
Oxford, 1992.
[24] NAYFEH A.H. Introduction to Perturbation
Techniques, Wiley, 1981.
[25] ODEN J., REDDY J. Variational Methods in Theoretical
Mechanics, Springer-Verlag, 1976.
[26] PEYTON JONES S.L. The implementation of functional
programming languages, Prentice Hall, 1987.
[27] PLACER J. Multiparadigm Research : A New Direction
In Language Design, ACM SIGPLAN Notices, 26(3),
1991, 9-17
[28] RAND R., ARMBRUSTER D. Perturbation Methods,
Bifurcation Theory and Computer Algebra, SpringerVerlag, 1987.
[29] RECHENBERG P. : Programming Languages as Thought
Models, Structured Programming, 11(3), 1990, 105-115
[30] REDDY J.N. Applied Functional Analysis and
Variational Methods in Engineering, McGraw-Hill,
1986.
[31] WOLFRAM S. The Mathematica Book, Wolfram Media
& Cambridge University Press, 1996.

8

