You can use the formula for geometric sum if the **modular multiplicative inverse** of the denominator of the geometric sum exists. Refer to @joffan 's answer.

Otherwise, you can solve this recursively. Note that:

1 + a + a^2 + a^3 + \dots + a^{2n+1} = (1 + a)(1 + (a^2) + (a^2)^2 + (a^2)^3 + \dots + (a^2)^n)

If we denote 1 + a + a^2 + a^3 + \dots + a^n as S(a,\ n), then

S(a, 2n + 1) = (1 + a) \cdot S(a^2, n)

S(a, 2n) = 1 + a \cdot (1 + a) \cdot S(a^2, n-1)

with base cases:

S(a, 0) = 1

S(a, 1) = 1 + a

Since at each step we are reducing the number of terms to half, the complexity is \mathcal{O}(\log{n}).

You just have to store the answer modulo M at every step.

Here is the Python

```
[1].
[1]: https://ideone.com/gWBvJl
```