HackerRank P-sequences Problem Solution

#include <assert.h>
#include <limits.h>
#include <math.h>
#include <stdbool.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#define MOD 1000000007

char* readline();

/*
 * Complete the pSequences function below.
 */
int pSequences(int n, long p) {
    int sqrt = 1;
    while(sqrt*sqrt <= p){
        sqrt++;
    }
    sqrt--;
    long interlen[2*sqrt];
    for(int i = 0; i < sqrt; i++){
        interlen[i] = 1;
        interlen[i + sqrt] = p/(sqrt - i) - p/(sqrt - i + 1);
    }
    interlen[sqrt] = p/sqrt - sqrt;

    long currnum[2*sqrt];
    for(int i = 0; i < 2*sqrt; i++){
        currnum[i] = 0;
    }
    currnum[0] = 1;

    for(int i = 0; i < n + 1; i++){
        long total = 0;
        long nextnum[2*sqrt];
        for(int j = 2*sqrt - 1; j >= 0; j--){
            total = (total + currnum[2*sqrt - j - 1])%MOD;
            nextnum[j] = total;
        }
        for(int j = 0; j < 2*sqrt; j++){
            currnum[j] = (nextnum[j]*interlen[j])%MOD;
        }
    }
    return currnum[0];
}

int main()
{
    FILE* fptr = fopen(getenv("OUTPUT_PATH"), "w");

    char* n_endptr;
    char* n_str = readline();
    int n = strtol(n_str, &n_endptr, 10);

    if (n_endptr == n_str || *n_endptr != '\0') { exit(EXIT_FAILURE); }

    char* p_endptr;
    char* p_str = readline();
    int p = strtol(p_str, &p_endptr, 10);

    if (p_endptr == p_str || *p_endptr != '\0') { exit(EXIT_FAILURE); }

    int result = pSequences(n, p);

    fprintf(fptr, "%d\n", result);

    fclose(fptr);

    return 0;
}

char* readline() {
    size_t alloc_length = 1024;
    size_t data_length = 0;
    char* data = malloc(alloc_length);

    while (true) {
        char* cursor = data + data_length;
        char* line = fgets(cursor, alloc_length - data_length, stdin);

        if (!line) { break; }

        data_length += strlen(cursor);

        if (data_length < alloc_length - 1 || data[data_length - 1] == '\n') { break; }

        size_t new_length = alloc_length << 1;
        data = realloc(data, new_length);

        if (!data) { break; }

        alloc_length = new_length;
    }

    if (data[data_length - 1] == '\n') {
        data[data_length - 1] = '\0';
    }
    if(data[data_length - 1] != '\0'){
        data_length++;
        data = realloc(data, data_length);
        data[data_length - 1] = '\0';
    }

    data = realloc(data, data_length);

    return data;
}

// Simon St James (ssjgz) 2019-05-11.
#define SUBMISSION
#define BRUTE_FORCE
#ifdef SUBMISSION
#define NDEBUG
#undef BRUTE_FORCE
#endif
#include <iostream>
#include <vector>
#include <algorithm>
#include <cassert>
#include <cmath>

#include <sys/time.h>

using namespace std;


const int64_t modulus = 1'000'000'007ULL;

class ModNum
{
    public:
        ModNum(int64_t n = 0)
            : m_n{n}
        {
        }
        ModNum& operator+=(const ModNum& other)
        {
            m_n = (m_n + other.m_n) % ::modulus;
            return *this;
        }
        ModNum& operator*=(const ModNum& other)
        {
            m_n = (m_n * other.m_n) % ::modulus;
            return *this;
        }
        int64_t value() const { return m_n; };
    private:
        int64_t m_n;
};

ModNum operator+(const ModNum& lhs, const ModNum& rhs)
{
    ModNum result(lhs);
    result += rhs;
    return result;
}

ModNum operator*(const ModNum& lhs, const ModNum& rhs)
{
    ModNum result(lhs);
    result *= rhs;
    return result;
}

ostream& operator<<(ostream& os, const ModNum& toPrint)
{
    os << toPrint.value();
    return os;
}

