Category Archives: Computer Science

A Probabilistic Model for COVID-19 Outbreak

As I am writing this blogpost, MIT is about to close until next fall and the undergrads are essentially getting evicted from the dorms.

Last two days were some of the most stressful days I have ever had. I have tons of work to do, but I cannot concentrate because I need to make sure my visa status will be okay, I need to pack my belongings and hopefully I’ll find a place to go.

COVID-19 is spreading rapidly across the US. This was expected. While MIT’s take on this is understandable, was it really necessary?

Let me try to answer that.

Let’s start with a small campus with \(N\) students where a small number of students \(\eta\) have been tested COVID-19 positive. I’ll also assume that the campus is in a state of lockdown — no student leaves the campus and no person from outside enters the campus.
Continue reading A Probabilistic Model for COVID-19 Outbreak

An Introduction to Variational Autoencoders

These people don’t exist! (source)

Variational Autoencoders (VAE) are really cool machine learning models that can generate new data. It means a VAE trained on thousands of human faces can new human faces as shown above!

Recently, two types of generative models have been popular in the machine learning community, namely, Generative Adversarial Networks (GAN) and VAEs. While GANs have had more success so far, a recent Deepmind paper showed that VAEs can yield results competitive to state of the art GAN models. Furthermore, VAE generated images retain more of the diversity of training dataset than GAN counterparts. Continue reading An Introduction to Variational Autoencoders

An Online Algorithm to Check for Bipartite Graphs

Bipartite Graphs

A graph \(G(V, E)\) is called bipartite if its vertices can be divided into two groups \(X\) and \(Y\) such that every edge connects one vertex from \(X\) and one vertex from \(Y\). The graph drawn below is bipartite.

Given a graph, can you determine if it is bipartite? Continue reading An Online Algorithm to Check for Bipartite Graphs

A Blind Robot Beside an Infinite Wall

Let’s think about the following problem:

Consider a wall that stretches to infinitely in both directions. There is a robot at position \(0\) and a door at position \(p\in\mathbb Z\) along the wall \((p\neq 0)\). The robot would like to get to the door, but it knows neither \(p\), nor the direction to the door. Furthermore, the robot cannot sense or see the door unless it stands right next to it. Give a deterministic algorithm that minimizes the number of steps the robot needs to take to get to the door.

This problem quite famous Continue reading A Blind Robot Beside an Infinite Wall

Multiplying Two Polynomials with Fast Fourier Transform

Polynomial multiplication is one of the most important problems in mathematics and computer science. It can be formally defined as follows:

You are given two polynomials roughly of equal size\[A(x)=a_0+a_1x+\dots+a_{n-1}x^{n-1}\]\[B(x)=b_0+b_1x+\dots+b_{n-1}x^{n-1}\]find a polynomial \[C(x)=c_0+c_1x+\dots+c_{2n-2}x^{2n-2}\] such that \(A(x)\cdot B(x)=C(x)\). At first, this may look like an easy problem, Continue reading Multiplying Two Polynomials with Fast Fourier Transform

Generating All Subsets of Size k in Python


Suppose you are given a set \(S\) with \(n\) elements and you need to generate every subset of size \(k\) from it. For instance, if \(S=\{3,4,5\}\) and \(k=2\), then the answer would be \(\{3,4\}\), \(\{3,5\}\), \(\{4,5\}\). So, how would you do that in Python?

First of all, this is a really simple exercise. Stuff like this often comes up when someone writes a moderately large script. I wanna talk about this one in particular because it has a combinatorial solution. Continue reading Generating All Subsets of Size k in Python

Lie Detectors and the Story of Halting Turing Machines


Is it possible to design a machine that detects lies? Sure, people have already built devices called polygraphs that monitor blood pressure, respiration, pulse etc. to determine if a person is giving false information. However, that is not same as “detecting lies” because polygraphs merely measure physical effects of telling lies. On the top of that, they are hella inaccurate. I am talking about REAL lie detectors, ie. machines that can instantly validate the statements themselves. Can we actually make such machines? (*insert new startup idea*) Continue reading Lie Detectors and the Story of Halting Turing Machines

ক্যালকুলেটর কীভাবে log-এর মান বের করে?

ব্যাপারটা একটু ক্যালকুলেটরের মত করে চিন্তা করা যাক। সে খুব সাদাসিধে একটি যন্ত্র। ধরা যাক, সে লগারিদম বোঝে না, কিন্তু বড়ো বড়ো সংখ্যার যোগ, বিয়োগ, গুণ, ভাগ বা বর্গমূল এক নিমিষে বের করে ফেলতে পারে। তাকে লগারিদম না শিখিয়ে আমার \(\log 2017\)-এর মান নির্ণয় করানো দরকার। বলো তো আমি এখন কী করতে পারি?

আমি হয়তো ক্যালকুলেটরকে বলতে পারি যে \(\ln (1+x)\)-এর টেইলর সিরিজ
\[\ln(1+x)=x-\frac{x^2}2+\frac{x^3}3-\cdots\;(|x| < 1)\]

