In my research, I've come across the following inhomogenous nonlinear ODE ($c \geq 0$ is an undetermined constant):
$$x_{ssss} - x_{ss} x = c$$
It has boundary conditions $$x_s(0) = x_{sss}(0) = x(1) = x_{ss}(1) = x_{sss}(1) = 0$$
I'm trying to understand if anyone has studied this equation, or if there are any analytical techniques I can use to better understand it. Thank you so much!
So I've been thinking about this one for a couple of days because from an asymptotic analysis perspective, it is both tractable and very nontrivial.
What follows is a way to analyse the selection problem for large $c$. It will involve WKB-type analysis and matched asymptotics. I've skipped a lot of detail, and no guarantees on all the constants and signs (the final result seems to line up with numerics though).
Start with the problem: $$ x''''(x) - xx''(s) = c, \qquad x'(0) = x'''(0) = 0, \qquad x(1) = x''(1) = x'''(1) = 0. $$ Assume that $c$ is large and positive, and scale $x = c^{1/2}X$, so that $$ \epsilon X'''' - XX'' = 1, \qquad \epsilon = c^{-1/2} $$ (the boundary conditions remain unchanged). This is in the more `natural' form for a singularly perturbed equation, with a small parameter multiplying the highest derivative term.
Now we consider an outer region (away from the boundary $s=1$), and an inner region (near $s=1$). In the outer region, the leading order (for small $\epsilon$) problem is $$ X_0X_0'' = -1. $$ This equation can be integrated to find a solution in implicit form $$ s = \mathrm{erf}\left[\sqrt{\log\left(-\frac{\sqrt{2/\pi}}{X_0}\right)}\right], \qquad -\sqrt{\frac{2}{\pi}} < X_0 < 0. $$ Here $X_0 = -\sqrt{2/\pi}$ corresponds to $s=0$ and $X_0 = 0$ corresponds to $s=1$. The constants are chosen so that both BCs at $s=0$ are satisfied, and $X_0(1) = 0$. Note that the other conditions at $s=1$ are not satisfied, indeed, the derivatives of $X_0$ will blow up at $s=1$, as if we use the asymptotic expansion of the error function we would find $$ X_0 \sim -\sqrt 2 (1-s)[-\log(1-s)]^{1/2}, \qquad s\to 1. $$ This is a symptom of the singular nature of the problem.
Now we scale into the inner problem near $s=1$ using: $$ X = \epsilon^{1/3}[-\log\epsilon]^{1/3} Y(t), \qquad s = 1-\epsilon^{1/3}[-\log\epsilon]^{-1/6} t. $$ While seemingly arbitrary, these scales are chosen so that a correct dominant balance is achieved in the inner problem, and matching to the outer behaviour is possible (the log terms arise because of the logarithmic behaviour of $X_0$). Under these scalings the inner problem becomes $$ Y'''' - YY'' = \ell, \qquad \ell = [-\log\epsilon]^{-1} \quad (\ll 1), \qquad Y(0) = Y''(0) = Y'''(0) = 0. $$ Expand $Y(t) = Y_0 +\ell Y_1(t) + \ldots$. To leading order, taking the boundary conditions into account, $Y_0$ is linear, and its coefficient may be found by matching to the outer problem, resulting in: $$ Y_0 = -\sqrt{\frac{2}{3}} t. $$ The problem for $Y_1$ (or more particularly, the problem for $Y_1''$) is now an inhomogeneous Airy equation: $$ Y_1'''' + \sqrt{\frac{2}{3}} t Y_1'' = 1, \qquad Y''(0) = Y'''(0) = 0. $$ (We only need $Y_1''$ for the selection problem). The solution using variation of parameters is: $$ Y_1'' = \left(-\frac{\pi}{a}\int_0^t \mathrm{Bi}(a\tau)\,\mathrm d\tau\right)\mathrm{Ai}(at) + \left(\frac{\pi}{a}\int_0^t \mathrm{Ai}(a\tau)\,\mathrm d\tau\right)\mathrm{Bi}(at), \qquad a = -(2/3)^{1/6} $$ Note because $a<0$ we are in the oscillatory regime of the Airy functions. This behaviour is crucial as it will match on to an oscillatory term in the outer problem, which will be what determines allowed values of $\epsilon$, as the oscillations must have the right frequency so that the BCs at $s=0$ are satisfied. To do this matching we need the $t\to\infty$ behaviour of the above solution, which is $$ Y_1'' \sim -\sqrt{\pi}\left(\frac{2}{3}\right)^{15/24} t^{-1/4} \cos\left[\left(\frac{2}{3}\right)^{5/4} t^{3/2} + \frac{\pi}{4}\right]. $$ The important piece of information here is the phase (term in square brackets), rather than the amplitude.
Now we need to return to the outer problem and consider the behaviour of a highly oscillatory term in the expansion by performing a WKB analysis, i.e. $$ X = X_0 + \ldots + \tilde X, \qquad \tilde X'' \sim \delta(\epsilon)U(s)\mathrm{e}^{\mathrm{i}\epsilon^{-1/2}\theta(s)} + c.c. $$ (c.c. short for complex conjugate). Here I'm thinking of $X''$ in particular to align with the inner problem, but it's the phase $\theta$ that is the important thing. Likewise, the magnitude $\delta(\epsilon)$, and the amplitude $U(s)$ are not terribly important. Substitute this expansion into the ODE and we find $$ -(\theta')^2 U + \mathrm{i}\epsilon^{1/2}[U\theta'' + 2U'\theta'] + X_0 U = O(\epsilon), $$ thus $$ \theta' = \pm\sqrt{-X_0}, \qquad U = (\theta')^{-1/2} = (-X_0)^{-1/4}. $$ Since we have $s$ in terms of $X_0$, it's easier to put the derivative in terms of $X_0$ and integrate to find $\theta$, which ultimately results in $$ \theta = \sqrt{\frac{\pi}{3}}\left(\frac{2}{\pi}\right)^{3/4} \mathrm{erf}c\left(\sqrt{\frac{3}{2}}\sqrt{\log\left(\frac{-\sqrt{2/\pi}}{X_0}\right)}\right) $$ (plus arbitrary constant representing phase shift, positive sign chosen arbitrarily). If we look at the behaviour of this $\theta$ as $s\to 1$ (equiv.~$X\to 0^-$) then $$ \theta \sim \frac{4}{3}(1-s)^{3/2} [-\log(1-s)]^{1/4}. $$ Matching with the inner problem then tells us what the phase shift must be. Ultimately, $$ \tilde X \sim \delta(s)U(s) \cos \left[\epsilon^{-1/2}\left(\theta|_{s=1} - \theta\right) + \frac{\pi}{4}\right] $$ Now to satisfy the condition at $s=0$, the phase must be an integer mulitple of $\pi$, leading to the selection condition $$ \epsilon^{-1/2} \sim \frac{(-n + 1/4)\pi}{\theta|_{s=1} - \theta|_{s=0}} = \frac{\pi^{5/4} \sqrt 3}{2^{3/4}}\left(n-\frac{1}{4}\right) \approx 4.31\left(n-\frac{1}{4}\right), \qquad n \in \mathbb N\to\infty. $$
This can then be related back to your values of $c$. It seems to match pretty well with some quick numerical checking (the first few values of $\epsilon$ I found numerically are $0.0937, 0.0175, 0.0071, 0.0038, \ldots$)