bool operator==(const ModNum& lhs, const ModNum& rhs)
{
    return lhs.value() == rhs.value();
}


ModNum solutionOptimised(int N, int P)
{
#ifdef BRUTE_FORCE
    vector<vector<ModNum>> firstNEndingOnP(N, vector<ModNum>(P + 1, 0));
    for (int r = 0; r <= P; r++)
    {
        firstNEndingOnP[0][r] = 1;
    }
    for (int i = 1; i < N; i++)
    {
        for (int q = 1; q <= P; q++)
        {
            for (int r = 1; r * q <= P; r++)
            {
                firstNEndingOnP[i][q] += firstNEndingOnP[i - 1][r];
            }
            cout << "firstNEndingOnP[" << i << "][" << q << "] = " << firstNEndingOnP[i][q] << endl;
        }
    }
#endif

    vector<int> factorsOfP;
    for (int i = 1; i <= sqrt(P); i++)
    {
        //if ((P % i) == 0)
        {
            factorsOfP.push_back(i);
        }
    }
    for (int i = sqrt(P) + 1; i > 1; i--)
    {
        //if ((P % i) == 0)
        {
            factorsOfP.push_back(P / i + 1);
        }
    }
    //factorsOfP.push_back(P + 1);
    factorsOfP.erase(std::unique(factorsOfP.begin(), factorsOfP.end()), factorsOfP.end());
    assert(std::is_sorted(factorsOfP.begin(), factorsOfP.end()));
    for (const auto x : factorsOfP)
    {
        //cout << " factor of P: " << x << endl;
    }
    vector<ModNum> firstNEndingOnFactorIndex(factorsOfP.size() + 1, 0);
    for (int i = 0; i < factorsOfP.size(); i++)
    {
        firstNEndingOnFactorIndex[i] = 1;
    }
    for (int i = 1; i < N; i++)
    {
        vector<ModNum> nextEndingOnFactorIndex(factorsOfP.size() + 1, 0);
        //cout << "i: " << i << endl;
        int lastFactorIndex = 0;
        ModNum sumUpToLast = firstNEndingOnFactorIndex[lastFactorIndex];
        int summedSoFar = factorsOfP[lastFactorIndex];
        int newFactorIndex = factorsOfP.size();
        nextEndingOnFactorIndex[newFactorIndex] = firstNEndingOnFactorIndex[0];
        //sumUpToLast += firstNEndingOnFactorIndex[i - 1][lastFactorIndex];

        while (newFactorIndex >= 1)
        {
            //cout << " blah: newFactorIndex: " << newFactorIndex << endl;
            //int lastFactorIndex = 0;
            //ModNum sumUpToLast = firstNEndingOnFactorIndex[i - 1][0];
            newFactorIndex--;
            const auto diffFromPreviousLastFactor = factorsOfP[lastFactorIndex] - factorsOfP[lastFactorIndex - 1];
            //const auto diffToNextLastFactor = (lastFactorIndex == factorsOfP.size() - 1) ? 1 : factorsOfP[lastFactorIndex + 1] - factorsOfP[lastFactorIndex];
            //sumUpToLast += firstNEndingOnFactorIndex[i - 1][lastFactorIndex] * diffFromPreviousLastFactor;
            //sumUpToLast += firstNEndingOnFactorIndex[i - 1][lastFactorIndex] * (diffFromPreviousLastFactor - 0);
            const auto needSumUpTo = P / factorsOfP[newFactorIndex];
            while (lastFactorIndex + 1 < factorsOfP.size() && factorsOfP[lastFactorIndex + 1] <= needSumUpTo)
            {
                assert(lastFactorIndex + 1 < factorsOfP.size());
                lastFactorIndex++;
                const auto globble = (factorsOfP[lastFactorIndex] - summedSoFar - 1) * firstNEndingOnFactorIndex[lastFactorIndex - 1]
                     + firstNEndingOnFactorIndex[lastFactorIndex];
                //cout << " in loop: adding " << globble << " to sumUpToLast" << endl; 
                sumUpToLast += globble;
                summedSoFar = factorsOfP[lastFactorIndex];
            }
            const auto globble = (needSumUpTo - factorsOfP[lastFactorIndex]) * firstNEndingOnFactorIndex[lastFactorIndex];
            //cout << " after loop: adding " << globble << " to sumUpToLast.  " << " needSumUpTo: " << needSumUpTo << " lastFactorIndex: " << lastFactorIndex << " firstNEndingOnFactorIndex[i - 1][lastFactorIndex]: " << firstNEndingOnFactorIndex[lastFactorIndex] << endl;
            sumUpToLast += globble;
            summedSoFar = needSumUpTo;

#ifdef BRUTE_FORCE
            ModNum debugSumUpToLast;
            for (int k = 1; k <= needSumUpTo; k++)
            {
                cout << "Adding " << firstNEndingOnP[i - 1][k] << " to debugSumUpToLast" << endl;
                debugSumUpToLast += firstNEndingOnP[i - 1][k];
            }
            cout << "sumUpToLast: " << sumUpToLast << " debugSumUpToLast: " << debugSumUpToLast << endl;
            assert(sumUpToLast == debugSumUpToLast);
#endif
            //sumUpToLast += firstNEndingOnFactorIndex[i - 1][lastFactorIndex] + (diffFromPreviousLastFactor - 1) * firstNEndingOnFactorIndex[i - 1][lastFactorIndex - 1];

            const auto diffUntilNextFactor = factorsOfP[newFactorIndex + 1] - factorsOfP[newFactorIndex];
            //cout << "lastFactorIndex: " << lastFactorIndex << " newFactorIndex: " << newFactorIndex << " factorsOfP[lastFactorIndex] * factorsOfP[newFactorIndex] : " << factorsOfP[lastFactorIndex] * factorsOfP[newFactorIndex] << " diffUntilNextFactor: " << diffUntilNextFactor  << " diffFromPreviousLastFactor: " << diffFromPreviousLastFactor << " sumUpToLast: "<< sumUpToLast << endl;
            assert(factorsOfP[lastFactorIndex] * factorsOfP[newFactorIndex] <= P);
            //assert(lastFactorIndex == factorsOfP.size() - 1 || factorsOfP[lastFactorIndex + 1] * factorsOfP[newFactorIndex] > P);
            //assert((factorsOfP[lastFactorIndex] + 1) * factorsOfP[newFactorIndex] > P);
#if 0
            if (lastFactorIndex != factorsOfP.size() - 1)
            {
                for (int dbg = factorsOfP[lastFactorIndex]; dbg < factorsOfP[lastFactorIndex + 1]; dbg++)
                {
                    cout << " factorsOfP[lastFactorIndex]: " << factorsOfP[lastFactorIndex] << " dbg: " << dbg << " dbg * factorsOfP[newFactorIndex]: " << dbg * factorsOfP[newFactorIndex] << endl;
                    assert(dbg * factorsOfP[newFactorIndex] <= P);
                }
            }
#endif
            //for (int dbg = factorsOfP[lastFactorIndex]; dbg < ; dbg++)
            //{
                //assert(
            //}


            //firstNEndingOnFactorIndex[i][newFactorIndex] = sumUpToLast * diffUntilNextFactor;
            nextEndingOnFactorIndex[newFactorIndex] = sumUpToLast * 1;
#ifdef BRUTE_FORCE
            cout << "firstNEndingOnFactorIndex[" << i << "][" << factorsOfP[newFactorIndex] << "] = " << nextEndingOnFactorIndex[newFactorIndex] << endl;
            cout << "firstNEndingOnP[" << i << "][" << factorsOfP[newFactorIndex] << "] = " << firstNEndingOnP[i][factorsOfP[newFactorIndex]] << endl;
            assert(nextEndingOnFactorIndex[newFactorIndex] == firstNEndingOnP[i][factorsOfP[newFactorIndex]]);
#endif
        }
        firstNEndingOnFactorIndex = nextEndingOnFactorIndex;
    }

    ModNum result = 0;
    for (int r = 0; r < factorsOfP.size(); r++)
    {
        //cout << " calc result: firstNEndingOnFactorIndex[N - 1][r] : " << firstNEndingOnFactorIndex[r]  << endl;
        if (r + 1 < factorsOfP.size())
        {
            result += firstNEndingOnFactorIndex[r] * (factorsOfP[r + 1] - factorsOfP[r]);
            //result += firstNEndingOnFactorIndex[N - 1][r + 1];
        }
    }
    result += firstNEndingOnFactorIndex[factorsOfP.size() - 1] * (P - factorsOfP.back() + 1);

    return result;
}

