How can I minimize the real part of the roots of this function involving both $x$ and $e^x$ terms?

The question

I have a function $$D(s) = s^2 + c s + k + K_d s e^{-s} + K_p e^{-s}.$$ The values of $$c$$ and $$k$$ are fixed, but I can choose $$K_d$$ and $$K_p$$. How do I choose these two values in order to minimize $$\max \{ \operatorname{Re}(s) : D(s) = 0 \}$$?

Motivation

I'm writing an autopilot program for a flight simulation video game. Given information about the current state of the aircraft, it attempts to calculate flight control inputs which will put it into a desired state.

I'm using PID controllers to accomplish this. However, it's difficult to make these work in practice. The main reason is that aircraft are inherently oscillatory in several ways, and poorly chosen PID parameters exacerbate these oscillations. Attempting to find PID parameters which control the aircraft effectively without producing oscillations is very tedious.

In order to try to calculate parameters automatically, I've decided to examine a "toy" control problem whose behavior is similar to the aircraft in the video game.

Problem

In the "toy" problem, the system is a damped harmonic oscillator. The oscillator is driven by a PD controller which attempts to drive the position of the oscillator to $$0$$. However, the controller suffers a delay of $$1$$ second.

The differential equation describing this system is

$$f''(t) = - c f'(t) - k f(t) - K_d f'(t - 1) - K_p f(t - 1).$$

Here, the $$c$$ term represents the damping force and the $$k$$ term represents the restoring force. The $$K_d$$ and $$K_p$$ terms form the driving force; the $$K_d$$ term attempts to slow the motion of the system, and the $$K_p$$ term attempts to push the position of the system towards $$0$$. The constants $$c$$ and $$k$$ cannot be changed, but we are able to select $$K_d$$ and $$K_p$$.

We can find the Laplace transform of $$f$$:

$$s^2 F(s) - s f(0) - f'(0) = -c (s F(s) - f(0)) - k F(s) - K_d (s e^{-s} F(s) - f(-1)) - K_p e^{-s} F(s)$$ $$s^2 F(s) + c s F(s) + k F(s) + K_d s e^{-s} F(s) + K_p e^{-s} F(s) = s f(0) + f'(0) + c f(0) + K_d f(-1)$$ $$F(s) = \frac{s f(0) + f'(0) + c f(0) + K_d f(-1)}{s^2 + c s + k + K_d s e^{-s} + K_p e^{-s}}$$

If I understand the Laplace transform correctly, the system converges whenever all of the poles of $$F(s)$$ have a negative real part; and it suffers from divergent oscillations whenever at least one of the poles of $$F(s)$$ has at least one positive real part.

So, the behavior is determined by the rightmost root of $$s^2 + c s + k + K_d s e^{-s} + K_p e^{-s}$$. If the real part of this root is negative, then the system will converge. Furthermore, the closer the real part is to negative infinity, the more quickly the system will converge. So, we want to choose $$K_d$$ and $$K_p$$ so as to make the real part of the rightmost root as small as possible.

Hence, the question at the top of this post.

My thoughts

The equation $$s^2 + c s + k + K_d s e^{-s} + K_p e^{-s} = 0$$ doesn't look like it admits an elementary solution. I could probably find its roots using some type of numerical search; is this the best way?

Even if I had a fast way to calculate the solutions to this equation, I'd then have to perform another search in order to find the one which minimizes the maximum real part. If I had to perform nested searches, then the whole process could get very slow.

Based on playing around with the function in graphing calculators, it looks like it usually has three roots near the origin (not necessarily distinct). Does this function always have exactly three roots near the origin when $$K_d$$ and $$K_p$$ are not both zero?

If a numerical search is the best way to go for both parts of the problem (locating the roots and minimizing them), maybe the best approach is going to be to use gradient descent in the outer loop to minimize the roots, and Newton's method in the inner loop to locate the roots.

$$\exp(-s)=\dfrac{\exp(-s/2)}{\exp(s/2)}\approx \dfrac{1-s/2+1/2!(s/2)^2-1/3!(s/2)^3+ ... + 1/m!(-s/2)^m}{1+s/2+1/2!(s/2)^2+1/3!(s/2)^3+ ...+1/m!(s/2)^m}.$$
You can gradually increase $$m$$ and try to find conditions for the polynomials.