#!/usr/bin/python

import os
import sys
import time

class interactive_progressbar():
    def __init__(self, *stages):
        '''
        '''
        self.stage_counter = 0
        self.n_stages = len(stages)
        self.stage_title = stages
        self.i = 0
        self.n = 0
        self.last_check = 0.0
        self.p_past = (0, 0.0)
        self.p_dpast = (0, 0.0)
        self.tf_a = (0, 0.0)
        self.tf_b = (0, 0.0)
        self.tf_c = (0, 0.0)
        self.start_time = 0.0
        self.width = int(os.getenv('COLUMNS', '80')) - 1
        self.prev_length = 0
        # Print initial progress bar.
        self.next_stage()
    
    def next_stage(self):
        '''
        '''
        self.last_check = 0.0
        self.p_past = (0, 0.0)
        self.p_dpast = (0, 0.0)
        self.tf_a = (0, 0.0)
        self.tf_b = (0, 0.0)
        self.tf_c = (0, 0.0)
        self.start_time = time.time()
        self.stage_counter += 1
        sys.stderr.write('\b' * self.prev_length)
        if self.stage_counter > self.n_stages:
            line = '{}/{}: {}'.format(
                self.n_stages, self.n_stages,
                self.stage_title[self.n_stages - 1]
            )
        else:
            line = '{}/{}: {}'.format(
                self.stage_counter, self.n_stages,
                self.stage_title[self.stage_counter - 1]
            )
        self.prev_length = len(line)
        sys.stderr.write(line)
        sys.stderr.flush()
        if self.stage_counter > self.n_stages:
            blank = self.width - len(line) - len(' [ Done! ]')
            sys.stderr.write(' [' + ' '*blank + ' Done! ]\n')
        
    
    def update(self, i, n):
        if time.time() - self.last_check > 1:
            self.last_check = time.time()
            # Calculate time left
            # uses a pessimistic approximation for O(n)=n^2
            time_past = time.time() - self.start_time
            d_time_past = (time_past - self.p_past[1]) / (i - self.p_past[0])
            dd_time_past = (d_time_past - self.p_dpast[1]) / (i - self.p_dpast[0])
            self.p_past = (i, time_past)
            self.p_dpast = (i, d_time_past)
            # t = f(i)
            # f(x) = a*x^3 + b*x^2  + c*x
            # f'(x) = 3*a*x^2 + 2*b*x + c
            # f''(x) = 6*a*x + 2*b
            a = (dd_time_past*i**2 - 2*d_time_past*i + 2*time_past) / (2*i**3)
            b = -(dd_time_past*i**2 - 3*d_time_past*i + 3*time_past) / (i**2)
            c = (dd_time_past*i**2 - 4*d_time_past*i + 9*time_past) / (2*i)
            if a < 0:
                a = 0
                b = d_time_past
                c = time_past - d_time_past*i
            if b < 0:
                b = 0
                c = time_past/i*n
            a = 0
            b = dd_time_past
            c = d_time_past - dd_time_past*i
            if i > 10:
                self.tf_a = self.tf_a[0] + 1, (self.tf_a[0]*self.tf_a[1] + a) / (self.tf_a[0] + 1)
                self.tf_b = self.tf_b[0] + 1, (self.tf_b[0]*self.tf_b[1] + b) / (self.tf_b[0] + 1)
                self.tf_c = self.tf_c[0] + 1, (self.tf_c[0]*self.tf_c[1] + c) / (self.tf_c[0] + 1)
            a, b, c = self.tf_a[1], self.tf_b[1], self.tf_c[1]
            time_left = int(round(a*n**3 + b*n**2 + c*n - time_past))
            progress = time_past / (a*n**3 + b*n**2 + c*n + 1)
            
            # Format time left
            time_left_t = time.gmtime(time_left)
            if time_left < 600:
                time_left_s = time.strftime('%M:%S', time_left_t)[1:]
            elif time_left < 3600:
                time_left_s = time.strftime('%M:%S', time_left_t)
            elif time_left < 86400:
                time_left_s = time.strftime('%H:%M:%S', time_left_t)
            else:
                time_left_s = str(time_left//86400) + '+'
                time_left_s += time.strftime('%H:%M', time_left_t)
            # Format progress bar
            stage_s = '{}/{}: {}'.format(
                self.stage_counter, self.n_stages,
                self.stage_title[self.stage_counter - 1]
            )
            crud = len(stage_s) + len(' [] ') + len(time_left_s)
            bar = '=' * int(round(progress * (self.width - crud)))
            blank = ' ' * int(round((1.0 - progress) * (self.width - crud)))
            # Print bar
            sys.stderr.write('\b' * self.prev_length)
            line = '{} [{}{}] {}'.format(stage_s, bar, blank, time_left_s)
            self.prev_length = len(line)
            sys.stderr.write(line)
            sys.stderr.flush()
            
           

class dummy_progressbar():
    def __init__(self, *args, **kwargs):
        pass
    def next_stage(self):
        pass
    def update(self, i, n):
        pass

def gcd(a, b):
    while b:
        a, b = b, a % b
    return a

def pascals_row(n, pb=None):
    '''    
    nCr(n, k) = pascals_row(n)[k]
    
    Save pascals_row(n) and use it as an optimized nCr function when
    `n` does not vary.
    '''
    prev = curr = [1]
    i = 0
    while i < n:
        i += 1
        curr = [1]
        j = 1
        while j < i:
            curr.append(prev[j-1] + prev[j])
            j += 1
        curr.append(1)
        prev = curr
        if pb is not None:
            pb.update(i + 1, n)
    return curr

def extend(numerator_a, denominator_a, numerator_b, denominator_b):
    '''
    na, nb, d = extend(na, da, nb, db)
    
    Eg.
    (3/4) + (2/3) = (9-8)/12
    (3/4) - (2/3) = (9-8)/12
    extend(3, 4, 2, 3) -> (9, 8, 12)
    '''
    # na1 = na0 * lcm(da0, db0)/da0
    # nb1 = nb0 * lcm(da0, db0)/db0
    # d1 = lcm(da0, db0)
    # lcm(x, y) = (x*y)/gcd(x, y)
    gcd_cache = gcd(denominator_a, denominator_b)
    new_numerator_a = numerator_a * (denominator_b / gcd_cache)
    new_numerator_b = numerator_b * (denominator_a / gcd_cache)
    new_denominator = denominator_a * (denominator_b / gcd_cache)
    # Dividing with gcd_cache before multiplying is a significant optimization.
    # 4 times faster for test_precission(3000, 2000, 100)
    return new_numerator_a, new_numerator_b, new_denominator
    
def integral(n, k, x, interactive):
    if interactive:
        pb = interactive_progressbar(
            "Pascal's triangle", "Numerator", "Denominator", "Finalization"
        )
    else:
        pb = dummy_progressbar()
    
    a, b, p_d = extend(k, n, 1, 2*x)
    p_low_n = max(0, a - b)
    if a + b > 1:
        p_high_n = p_d
    else:
        p_high_n = a + b
    ncr_row = pascals_row(n - k, pb)
    pb.next_stage()
    # Calculate the main numerator.
    # sum(i=0, n-k, nCr(n-k, i) * (-1)^i * (high^(n-i+1)-low^(n-i+1))/(n-i+1))
    numerator_sum_n, numerator_sum_d = 0, 1
    p_low_exp_n = p_low_n**(k+1)
    p_high_exp_n = p_high_n**(k+1)
    p_exp_d =  p_d**(k+1)
    i = n - k
    while i >= 0:
        # Calculate i'th term in the numerator sum.
        temp_n = ncr_row[i] * (p_high_exp_n - p_low_exp_n)
        temp_d = (n - i + 1) * p_exp_d
        # Add to main numerator
        numerator_sum_n, diff, numerator_sum_d = extend(
            numerator_sum_n, numerator_sum_d,
            temp_n, temp_d
        )
        if i % 2: # (-1)^i
            numerator_sum_n -= diff
        else:
            numerator_sum_n += diff
        # Next p_low^(n-i+1) nad p_high^(n-i+1)
        p_low_exp_n *= p_low_n
        p_high_exp_n *= p_high_n
        p_exp_d *= p_d
        # Loop
        pb.update(n-k-i+1, n-k+1)
        i -= 1
    pb.next_stage()
    # Calculate the main denominator.
    # sum(i=0, n-k, nCr(n-k, i) * (-1)^i / (n-i+1))
    denominator_sum_n, denominator_sum_d = 0, 1
    i = 0
    while i <= n - k:
        temp_n = ncr_row[i]
        temp_d = n - i + 1
        # Add to main denominator
        denominator_sum_n, diff, denominator_sum_d = extend(
            denominator_sum_n, denominator_sum_d,
            temp_n, temp_d
        )
        if i % 2: # (-1)^i
            denominator_sum_n -= diff
        else:
            denominator_sum_n += diff
        # Loop
        pb.update(i+1, n-k+1)
        i += 1
    pb.next_stage()
    # Compute integral and target.
    integral_n = numerator_sum_n * denominator_sum_d
    integral_d = denominator_sum_n * numerator_sum_d
    # Compare
    # NOTE:
    # a/b > c/d  | *b *d
    # ad > bc is wrong if d xor b is negative
    # If d xor b is negative, ad < bc
    if integral_d < 0:
        integral_n, integral_d = -integral_n, -integral_d
    # target_d is x and x is not negative, no test needed
    pb.next_stage()
    return integral_n, integral_d

def test_precission(n, k, x, interactive=True):
    integral_n, integral_d = integral(n, k, x, interactive)
    target_n = x - 1
    target_d = x
    return (integral_n * target_d) > (target_n * integral_d)