Count the Number of Infection Sequences
Source: leetcode
Difficulty: Hard
Tags: array, math, combinatorics

You are given an integer `n` and an array `sick` sorted in increasing order, representing positions of infected people in a line of `n` people.


At each step, **one** uninfected person **adjacent** to an infected person gets infected. This process continues until everyone is infected.


An **infection sequence** is the order in which uninfected people become infected, excluding those initially infected.


Return the number of different infection sequences possible, modulo `109+7`.


 


**Example 1:**



**Input:** n \= 5, sick \= \[0,4]


**Output:** 4


**Explanation:**


There is a total of 6 different sequences overall.


- Valid infection sequences are `[1,2,3]`, `[1,3,2]`, `[3,2,1]` and `[3,1,2]`.
- `[2,3,1]` and `[2,1,3]` are not valid infection sequences because the person at index 2 cannot be infected at the first step.



**Example 2:**



**Input:** n \= 4, sick \= \[1]


**Output:** 3


**Explanation:**


There is a total of 6 different sequences overall.


- Valid infection sequences are `[0,2,3]`, `[2,0,3]` and `[2,3,0]`.
- `[3,2,0]`, `[3,0,2]`, and `[0,3,2]` are not valid infection sequences because the infection starts at the person at index 1, then the order of infection is 2, then 3, and hence 3 cannot be infected earlier than 2\.



 


**Constraints:**


- `2 <= n <= 105`
- `1 <= sick.length <= n - 1`
- `0 <= sick[i] <= n - 1`
- `sick` is sorted in increasing order.




C++ solution:
constexpr int m = 100000;
long long fact[m+1], inv_fact[m+1], mod = 1000000007;
bool init = false;
class Solution {
    // Time complexity O(slogm); Memory complexity O(m);
public:
    int numberOfSequence(int n, vector<int>& sick) {
        // (sum(segments))!/(multiply(for each segments!)) * 2 ^ (sum(segments)-len(segments));
        if (!init) {
            init = fact[0] = 1;
            // we only need modular factorials and inverse factorials up to mx inclusive:
            for (int i=1;i<=m;i++) 
                fact[i] = i * fact[i-1] % mod;
            // FLT: a^(p-1)%p = 1
            // a * a^(p-2)%p = a^-1
            // thus:
            // a * a^(p-2) = a^(p-1) * a^-1 = 1 * a^(-1) (mod p)
            // a^-1 = a^(p-2)
            inv_fact[m] = modPow(fact[m], mod - 2, mod);
            for (int i=m-1;i>=0;i--) {
                inv_fact[i] = inv_fact[i+1] * (i+1) % mod;
            }
        }
        long long res = 1;
        res = (res * inv_fact[sick[0]]) % mod;
        res = (res * inv_fact[n-sick.back()-1]) % mod;
        for (int i = 1;i<sick.size();i++) {
            int dist = sick[i] - sick[i-1]-1;
            if (dist == 0)
                continue;
            res = (res * modPow(2, dist-1, mod)) % mod;
            res = (res * inv_fact[dist]) % mod;
        }
        return (res * fact[n-sick.size()]) % mod;
    }
private:
    int modPow(int x, int y, int mod) {
        if (y == 0)
            return 1;
        long p = modPow(x,y/2,mod)%mod;
        p = (p*p)% mod;
        return y%2 ? (p*x) % mod : p;
    }
    long long binom_coeff(int n, int k) {
        return max(1LL, fact[n] * inv_fact[k] % mod * inv_fact[n-k] % mod);
    }
};

Go solution:
const MOD = int(1e9+7)

// (k1 + k2 + k3)! / (k1! * k2! * k3!)
func multinomial(k []int) int {
    n := 0
    denominator := 1
    for i := range k {
        n += k[i]
        n %= MOD
        denominator *= factorial(k[i])
        denominator %= MOD
    }
    return factorial(n) * modInverse(denominator) % MOD
}

var factorialMemo = []int{1}
func factorial(num int) int {
    for i := len(factorialMemo); i <= num; i++ {
        factorialMemo = append(factorialMemo, (factorialMemo[i-1] * i) % MOD)
    }
    return factorialMemo[num]
}

func modInverse(a int) int {
    d, modInv, _ := gcd(a, MOD)
    if d != 1 {
        panic(fmt.Sprintf("gcd(%d, %d)=%d, must be 1", a, MOD, d))
    }
    return (modInv + MOD) % MOD
}

