Solving another non-trivial recurrence relation

Question

This is a generalization of question Yet another non-trivial recurrence relation to solve. . Let $E \in {\mathbb R}$ and $\lambda^{C} \ge 0$, $\lambda^{M} \ge 0$ and $\Lambda \ge 0$.In addition to that take $q \in (0,1)$.Without loss of generality we can set $\lambda^{C}=1$. Now, we consider a following recurrence relation: \begin{eqnarray} E \cdot X_{E}(n) &=& X_{E}(n+1) \cdot \left( \lambda^{M} +(n+1) \cdot \lambda^{C}\right) + \\ &&X_{E}(n) \cdot \left(-\lambda^{M} 1_{n \ge 1} - n \cdot \lambda^{C} - \Lambda \right) +\\ && \Lambda \cdot \sum\limits_{i=1}^n X_{E} (n-i) \cdot q \cdot (1-q)^{i-1} \quad (i) \end{eqnarray} As considered in 1 the equation above describes a queuing system with orders conforming to three Poisson processes (cancellations, limit orders and market orders) and order sizes conforming to a negative binomial distribution with parameter $q$.

Now using the technique of Z-transforms we found a solution to the recurrence above. It reads: \begin{eqnarray} &&X_E(n) = (1-q)^n \cdot \frac{(1+\lambda^{M} + \frac{\Lambda}{1-q})^{(n)}}{(1+\lambda^{M})^{(n)}} \cdot {\mathfrak F}_{E}(n) \quad (ii) \end{eqnarray} where if $E=0$ then: \begin{equation} {\mathfrak F}_{E=0}(n) := 1_{n=0} + \frac{\Lambda}{\Lambda + (n+\lambda^{M})(1-q)} 1_{n>0} \end{equation} and if $E=-i$ where $i \in {\mathbb N}_+$ then \begin{eqnarray} {\mathfrak F}_{E=-i}(n) := 1_{n=0} + \frac{1_{n\ge 1}}{(i-1)!} \cdot \sum_{l=0}^n \frac{{\mathfrak P}^{(l,i)}(n)}{l! (1-q)^{n-l+i} (\frac{\Lambda}{1-q} + \lambda^{M} +l - i+1)^{(n-l+i)}} \end{eqnarray} where \begin{eqnarray} &&{\mathfrak P}^{(l,i)}(n) := \sum\limits_{p=l}^n \frac{p!(p-l+1)^{(i-1)} (-q)^{p-l}}{(n-p)!} \left[1_{p=n}+1_{p<n} (p-l)^{(n-p-1)} (-1+p-l-p q+n q) \right] \cdot (1+\lambda^{M}+p)^{(n-p)} \end{eqnarray} The quantities ${\mathfrak P}^{(l,i)}(n)$ are polynomials of order $(n-l)$ in the variable $\lambda^{M}$.