সে ডানপাশের ধারাটির প্রথম কিছু পদের যোগফল নিলে \(\ln (1+x)\)-এর মোটামুটি নিখুঁত মান পাবে। এরপর সহজেই যেকোন ধ্বনাত্নক সংখ্যার \(\log\) নির্ণয় করা সম্ভব, কেননা

\[\log 2017 = \frac{\ln 2017}{\ln 10}=\frac{\ln\left(1-\frac{2016}{2017}\right)}{\ln\left(1-\frac{9}{10}\right)}\]

কিন্তু আমার এভাবে করাটা পছন্দ হল না। কারণ কাজটা বড্ড রসকষহীন।

আরও মজার কিছু ভাবা যাক।

খুব সহজেই দেখা যায় যে \[10^3 < 2017 < 10^4 \Longrightarrow 3 <\log 2017 < 4\] কিন্তু দশমিকের পরের অঙ্কটি কত হবে? \(3.1\) থেকে \(3.9\) পর্যন্ত \(10\) এর পাওয়ার নিয়ে আমি দেখতে পারি কোনটা \(2017\) এর সবচেয়ে কাছাকাছি। এরপর দশমিকের পরের দ্বিতীয় অঙ্কের জন্যও একই কাজ করতে পারি। এভাবে একটি একটি করে অঙ্ক আমি নির্ণয় করে যেতে পারি। কিন্তু প্রতিটি অঙ্কের জন্য এতবার করে \(10\) এর পাওয়ার নেওয়া আসলে অপচয়। একই কাজ বাইনারিতে করলে কেমন হয়? সেক্ষেত্রে প্রতিটি অঙ্ক \(0\) অথবা \(1\).  হিসেব করে দেখা যায় যে \(2017 < 10^{3}\times 10^{\frac{1}{2}}\approx 3162\). অতএব, বাইনারিতে দশমিকের পরের অঙ্কটি  \(0\). এর পরের অঙ্কটি? দেখা যায় যে \(2017 > 10^3\times 10^{\frac{1}{4}} \approx 1778\). অতএব, দশমিকের পরের দ্বিতীয় অঙ্কটি \(1\). বাইনারিতেও আমি একই কাজ পরবর্তী অঙ্কগুলোর জন্যে করে যেতে পারি। এভাবে \(\log\)-এর মান যেকোন ঘর পর্যন্ত নিখুঁতভাবে নির্ণয় করা সম্ভব।

এবার তাহলে ঝটপট কোড লিখে ফেলা যাক!

from math import sqrt
def log(base, num, precision = 40): 
    n, e = float(base), 1.0 
    # Find a power of base greater than num 
    while n <= num: 
        n, e = n**2, e * 2 
    # that's initial guess for log(base, num) 
    approx, log_approx = n, e 
    for i in range(precision): 
        # improve guess 
        n, e = sqrt(n), e / 2 
        if approx / n >= num: 
            approx = approx / n 
            log_approx -= e 
    return log_approx

পুনশ্চ ১: তুমি কি ধরতে পেরেছ উপরের কোডটা যে আসলে বাইনারি সার্চ?

পুনশ্চ ২: সব ক্যালকুলেটর দ্বিতীয় পদ্ধতি অনুসরণ করে না। 

ইন্টিজার ফ্যাক্টোরাইজেশন: পর্ব ১ পোলার্ডের রো (rho) এলগোরিদম

ইন্টিজার ফ্যাক্টোরাইজেশন পর্ব ০

গত পর্বে আমরা গায়ের জোরে কিছু সংখ্যার উৎপাদক বের করার চেষ্টা করেছিলাম। দুঃখজনকভাবে আমাদের বাচ্চা প্রোগ্রামটি \(20\) ডিজিটের একটি সংখ্যাকে উৎপাদকে বিশ্লেষণ করতে গিয়ে ক্লান্ত হয়ে গিয়েছিল। \(50\) বা \(100\) ডিজিটের কোন সংখ্যাকে উৎপাদকে বিশ্লেষণ করতে হলে সেই প্রোগ্রামের কী অবস্থা হবে ভাব একবার!

মানুষ যখন দেখল বড় বড় সংখ্যার জন্য \(2\) থেকে \(N-1\) [কিংবা \(\sqrt N\)] পর্যন্ত সব সংখ্যা দিয়ে \(N\)-কে এমন ভাগ করা অনেক সময়সাপেক্ষ, তখন তারা ট্রিক খুঁজতে লাগল কী করে আরও কম সংখ্যক সংখ্যা দিয়ে ভাগ করে \(N\)-এর উৎপাদক বের করা যায়। এমন অনেকগুলো ট্রিক রয়েছে, যার মাঝে একটি হল জন্মদিন সমস্যা(পরিশিষ্টে দেখ)। একে কাজে লাগিয়ে জন্ম হয়েছে একটি চালাকচতুর এলগোরিদমের।

পোলার্ডের \(\rho\) এলগোরিদম

মনে কর, আমাদেরকে \(N=55\)-কে উৎপাদকে বিশ্লেষণ করতে হবে। (আমরা যদিও জানি যে \(55\) এর একটি উৎপাদক \(p=5\), তবুও আমরা সেটি এলগোরিদম বের করব।)

