second order ODE with Heaviside coefficient I am wondering if there is an analytical solution for such a boundary value problem,
$$
\frac{\partial ^{2}c}{\partial x^2}+\frac{\partial c}{\partial x} -(u_{c1}(x)-u_{c2}(x))c =0
$$
where $u_{c1}(x)$ and $u_{c2}(x)$ are Heaviside function with the jump points at $x=c_1$ and $x=c_2$, respectively. The boundary condition is,
$$
c(x=0)=c_0
$$
$$
c(x\rightarrow\infty)=0
$$
Any comments or suggestions are appreciated! 

Let us follow the steps of this answer. For the sake of simplicity, let us assume $c_1<c_2$ and $c_0=0$. Using the definition of the Heaviside function, we can split the equation as follows:
$$
\begin{cases}
c''+c'=0, & 0\leq x<c_1\\
c''+c'-c=0, & c_1\leq x<c_2\\
c''+c'=0 &  x\geq c_2.
\end{cases}$$
Solving the first equation and applying the first boundary condition, $c(0)=0$, we get
\begin{equation}\label{1}
c(x)=\alpha_1(1-e^{-x}),\ 0\leq x<c_1.\tag{1}
\end{equation}
In the same way, solving the $3^{rd}$ equation, and applying the second boundary condition, we have
\begin{equation}\label{2}
c(x)=\beta_3e^{-x}, \ x\geq c_2.\tag{2}
\end{equation}
Now, we know that the second equation has the following general solution:
\begin{equation}\label{3}
c(x)=\alpha_2e^{\lambda_1x}+\beta_2e^{\lambda_2x}, \ c_1<x\leq c_2,\tag{3}
\end{equation}
where $\lambda_1=-\frac{1}{2} \left(1+\sqrt{5}\right)$ and $\lambda_2=\frac{1}{2} \left(-1+\sqrt{5}\right)$. If we seek for a continuous solution, we must evaluate the constant $c_1$ into \eqref{1} and \eqref{3} and then equals the result. The same for $c_2$. In carrying out these details, we obtain the following linear system:
\begin{equation}
\begin{cases}
\alpha_2e^{\lambda_1c_1}+\beta_2e^{\lambda_2c_1}=\alpha_1(1-e^{-c_1})\\
\alpha_2e^{\lambda_1c_2}+\beta_2e^{\lambda_2c_2}=\beta_3e^{-c_2}
\end{cases}
\end{equation}
Solving this system, we will find $\alpha_2$ and $\beta_2$ as function of $\alpha_1$ and $\beta_3$, as follows:
\begin{equation}\label{4}\alpha_2=-\frac{e^{-c_1-c_2} \left(-\beta_3 e^{c_1 +\lambda_2c_1}+\alpha_1 e^{c_1+c_2 \text{$\lambda_2$}+c_2}-\text{$\alpha_1 $} e^{c_2 \text{$\lambda_2 $}+c_2}\right)}{e^{c_1 \text{$\lambda_2 $}+c_2 \text{$\lambda_1 $}}-e^{c_1 \text{$\lambda_1 $}+c_2 \text{$\lambda_2 $}}}\tag{4}\end{equation}
and
\begin{equation}\label{5}\beta_2=\frac{\text{$\alpha_1 $} \left(1-e^{-c_1}\right) e^{c_2 \text{$\lambda_1 $}}-\text{$\beta_3 $} e^{c_1 \text{$\lambda_1 $}-c_2}}{e^{c_1 \text{$\lambda_2 $}+c_2 \text{$\lambda_1 $}}-e^{c_1 \text{$\lambda_1 $}+c_2 \text{$\lambda_2 $}}}.\tag{5}\end{equation}
So, the general solution to this problem will be given by
$$c(x)=\begin{cases}
\alpha_1(1-e^{-x}),&  0\leq x<c_1,\\
\alpha_2e^{\lambda_1x}+\beta_2e^{\lambda_2x},&  \ c_1<x\leq c_2,\\
\beta_3e^{-x},&  \ x\geq c_2,
\end{cases} $$
where $\alpha_1$ and $\beta_3$ are arbitrary constants and $\alpha_2$ and $\beta_2$ are constants depending on $\alpha_2$ and $\beta_2$  given in \eqref{4} and \eqref{5}.
For instance, if we set $c_1=c_2=\alpha_1=\beta_3=1$, we will obtain the following graph:

