Mathematica code for Gram Schmidt Orthogonalization for a general Inner Product (in particular, weighted Bergman spaces)

Question

I am trying to write a Mathematica code to find an orthonormal basis of the weighted Bergman space defined by: $$A^2(\varOmega) = \Big\{f \in \operatorname{Hol}(\varOmega) \; \Big| \; \|f\|= \textstyle\sqrt{\int_{\varOmega}|f(z)|^2 \, (1-|z|^2) \, \mathrm{d}A(z)} < \infty\Big\}$$

with the normalized Lebesgue measure on region $\Omega$. The inner product is like the usual $L^2$ inner product, $\langle f,g\rangle = \int_{\varOmega}\overline{f(z)} \, g(z) \, (1-|z|^2) \, \mathrm{d}A(z) \} $

I'm looking at $\Omega$ of the form:

Ω = RegionDifference[Disk[{0, 0}, 1], Disk[{1, 0}, b]]

which looks like this:

Can someone help me out with the code? I'm quite new to Mathematica and don't know much about it. I'm trying to see what orthogonalization of the set $\{1,z,z^2,\dots\}$ looks like.

Do you already have a function for evaluating your inner product? — J. M.'s missing motivation, Jul 13 '22 at 17:34
This is a subset of the complex plane but the integral is interpreted as an area integral over $R^2$. — Aditya Ghosh, Jul 14 '22 at 08:28

score 8 · Accepted Answer · edited Jul 14 '22 at 08:44

8

If floating point approximations suffice, then this leads quite quickly to results:

b = 1/10;
Ω = DiscretizeRegion[RegionDifference[Disk[{0, 0}, 1], Disk[{1, 0}, b]],
                     MaxCellMeasure -> (1 -> 0.01)]; 
innerprod[f_, g_] :=
Integrate[f Conjugate[g] (1 - (x + I y) Conjugate[(x + I y)]), {x, y} ∈ Ω]; 
basis = (x + I y)^Range[0, 5]; 
orthobasis = Orthogonalize[basis, innerprod];

If you want exact computations, then you can try to leave away DiscretizeRegion. But be prepared to wait very long time.

edited Jul 14 '22 at 08:44

J. M.'s missing motivation

124,525
11
401
574

answered Jul 13 '22 at 19:25

Henrik Schumacher

106,770
7
179
309

Thank you for this! It works although it is quite slow. Just a couple of follow-up questions. Is there a way to extract the coefficients and tabulate them? – Aditya Ghosh Jul 14 '22 at 08:26
1

have a look at CoefficientRules. – Henrik Schumacher Jul 14 '22 at 08:29

score 4 · Answer 2 · answered Jul 14 '22 at 10:01

In case you need/want symbolic results, you have to do a bit more work to get results in a reasonable amount of time.

I'll use the same parameters used in Henrik's answer:

With[{b = 1/10}, Ω = RegionDifference[Disk[{0, 0}, 1], Disk[{1, 0}, b]]]

The trick is to recognize that computing the inner product $\langle f,g\rangle$ when both functions $f,g$ are polynomials entails a lot of repetitive work, due to the evaluation of expressions of the form $\langle z^m,z^n\rangle$. It thus makes sense to evaluate and memoize these subexpressions so that they get used in the inner product computation:

iProdPow[m_, n_] := iProdPow[m, n] = 
   Simplify[Integrate[ComplexExpand[(x + I y)^m Conjugate[(x + I y)^n]
                                    (1 - (x + I y) Conjugate[x + I y]),
                                    TargetFunctions -> {Re, Im}],
                      {x, y} ∈ Ω]]
iProd[f_, g_, z_] /; PolynomialQ[f, z] && PolynomialQ[g, z] := 
   With[{fe = CoefficientList[f, z], ge = CoefficientList[g, z]}, 
        Simplify[Sum[fe[[j + 1]] ge[[k + 1]] iProdPow[j, k],
                     {j, 0, Length[fe] - 1}, {k, 0, Length[ge] - 1}]]]

I've found Orthogonalize[] to be a bit slow for this application, so I decided to ~~steal~~ borrow and modify acl's Gram-Schmidt code instead:

With[{n = 3},
     res = Block[{u, v},
                 v = z^Range[0, n]; u = ConstantArray[0, n + 1]; u[[1]] = v[[1]];
                 Do[u[[i]] = 
                    Collect[v[[i]] - 
                             Sum[u[[j]] iProd[v[[i]], u[[j]], z]/
                                 iProd[u[[j]], u[[j]], z], {j, 1, i - 1}],
                            z, Simplify],
                    {i, 2, n + 1}];
                 u/Sqrt[iProd[#, #, z] & /@ u]]];

The resulting expressions are quite complicated:

LeafCount /@ res
   {33, 158, 1644, 28985}

but evaluating N[] and comparing them with Henrik's results (and squinting a little bit) shows that the results from this method are consistent with his.

I ran your and Henrik's code for n = 10 case. Henrik's code took 57s to run while yours took 120s (I changed your code to NIntegrate as I only wanted approx values). I'm not sure why your code takes longer. Seems like yours would save time. For n = 20, the program doesn't terminate (in a reasonable time) for either code. — Aditya Ghosh, Jul 14 '22 at 20:55

Mathematica code for Gram Schmidt Orthogonalization for a general Inner Product (in particular, weighted Bergman spaces)

2 Answers2