আমরা করব কি, র‍্যান্ডম পূর্ণসংখ্যার একটি অসীম ধারা (infinite sequence) \(a_0, a_1, a_2, a_3, a_4,\ldots\) [সংক্ষেপে \(\{a_i\}_{i = 0}^{\infty}\)] থেকে খেয়াল খুশিমত কিছু পদ, যেমন: \(a_0, a_3, a_7, a_{100}\), তুলে নেব। জন্মদিন সমস্যা অনুযায়ী, এই চারটি পদের মাঝে \(\pmod {5}\)-এ কনগ্রুয়েন্ট অন্তত দুটি পদ পাওয়ার সম্ভাবনা\[1-e^{-\frac{4^2}{2\times 5}}\approx 80\%\]এমন দুটি পদ \(a_m\) এবং \(a_n\) যদি সত্যি সত্যি থাকে, তাহলে \[a_m\equiv a_n\pmod{5}\]\[\implies a_m-a_n\equiv 0\pmod{5}\]\[\implies d = \gcd(55, a_m-a_n)\ge 5\] অর্থাৎ, আমরা \(55\)-এর একটি উৎপাদক \(d\) নির্ণয় করে ফেলব। অবশ্য কেউ প্রশ্ন করতেই পারে \(d = 55\) হলে? হ্যাঁ, সেটা হতে পারে। কিন্তু চারটি সংখ্যার মাঝে \(\pmod{55}\)-এ কনগ্রুয়েন্ট দুটি সংখ্যা পাওয়ার সম্ভাবনা মাত্র\[1-e^{-\frac{4^2}{2\times 55}}\approx 13.5\%\]

অতএব, আমরা যদি \(\{a_i\}_{i = 0}^{\infty}\) থেকে প্রতিবার যেকোন দুইটি পদ \(a_i\) এবং \(a_j\) নিয়ে \(\gcd(N, a_i-a_j)\) নির্ণয় করে যেতে থাকি, তাহলে খুবই ভাল সম্ভাবনা আছে যে আমরা \(a_m-a_n\) পেয়ে যাব এবং \(55\)-এর একটি প্রকৃত উৎপাদক \(d\) নির্ণয় করে ফেলব।

তবে ঝামেলা হচ্ছে খাঁটি র‍্যান্ডম সংখ্যার অসীম ধারা ব্যবহার করলে সেটির অনেকগুলো পদ আমাদেরকে মেমরিতে রাখতে হবে, ফলে \(50\) বা \(100\) ডিজিটের কোন সংখ্যার জন্য কোটি কোটি পেটাবাইট মেমরির প্রয়োজন হবে। তাই আমরা এমন কোন সংখ্যার ধারা নেব যার মান প্রকৃতপক্ষে ডিটারমিনিস্টিক, অর্থাৎ সূত্র মেনে চলে, কিন্তু \(\pmod p\)-তে মোটামুটি র‍্যান্ডমলি আসে। (এমন সংখ্যাকে বলে সুডোর‍্যান্ডম সংখ্যা)

উদাহরণস্বরূপ, অসীম ধারা হিসেবে \(a_0=1\), এবং \(a_i=a_{i-1}^2+1\), অর্থাৎ \(1, 2, 5, 26, 677,\ldots\) ধারাটি নেওয়া নেব এবং প্রতিটি পদের জন্য \(\gcd(N, a_i-a_0)\) নির্ণয় করব। লক্ষ্য কর যে, আমাদের \(a_i\)-এর প্রকৃত মানের পরিবর্তে \(a_i\pmod N\) জানলেই চলে।

def rhonaive(N):
    a_0 = d = 1
    a_i = a_0
    while d == 1:
        a_i = (a_i * a_i + 1) % N
        d = gcd(a_i - a_0, N)
    return d

অনেকসময় দেখা যেতে পারে আমাদের নেওয়া সুডোর‍্যান্ডম সংখ্যার ধারাটি “যথেষ্ট র‍্যান্ডম” নয়। তেমন ক্ষেত্রে আমাদের \(a_i = a_{i-1}^2+1\) এর পরিবর্তে আরেকটু জটিল কোন সম্পর্ক ব্যবহার করতে হবে। তাই, সাধারণভাবে(generally) আমরা \(a_i=f(a_{i-1})\) ধরতে পারি যেখানে \(f(x)\) হচ্ছে একটি সুবিধাজনক বহুপদী। যেমন- উপরের প্রোগ্রামে \(f(x)=x^2+1\).

লক্ষ্য কর যে, \(\{a_i\}_{i = 0}^{\infty}\) ধারার প্রতিটি পদ শুধুমাত্র আগের পদের উপর নির্ভরশীল। তাই \(a_m \equiv a_n\pmod p\) একবার পাওয়া গেলে ধারাতে পূর্ববর্তী পদগুলো সাইক্লিকভাবে ঘুরে ফিরে আসতে থাকবে। অর্থাৎ,\(a_{m+i}\equiv a_{n+i}\pmod p\).