func gcd(a int, b int) (int, int, int) {
    if a == 0 {
        return b, 0, 1
    }
    d, x, y := gcd(b % a, a)
    return d, y - (b / a) * x, x
}

func pow(a int, p int) int {
    res := 1
    for p > 0 {
        if p % 2 == 1 {
            res *= a
            res %= MOD
        }
        a *= a
        a %= MOD
        p /= 2
    }
    return res
}

func numberOfSequence(n int, sick []int) int {
    options := 1
    lenghts := []int{sick[0], n - sick[len(sick)-1] - 1} // non-infected leftmost and rightmost lengths
    for i := 1; i < len(sick); i++ { // all lenghts of non-infected between infected
        l := sick[i] - sick[i-1] - 1
        if l > 0 {
            lenghts = append(lenghts, l)

            // neither leftmost nor rightmost
            // can be infected from both sides with 2^(len-1) options
            options *= pow(2, l-1)
            options %= MOD
        }
    }

    return options * multinomial(lenghts) % MOD
}

Java solution:
class Solution {
    private int MOD = (int) Math.pow(10, 9) + 7;

    public int numberOfSequence(int n, int[] sick) {
        List<Integer> segments = new ArrayList<>();
        int startSegmentLen = sick[0];
        int endSegmentLen = n - sick[sick.length - 1] - 1;
        // System.out.println("start, end, segments: " + startSegmentLen + "," + endSegmentLen + ",");

        for (int i = 1; i < sick.length; i++) {
            if (sick[i] - sick[i-1] - 1 > 0) {
                segments.add(sick[i] - sick[i-1] - 1);
            }
            // System.out.println(segments.get(segments.size() - 1));
        }

        int[] fac = prepareFactorial(n);
        int[] pow2 = preparePow2(n);

        int totalNotSick = n - sick.length;

        int ans = fac[totalNotSick];

        // order within segment does not matter outside of segment;
        // opening and closing segments may be unbounded so possible sequence is 1
        ans = multiply(ans, inverse(fac[startSegmentLen])); 
        ans = multiply(ans, inverse(fac[endSegmentLen]));

        for (int segment : segments) {
            ans = multiply(ans, pow2[segment - 1]); // ordering within segment to get sick
            ans = multiply(ans, inverse(fac[segment])); // order within segment does not matter outside of segment
        }

        return ans;
    }

    private int inverse(int n) {
        return exp(n, MOD - 2);
    }

    private int exp(int base, int exp) {
        int res = 1;
        while (exp > 0) {
            if (exp % 2 == 1) {
                res = multiply(res, base);
            }
            base = multiply(base, base);
            exp /= 2;
        }
        return res;
    }

    private int multiply(int a, int b) {
        return (int) ((1L * a * b) % MOD);
    }

    private int[] preparePow2(int n) {
        int[] pow2 = new int[n + 1];
        pow2[0] = 1;
        for (int i = 1; i <= n; i++) {
            pow2[i] = multiply(pow2[i - 1], 2);
        }
        return pow2;
    }

    private int[] prepareFactorial(int n) {
        int[] fac = new int[n + 1];
        fac[0] = 1;
        for (int i = 1; i <= n; i++) {
            fac[i] = multiply(fac[i - 1], i);
        }
        return fac;        
    }
}

JavaScript solution:
//#region Modulo
class Modulo {
  /**
   * @param {number} modulo
   */
  constructor(modulo) {
    /** @type {number} @readonly */
    this.modulo = modulo;
  }

  /**
   * @param  {...number} numbers
   */
  add(...numbers) {
    let result = 0;
    for (let number of numbers) {
      result = (result + (number % this.modulo)) % this.modulo;
    }

    if (result < 0) result += this.modulo;
    return result;
  }

  /**
   * @private
   * @param {number} a
   * @param {number} b
   * @returns {number}
   */
  _quickMul(a, b) {
    a = ((a % this.modulo) + this.modulo) % this.modulo;
    b = ((b % this.modulo) + this.modulo) % this.modulo;
    if (a === 0 || b === 0) return 0;
    if (Math.log2(a) + Math.log2(b) < 50) return (a * b) % this.modulo;

    let result = 0;
    while (b) {
      while (b % 2 === 0) {
        a = (a * 2) % this.modulo;
        b /= 2;
      }

      if (b % 2 !== 0) {
        result = (result + a) % this.modulo;
        b--;
      }
    }

    return result;
  }

