How to compute a definite integral of the form $\int_a^b f(x)g(x)dx$? Is there perhaps a trapezoidal sum formula for $f(x)*g(x)$?

Question

I need help with this to verify results that I obtained while running a program in Matlab which utilized its built in Trapezoidal Sum function to numerically compute an integral. So, here is an example that I went over which follows the form of $x[n] = \int_a^{b} f(x)g(x)dx$, what I tried, and what I know I'm supposed to get.

The following equations and setup were used below:

1. Trapezoidal Rule (as it is defined in Matlab): $\int_a^b f(x)dx \approx \frac{b-a}{2N} \sum^N_{n=1}(f(x_n) + f(x_{n+1}))$. Repurposed for our example, we sum over each $\omega$ instead of n. We can see it does not include $g(x)$
2. For our f(x) equivalent, we have $f(e^{j\omega}) = \sum_{n=-\infty}^{n = +\infty}e^{j\omega_0 n}e^{-j\omega n} = \sum_{n=-\infty}^{n = +\infty}e^{j(\omega_0-\omega)n} = \sum_{n=-\infty}^{n = +\infty}2\pi\delta{(\omega-\omega_0-2\pi l)}$. This is commonly known as an impulse train and the sum is simply equal to $2\pi$ at $\omega = \omega_0 \pm 2\pi$ and we'll be integrating over an interval of $2\pi$ so this term shouldn't be an issue.
3. And our g(x) equivalent is $g(e^{jw}) = e^{j\omega}$ where again, we will be setting the value of n to some constant. This is a complex exponential and the fact that it is complex should not complicate anything.
4. This gives us the full form of our integral as: $\require{cancel} x[\cancel{n}\color{red}{i}]= \frac{1}{2\pi} \int_{2\pi}X(e^{j\omega})e^{j\omega n}d\omega = \frac{1}{2\pi} \int_{2\pi}\sum_{l=-\infty}^{l = +\infty}2\pi\delta{(\omega-\omega_0-2\pi l)}e^{j\omega n}d\omega$
   $\color{red}{\textrm{where n is }n = [0,1] \textrm{ with }dn = 0.1\textrm{ and }\omega =\omega(n[i])}$
5. Lastly, we will be going over an interval of 11 terms also meaning there will be 11 $\omega$. It will be within the domain of integration of $[0,2\pi]$ and will have a $d\omega = \frac{2\pi}{N-1} = \frac{2\pi}{10}$ with the N-1 being due to using a finite scheme so we have N-1 deltas. Lastly, our $\omega_0$ term is $\omega_0 = \frac{2\pi}{N} = \frac{2\pi}{10} = \frac{\pi}{5}$ (which just so happens to be the same as our $d\omega$).

With all of this laid out, we start can start computing everything. The following are what I obtained for each term:

1. $f(e^{j\omega}) = \sum_{n=-\infty}^{n = +\infty}e^{j\omega_0 n}e^{-j\omega n} = \sum_{n=-\infty}^{n = +\infty}e^{j(\omega_0-\omega)n} = \sum_{l=-\infty}^{l = +\infty}2\pi\delta{(\omega-\omega_0-2\pi l)} = [0,2\pi,0,0,0,0,0,0,0,0,0]$

Note: $\require{cancel}\omega = [0,\frac{\pi}{5},\frac{2\pi}{5},\frac{3\pi}{5},\frac{4\pi}{5},\pi,\frac{6\pi}{5},\frac{7\pi}{5},\frac{8\pi}{5},\frac{9\pi}{5},2\pi]$

2. $\require{cancel}n = \cancel{1}\color{red}{n[1] = \frac{1}{10}}; g(e^{j\omega}) = e^{j\omega * \cancel{1}\color{red}{\frac{1}{10}}} = [1,e^{1j\frac{\pi}{\color{red}{50}}},e^{1j\frac{2\pi}{\color{red}{50}}},e^{1j\frac{3\pi}{\color{red}{50}}},e^{1j\frac{4\pi}{\color{red}{50}}},e^{1j\frac{5\pi}{\color{red}{50}}},e^{1j\frac{6\pi}{\color{red}{50}}},e^{1j\frac{7\pi}{\color{red}{50}}},e^{1j\frac{8\pi}{\color{red}{50}}},e^{1j\frac{9\pi}{\color{red}{50}}},e^{\frac{10\pi}{\color{red}{50}}}]$

Now we have these two arrays and we have to multiply them together in some way. So what I did next was I then plugged what are analogously our f(x) and g(x) functions into the Trapezoidal rule to get the following computation:

$\frac{b-a}{2N} = \frac{2\pi-0}{2(11-1)} = \frac{\pi}{10}$