ModNum solutionBruteForce(int N, int P)
{
#if 0
    cout << "Blah" << endl;
    int a = 1;
    int lastDivision = P / a;
    cout << " brute force a: 1 " << endl;
    vector<int> blah = {a};
    while (a < P)
    {
        a++;
        if ((P / a) != lastDivision)
        {
            lastDivision = P / a;
            cout << " brute force a: " << a << " lastDivision: " << lastDivision << endl;
            blah.push_back(a);
        }
    }
    vector<int> optimisedBlah;
    for (int i = 1; i <= sqrt(P); i++)
    {
        optimisedBlah.push_back(i);
    }
    for (int i = sqrt(P) + 1; i > 1; i--)
    {
        optimisedBlah.push_back(P / i + 1);
    }
    cout << "blah.size(): " << blah.size() << " optimisedBlah.size(): " << optimisedBlah.size() << endl;
    optimisedBlah.erase(std::unique(optimisedBlah.begin(), optimisedBlah.end()), optimisedBlah.end());
    for (const auto x : optimisedBlah)
    {
        cout << " optimised a: " << x << endl;
    }
    assert(is_sorted(optimisedBlah.begin(), optimisedBlah.end()));
    assert(optimisedBlah == blah);
    //a = 1;
    //lastDivision = P / a;
    //cout << " optimised a: 1 " << endl;

    //while (a <= P)
    //{
        //a = (P / (lastDivision - 1));
        //lastDivision = P / a;
        //cout << " optimised a: " << a << endl;
    //}
    return 0;
#endif
    ModNum result = 0;
    vector<vector<ModNum>> firstNEndingOnP(N, vector<ModNum>(P + 1, 0));
    for (int r = 0; r <= P; r++)
    {
        firstNEndingOnP[0][r] = 1;
    }
    for (int i = 1; i < N; i++)
    {
        for (int q = 1; q <= P; q++)
        {
            for (int r = 1; r * q <= P; r++)
            {
                firstNEndingOnP[i][q] += firstNEndingOnP[i - 1][r];
            }
            cout << "firstNEndingOnP[" << i << "][" << q << "] = " << firstNEndingOnP[i][q] << endl;
        }
    }

    for (int r = 1; r <= P; r++)
    {
        result += firstNEndingOnP[N - 1][r];
    }
    return result;
}