In[427]:= 
(*The recurrence relations (equations (7.24) (7.25) from Abergel \
"Limit Order Book"*)
(*See also \
https://math.stackexchange.com/questions/3550307/solving-another-non-\
trivial-recurrence-relation *)
lM =.; lC = 1; L =.; q =.; ii = 
 RandomInteger[{1, 10}]; EE = -ii; M = 10;
{t0, PP} = 
  Timing[Table[
    Sum[(p)!  /(n - p)! Pochhammer[p - l + 1, ii - 1] (-q)^(p - l)
       Which[p == n, 1, True, 
       Pochhammer[p - l, 
         n - p - 1] (-1 + p - l - p q + n q)] Pochhammer[1 + lM + p, 
       n - p], {p, l, n}], {n, 0, M}, {l, 0, n}]];
Print["Polynomials precomputed in t0=", t0, " secs."];


XX = Table[((1 - q)^n Pochhammer[1 + lM + L/(1 - q), n])/
    Pochhammer[1 + lM, 
     n] (Sum[PP[[1 + n, 
       1 + l]]/((l)! (1 - q)^(n - l + ii) Pochhammer[
        L/(1 - q) + lM + l - ii + 1, n - l + ii] ), {l, 0, n}] ) , {n,
     0, M}];
Table[-EE XX[[1 + n]] + XX[[2 + n]] (lM + (n + 1) lC) - 
   XX[[1 + n]] (If[n == 0, 0, lM] + n lC + L) + 
   L Sum[ XX[[1 + n - i]] q (1 - q)^(i - 1), {i, 1, n}], {n, 0, 
   Length[XX] - 2}] // FullSimplify



During evaluation of In[427]:= Polynomials precomputed in t0=0. secs.

Out[431]= {0, 0, 0, 0, 0, 0, 0, 0, 0, 0}

For $E=0$ the solution is a strictly positive waveform with one mode only whereas when $E= - i \cdot \lambda^{C}$ where $i \in {\mathbb N}_+$ then we have : \begin{eqnarray} \sum\limits_{n=0}^\infty X_{E=-i \cdot \lambda^{C}} (n) \simeq \delta_{E,0} \quad (iii) \end{eqnarray}

Below I plot the solution in question for $\lambda^{M} = 0.1$, $\Lambda=7.1$, $q=0.75$ and $i=0$ (Left figure) and then $i=1,2,3$ (Right figure in Blue, Purple and Cyan respectively).

In[456]:= (*This is plotting the stationary state (EE=0) along with \
the transient modes (EE=-1,-2,-3,..)*)
SetOptions[ListPlot, ImageSize -> 800, 
  LabelStyle -> {15, FontFamily -> "Arial"}, 
  BaseStyle -> {15, FontFamily -> "Bold"}];
{lM, L, q} = {1/10, 71/10, 3/4}; EE = -2; M = 50;
F21 = q Hypergeometric2F1[1, 1 + lM - L/(q - 1), 1 + lM, 1 - q];
PreFct = 1/F21 ;
{t0, PP} = 
  Timing[Table[
    Sum[(p)!  /(n - p)! Pochhammer[p - l + 1, ii - 1] (-q)^(p - l)
       Which[p == n, 1, True, 
       Pochhammer[p - l, 
         n - p - 1] (-1 + p - l - p q + n q)] Pochhammer[1 + lM + p, 
       n - p], {p, l, n}], {ii, 1, 3}, {n, 0, M}, {l, 0, n}]];
myModes = 
  PreFct Table[((1 - q)^n Pochhammer[1 + lM + L/(1 - q), n])/
     Pochhammer[1 + lM, 
      n] (If[ii == 0, If[n == 0, 1, L/(L + (n + lM) (1 - q))], 
       Sum[PP[[0 + ii, 1 + n, 
          1 + l]]/((l)! (1 - q)^(n - l + ii) Pochhammer[
          L/(1 - q) + lM + l - ii + 1, n - l + ii] ), {l, 0, 
         n}]] ) , {ii, 0, 3}, {n, 0, M}];
pl = ListPlot[myModes[[1]], PlotRange -> All, 
   PlotMarkers -> Automatic, Joined :> True, ImageSize -> 800];
pl1 = ListPlot[myModes[[Range[2, Length[myModes]]]], PlotRange -> All,
    PlotMarkers -> Automatic, Joined :> True, ImageSize -> 800];
pl = GraphicsGrid[{{pl, pl1}}];
(Total[#] & /@ myModes) // N
Export["myModes.jpg", pl, "JPEG"];
Import["myModes.jpg"]


Out[465]= {1., 1.26251*10^-7, -4.88867*10^-9, 3.9898*10^-10}

Now my question is twofold. Firstly, can we simplify expression $(ii)$ so that it can be evaluated faster on the computer. Secondly, I said before if $E/\lambda^{C}$ is a negative integer the waveforms, under a proper normalization, are orthogonal to each other. They are a sort of Sturm-Liouville system in discrete time. Can we actually prove that orthogonality condition rigorously ? In other words we want to prove identity $(iii)$ and also we want to find a weighting function $\left( \omega(n) \right)_{n=0}^\infty$ such that : \begin{eqnarray} \sum\limits_{n=0}^\infty X_{-i \cdot \lambda^{C}}(n) X_{-j \cdot \lambda^{C}}(n) \cdot \omega(n) = \delta_{i,j} \end{eqnarray}

1 "The Order Book as a Queueing System" in: F Abergel et al, Limit Order Books, Physics of Society: Econophysics and Sociophysics, Cambridge University Press 2016