তবে এই সাইকেল \(a_0\) থেকে শুরু হতে হবে এমন কোন কথা নেই। যেমন- \(p=29\) এবং \(f(x)=x^2+1\) হলে, \(\pmod p\)-তে \(\{a_i\}_{i = 0}^{\infty}\) ধারাটি হয়

\(i\) \(0\) \(1\) \(2\) \(3\) \(4\) \(5\) \(6\) \(7\) \(8\) \(9\) \(10\) \(11\)
\(a_i\) \(1\) \(2\) \(5\) \(26\) \(10\) \(14\) \(23\) \(8\) \(7\) \(21\) \(\color{red}7\) \(\color{red}{21}\)

দেখতেই পাচ্ছ, \(a_0\)-র কনগ্রুয়েন্ট কোন পদ দ্বিতীয়বার পাওয়া যায়নি। তাই, প্রোগ্রাম ১-এ সবসময় \(a_0\) থেকে সাইকেল শুরু হবে ধরে নেওয়াটা ভুল ছিল।

এমন ধারাতে সাইকেল (অর্থাৎ \(\pmod p\)-তে কনগ্রুয়েন্ট দুটি পদ) নির্ণয় করতে আরেকটু বুদ্ধি খাটাতে হয়। আমি এখানে দুটি এলগোরিদমের কথা বলে দিচ্ছি। প্রথমটি ফ্লয়েডের এলগোরিদম, দ্বিতীয়টি ব্রেন্টের এলগোরিদম।

ফ্লয়েডের এলগোরিদম

মনে কর \(a_i\) ধারার সাইকেলের দৈর্ঘ্য \(l\). \(n\) যদি যথেষ্ট বড় হয়, তবে,\[a_n\equiv a_{n+l}\equiv a_{n+2l}\equiv\cdots\equiv a_{n+kl}\pmod p\]
আমরা যদি যথেষ্ট বড় কোন \(k\)-র জন্য \(n=kl\) নেই, তবে,\[a_n \equiv a_{n+kl}=a_{2kl}=a_{2n}\pmod p\]
অর্থাৎ, \(a_i\) ধারায় এমন অসীম সংখ্যক পদ \(a_n\) আছে যাদের জন্য \(a_n\equiv a_{2n}\pmod p\). সুতরাং \(\gcd(N, a_{2i}-a_i)\) নির্ণয় করে যেতে থাকলেই একটা সময় \(\pmod p\)-তে কনগ্রুয়েন্ট দুটি পদ পাওয়া যাবে। ফ্লয়েডের এলগোরিদম এই ধারণাটি ব্যবহার করে।

ফ্লয়েডের এলগোরিদমে একটি খরগোশ আর একটি কচ্ছপের মাঝে দৌড় প্রতিযোগিতার চলে। কচ্ছপটি ধারার একটি একটি পদে \(a_0, a_1, a_2,a_3,\ldots\) পর্যায়ক্রমে যেতে থাকে। একইসময়ে, খরগোশটি ধারার \(a_0, a_2, a_4, a_6,\ldots\) পদগুলোতে লাফিয়ে লাফিয়ে যেতে থাকে। আমরা খরগোশের অবস্থান এবং কচ্ছপের অবস্থানের পার্থক্যের সাথে \(N\)-এর গসাগু নির্ণয় করতে থাকি। যখনই গসাগু \(1\) এর চেয়ে বড় হয়ে যায়, এলগোরিদমটি টার্মিনেট করে।

def rhofloyd(N):
    #both hare and tortoise at a_0
    tortoise = hare = d = 1
    while d == 1:
        tortoise = f(tortoise) % N #tortoise steps once
        hare = f(f(hare)) % N #hare steps twice
        d = gcd(hare - tortoise, N)
    return d

ব্রেন্টের এলগোরিদম

ব্রেন্টের এলগোরিদমে অনেকগুলো রাউন্ডে খরগোশ আর কচ্ছপের দৌড় প্রতিযোগিতা চলে। প্রতিযোগিতার \(i\)-তম রাউন্ডে খরগোশ আর কচ্ছপ উভয়েই \(a_{2^i}\)-তম পদে অবস্থান করে। প্রতিযোগিতা শুরুর পর খরগোশ তার জায়গায় আরাম করে ঘুম দেয়, আর কচ্ছপটি টুক টুক করে একটি একটি পদ সামনে এগুতে থাকে। আমরা প্রতিবারই খরগোশের অবস্থান এবং কচ্ছপের অবস্থানের পার্থক্যের সাথে \(N\)-এর গসাগু \(d\) নির্ণয় করি। যখন কচ্ছপ খরগোশের চেয়ে \(2^i\) পদ সামনে এগিয়ে যায়, খরগোশের তড়াক করে জেগে ওঠে এবং সে এক লাফে কচ্ছপের অবস্থানে চলে আসে। এরপর প্রতিযোগিতার \(i+1\)-তম রাউন্ড শুরু হয়।

আগের মতই, \(d> 1\) হয়ে গেলে এলগোরিদমটি টার্মিনেট করবে।