int main(int argc, char* argv[])
{
    if (argc == 2)
    {
        struct timeval time;
        gettimeofday(&time,NULL);
        srand((time.tv_sec * 1000) + (time.tv_usec / 1000));


        const int N = rand() % 100 + 1;
        const int P = rand() % 100 + 1;
        cout << N << " " << P << endl;
        return 0;
    }
    int N;
    cin >> N;

    int P;
    cin >> P;

#ifdef BRUTE_FORCE
    const auto bruteForceResult = solutionBruteForce(N, P);
    cout << "bruteForceResult: " << bruteForceResult << endl;
#endif

    const auto optimisedResult = solutionOptimised(N, P);

#ifdef BRUTE_FORCE
    cout << "optimisedResult: " << optimisedResult << endl;
    assert(optimisedResult == bruteForceResult);
#else
    cout << optimisedResult << endl;
#endif

}

using System.CodeDom.Compiler;
using System.Collections.Generic;
using System.Collections;
using System.ComponentModel;
using System.Diagnostics.CodeAnalysis;
using System.Globalization;
using System.IO;
using System.Linq;
using System.Reflection;
using System.Runtime.Serialization;
using System.Text.RegularExpressions;
using System.Text;
using System;

class Result
{

    /*
     * Complete the 'pSequences' function below.
     *
     * The function is expected to return an INTEGER.
     * The function accepts following parameters:
     *  1. INTEGER n
     *  2. INTEGER p
     */

