I'm working modeling the behavior of periodic variable stars and I have a question about reducing the expression of a parameter involved in this analysis.

Let $f(t)$ be a Fourier series define as:

$$f(t)= \langle m \rangle + \displaystyle{\sum_{k=1}^{N} A_{k}\sin \left(\frac{2\pi kt}{P} + \phi_{k}\right)}.$$

The variable $t$ is real and defined in the so-called **phase space**, where $0\leq t \leq 1$. These series describe a periodic process with period $P$, $A_i$ and $\phi_i$ the individual Amplitude and phases for all the $N$ harmonics; and $\langle m \rangle$ the average value of the time-series that contains the process. All the shown variables are known. I want to figure out the length of this Fourier series using the definition:

$${\displaystyle s=\int _{a}^{b}{\sqrt {1+\left(f'(x)\right)^{2}}} dx.}$$

If we define $x_k:= \frac{2\pi k}{P}$, the derivative of the Fourier series $f'(t)$ will be:

$$ f'(t) = \displaystyle{ \sum_{k=1}^{N} A_{k} x_{k} \cos(x_{k} t + \phi_{k}).}$$

Then: \begin{split} (f'(t))^{2}=\displaystyle{\left(\sum_{k=1}^{N} A_{k} x_{k} \cos(x_{k} t + \phi_{k})\right)^{2}}\\ = \left(\sum_{i=1}^{N} A_{i} x_{i} \cos(x_{i} t + \phi_{i})\right) \times \left(\sum_{j=1}^{N} A_{j} x_{j} \cos(x_{j} t + \phi_{j})\right) \\ = \sum_{i,j=1}^{N} a_i a_j \cos(x_{i} t + \phi_{i}) \cos(x_{j} t + \phi_{j}), \end{split}

where we also define $a_k := A_k x_k$. Finally our lenght will be:

$${\displaystyle s=\int _{a}^{b}{\sqrt {1+ \sum_{i,j=1}^{N} a_i a_j \cos(x_{i} t + \phi_{i}) \cos(x_{j} t + \phi_{j}) }}\ \ dt.}\ \ (1)$$

The integration limits are $a=0$ and $b=1$, given we live in this convenient phase space. Now comes the big question:

Did you have some idea or trick to attack and solve this integral? I tried with some trigonometric manipulations, but I got not great progress. I'm trying to obtain the less-expensive function to compute and so far, is the expression (1).

Thanks for reading this and every kind of feedback will be awesome :)