লক্ষ্য কর যে, \(i\) এর মান বাড়ার সাথে সাথে \(2^{i}\) এবং \(2^{i+1}\)-এর মধ্যবর্তী পার্থক্যও বাড়তে থাকে। একটা সময় এই পার্থক্য \(\{a_i\}_{i = 0}^\infty\) ধারার সাইকেলের দৈর্ঘ্য \(l\) থেকে বড় হয়ে যায়। তখন কচ্ছপের কোন এক অবস্থানের জন্য খরগোশ ও কচ্ছপের অবস্থানের পার্থক্য হবে ঠিক \(l\)টি পদ। এই অবস্থাতেই আসলে এলগোরিদমটি টার্মিনেট করবে।

def rhobrent(N):
    #both hare and tortoise are at a_(2^0)
    #power is the nearest 2^i behind tortoise
    tortoise = hare = f(1)
    power, lead, d = 1, 0, 1
    while d == 1:
        if lead == power: #tortoise reaches 2^(i+1)-th term
            hare = tortoise #hare jumps to tortoise's position
            power *= 2 #update power
            lead = 0 #reset lead, since both hare and tortoise 
                     #are in the same position now
        tortoise = f(tortoise) % N #tortoise steps once
        lead += 1 #tortoise leads the race by one more step
        d = gcd(hare - tortoise, N)
    return d

ব্রেন্টের এনালাইসিস অনুসারে, ব্রেন্টের এলগোরিদম ফ্লয়েডের এলগোরিদম থেকে \(24\%\) দ্রুততর।

পূর্ণাঙ্গ প্রোগ্রাম

ব্রেন্টের এলগোরিদম দিয়ে যেহেতু সাইকেল শনাক্তকরণ সবচেয়ে দ্রুত করা যায়, তাই সেটি দিয়েই নিচে একটি পূর্ণাঙ্গ প্রোগ্রাম দেওয়া হল।

from fractions import gcd
from time import time

def f(x):
    return (x**2 + x + 1) % N

def rhobrent(N):
    tortoise = hare = f(1)
    power, lead, d = 1, 0, 1
    while abs(d) == 1:
        if lead == power:
            hare, lead = tortoise, 0
            power *= 2
        tortoise = f(tortoise)
        lead += 1
        d = gcd(tortoise - hare, N)
    return d

N = int(raw_input())
exec_time = time()
d = 1 if(-2 < N < 2) else rhobrent(N)
exec_time = time() - exec_time
print d, '*', N / d
print "Execution time: %.3f s" % exec_time
#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include <iso646.h> //and

#define MAX 2147483647

long long f(long long x);
long long gcd(long long a, long long b);
long long rhobrent(long long N);

long long N;

