Analytically solving of two coupled second order nonlinear ODEs

I have two coupled non linear second order differential equation.

-y1''[x] == + Exp[k1 (y1[x] - y2[x])] - Exp[-k2 (y1[x] - y2[x])],
-y2''[x] == - Exp[k1 (y1[x] - y2[x])] + Exp[-k2 (y1[x] - y2[x])]


I have been trying to solve it for past two weeks. Of course, I could get the numeric solution, but I am not able to get analytical solution. So, I was wondering if someone could help me how to solve this analytically so that I can get solution in terms of y1 and y2. Any solution to this problem would be of interest.

I am completely new to Mathematica, so any example worksheets on this sort of equation would be gratefully appreciated!

• Your code has syntax issues you have to use Exp instead of exp. Also, function arguments should be inside square brackets.[ ... ] But that apart, I don't think DSolve can solve non-linear coupled ODEs. – Lotus Sep 6 '18 at 6:14
• Analytical solving seems impossible. – Αλέξανδρος Ζεγγ Sep 6 '18 at 6:20
• However, I find that $y_1 + y_2 = a_+ x + b_+, y_1 - y_2 \sim -a_- \sin(b_- x)$, the 2nd one being a fitted approximation, where $a_\pm, b_\pm$ are constants to be determined. – Αλέξανδρος Ζεγγ Sep 6 '18 at 7:02
• @ Αλέξανδρος Ζεγγ Is your second approximation the solution of the linearized ode assuming y1-y2 to be small? If so I would expect y1-y2->0 for x-> \[Infinity] , which is only possible for b_==0! – Ulrich Neumann Sep 6 '18 at 7:46
• @UlrichNeumann No, I didn't make such an assumption. Actually, I used NDSolve with some values set for $k_1, k_2$ as well as some initial values for $y_1(0), y_1'(0), y_2(0), y_2'(0)$ and I found that then $y_1 - y_2$ looked like something as a sine. – Αλέξανδρος Ζεγγ Sep 6 '18 at 9:43

Although, as noted in comments above, no symbolic solution appears to exist for general parameters, some progress can be made. Define

eq1 = y1''[x] + Exp[k1 (y1[x] - y2[x])] - Exp[-k2 (y1[x] - y2[x])];
eq2 = y2''[x] - Exp[k1 (y1[x] - y2[x])] + Exp[-k2 (y1[x] - y2[x])];


Then, also as noted in a comment,

eq1 + eq2 == 0
(* y1''[x] + y2''[x] == 0 *)


suggesting the substitution

eqw = Simplify[eq1 /. {y2 -> Function[x, -w[x] + c1 + c2 x],
y1 -> Function[x, w[x] + c1 + c2 x]}]
(* E^(2 k1 w[x]) - E^(-2 k2 w[x]) + w''[x] *)


which DSolve can solve, up to a point, yielding an implicit solution in terms of an integral.

DSolve[eqw == 0, w, x] // First
(* Integrate[1/Sqrt[2*(-E^(2*k1*K[1])/(2*k1) - 1/(2*E^(2*k2*K[1])*k2)) + C[1]],
{K[1], 1, w[x]}]^2 == (x + C[2])^2 *)


Explicit solutions exist in at least two cases,

FullSimplify[eqw /. k2 -> k1];
Flatten@DSolve[% == 0, w[x], x];
(* {w[x] -> -((I JacobiAmplitude[I Sqrt[k1] Sqrt[-2 + k1 C[1]] (x + C[2]),
-(4/(-2 + k1 C[1]))])/k1),
w[x] -> (I JacobiAmplitude[I Sqrt[k1] x Sqrt[-2 + k1 C[1]] +
I Sqrt[k1] Sqrt[-2 + k1 C[1]] C[2], -(4/(-2 + k1 C[1]))])/k1} *)

FullSimplify[eqw /. k2 -> -k1];
Flatten@DSolve[% == 0, w[x], x]
(* {w[x] -> C[1] + x C[2]} *)