Let $a, b, x, y \in \mathbf{R}$, satisfy the system of equations
$$
\left\{\begin{array}{l}
a x+b y=3, \\
a x^{2}+b y^{2}=7, \\
a x^{3}+b y^{3}=16, \\
a x^{4}+b y^{4}=42 .
\end{array}\right.
$$

Find the value of $a x^{5}+b y^{5}$.
$(1990$, American Mathematical Invitational)
Solution 1:
From
$a x^{3}+b y^{3}$
we get
$$
=\left(a x^{2}+b y^{2}\right)(x+y)-(a x+b y) x y,
$$
thus
$$
\begin{array}{l}
16=7(x+y)-3 x y . \\
\text { From } a x^{4}+b y^{4} \\
=\left(a x^{3}+b y^{3}\right)(x+y)-\left(a x^{2}+b y^{2}\right) x y,
\end{array}
$$

From $a x^{4}+b y^{4}$
we get
$$
42=16(x+y)-7 x y \text {. }
$$

From equations (2) and (3), we solve to get
$$
\begin{array}{l}
x+y=-14, x y=-38 . \\
\text { Therefore, } a x^{5}+b y^{5} \\
=\left(a x^{4}+b y^{4}\right)(x+y)-\left(a x^{3}+b y^{3}\right) x y \\
=42 \times(-14)-16 \times(-38)=20 .
\end{array}
$$

Solution 2: This problem can be handled using the perspective of a second-order linear recurrence sequence. The general term of a second-order linear recurrence sequence is
$$
a_{n}=a x^{n}+b y^{n} \text {. }
$$

Conversely, $\left\{a x^{n}+b y^{n}\right\}$ is a second-order recurrence sequence. The recurrence relation is
$$
a_{n+1}=c a_{n}+d a_{n-1},
$$

where $c=x+y, d=-x y$.
Thus,
$$
\left\{\begin{array}{l}
16=7 c+3 d, \\
42=16 c+7 d, \\
a_{5}=42 c+16 d .
\end{array}\right.
$$

Considering this as a system of equations in $1, c, d$,
$$
\left\{\begin{array}{l}
16-7 c-3 d=0, \\
42-16 c-7 d=0, \\
a_{5}-42 c-16 d=0 .
\end{array}\right.
$$

Viewing $(1, c, d)$ as a non-zero solution to the system of equations (4), the determinant of the coefficient matrix is 0, i.e.,
$$
\left|\begin{array}{lll}
16 & -7 & -3 \\
42 & -16 & -7 \\
a_{5} & -42 & -16
\end{array}\right|=0 \text {. }
$$

Solving this, we get $a_{5}=20$.
This is a clever solution, but some students might ask if there is a direct method, using "brute force" to solve for $a, b, x, y$ from the system of equations (1) and then substitute back into $a x^{5}+b y^{5}$. Of course, this is a complex approach for ordinary people, but for the Indian mathematician Ramanujan (1887-1920), it would be a piece of cake.
Ramanujan proposed and solved the following tenth-order system of equations:
$$
\begin{array}{l}
x+y+z+u+v=2, \\
p x+q y+r z+s u+t v=3, \\
p^{2} x+q^{2} y+r^{2} z+s^{2} u+t^{2} v=16, \\
p^{3} x+q^{3} y+r^{3} z+s^{3} u+t^{3} v=31, \\
p^{4} x+q^{4} y+r^{4} z+s^{4} u+t^{4} y=103, \\
p^{5} x+q^{5} y+r^{5} z+s^{5} u+t^{5} v=235, \\
p^{6} x+q^{6} y+r^{6} z+s^{6} u+t^{6} v=674, \\
p^{7} x+q^{7} y+r^{7} z+s^{7} u+t^{7} v=1669, \\
p^{8} x+q^{8} y+r^{8} z+s^{8} u+t^{8} v=4526, \\
p^{9} x+q^{9} y+r^{9} z+s^{9} u+t^{9} v=11595 .
\end{array}
$$

Ramanujan first considered the general system of equations
$$
x_{1}+x_{2}+\cdots+x_{n}=a_{1} \text {, }
$$
$$
\begin{array}{l}
x_{1} y_{1}+x_{2} y_{2}+\cdots+x_{n} y_{n}=a_{2}, \\
x_{1} y_{1}^{2}+x_{2} y_{2}^{2}+\cdots+x_{n} y_{n}^{2}=a_{3}, \\
\cdots \cdots \cdots \\
x_{1} y_{1}^{2 n-1}+x_{2} y_{2}^{2 n-1}+\cdots+x_{n} y_{n}^{2 n-1}=a_{2 n} . \\
\text { Let } F(\theta)=\frac{x_{1}}{1-\theta y_{1}}+\frac{x_{2}}{1-\theta y_{2}}+\cdots+\frac{x_{n}}{1-\theta y_{n}} . \\
\text { But } \frac{x_{i}}{1-\theta y_{i}}=x_{i}\left(1+\theta y_{i}+\theta^{2} y_{i}^{2}+\theta^{3} y_{i}^{3}+\cdots\right)
\end{array}
$$
$(i=1,2, \cdots, n)$, hence
$$
\begin{array}{l}
F(\theta)=\sum_{i=1}^{n} x_{i} \cdot \sum_{k=0}^{\infty}\left(\theta y_{i}\right)^{k} \\
=\sum_{k=0}^{\infty}\left(\sum_{i=1}^{n} x_{i} y_{i}^{k}\right) \theta^{k}=\sum_{k=0}^{\infty} a_{k+1} \theta^{k}
\end{array}
$$

Converting it to a rational function with a common denominator, we get
$F(\theta)=\frac{A_{1}+A_{2} \theta+A_{3} \theta^{2}+\cdots+A_{n} \theta^{n-1}}{1+B_{1} \theta+B_{2} \theta^{2}+\cdots+B_{n} \theta^{n}}$.
Thus, $\sum_{k=0}^{\infty} a_{k+1} \theta^{k} \cdot \sum_{s=0}^{n} B_{s} \theta^{s}=\sum_{i=1}^{n} A_{t} \theta^{t-1}\left(B_{0}=1\right)$.
Hence, $A_{t}=\sum_{k=0}^{t-1} a_{k+1} B_{t-1-k}(t=1,2, \cdots, n)$, $0=\sum_{s=0}^{n} B_{s} a_{n+l-s}(t=1,2, \cdots, n)$.
Since $a_{1}, a_{2}, \cdots, a_{n}, a_{n+1}, \cdots, a_{2 n}$ are known,
we can solve for $B_{1}, B_{2}, \cdots, B_{n}$ from the last $n$ equations, and then substitute into the first $n$ equations to solve for $A_{1}, A_{2}, \cdots, A_{n}$. Knowing $A_{i} 、 B_{i}(i=1,2, \cdots, n)$, we can construct the rational function $F(\theta)$, and then expand it into partial fractions. Thus, we get
$$
F(\theta)=\frac{p_{1}}{1-q_{1} \theta}+\frac{p_{2}}{1-q_{2} \theta}+\cdots+\frac{p_{n}}{1-q_{n} \theta} .
$$

Clearly, $x_{i}=p_{i}, y_{i}=q_{i}(i=1,2, \cdots, n)$.
This is the solution to the general system of equations.
For the considered case, we have
$F(\theta)=\frac{2+\theta+3 \theta^{2}+2 \theta^{3}+\theta^{4}}{1-\theta-5 \theta^{2}+\theta^{3}+3 \theta^{4}-\theta^{5}}$.
Expanding it into partial fractions, we get the values of the unknowns:
$x=-\frac{3}{5}, p=-1$;
$y=\frac{18+\sqrt{5}}{10}, q=\frac{3+\sqrt{5}}{2} ;$
$z=\frac{18-\sqrt{5}}{10}, r=\frac{3-\sqrt{5}}{2}$;
$u=-\frac{8+\sqrt{5}}{2 \sqrt{5}}, s=\frac{\sqrt{5}-1}{2}$;
$v=\frac{8-\sqrt{5}}{2 \sqrt{5}}, t=-\frac{\sqrt{5}+1}{2}$.
Readers can apply this method to the problem mentioned at the beginning.