Update: this answer appears to be flawed. See the comment from Kristoffer Arnsfelt Hansen.
I don't know how to solve your problem, but here is a technique to solve a simpler version of your problem: namely, given an edge $e$, test whether there exists any simple path from $s$ to $t$ that includes edge $e$. (This corresponds to the special case of your problem where $l=\infty$.)
You can solve this simpler problem using "max-flow with lower bounds" as a subroutine. In the standard max-flow problem, the capacity of each edge gives us an upper bound on the amount of flow going through that edge, and we require that the amount of flow on edge be lower-bounded by 0. In "max-flow with lower bounds", we are allowed to specify both a lower bound and an upper bound on the amount of flow through that edge. It is known that "max-flow with lower bounds" can be solved in polynomial time.
Now, suppose we have an edge $e \in E$, and we want to test whether there exists a simple path from $s$ to $t$ that includes edge $e$. We're going to set up a max-flow with lower bounds problem. In particular, take graph $G$ and add a new node $s_0$ with edge $s_0 \to s$ and a new node $t_1$ with edge $t \to t_1$. Make the capacity (upper-bound) on each edge 1. The lower-bound on all edges will be 0, except that the lower-bound on edge $e$ is 1. Now check whether there exists a feasible flow from $s$ to $t$ that satisfies all the bounds (this test can be done in polynomial time, as mentioned above). If there is no flow, then there is no simple path from $s$ to $t$. If there is such a flow, then tracing out that flow yields a simple path from $s$ to $t$ that includes edge $e$, so there does exist such a simple path.
How do we solve a "max-flow with lower bounds" problem? In this case, only one edge has a non-zero lower bound. Therefore, we can use a standard approach to network flow, where at each point we choose an augmenting path by computing shortest paths in the residual graph -- except that here we ask (roughly) that one of the augmenting paths includes the edge $e$.
I learned this idea from the following paper:
- Finding a Simple Path with Multiple Must-include Nodes. Hars Vardhan, Shreejith Billenahalli, Wanjun Huang, Miguel Razo, Arularasi Sivasankaran, Limin Tang,Paolo Monti, Marco Tacca, and Andrea Fumagalli. The University of Texas at Dallas, Technical Report UTD/EE/2/2009. June 2009.
(Make sure to read the techreport version, not the published version. This idea is found in the second paragraph of the introduction.)
Unfortunately, I don't know how to extend this technique to solve your original problem, with your upper-bound $l$ on the length of the simple path.
Alternatively, we could solve your problem in a straightforward way using integer linear programming (ILP). In practice, ILP solvers are pretty good on many problems. However, their worst-case running time is still exponential, so this is not going to give an algorithm with polynomial worst-case running time. Let me know if you want me to elaborate on how to formulate this using ILP.