Consider the Hopf fibrations $S^1\to S^{2n+1}\to CP^n$ and $S^3\to S^{4n+3}\to HP^n$. These are Riemannian submersions with totally geodesic fibers. Consider now their canonical variations (the so-called Berger spheres), i.e., for each $t>0$, let $S^{2n+1}_t$ (resp $S^{4n+3}_t$) be the Riemannian manifold $S^{2n+1}$ (resp $S^{4n+3}$) endowed with the metric obtained by scaling the round metric by a factor of $t^2$ in the direction of the fibers $S^1$ (resp $S^3$). For each $t>0$, they are the total space of a Riemannian submersion with totally geodesic fibers isometric to $tS^1$ (resp $tS^3$). In this context, the Laplacian $\Delta_t$ of the canonical variation can be related to the original Laplacian $\Delta_1$ in terms of the vertical Laplacian $\Delta_v$, by the formula $$\Delta_t=\Delta_1+(\tfrac{1}{t^2}-1)\Delta_v,$$ see Berard-Bergery and Bourguignon [Illinois Journ of Math, 1982]. The operator $\Delta_v$ acts on functions of the total space and is defined by $$(\Delta_v f)(p)=\Delta_{F_p}(f\vert_{F_p})$$ where $\Delta_{F_p}$ is the Laplacian of the fiber $F_p$.
Since the fibers are totally geodesic, the spectrum of $\Delta_v$ coincides with the spectrum of the Laplacian of the fiber, and all the above operators commute. In particular, $L^2$ of the total space admits a decomposition in simultaneous eigenspaces of $\Delta_t$ and $\Delta_v$. As a consequence, we have the following inclusion of spectra: $$Spec(\Delta_t) \subset Spec(\Delta_1)+(\tfrac{1}{t^2}-1) Spec(\Delta_F).$$
I am interested in computing $\lambda_1(t)$, the first non-zero eigenvalue of $\Delta_t$, for all $t$. By the above, for all $t>0$, we know $\lambda_1(t)$ is of the form $\mu_k+(\tfrac{1}{t^2}-1)\phi_j$ where $\mu_k\in Spec(\Delta_1)$ is an eigenvalue of the Laplacian on the original total space and $\phi_j\in Spec(\Delta_F)$ is an eigenvalue of the Laplacian of the fiber. The problem is that not all combinations of $\mu_k$'s and $\phi_j$'s give an eigenvalue of $\Delta_t$ and the first non-zero eigenvalue of $\Delta_t$ might be given by different combinations of $\mu_k$'s and $\phi_j$'s for each $t$. Recall that since total space and fibers are spheres, both $\mu_k$'s and $\phi_j$'s are easy to compute, namely the $k$th eigenvalue of the $m$-sphere is $k(k+m-1)$.
In the case of the first family with 1-dim fibers, this computation follows from a paper of Tanno [Tohoku Math Journ, 1979]. The trick is to look at a trajectory of the vector tangent to the $S^1$ fiber (a great circle) and solve the eigenvalue equation there (which becomes an ODE). Using this he proves that the only combinations of $\mu_k$'s and $\phi_j$'s permitted are when $0\leq j\leq k$ and $k-j$ is even. He then also finds explicit eigenfunctions when $k=j$; $j=1$ and $k$ odd; $j=2$ or $j=0$ and $k$ even. With these, it is very easy to compute $\lambda_1(t)$ for this family where the fibers are 1-dim.
However, I do not know any way of extending the result to the case of 3-dim fibers, or also the case of 7-dim fibers, $S^7\to S^{15}\to S^8$. It seems that Tanno's technique uses strongly the fact that the fibers are 1-dim. Any suggestion on how to compute this other first eigenvalues would be greatly appreciated.