Orthonormalise the cosine sine functions w.r.t a non uniform measure

Python usage
1

I wouls like to know if there is any means in OpenTURNS to build on orthonormal sequence of functions from the usual cosine,sine functions when the scalar product is not w.r.t the uniform measure on [-\pi,\pi] (for instance a truncated gaussian).

2

Hello,

We have the Fourier basis that you can use with Functional Chaos algos for example, maybe thats what you need:
https://openturns.github.io/openturns/latest/user_manual/_generated/openturns.FourierSeriesFactory.html

3

Hi,

You probably ask for an orthonormalization algorithm able to produce an orthonormal basis wrt a given measure but starting from the Fourier basis. There is no such algorithm in OpenTURNS and no way to build a user-defined orthogonal basis, so you are currently stuck. We used to have several such algorithms in the past (ChebychevAlgorithm, ModifiedGramSchmidtAlgorithm) but they were dedicated to the construction of polynomial bases given a measure, and have been superseeded by AdaptiveStieltjesAlgorithm, which was both more robust (by working directly on the 3 terms reccurence coefficients) and more numerically stable.

What I can do is first to implement a FiniteOrthogonalBasis class, allowing you to create an orthonormal basis from a collection of functions. This way you can build your basis outside of OT and use it in the algorithms. Then, I can implement a modified Gram Schmidt algorithm able to build an orthonormal basis given a Basis and a Distribution. But you will have to wait for a release or two…

Cheers

Régis

4

Hi,

I wrote this little piece of code to illustrate the procedure:

import openturns as ot
import math as m

# This basis is orthonormal wrt U(-pi, pi)
fourier = ot.FourierSeriesFactory()
# This is in order to build Function and not UniVariateFunction
basisF = ot.OrthogonalProductFunctionFactory([fourier])
size = 10

# Target measure
measure = ot.TruncatedNormal(0.0, 1.0, -1.5, 2.5)

# Compute int_R f(t)g(t)m(t)dt
# The dot product is computed using a Gauss integration wrt the target measure
Nnodes = 32
nodes, weights = ot.GaussProductExperiment(measure, [Nnodes]).generateWithWeights()

def dotF(f, g):
    return sum([f(x)[0]*g(x)[0]*w for x, w in zip(nodes, weights)])

# Here we verify that the Fourier basis is not orthonormal wrt the target measure
gram = ot.SymmetricMatrix(size)
f = [basisF.build(i) for i in range(size)]
for i in range(size):
    for j in range(i+1):
        gram[i,j] = dotF(f[i], f[j])

# Different from the identity matrix!
print(gram)

# Here we build the orthonormal basis associated to the Fourier basis and the target measure
# We use the modified Gram-Schmidt algorithm, much more stable than
# the crude Gram-Schmidt algorithm, see e.g https://laurenthoeltgen.name/post/gram-schmidt/
newF = [ot.Function(fi) for fi in f]
for j in range(size):
    newF[j] = ot.LinearCombinationFunction([newF[j]], [1.0 / m.sqrt(dotF(newF[j], newF[j]))])
    for i in range(j+1, size):
        newF[i] = ot.LinearCombinationFunction([newF[i], newF[j]], [1.0, -dotF(newF[i], newF[j])])

# Here we verify that the new basis is orthonormal wrt the target measure
for i in range(size):
    for j in range(i+1):
        gram[i,j] = dotF(newF[i], newF[j])

print(gram)

# Draw these functions
a = measure.getRange().getLowerBound()[0]
b = measure.getRange().getUpperBound()[0]

g = ot.Graph(newF[0].draw(a, b, 1024))
for i in range(1, size):
    g.add(newF[i].draw(a, b, 1024))
ot.Show(g)

# And print them... not really readable!
for i in range(size):
    print(newF[i])

You can see that it works pretty well, even if the expressions of the resulting functions are not easy to read!

Cheers

Régis