How to model packet delay variation

Lately our lab are working on simulating packet delay variation in NS-3-DCE. My adviser was given a project to simulate a cellular network(base station to core network) throughput given delay and packet delay variation (PDV). 
------------------------------------------------------------------------------------------------------------------------------------------------------

Mathematical Modeling

We model transmission delay between node/hop as exponential random variable and processing time from each node/hop as a fixed value (t(px)) where x represent node index.

From our assumption we model total end to end delay as total transmission delay plus total processing time as shown in the equation 1 below.

Equation 1

Furthermore we denote d as the total processing time. Total transmission delay which is sum of exponential random variable can be model as Erlang random variable. We construct Erlang random varibale from Gamma Distribution with integer shape parameter (k). We use Gamma distribution because it is available in NS-3.  Gamma Distribution is characterized by shape parameter (k) and scaling variable (θ). It has mean/average/expectation equal to  kθ and variant equal to k*θ^2.

Equation 2

Based on equation 2. we can model average delay as the expectation (average) of end to end delay.

Equation 3

Next we model packet delay variation (PDV). PDV is the difference between end to end delay of two consecutive packet.

In the illustration it is t(n)-t(n-1). Then to make sure there is no reordering if we find that current packet arrives earlier than previous packet we will set current packet arrival time equal to previous packet.

So PVD is the expectation (average) of  t(n)-t(n-1).

Equation 4

Where both X(n) and X(n-1) follow Gamma distribution (Γ(k,θ)). Then we can write it as  

Equation 5

We can take the scaling parameter out from the expectation, so we get expectation of the difference between two unit scale gamma distribution multiply by scaling parameter θ.

Equation 6

We denote expectation of the difference between two unit scale gamma distribution as Pk. We calculate Pk using matlab. First, we select value of k then we generate two unit scale gamma random number with shape parameter k. then we find the absolute difference between those value then try it again several times and find the average. We calculated Pk for several value of k. This value serve as look up table.

Finally our goal is to create average delay and PDV to be the same as what is given form the user.

So for each average delay and PDV we need to set simulation parameter d and Gamma distribution parameter k and θ. To make this happen we set the value of k and given k we find average delay and θ. From Equation 6 we can find θ, and given θ value we can use Equation 3 to find d.