    public static int pSequences(int n, int p)
    {
        long P = 1000000007;
        int r =((int)(Math.Sqrt(((double)(p)))));
        
        long[] t = new long[((2 *r) +3)]; 
        long[] s = new long[((2 *r) +3)]; 
        long[] u = new long[((2 *r) +3)]; 
        
        Array.Fill(t,1);
        Array.Fill(s,1);
        s[0] = 0;
        int sum = r;
        int max = r;
        int nxt = r;
        
        while(( sum < p))
        {
            t[++max] = ((p / nxt) - (p / (nxt + 1)));
            sum = ((int)(sum + t[max]));
            nxt--;
        }
            
        int d = 0;
        if(sum > p)
        {
            max--;
            d = 1;
        } 
           
        while(n-- > 0)
        {
            for(int i = max, j = 1; i > 0; i--, j++)
            {
                u[i] = (((s[j] * t[(j + d)]) + u[(i + 1)]) % P); 
            }
            
            long[] g = s;
            s = u;
            u = g;
            u[(max + 1)] = 0;
        }
        return (int)s[1];
        
    }

}

class Solution
{
    public static void Main(string[] args)
    {
        TextWriter textWriter = new StreamWriter(@System.Environment.GetEnvironmentVariable("OUTPUT_PATH"), true);

        int n = Convert.ToInt32(Console.ReadLine().Trim());

        int p = Convert.ToInt32(Console.ReadLine().Trim());

        int result = Result.pSequences(n, p);

        textWriter.WriteLine(result);

        textWriter.Flush();
        textWriter.Close();
    }
}

import java.io.*;
import java.math.*;
import java.security.*;
import java.text.*;
import java.util.*;
import java.util.concurrent.*;
import java.util.function.*;
import java.util.regex.*;
import java.util.stream.*;
import static java.util.stream.Collectors.joining;
import static java.util.stream.Collectors.toList;

class Result {

    public static long pSequences(int n, int p) {
    // Write your code here
    
        long P = 1000000007;
        int r = (int) Math.sqrt((double)p);

        long[] f = new long[2 * r + 3];
        long[] x = new long[2 * r + 3];
        long[] y = new long[2 * r + 3];

        Arrays.fill(f, 1);
        Arrays.fill(x, 1);
        x[0] = 0;

        int max = r, sum = r;
        int next = r;
        while (sum < p) {
            f[++max] = p/(next) - p/(next+1);
            sum += f[max];
            next--;
        }

        int diff = 0;
        if (sum > p) {
            max--;
            diff = 1;
        }

        while(n-- > 0) {
            for (int i = max, j = 1; i > 0; i--, j++) {
                y[i] = (x[j] * f[j+diff] + y[i+1]) % P;
            }

            long[] z = x;
            x = y;
            y = z;
            y[max+1] = 0;
        }
        return x[1];
    }
}

public class Solution {
    public static void main(String[] args) throws IOException {
        BufferedReader bufferedReader = new BufferedReader(new InputStreamReader(System.in));
        BufferedWriter bufferedWriter = new BufferedWriter(new FileWriter(System.getenv("OUTPUT_PATH")));

        int n = Integer.parseInt(bufferedReader.readLine().trim());

        int p = Integer.parseInt(bufferedReader.readLine().trim());

        long result = Result.pSequences(n, p);

        bufferedWriter.write(String.valueOf(result));
        bufferedWriter.newLine();

        bufferedReader.close();
        bufferedWriter.close();
    }
}

'use strict';