int main(void){
    long long d;
    scanf("%lld", &N);
    clock_t exec_time = clock();
    d = (N < 2 and N > -2) ? 1 : rhobrent(N);
    exec_time = clock() - exec_time;
    double t = (double) exec_time / CLOCKS_PER_SEC;
    printf("%lld * %lld\n", d, N/d);
    printf("Execution time: %.3lf s\n", t);
    return 0;

long long f(long long x){
    if(x < MAX) 
        return (x * x + x + 1) % N;
    long long y = x, output = x + 1;
    while(y > 0){
        if(y % 2 == 1)
            output += x;
        x *= 2;
        y /= 2;
        output %= N;
        x %= N;
    return output;

long long gcd(long long a, long long b){
    long long r;
    while(b != 0){
        r =  a % b;
        a = b;
        b = r;
    return (a > 0) ? a : -a;

long long rhobrent(long long N){
    long long tortoise = f(1), hare = f(1), power = 1, lead = 0, d = 1;
    while(d == 1){
        if(lead == power){
            hare = tortoise;
            power *= 2;
            lead = 0;
        tortoise = f(tortoise);
        d = gcd(N, tortoise - hare);
    return d;


  1. উপরের কোডে \(f(x)=x^2+x+1\) এবং ফাংশানের ভেতরেই \(N\) দিয়ে মডুলাসের কাজ শেষ করা হয়েছে।
  2. প্রায় সকল আধুনিক কম্পাইলারে long long int-এর ব্যাবধী \([-2^{63}, 2^{63})\). তাই সি প্রোগ্রামে \(f(x)\) ফাংশনটি ঠিকঠাকমত কাজ করার জন্য \[x^2+x+1\leq 2^{63}-1\Longrightarrow x<2^{31.5}-1\] হতে হবে। এজন্য \(x \ge 2^{31}\) -এর জন্য \(y = x\) নামের একটি নতুন ভেরিয়েবল নিয়ে এর প্রতিটি বিটের সাথে আলাদাভাবে \(x\)-এর গুণফল বের করে সেই সমষ্টির মডুলাস নেওয়া হয়েছে।
  3. \(\pmod N\) এবং \(\pmod p\)-তে \(N\)-এর সাইকেল একই হলে উপরের প্রোগ্রাম \(N\) কৃত্রিম হওয়া সত্বেও এর উৎপাদক খুঁজে পাবে না। এরকম হলে \(f(x)\) পলিনোমিয়ালটি পরিবর্তন করে আবার চেষ্টা করতে পার।


মনে কর, \(\{a_i\}_{i=0}^{\infty}\) একটি সত্যিকারে র‍্যান্ডম ধারা। জন্মদিনের সমস্যা অনুযায়ী, \(\pmod p\)-তে কনগ্রুয়েন্ট দুটি পদ পাওয়ার সম্ভাবনা \(99\%\) হতে হলে \(\sqrt{-2\times\ln(1-0.99)}\sqrt{p}\)-টি পদ নিতে হবে। অতএব, রো এলগোরিদমের কমপ্লেক্সিটি \(O(\sqrt p)\leq O\left(N^{\frac 1 4}\right)\). [অর্থাৎ \(N\) সংখ্যাটিতে \(k\)টি ডিজিট থাকলে রো এলগোরিদমের রানটাইম \(O\left(10^{\frac k 4}\right)\)]
গতপর্বের চেয়ে আরেকটু ভালো বাউন্ড, তাই না? পরবর্তী পর্বে আমরা জন পোলার্ডেরই আরেকটি এলগোরিদম দেখব, যেটির নাম \(p-1\) এলগোরিদম


জন্মদিন সমস্যা

তোমার ফেসবুক ফ্রেন্ডলিস্টে নূন্যতম কতজন থাকলে একথা নিশ্চিত করে বলা যাবে যে অন্তত দুইজন ফ্রেন্ডের জন্মদিন একই তারিখে হওয়ার সম্ভাবনা \(99\%\)?

মনে কর, তোমার ফ্রেন্ডলিস্টে \(m\) জন ফ্রেন্ড আছে, এবং কোন ফ্রেন্ডের জন্মদিন ২৯ ফেব্রুয়ারি নয়। কোন একজন ফ্রেন্ডের জন্মদিন ৬ অক্টোবর হলে আর কোন ফ্রেন্ডের জন্মদিন ৬ অক্টোবর না হওয়ার সম্ভাবনা \(1-\frac {1}{365}\). আরেক ফ্রেন্ডের জন্মদিন মনে কর ৪ আগস্ট। সুতরাং, তৃতীয় কোন ফ্রেন্ডের জন্মদিন ৬ অক্টোবর বা ৪ আগস্ট না হওয়ার সম্ভাবনা \(1-\frac {2}{365}\). অতএব, \(m\) জন ফ্রেন্ডের ভিন্ন ভিন্ন দিনে জন্মদিন পড়ার সম্ভাবনা \[\left(1-\frac {1}{365}\right)\cdots\left(1-\frac {m-1}{365}\right)\]
তাহলে, অন্তত দুজন ফ্রেন্ডের জন্মদিন একই দিনে পড়ার সম্ভাবনা \[1 – \left(1-\frac {1}{365}\right)\cdots\left(1-\frac {m-1}{365}\right)\]
\(x\)-এর মান \(0\)-এর কাছাকাছি হলে আমরা টেইলরের সিরিজ থেকে \(1-x\approx e^{-x}\) লিখতে পারি। এভাবে উপরের সম্ভাবনাটিকে সরল করে লিখা যায় \[1-e^{-\frac {1}{365}}\cdot e^{-\frac {2}{365}}\cdots e^{-\frac {m-1}{365}}\]\[\approx 1-e^{-\frac{m^2}{2\times 365}}\]এখন \(1-e^{-\frac{m^2}{2\times 365}}\ge 0.99\) হতে হলে \[m\ge \sqrt{-2\times 365\ln(1-0.99)}\approx 58\] হতে হবে।

তোমার ফ্রেন্ডলিস্টে মাত্র \(58\) জন মানুষ থাকলেই তুমি \(99\%\) গ্যারান্টি দিয়ে বলতে পারবে যে অন্তত দুইজনের জন্মদিন একই তারিখে! এই বিখ্যাত সমস্যাটির নাম জন্মদিন সমস্যা

লক্ষ্য কর জন্মদিনের সমস্যার সমাধান থেকে এটা বলা যায় যে, \(0\) থেকে \(N-1\) এর মাঝে \(m\) সংখ্যা নিলে কোন একটি সংখ্যা অন্তত দুইবার নেওয়ার সম্ভাবনা \(1-e^{-\frac{m^2}{2N}}\).

ইন্টিজার ফ্যাক্টোরাইজেশন পর্ব ০

ইন্টিজার ফ্যাক্টোরাইজেশন: পর্ব ১ পোলার্ডের রো (rho) এলগোরিদম

আমরা কেন ইন্টিজার ফ্যাক্টোরাইজেশন করি

ইন্টেজার ফ্যাক্টোরাইজেশন অর্থ হচ্ছে একটা পূর্ণসংখ্যাকে উৎপাদকে বিশ্লেষণ করা। যেমন \(36 = 4\times 9\). এখানে \(4\) আর \(36\) উভয়েই \(36\)-এর উৎপাদক।

ছোটবেলায় আমাদের সবাই ইন্টেজার ফ্যাক্টোরাইজেশন করতাম। কেন জানো? পরীক্ষায় নম্বর পাওয়ার জন্য! প্রশ্নে একটা সংখ্যা দেওয়া থাকত আর আমরা দুই থেকে শুরু করে একটা একটা সংখ্যা দিয়ে প্রশ্নে দেওয়া সংখ্যাটাকে ভাগ দিতাম। যখন ভাগ যেত, তখন বোঝা যেত যে একটা উৎপাদক (বা গুণনীয়ক) পেয়ে গেছি। খুবই রুটিনমাফিক একটা কাজ, তাই না?

তবে বাচ্চাদের পরীক্ষার প্রশ্নে আসা ছাড়াও ইন্টিজার ফ্যাক্টোরাইজেশনের অন্যরকম একটা গুরুত্ব আছে। যেমন- অনলাইনে তোমার, আমার, সকলের ব্যক্তিগত তথ্য পাজি লোকদের হাত থেকে সুরক্ষিত রাখতে গোপন সংকেতের আড়ালে রাখা হয়। এরকম অনেক গোপন সংকেত বানাতে ও ভাঙতে দুইশ, তিনশ ডিজিটের বিশাল বিশাল সংখ্যাকে উৎপাদকে বিশ্লেষণ করতে হয়, কিংবা প্রাইম খুঁজে বের করতে হয়। এসব হিসেব এত বড় হয় যে খাতা কলমে করতে গেলে আমাদের চুল পেকে যাবে, পাতার পর পাতা ফুরিয়ে যাবে। সেজন্য, আমরা তখন কম্পিউটারের সাহায্য নেই।

কিন্তু সমস্যা কি জানো? কম্পিউটার খুব বোকা। যদি আমি বলি, “কম্পিউটার, আমার একটা হিসেবে \(731\)-এর উৎপাদকে বিশ্লেষণ লাগবে। তুমি চট করে সেটা বের করে দাও তো দেখি!”, সে হা করে বসে থাকবে। কারণ সে তো আর জানে না কীভাবে উৎপাদকে বিশ্লেষণ করতে হয়!

কম্পিউটারকে দিয়ে উৎপাদকে বিশ্লেষণ করানোর জন্য তাকে তোমার প্রাইমারির বাচ্চাদের মত একটু একটু করে বুঝিয়ে বলতে হবে, “কম্পিউটার, \(2\) থেকে \(730\) পর্যন্ত প্রতিটি সংখ্যা দিয়ে \(731\)-কে ভাগ দাও। যদি কোনটা দিয়ে নিঃশেষে ভাগ যায়, তবে সেই সংখ্যা এবং ভাগফল আমাকে দেখাও।”

def factor(N):
    for d in range(2, abs(N) + 1):
        if N % d == 0: return d

N = int(raw_input())
d = 1 if(2 > N > -2) else factor(N)
print d, '*', N / d

কম্পিউটার যেহেতু বাংলা পারে না, তাই সে বোঝে এমন একটি ভাষা পাইথনে আমরা আমাদের কথাগুলো লিখেছি। এই কোডে একটি পূর্ণসংখ্যা \(N\)-কে ইনপুট হিসেবে দেওয়া হবে। আর কম্পিউটার সাথে সাথে \(N\)-কে উৎপাদকে বিশ্লেষণ করে খুশিমনে আমাদের জানিয়ে দেবে। যত বড় সংখ্যাই হোক না কেন, আমাদের ছেলেবেলায় শেখা নিয়ম দিয়ে তার উৎপাদক বের করা সম্ভব।

উপরের প্রোগ্রামে \(731\) ইনপুট দিলে কম্পিউটার আমাদের বলে দেবে যে \(731 = 17\times 43\).

এবার একটা কাজ কর। উপরের প্রোগ্রামটি চালিয়ে \(N = 18446744073709551619\) ইনপুট দাও। কী হল?

প্রোগ্রামটি আর থামছে না।

এর কারণ হচ্ছে উপরের প্রোগ্রামে কম্পিউটারকে \(2\) থেকে \(N-1\) পর্যন্ত সকল সংখ্যা দিয়েই \(N\)-কে ভাগ (প্রকৃতপক্ষে মডুলাস নির্ণয়) করতে হচ্ছে। কম্পিউটারের প্রতিটা ভাগ করতে কিছুটা সময় লাগে। তাই \(N\) যত বড় হবে, কম্পিউটারকে তত বেশি সংখ্যা দিয়ে ভাগ করতে হবে, এবং প্রোগ্রামের রান করতে সময় (রানটাইম) বেশি লাগবে। রানটাইম \(N\)-এর সাথে পাল্লা দিয়ে কত দ্রুত বাড়বে সেটার একটা ধারণাও আমরা পেতে পারি।

ধর \(N\) সংখ্যাটিতে \(k\) সংখ্যক ডিজিট আছে। আমরা আসন্নভাবে \(N\approx 10^k\) ধরতে পারি। যদি \(N\) একটি মৌলিক সংখ্যা হয়, তবে প্রোগ্রাম ১-এ \(N\)-কে \(2\) থেকে \(N-1\) পর্যন্ত সকল সংখ্যা দিয়েই ভাগ করা হবে। যদি প্রতিবার ভাগ করতে সর্বোচ্চ \(c\) পরিমাণ সময় লাগে, তবে উপরের প্রোগ্রামটি রান করতে সর্বোচ্চ \(Nc=10^{k}c\) সময় লাগবে। এর অর্থ হচ্ছে, ডিজিটের সংখ্যা বাড়ার সাথে সাথে উপরের প্রোগ্রামের রানটাইমও সূচকীয়ভাবে বৃদ্ধি পায়। [বা, বিগ O নোটেশনে উপরের প্রোগ্রামের রানটাইম \(O(10^{k})\)]

প্রোগ্রামিং-এ রানটাইম সূচকীয়ভাবে বাড়া খুবই দুঃখের একটি কথা। কারণ ইনপুটের একটুখানি পরিবর্তনে প্রোগ্রামের রানটাইমে রাতদিন পরিবর্তন হয়ে যেতে পারে।

মনে কর তোমার কাছে যদি কোন মহাক্ষমতাধর সুপার কম্পিউটার আছে যে \(200\) ডিজিটের একটা সংখ্যার উৎপাদকে বিশ্লেষণ প্রোগ্রাম ১ দিয়ে \(1\) সেকেন্ডে বের করতে পারে। তাহলে \(210\) ডিজিটের একটা সংখ্যার উৎপাদকে বিশ্লেষণ করতে কম্পিউটারটির সময় লাগবে \(10^{210-200}=10^{10}\) সেকেন্ড \(\approx 317\) বছর!

আমদের জীবনটা খুব ছোট, জানো? একটা প্রোগ্রামের আউটপুট দেখে যাওয়ার জন্য আমাদের \(317\) বছর বেঁচে থাকা সম্ভব না। তাই খুব ফাস্ট কোন এলগোরিদম দরকার।

প্রোগ্রাম ১-কেই অনেকখানি ফাস্ট বানানো যায় যদি আমরা \(N-1\) এর পরিবর্তে \(\sqrt{N}\) পর্যন্ত প্রতিটা সংখ্যা দিয়ে ভাগ দেই। কেননা প্রতিটি কৃত্রিম সংখ্যার \(\sqrt N\)-এর মাঝেই অন্তত একটা উৎপাদক থাকে। আর যদি উৎপাদক না পাওয়া যায়, তবে সংখ্যাটি প্রাইম। এই পরিবর্তনের ফলে আমাদের এলগোরিদমের রানটাইম হবে \(10^{\frac k 2}c\) [বা \(O\left(10^{\frac k 2}\right)\)]. কিন্তু এটাও যে সূচকীয় হারে বাড়ে!

এর থেকে ফাস্ট কি কোন এলগোরিদম নেই? উত্তরটা হচ্ছে, হ্যাঁ আছে। তবে সেটি জানতে হলে তোমার মডুলার এরিথমেটিক সম্পর্কে হালকা পাতলা ধারণা লাগবে।

যদি \(a\) এবং \(b\) কে \(n\) দিয়ে ভাগ করলে একই ভাগশেষ থাকে, তাহলে আমরা লিখি \(a \equiv b\pmod{n}\). এটাকে পড়া হয় \(a\) is congruent to \(b\) modulo \(n\). যেমন: \(57\equiv 112\pmod{55}\).

লক্ষ্য কর যে, যদি \(n\) দিয়ে \(a\)-কে ভাগ করলে \(r\) ভাগশেষ থাকে, তবে \(a\equiv r\pmod{n}\) লিখা যায়। যেমন: \(57\equiv 2\pmod{55}\).

যেকোন পূর্ণসংখ্যাকে যদি \(55\) দিয়ে ভাগ কর, তবে তুমি \(0\) থেকে \(54\) এর মাঝে কোন একটি ভাগশেষ পাবে। তাহলে, তুমি যদি অনেকগুলো র‍্যান্ডম পূর্ণসংখ্যার নাও, এদের মাঝে \(\pmod {55}\)-এ কনগ্রুয়েন্ট, অর্থাৎ \(55\) দিয়ে ভাগ করলে একই ভাগশেষ থাকে, এমন দুটি পূর্ণসংখ্যা পাওয়ার কিছুটা সম্ভাবনা থাকে, তাই না? যত বেশি সংখ্যা নেবে, এই সম্ভাবনাও তত বাড়বে।

আবার, একই সংখ্যাগুলির মাঝে \(\pmod {5}\)-এ কনগ্রুয়েন্ট দুটি সংখ্যা পাওয়ার সম্ভাবনাও থাকে। কিন্তু এই সম্ভাবনা \(\pmod{55}\)-এ কনগ্রুয়েন্ট হওয়ার সম্ভাবনার চেয়ে অনেক বেশি। কেননা, \(\pmod{5}\)-এ তো খালি \(0\) থেকে \(4\) পর্যন্ত ভাগশেষ আসতে পারে।

এই ধারণাকে ব্যবহার করে উৎপাদকে বিশ্লেষণের একটি চমৎকার এলগোরিদম রয়েছে, যার নাম পোলার্ডের রো(\(\rho\)) এলগোরিদম। পরবর্তী পর্বে এটি নিয়ে আলোচনা করব।