# Double Integral for Gaussian Quadrature

I am trying to apply the Gaussian-Legendre Quadrature rule to a Double Integral, namely

$$\int^1_0 \int^1_0 \text{sin}(x^2+y^2)dxdy$$

I have done the following:

Define $$\phi_n(x)$$ is the Legendre polynomial of degree $$n$$. Using the fact that \begin{align*} \int^b_{a} f(x)dx &= \frac{b-a}{2}\int^1_{-1}f\left(\frac{b-a}{2}x+\frac{a+b}{2} \right)dx\\ &= \frac{b-a}{2}\sum^{n}_{i=1}w_if\left(\frac{b-a}{2}x_i+\frac{a+b}{2}\right) \end{align*}

with

$$w_i = \frac{2}{(1-x_i^2)[\phi_n '(x_i)]^2}$$

where each $$x_i$$ is a root of the polynomial $$\phi_n(x)$$, I have applied this logic to the double integral case and obtained

$$\int^1_0\int^1_0 F(x,y)dxdy = \frac{1}{4}\sum^n_{i=1}\sum^n_{j=1}w_i w_jF\left(\frac{1}{2}x_i+\frac{1}{2},\frac{1}{2}y_i+\frac{1}{2}\right)$$

with

$$w_i = \frac{2}{(1-x_i^2)[\phi_n '(x_i)]^2} \ \ \text{ and } \ \ w_j = \frac{2}{(1-y_j^2)[\phi_n '(y_j)]^2}$$

However, this logic has not worked under implementation (using C++ code).

Is the logic here sound? Or have I oversimplied the issue?

Note, I looked throughly online and on this website for an and couldn't find anything conclusive - if there is a link, I would appreciate it if you could direct me over.

• Of course there's no need to do this one as a double integral, if you use the expansion $\sin(x^2+y^2) = \sin(x^2)\cos(y^2)+\cos(x^2)\sin(y^2)$. – Robert Israel Apr 15 at 2:43
• Oh yes, I am aware - however the assignment specifically requested that we make use of double integrals. Thank you though – Naji Apr 15 at 2:47

Have you checked your quadrature method on $$1$$-variable integrals? How does it do on $$\int_0^1 \int_0^1 \sin(x^2)\; dx\; dy$$ and $$\int_0^1 \int_0^1 \cos(y^2)\; dx\; dy$$?
• The approximation does not seem to be converging to $\frac{1}{6\sqrt{2}}$ as required, it seems to rather converge (slowly) to $0$ as $n$ increases. Following your advice, I evaluated both double integrals that you posted and they don't seem giving the correct approximation. – Naji Apr 15 at 2:56
• Seems like some problem with your programming then. BTW, the answer should not be $1/(6 \sqrt{2})$, it should be $\pi FresnelS(\sqrt{2/\pi}) FresnelC(\sqrt{2/\pi})$ or approximately $0.5612903983$. – Robert Israel Apr 15 at 4:14