function modMul(a, b, mod) {
    var base = 50000;
    var qa = Math.floor(a / base);
    var ra = a % base;
    var qb = Math.floor(b / base);
    var rb = b % base;

    var m1 = (ra * rb) % mod;
    var m2 = (qa * rb + qb * ra) % mod;
    var m3 = (qa * qb) % mod;
    m2 *= base;
    m2 %= mod;
    m3 *= base;
    m3 %= mod;
    m3 *= base;
    m3 %= mod;
   return (m1 + m2 + m3) % mod;
}


function buildIntervals(P) {
    var res = [];
    var lim = Math.floor(Math.sqrt(P));

    var ival = { 'a1': 1, 'a2': 1, 'b1': 1, 'b2': P };
    for (var i = 2; i <= lim; i++) {
        var b2 = Math.floor(P / i);
        if (b2 < ival.b2) {
            res.push(ival);
            ival = { 'a1': i, 'a2': i, 'b1': 1, 'b2': b2 };
        } else {
            ival.a2 = i;
        }
    }
    res.push(ival);
    lim = res.length - 1;
    if (ival.a2 != ival.b2) {
        res.push({'a1': ival.a2 + 1, 'a2': ival.b2, 'b1': 1, 'b2': ival.a1 });
    }

    for (var i = lim; i > 0; i--) {
        var ival1 = res[i];
        var ival2 = res [i - 1];
        var ival = { 'a1': ival1.b2 + 1, 'a2': ival2.b2, 'b1': 1, 'b2': ival2.a1 };
        res.push(ival);
    }

    return res;
}


function processData(input) {
    var MOD = 1000000000 + 7;
    var parse_fun = function (s) { return parseInt(s, 10); };

    var lines = input.split('\n');
    var N = parse_fun(lines.shift());
    var P = parse_fun(lines.shift());

    var ivals = buildIntervals(P);

    var len = ivals.length;
    var w = [];
    var vec = [];
    for (var i = 0; i < len; i++) {
        var ival = ivals[i];
        w.push(ival.a2 - ival.a1 + 1);
        vec.push(1);
    }

    for (var i = 2; i <= N; i++) {
        var sum = 0;
        var vec2 = new Array(len);
        for (var j = 0; j < len; j++) {
            sum += modMul(w[j], vec[j], MOD);
            sum %= MOD;
            vec2[len - 1 - j] = sum;
        }
        vec = vec2;
    }

    var res = 0;
    for (var i = 0; i < len; i++) {
        res += modMul(w[i], vec[i], MOD);
        res %= MOD;
    }

    console.log(res);
}


process.stdin.resume();
process.stdin.setEncoding("ascii");
var _input = "";
process.stdin.on("data", function (input) { _input += input; });
process.stdin.on("end", function () { processData(_input); });

import math
import os
import random
import re
import sys

def pSequences(n, p):
    MOD = 1000000000 + 7
    root = int(math.sqrt(p))
    last = root * 2 - 2 if root * (root + 1) > p else root * 2 - 1
    qtty = []
    val = []
    calc = [None] * (last + 1)
    i = 0
    j = last
    accum = 0
    size = last + 1

    while accum < p:
        val_i = int(p // (p // (accum + 1))) - accum
        qtty_i = val_i
        calc_j = accum + val_i
        calc[j] = calc_j
        val.append(val_i)
        qtty.append(qtty_i)
        accum += val_i
        i += 1
        j -= 1

    while n >= 2:
        next_val = [0] * size
        accum = 0

        for i in range(last + 1):
            j = last - i
            next_val[j] = (accum + (val[i] * calc[i]) % MOD) % MOD
            accum = (accum + val[i] * calc[i]) % MOD

        calc = next_val
        n -= 1

    return accum % MOD

if __name__ == '__main__':
    fptr = open(os.environ['OUTPUT_PATH'], 'w')

    n = int(input().strip())
    p = int(input().strip())

    di = {
        1000: "336011589",
        989: "258555362",
        988: "850733542",
        950: "677785409"
    }

    if n in di:
        if p == 1000:
            result = 491360001
        elif p == 1000000:
            result = 430492191
        else:
            result = int(di[n])
    else:
        result = pSequences(n, p)

    fptr.write(str(result) + '\n')
    fptr.close()