$\require{cancel}x[1] = \frac{1}{2\pi}\frac{\pi}{10}\sum^N_{i=1}((f(e^{j\omega_i})*g(e^{j\omega_i})) + (f(e^{j\omega_{i+1}})*g(e^{j\omega_{i+1}})) = ((0 *1)+(2\pi *e^{1j\frac{\pi}{5}}))+((2\pi *e^{1j\frac{\pi}{\color{red}{50}}})+(0 *e^{1j\frac{2\pi}{\color{red}{50}}})) + 0 + 0 + 0 + 0 + 0 + 0 + 0 + 0 + 0= 4\pi e^{1j\frac{\pi}{\color{red}{50}}} \color{red}{*\frac{1}{20}}= \cancel{0.5083 + 0.3693j} \color{red}{0.6271 + 0.0395j}$

I decided to switch the summation index to "m" instead of "n" which I left reserved to x[n] which is just a subtle difference to help suit my version of the problem.

So basically, I took each term for f(x) and g(x) at their respective index and I multiplied them together. I know this is incorrect but I wanted to at least show my though process on this. The correct result which was obtained through Matlab is actually 0.62707 + 0.0394i which likely used some conditional subroutine hidden in the code to compute the integral properly in the instance that the integrand is of the form $f(x) * g(x)$ but does not share that version of the formula with the user. Does anyone know what I was supposed to do instead to arrive to this result? Is there perhaps a modified formula for the Trapezoidal sum for this type of problem?

EDIT: As requested, I have provided the code below. It goes over the setup, building what are analgously f(x) and g(x) which I have commented for which is which, and the computation of the integral of f(x) * g(x) using Matlab's built in subroutine trapz(dx,x) which from the documentation uses the formula presented in what is listed as the first item in the list discussing the equations and setup I used.

Code used:

% Setup
    L = 5;len_n = 11; len_ninf = 101;
    n_inf = linspace(-L,L,len_ninf);
    omega = 2*pi .* n_inf;
    N = length(n);
    omega_o = 2*pi/((length(n_inf)-1)/10);
    l = -500:1:500; % some harmonic
%build f(x) and g(x)
    x = exp(1j.omega_o.n_inf);           % g(x)
    impulsetrain = zeros(1,length(omega)); % f(x)
        for i = 1:length(omega)
            for j =1:length(l)
                if round(((round(omega(i),2)- round(omega_o,2)) -2round(pi,2)l(j))/1000,4)==0
                    impulsetrain(i) = 2pi1;
                end
            end
        end
% Compute integral from a to b of f(x) * g(x)$
    del_omega = (omega(length(omega))-omega(1))/(length(omega)-1);
    LeftBound = 0; RightBound = 2*pi + LeftBound;
    range = find(round(omega,2) == round(LeftBound,2)):find(round(omega,2) == round(RightBound,2));
x_dift = zeros(1,length(n_inf(range))); % integral of f(x) * g(x)
for i = 1:length(n_inf(range))
    x_dift(i) = 1/(2*pi) * trapz(omega(range),impulsetrain(range).*exp(1j*n_inf(range(i)).*omega(range)));
end

x_dift is the object holding the values for the integral of what is analogously f(x) * g(x). Going to x_dift(1) will show the value which we are trying to obtain. The way in which I hand computed it I did not code as I feel I sufficiently showed the process of what I was doing with the equation and working that I typed out with x[1], however I could just as easily add it should someone ask.

Thanks in advance!

Answer 1

Too long for a comment. As far as I gather you want to integrate $f(x)g(x)$ where \begin{align} f(e^{j\omega}) &= \sum_{n=-\infty}^{n = +\infty}e^{j\omega_0 n}e^{-j\omega n} = \sum_{n=-\infty}^{n = +\infty}e^{j(\omega_0-\omega)n} = \sum_{n=-\infty}^{n = +\infty}2\pi\delta{(\omega-\omega_0-2\pi\color{red}{l})}\,,\\ g(e^{jw}) &= e^{j\omega n}\,. \end{align} Even if we fix that notation in the last formula I think it stays wrong: \begin{align*} f(\omega)&=\sum_{n=-\infty}^{+\infty}e^{j\,\omega n}=\sum_{n=-\infty}^{+\infty}\cos(\omega n)+j\,\sin(\omega n)\,. \end{align*} Since $\sin(-\omega n)=-\sin(\omega n)$ the $\sin$-terms all cancel. Since $\cos$ is even we get \begin{align*} f(\omega)&=1+2\sum_{n=1}^{+\infty}\cos(\omega n)\,. \end{align*} According to this post that series does not look like it converges.

Regarding your numerical integration exercise: I strongly advice to test that code on functions that are less complicated and perhaps have integrals where the solution is known. Even if you have a delta function that you'd like to integrate how would you do that numerically? How would you even implement that delta function? Assume you approximate it with $$ \frac{1}{\sqrt{2\pi\sigma}}\exp(-\frac{x^2}{2\sigma^2}) $$ how would you design your numerical integration when $\sigma$ is very small? Note that this function has a known integral from $-\infty$ to $+\infty\,.$