Solving a system of integral equations of the form $\int_0^r K(r,s)f(s)\mathrm{d}s + \int_r^1 Q(r,s)g(s) \mathrm{d}s = 1$ with $s,r\in[0,1]$

Question

Consider the following system of integral equations for the unknown functions $A(s)$ and $B(s)$, \begin{align} \int_0^r \left( \phi_1 \left(\tfrac{s}{r} \right) A(s) + \phi_2\left(\tfrac{s}{r} \right) B(s) \right) \mathrm{d} s + \int_r^1 \left( \phi_3\left(\tfrac{r}{s} \right) A(s) + \phi_4\left(\tfrac{r}{s} \right) B(s) \right) \mathrm{d} s = \frac{1}{2} \, , \tag{1} \\ \int_0^r \left( \psi_1\left(\tfrac{s}{r} \right) A(s) + \psi_2\left(\tfrac{s}{r} \right) B(s) \right) \mathrm{d} s + \int_r^1 \left( \psi_3\left(\tfrac{r}{s} \right) A(s) + \psi_4\left(\tfrac{r}{s} \right) B(s) \right) \mathrm{d} s = 0 \, , \tag{2} \end{align} where $s, r \in [0,1]$. Here, the kernel functions associated with Eq. (1) are given by \begin{align} \phi_1 (x) &= x K (x) \, , \\ \phi_2 (x) &= \frac{2E(x) - \left( 2-x^2\right)K(x)}{x} \, , \\ \phi_3 (x) &= K(x) \, , \\ \phi_4 (x) &= 2E(x) - K(x) \, . \end{align} In addition, the kernel functions associated with Eq. (2) are given by

\begin{align} \psi_1 (x) &= x \phi_4(x) \, , \\ \psi_2 (x) &= \frac{\left( x^2+2\right) K(x) - 2 \left( x^2+1\right) E(x)}{3x} \, , \\ \psi_3 (x) &= \frac{\phi_2(x)}{x} \, , \\ \psi_4 (x) &= \frac{\psi_2 (x)}{x} \, . \end{align}

Here, $K$ and $E$ are the complete elliptic integrals of the first and second kind, respectively.

The convention used here is that of Maple such that $K(x) = (\pi/2) \left( 1 + x^2/4+ \dots \right)$ and $E(x) = (\pi/2) \left( 1 - x^2/4+ \dots \right)$ for $x \ll 1$.

By setting $r = 0$ in Eq. (1) so that only the second integral remains, we obtain \begin{equation} \int_0^1 \left( A(s) + B(s) \right) \mathrm{d} s = \frac{1}{\pi} \, . \end{equation}

Likewise, by dividing both terms of Eq. (2) by $r^2$ and taking the limit when $r \to 0$ we obtain \begin{equation} \int_0^1 \left( A(s) - 3B(s) \right) \frac{\mathrm{d}s }{s^2} = 0 \, . \end{equation}

I attempted to solve the given system numerically using finite differences. The functions exhibit rapid decay as the variable $s$ increases. Notably, at $s=0$, both $A(s)$ and $B(s)$ assume significantly large values. These numerical outcomes were derived using 1, 000 discretization points.

I am now interested in exploring the possibility of obtaining closed-form analytical expressions for $A(s)$ and $B(s)$ through an analytical approach. I would greatly appreciate any assistance in this endeavor.

Thank you in advance!