-rw-r--r-- 3102 cryptattacktester-20231020/queue_prob.cpp raw
#include <cassert>
#include <vector>
#include "queue_prob.h"

using namespace std;

// generator runs G times
// each time producing an event with probability p
// which is then inserted into a queue
// output: chances of 0, 1, 2, ..., Q-1 entries in the queue
// i.e., coefficients of x^0,...,x^(Q-1) in (1-p+px)^G
vector<bigfloat> queue_distribution(bigint G,bigfloat p,bigint Q)
{
  assert(Q >= 0);
  assert(G >= 0);

  vector<bigfloat> phi(Q,0);
  if (Q >= 1) phi.at(0) = 1;

  if (G == 0) return phi;

  bigfloat p1m = 1-p;

  // XXX: can be faster to use the bigfloat version below

  vector<bigfloat> base(Q,0);
  if (Q >= 1) base.at(0) = p1m;
  if (Q >= 2) base.at(1) = p;

  for (bigint pos = nbits(G)-1;pos >= 0;--pos) {
    vector<bigfloat> newphi(Q);
    for (bigint e = 0;e < Q;++e) newphi.at(e) = 0;
    for (bigint e = 0;e < Q;++e)
      for (bigint d = 0;d+e < Q;++d)
        newphi.at(d+e) += phi.at(d)*phi.at(e);

    if (G.bit(pos)) {
      for (bigint e = 0;e < Q;++e) phi.at(e) = 0;
      for (bigint e = 0;e < Q;++e)
        for (bigint d = 0;d+e < Q;++d)
          phi.at(d+e) += newphi.at(d)*base.at(e);
    } else {
      for (bigint e = 0;e < Q;++e) phi.at(e) = newphi.at(e);
    }
  }

  return phi;
}

vector<bigfloat> queue_distribution(bigfloat G,bigfloat p,bigint Q)
{
  assert(Q >= 0);

  vector<bigfloat> phi(Q,0);
  if (Q >= 1) phi.at(0) = 1;

  if (G.definitely_zero()) return phi;

  bigfloat p1m = 1-p;

  vector<bigfloat> ppow;
  ppow.push_back(1);
  for (bigint e = 0;e < Q;++e) ppow.push_back(ppow.at(e)*p);

  vector<bigfloat> p1mpow;
  p1mpow.push_back(1);
  for (bigint e = 0;e < Q;++e) p1mpow.push_back(p1mpow.at(e)*p1m);
  bigfloat p1mpowGQ1 = pow(p1m,G-(Q-1));

  for (bigint e = 0;e < Q;++e) {
    // phi_e = binomial(G,e) p^e (1-p)^(G-e)
    phi.at(e) = binomial_bigfloat(G,e)*ppow.at(e)*p1mpow.at(Q-1-e)*p1mpowGQ1;
  }

  return phi;
}

// generator runs G times
// each time producing an event with probability p
// which is then inserted into a queue
// but queue is limited to length Q
// output: average number of events in queue
bigfloat queue_average(bigint G,bigfloat p,bigint Q)
{
  if (G <= 0) return 0;
  if (G <= Q) return p*bigfloat(G);

  // define phi = (1-p+px)^G
  // then phi_e = chance of exactly e events being produced
  // so on average queue keeps X events
  // where X = sum_{e} min{e,Q} phi_e
  //
  // have sum_e Q phi_e = Q phi(1) = Q
  // so sum_{e} (Q-min{e,Q}) phi_e = Q-X
  // and this sum has the advantage of being supported on e<Q

  auto phi = queue_distribution(G,p,Q);

  bigfloat result = 0;
  for (bigint e = 0;e < Q;++e)
    result += (Q-e)*phi.at(e);

  result = Q-result;
  result = bigfloat_guarantee_nonnegative(result);
  return result;
}

bigfloat queue_average(bigfloat G,bigfloat p,bigint Q)
{
  if (!G.definitely_positive()) return 0;
  if (!bigfloat(Q).definitely_less(G)) return p*G;

  auto phi = queue_distribution(G,p,Q);

  bigfloat result = 0;
  for (bigint e = 0;e < Q;++e)
    result += (Q-e)*phi.at(e);

  result = Q-result;
  result = bigfloat_guarantee_nonnegative(result);
  return result;
}