I am trying to solve this PDE:
$$ \frac{\partial c}{\partial t} = D \frac{\partial^2 c}{\partial x^2} $$
subject to the initial and boundary conditions for $c(t,x)$ $$ c(0,x) = c_{\infty} \\ c(t,\infty) = c_{\infty} \\ \frac{\partial c}{\partial x} = K \left[u(t) - u(t - t_i) \right] $$
where $u(t)$ is the unit step function: $$ u(t) = \begin{cases} 1, & \text{for $t \ge 0$}\\ 0, & \text{for $t < 0$} \end{cases} $$
I found this problem in a paper: Verification and Analysis of Transference Number Measurements by the Galvanostatic Polarization Method by Hooman Hafezi and John Newman.
They suggest that the solution takes the form:
$$ c(t,x) = c_{\infty} - K\sqrt{4D t}\ i\ \mathrm{erfc}{\left(\frac{x}{\sqrt{4D t}}\right)} + K u(t-t_i) \Bigg\{\sqrt{\frac{4D(t-t_i)}{\pi}} \exp{\left(\frac{-x^2}{\sqrt{4D(t-t_i)}}\right)} - x\ \mathrm{erfc}{\left(\frac{x}{\sqrt{4D(t-t_i)}}\right)}\Bigg\} $$
but I am having difficulty in reproducing that solution, and perhaps am just naive, but I suspect there is a typo because I don't think it makes sense for there to be an imaginary part to the solution ($ i \ \mathrm{erfc}$)
After a few days of working on this problem, I am pretty confident I have found a solution and I will share it in case anyone else is interested.
If we change the boundary condition at $x=0$ slightly, to be a constant flux without a step change, i.e.:
$$ \frac{\partial c}{\partial t} = D \frac{\partial^2 c}{\partial x^2} $$
$$ c(0,x) = c_{\infty} \\ c(t,\infty) = c_{\infty} \\ \frac{\partial c}{\partial x}(x=0) = K $$
The problem is relatively straight-forward to solve using Laplace transforms.
$$ \mathcal{L} \left[ \frac{\partial c}{\partial t} \right] = \mathcal{L} \left[ D\frac{\partial^2 c}{\partial x^2} \right] \\ s \Theta(s,x) - c(0, x) = D \frac{\partial^2}{\partial x^2} \Theta(s,x) \\ D \frac{\partial^2 \Theta}{\partial x^2} - s \Theta = - c_{\infty} $$
The form of the solution is then $$ \Theta (s,x) = \frac{c_{\infty}}{s} + c_1 \exp{\left( \sqrt{\frac{s}{D}} x \right)} + c_2 \exp{\left( -\sqrt{\frac{s}{D}} x \right)} $$
$c_1 = 0$ is infered from initial value theorem, or from the boundary condition at $x = \infty$
$$ \Theta (s,x) = \frac{c_{\infty}}{s} + c_2 \exp{\left( -\sqrt{\frac{s}{D}} x \right)} $$
Using the boundary condition at $x=0$, we can solve for $c_2$
$$ \mathcal{L} \left[ \frac{\partial c}{\partial x} \right] = \mathcal{L} \left[ K \right] \\ \frac{\partial \Theta}{\partial x} = \frac{K}{s} \\ -c_2 \sqrt{\frac{s}{D}} \exp{\left( -\sqrt{\frac{s}{D}} (0) \right)} = \frac{K}{s} \\ c_2 = -\frac{K \sqrt{D}}{s \sqrt{s}} $$
Re-writing $\Theta(s,x)$ with the values for the constants: $$ \Theta(s,x) = \frac{c_{\infty}}{s} - \frac{K \sqrt{D}}{s \sqrt{s}} \exp{\left( -\sqrt{\frac{s}{D}} x \right)} $$
Taking the inverse Laplace transform: $\mathcal{L}^{-1} \left[\Theta(s,x) \right]$, gives us the concentration as a function of time and position.
$$ c(t,x) = c_{\infty} - K \left[ \sqrt{\frac{4Dt}{\pi}} \exp{\left(-\frac{x^2}{4Dt} \right)} - x\ \mathrm{erfc}{ \left( \frac{x}{\sqrt{4Dt}} \right)}\right] $$
Using the boundary condition at $x=0$ posed in the question: $$ \frac{\partial c}{\partial x} = K \left[u(t) - u(t - t_i) \right] $$
Then the transformed boundary condition becomes $$ \mathcal{L} \left[\frac{\partial c}{\partial x} \right] = \mathcal{L} \big[ K \left[ u(t) - u(t - t_i)\right] \big] \\ \frac{\partial \Theta}{\partial x} = \frac{K}{s} \left[1 - e^{-t_i s} \right] \\ $$
which is nearly identical to the original transform, but there is an additional term:
$-e^{-t_i s}$
It remains that $c_1 = 0$, while $c_2$ changes slightly:
$$ c_2 = -\frac{K \sqrt{D}}{s \sqrt{s}} \left[1 - \exp{(-t_i s)} \right] $$
$$ \Theta(s,x) = \frac{c_{\infty}}{s} - \left[1 - e^{-t_i s} \right] \frac{K}{s} \sqrt{\frac{D}{s}} \exp{\left(-\sqrt{\frac{s}{D}}x\right)} $$
Using line 27 from this Laplace transform table $$ \mathcal{L}\left[u(t - c) f(t-c)\right] = e^{-cs} F(s) $$ where $$ F(s) = \mathcal{L}\left[f(t)\right] $$
the inverse Laplace transform of $\Theta(s,x)$ becomes: $$ c(t,x) = c_{\infty} - K \left\{\left[ \sqrt{\frac{4Dt}{\pi}} \exp{\left(\frac{-x^2}{4Dt}\right)} - x\ \mathrm{erfc}{\left( \frac{x}{\sqrt{4Dt}} \right)} \right] \\ - u(t-t_i) \left[ \sqrt{\frac{4D(t-t_i)}{\pi}} \exp{\left(\frac{-x^2}{4D(t-t_i)}\right)} - x\ \mathrm{erfc}{\left( \frac{x}{\sqrt{4D(t-t_i)}} \right)} \right] \right\} $$
Here is a plot of these profiles with $t < t_i$ and $t > t_i$ ($t_i = 10$, $K = -1$, $D = 1$):
What we can see is that when $t<t_i$, at $x=0$, $\Delta c$ varies linearly with respect to $\sqrt{t}$ ("The Measurement of a Complete Set of Transport Properties for a Concentrated Solid Polymer Electrolyte Solution" by Yanping Ma, et al.); and when $t>t_i$, $\Delta c$ (at $x=0$) varies linearly with respect to $\tau$, where $\tau = \frac{\sqrt{t_i}}{\sqrt{t} + \sqrt{t - t_i}}$ (as suggested in "Verification and Analysis of Transference Number Measurements" by the Galvanostatic Polarization Method by Hooman Hafezi and John Newman).
I may be naive, but it appears the form of $c(t,x)$ derived in this post is significantly different than the form derived by Hooman Hafezi and John Newman, albeit they give identical results at $x=0$.