5
$\begingroup$

For a function $f\colon [-1,1]\to\mathbb{R}$, the Gagliardo seminorm of $f$ is defined to be

$$ |f| = \int_{-1}^1\int_{-1}^1 \frac{(f(x)-f(y))^2}{(x-y)^2}\, \mathrm{d}x\, \mathrm{d} y. $$

Given $(x,x^2,\ldots,x^n)$ for some $n$, I am interested in finding an orthonormal basis of these functions with respect to this seminorm. (We have to start at $x^1$ and not $x^0$ since the seminorm of the latter is zero).

The function Orthogonalize is normally how I do something like this. Indeed, in the help entry for Orthogonalize under Generalizations and Extensions, there is something similar to what I want to do. I can't figure out the notation/syntax for this norm.

What do I need to do?

$\endgroup$
  • 4
    $\begingroup$ To orthogonalize a set of vectors/functions you need an inner product, not just a norm. $\endgroup$ – Marius Ladegård Meyer Apr 12 '16 at 13:05
  • 4
    $\begingroup$ @MariusLadegårdMeyer Yeah, I ran into that too. We could use the polarization identity to construct the inner product from the norm, (hopefully). Or I was thinking that we could use $\langle f,g \rangle = \int_{-1}^1\int_{-1}^1 (f(x)-f(y))(g(x)-g(y)) / (x-y)^2\, \mathrm{d}x\, \mathrm{d} y $ for the inner product. $\endgroup$ – fred Apr 12 '16 at 13:10
9
$\begingroup$

Perhaps this?

ip[f_, g_, v_] := Module[{x, y, int},
  int = (((f /. v -> x) - (f /. v -> y)) ((g /. v -> x) - (g /. v -> y)))/(x - y)^2;
  Integrate[int, {x, -1, 1}, {y, -1, 1}]
  ];

Orthogonalize[x^Range[4], ip[##, x] &]
(* 
  {x/2, 1/2 Sqrt[3/2] x^2, 3/2 Sqrt[5/13] (-((2 x)/3) + x^3), 
   15/8 Sqrt[21/31] (-((14 x^2)/15) + x^4)}
*)

The inner product ip assumes that f and g are real-valued functions of the parameter v.


Update: For the sake of speed on polynomials

If you're dealing only with polynomials, one can get a speed-up by precomputing the integrals of monomials.

Block[{i, j, k},   (* protect i,j,k since they will be evaluated during the definition *)
  Assuming[k ∈ Integers && k >= 0,
   ip`power[i_Integer?EvenQ] =        (* integral of even powers *)
    Integrate[x^(2 k) , {x, -1, 1}] /. k -> i/2;
   ip`moment[i_Integer?EvenQ, j_Integer?EvenQ] =
    ip`power[i] ip`power[j];          (* the integral of x^i y^j is a product *)
   ]];
ip`moment[_Integer, _Integer] := 0;   (* odd powers integrate to zero *)

ipfast[f_, g_, v_] := Module[{x, y, int, coeff, mom},
  int = (((f /. v -> x) - (f /. v -> y)) ((g /. v -> x) - (g /. v -> y)))/(x - y)^2;
  mom = SparseArray@CoefficientList[Factor@int, {x, y}];  (* get coefficient list *)
  coeff = mom["NonzeroValues"];                           (* get coeffs of monomials *)
  mom = mom["NonzeroPositions"] - 1;                      (* get powers of monomials *)
  coeff.(ip`moment @@@ mom)                               (* integrate and add results *)
  ]

b = x^Range[5];
res = Orthogonalize[b, ip[##, x] &]; // AbsoluteTiming // First
res2 = Orthogonalize[b, ipfast[##, x] &]; // AbsoluteTiming // First
(*
  7.20899
  0.021552
*)

res - res2 // Simplify
(*  {0, 0, 0, 0, 0}  *)
|improve this answer|||||
$\endgroup$
5
$\begingroup$

I seem to be getting different results depending on whether I use OP's suggested inner product in the comments, or if I use the polarization identity.

Before everything else, however, here is a routine by Velvel Kahan for evaluating the divided difference of a polynomial, based on Horner's method:

polynomialDividedDifference[poly_, {x_, a_, b_}] /; PolynomialQ[poly, x] :=
          Module[{d = 0, y = 0},
                 Do[y = b y + Coefficient[poly, x, k]; d = a d + y,
                    {k, Max[0, Exponent[poly, x]], 1, -1}];
                 Expand[d]]

With this, here is the result of using the inner product

$$\langle f,g\rangle=\int_{-1}^1\int_{-1}^1\frac{(f(u)-f(v))(g(u)-g(v))}{(u-v)^2}\,\mathrm du\mathrm dv$$

Orthogonalize[x^Range[4], 
              Integrate[polynomialDividedDifference[#1, {x, u, v}]
                        polynomialDividedDifference[#2, {x, u, v}],
                        {u, -1, 1}, {v, -1, 1}] &]
   {x/2, 1/2 Sqrt[3/2] x^2, 3/2 Sqrt[5/13] (-((2 x)/3) + x^3), 
    15/8 Sqrt[21/31] (-((14 x^2)/15) + x^4)}

(same as in Michael's answer)

Contrast this with the result of using the inner product

$$\langle f,g\rangle=\frac14\left(\left(\int_{-1}^1\int_{-1}^1\frac{(f(u)-f(v)+g(u)-g(v))^2}{(u-v)^2}\,\mathrm du\mathrm dv\right)^2-\left(\int_{-1}^1\int_{-1}^1\frac{(f(u)-f(v)-(g(u)-g(v)))^2}{(u-v)^2}\,\mathrm du\mathrm dv\right)^2\right)$$

based on the polarization identity:

Orthogonalize[x^Range[4],
              (Integrate[polynomialDividedDifference[#1 + #2, {x, u, v}]^2,
                         {u, -1, 1}, {v, -1, 1}]^2 - 
               Integrate[polynomialDividedDifference[#1 - #2, {x, u, v}]^2,
                         {u, -1, 1}, {v, -1, 1}]^2)/4 &]
   {x/8, (3 x^2)/16, (1036800 (-((719 x)/2880) + x^3))/3838561,
    (17920000 (-((3177 x^2)/6400) + x^4))/33294461}

The OP will have to choose which basis is most appropriate for his situation.

|improve this answer|||||
$\endgroup$
  • $\begingroup$ Hmmm, perhaps things are different because we only have a semi-norm instead of a norm. I will have to think about this. $\endgroup$ – fred Apr 12 '16 at 17:41
  • $\begingroup$ Right, that's why I said that you need to do further investigations on which inner product is the right one for your application. $\endgroup$ – J. M.'s technical difficulties Apr 12 '16 at 17:43

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.