Error in Creating Orthogonal Polynomials

Question

I'm trying to create my own set of polynomials orthogonal to weight $w(x)=x^{14}$ on $[-a,a]$. My code:

N1=10;
Int0[n_,x_]=Integrate[x*(DUO[n])^2*x^14,x];
Int1[n_,x_]=Integrate[(DUO[n])^2*x^14,x];

DUO[0]=1;
DUO[1]=(x-Divide[Int0[0,N1]-Int0[0,-N1],Int1[0,N1]-Int1[0,-N1]])*DUO[0];

DUO[n_]:=DUO[n]=((x-Divide[Int0[n-1,N1]-Int0[n-1,-N1],Int1[n-1,N1]-Int1[n-1,-N1]])*DUO[n-1])-(Divide[Int1[n-1,N1]-Int1[n-1,-N1],Int1[n-2,N1]-Int1[n-2,-N1]]*DUO[n-2]);
DUO[2]

However, when evaluating DUO[2], I keep on getting random and non-sensical errors. I'm following the Gram-Schmidt Process, and I checked it over several times, and I can't find anything wrong in my implementation (other than it being horribly inefficient :) ). So could someone help me find what's wrong here?

Answer 1

The problem is that Int0 and Int1 are defined with = (Set) rather than := (SetDelayed), so the RHS evaluates immediately. This is problematic because DUO[n] does not yet explicitly depend on x. For example, your definition for Int0 becomes

Int0[n_, x_] = 1/16 x^16 DUO[n]^2

(to see this, run ?Int0).

Also we should directly evaluate the definite integral rather than take the difference between indefinite integrals. We define our inner product:

ip[f_, g_] := Integrate[f g x^14, {x, -a, a}];
ip[f_] := ip[f, f];

Then we implement (18), (19), and (24) in the linked mathworld page:

p[0] = 1;
p[1] = (x - ip[x p[0], p[0]] / ip[p[0]]) * p[0];
p[n_] := p[n] =
  Subtract[
    (x - ip[x p[n - 1], p[n - 1]] / ip[p[n - 1]]) * p[n - 1],
    ip[p[n - 1]] / ip[p[n - 2]] * p[n - 2]
  ];

We get:

Table[{n, p[n]}, {n, 0, 4}] // TableForm

$$ \begin{array}{cc} 0 & 1 \\ 1 & x \\ 2 & x^2-\frac{15 a^2}{17} \\ 3 & x \left(x^2-\frac{15 a^2}{17}\right)-\frac{4 a^2 x}{323} \\ 4 & x \left(x \left(x^2-\frac{15 a^2}{17}\right)-\frac{4 a^2 x}{323}\right)-\frac{289}{399} a^2 \left(x^2-\frac{15 a^2}{17}\right) \\ \end{array}$$

Check orthogonality:

Table[ip[p[m], p[n]], {m, 0, 4}, {n, 0, 4}] // TableForm

$$\begin{array}{ccccc} \frac{2 a^{15}}{15} & 0 & 0 & 0 & 0 \\ 0 & \frac{2 a^{17}}{17} & 0 & 0 & 0 \\ 0 & 0 & \frac{8 a^{19}}{5491} & 0 & 0 \\ 0 & 0 & 0 & \frac{8 a^{21}}{7581} & 0 \\ 0 & 0 & 0 & 0 & \frac{128 a^{23}}{3661623} \\ \end{array}$$