  /**
   * @param  {...number} numbers
   */
  mul(...numbers) {
    let result = 1;
    for (let number of numbers) {
      result = this._quickMul(result, number);
      if (result === 0) return 0;
    }

    if (result < 0) result += this.modulo;
    return result;
  }

  /**
   * @param {number} a
   * @param {number} b
   * @returns {number}
   */
  pow(a, b) {
    a = ((a % this.modulo) + this.modulo) % this.modulo;
    if (a === 0) return 0;

    let result = 1;
    while (b) {
      while (b % 2 === 0) {
        a = this._quickMul(a, a);
        b /= 2;
      }

      if (b % 2 !== 0) {
        result = this._quickMul(result, a);
        b--;
      }
    }

    return result;
  }
}
//#endregion

//#region Division
/** @private @type {undefined | number} */
Modulo.prototype._phi = undefined;

/**
 * @returns {number}
 */
Modulo.prototype.getPhi = function () {
  if (this._phi !== undefined) return this._phi;

  let result = this.modulo;
  let temp = this.modulo;
  for (let i = 2; i <= Math.sqrt(temp); i++) {
    if (temp % i === 0) {
      result /= i;
      result *= i - 1;
    }
    while (temp % i === 0) temp /= i;
  }
  if (temp > 1) {
    result /= temp;
    result *= temp - 1;
  }

  this._phi = result;
  return result;
};

/**
 * @param {number} a
 * @returns {number}
 */
Modulo.prototype.getInverse = function (a) {
  return this.pow(a, this.getPhi() - 1);
};

/**
 * @param {number} a
 * @param {number} b
 * @returns {number}
 */
Modulo.prototype.div = function (a, b) {
  return this._quickMul(a, this.getInverse(b));
};
//#endregion

let mod = new Modulo(1000000007);

let FACTS = new Uint32Array(100001);
FACTS[0] = 1;
for (let i = 1; i <= 1e5; ++i) FACTS[i] = mod.mul(FACTS[i - 1], i);

/**
 * @param {number} n
 * @param {number[]} sick
 * @return {number}
 */
var numberOfSequence = function (n, sick) {
  let m = sick.length;
  let e = FACTS[n - m];
  let d = mod.mul(FACTS[sick[0]], FACTS[n - 1 - sick[m - 1]]);

  for (let i = 1; i < m; ++i) {
    let dis = sick[i] - sick[i - 1] - 1;
    if (dis <= 1) continue;
    e = mod.mul(e, mod.pow(2, dis - 1));
    d = mod.mul(d, FACTS[dis]);
  }

  return mod.div(e, d);
};

Python solution:
MOD = 10 ** 9 + 7
m2 = {}
m1 = {}

def fac(n):
    if n not in m1: 
        r = 1 if n < 2 else n * fac(n - 1) % MOD
        m1[n] = r
    return m1[n]

def ifac(n):
    if n not in m2:
        r = pow(fac(n), -1, MOD)
        m2[n] = r
    return m2[n]

class Solution:
    def numberOfSequence(self, n, s):
        r = fac(n - len(s)) * ifac(s[0]) * ifac(n - s[-1] - 1) % MOD
        for g in [b - a - 1 for a, b in pairwise(s) if b - a > 1]:
            r *= ifac(g) * pow(2, g - 1, MOD)
            r %= MOD
        return r

Ruby solution:
def number_of_sequence(n, sick)
  mod = 10**9 + 7
  unless defined?(@fact)
    max_m = 10**5
    @fact = Array.new(max_m + 1, 1)
    (1..max_m).each { |i| @fact[i] = @fact[i-1] * i % mod }
    @inv_fact = Array.new(max_m + 1, 1)
    @inv_fact[max_m] = @fact[max_m].pow(mod - 2, mod)
    (max_m - 1).downto(0) { |i| @inv_fact[i] = @inv_fact[i + 1] * (i + 1) % mod }
  end

  total = n - sick.size
  return 0 if total < 0

  left_gap = sick[0]
  right_gap = n - sick[-1] - 1

  res = 1
  res = res * @inv_fact[left_gap] % mod
  res = res * @inv_fact[right_gap] % mod

  (1...sick.size).each do |i|
    a = sick[i-1]
    b = sick[i]
    gap = b - a - 1
    next if gap <= 0

    res = res * 2.pow(gap - 1, mod) % mod
    res = res * @inv_fact[gap] % mod
  end

  res = res * @fact[total] % mod
  res
end