Randomized অ্যালগরিদম — Las Vegas আর Monte Carlo

১ · এলোমেলোতা কেন সাহায্য করে — প্রতিপক্ষ, লটারি আর খাতা যাচাই

এ পর্যন্ত প্রায় সব অ্যালগরিদম ছিল deterministic: একই ইনপুটে প্রতিবার হুবহু একই ধাপ। এর একটি দুর্বলতা আছে। কেউ যদি আপনার কোড জানেন, তিনি ঠিক সেই ইনপুট বানাতে পারেন, যেখানে আপনার অ্যালগরিদম সবচেয়ে ধীর। পর্ব ২২-এ দেখেছেন, সবসময় প্রথম উপাদানকে pivot নিলে quick sort সাজানো অ্যারেতে O(n²) হয়ে যায়। পর্ব ৩২-এ দেখেছেন hash flooding: একই ঘরে পড়ে এমন হাজার হাজার key পাঠিয়ে সার্ভার থামিয়ে দেওয়া। এমন ইচ্ছাকৃত খারাপ ইনপুটকে বলে adversarial input.

Randomized অ্যালগরিদম চলার সময় নিজে এলোমেলো সংখ্যা নেয় — যেন মুদ্রা ছোড়ে। প্রতিপক্ষ ইনপুট বেছে নেন আগে, কিন্তু মুদ্রা পড়ে পরে। তাই তিনি জানেন না কোন pivot বা কোন hash function বাছা হবে। আমরা যে নিশ্চয়তা দিই, তা হয় এমন: “প্রতিটি ইনপুটের জন্য, আমাদের মুদ্রার ওপর গড়ে, খরচ এত” — এটিই expected running time. পর্ব ৫-এর average case-এর সাথে পার্থক্যটি মনে রাখুন: সেখানে ধরে নিই ইনপুট এলোমেলো; এখানে ইনপুট নিয়ে কোনো অনুমান নেই।

এলোমেলোতা আরও দুটি উপহার দেয়। সরলতা: treap বা skip list (পর্ব ৩৮) red-black tree-র চেয়ে অনেক ছোট কোডে একই প্রত্যাশিত গতি দেয়। গতি: কিছু সমস্যায় জানা সবচেয়ে ভালো deterministic পদ্ধতির চেয়ে randomized পদ্ধতি অনেক দ্রুত বা সরল — যেমন matrix গুণফল যাচাই বা বড় সংখ্যার primality test.

Las Vegas আর Monte Carlo — দুটোই জুয়ার শহরের নাম; নামকরণের গল্প ভাবনার প্রশ্ন ২-এ। মনে রাখার সহজ উপায়: Las Vegas-এ আপনি জেতা পর্যন্ত খেলেন (সময় অনিশ্চিত, ফল নিশ্চিত); Monte Carlo-তে নির্দিষ্ট দান খেলে উঠে পড়েন (সময় নিশ্চিত, ফল অনিশ্চিত)।

২ · সম্ভাবনার যন্ত্রপাতি — expectation, indicator আর দুটি সীমা

এই পর্বের প্রতিটি প্রমাণ পাঁচ-ছয়টি সরল হাতিয়ারে দাঁড়িয়ে। প্রতিটি একবার পরিষ্কার করে নিই।

Random variable আর expectation. Random variable X হলো এলোমেলো পরীক্ষার ফল থেকে একটি সংখ্যা — যেমন “quicksort-এ মোট তুলনা”. তার expected value E[X] = Σ x · P(X = x), অর্থাৎ বহুবার চালালে গড়ে যা আসে। একটি ছক্কার E = (1 + 2 + … + 6)/6 = 3.5.

Indicator variable. কোনো ঘটনা A ঘটলে I_A = 1, না ঘটলে 0 — এটি indicator random variable. এর মূল গুণ: E[I_A] = P(A). কোনো জিনিস “গোনা” মানে অনেকগুলো indicator যোগ করা।

Linearity of expectation. যেকোনো random variable-এর জন্য E[X + Y] = E[X] + E[Y] — তারা স্বাধীন হোক বা না হোক. এটিই linearity of expectation, এই পর্বের সবচেয়ে শক্তিশালী হাতিয়ার। জটিল জিনিসের গড় বের করতে সেটিকে indicator-এর যোগফল হিসেবে লিখুন, প্রতিটির সম্ভাবনা বের করুন, যোগ করুন — ঘটনাগুলোর পারস্পরিক সম্পর্ক নিয়ে মাথা ঘামাতে হয় না।

স্বাধীনতা. দুটি ঘটনা independent হলে P(A ও B) = P(A) · P(B). আলাদা আলাদা মুদ্রা নিক্ষেপ স্বাধীন — তাই k বার স্বাধীনভাবে চালানো পরীক্ষার প্রতিটিতে ব্যর্থ হওয়ার সম্ভাবনা গুণফল। Amplification-এর পুরো ভিত্তি এটি। সাবধান: linearity-তে স্বাধীনতা লাগে না, কিন্তু গুণফলের নিয়মে লাগে।

Markov-এর অসমতা. X ≥ 0 হলে যেকোনো a > 0-এর জন্য P(X ≥ a) ≤ E[X]/a. এটিই Markov's inequality. যুক্তি এক লাইনের: X ≥ a · I_X≥a, দুই পাশে expectation নিন। ব্যবহার: একটি Las Vegas অ্যালগরিদমের প্রত্যাশিত সময় T হলে, 2T সময়ের বেশি লাগার সম্ভাবনা ≤ 1/2. তাই “2T ধাপ পর থামাও, না পারলে ‘জানি না’ বলো” — Las Vegas থেকে Monte Carlo বানানো যায়।

Union bound. P(A₁ বা A₂ বা … বা A_m) ≤ P(A₁) + … + P(A_m) — ঘটনাগুলো স্বাধীন না হলেও। এটিই union bound: “অনেকগুলো খারাপ ঘটনার যেকোনো একটি” ঘটার সম্ভাবনা বের করতে। পর্ব ৩৮-এ skip list-এর উচ্চতা এভাবে বেঁধেছিলেন।

Chernoff bound — শুধু নামটুকু. অনেকগুলো স্বাধীন 0/1 variable-এর যোগফল তার গড় থেকে বেশি দূরে যাওয়ার সম্ভাবনা সূচকীয় হারে কমে: গড় μ হলে, 0 < δ ≤ 1-এর জন্য P(X ≥ (1 + δ)μ) ≤ e^−δ²μ/3. এটি Chernoff bound. Markov শুধু গড় জানে, তাই দুর্বল; Chernoff স্বাধীনতা কাজে লাগায়, তাই অনেক শক্ত। “৫০টি স্বাধীন চালানোর সংখ্যাগরিষ্ঠ ভোট প্রায় কখনো ভুল নয়” — এমন দাবি Chernoff দিয়েই প্রমাণ হয়।

৩ · Las Vegas — randomized quicksort-এর প্রমাণ, quickselect আর treap

পর্ব ২২-এ randomized quicksort দেখেছেন: প্রতিটি ধাপে pivot এলোমেলোভাবে বাছা। এখন indicator দিয়ে তার প্রত্যাশিত তুলনা-সংখ্যা প্রমাণ করি — প্রমাণটি বিশ্ববিদ্যালয়, চাকরি আর ইন্টারভিউর পরীক্ষায় খুব প্রচলিত।

সাজানো ক্রমে উপাদানগুলোর নাম দিন z₁ < z₂ < … < z_n (সব আলাদা ধরি)। i < j হলে X_ij = 1, যদি z_i আর z_j কখনো সরাসরি তুলনা হয়। দুটি উপাদান সর্বোচ্চ একবার তুলনা হয় (তুলনা হয় শুধু pivot-এর সাথে, আর pivot পরে আর কোনো ভাগে থাকে না)। তাই মোট তুলনা X = Σ_i<j X_ij.

মূল প্রশ্ন: P(z_i আর z_j তুলনা হয়) কত? সেট Z_ij = {z_i, z_i+1, …, z_j} দেখুন, এতে j − i + 1 টি উপাদান। এই সেট থেকে প্রথম যে উপাদান pivot হয়, সে-ই সিদ্ধান্ত নেয়। (ক) সে z_i বা z_j হলে দুজনের তুলনা হয়। (খ) সে মাঝের কোনো z_k হলে z_i যায় বাঁ ভাগে, z_j ডান ভাগে — আর কখনো দেখা হয় না। সেট থেকে কেউ pivot হওয়ার আগ পর্যন্ত পুরো সেট একই ভাগে থাকে, আর pivot এলোমেলো, তাই সেটের প্রত্যেকের “প্রথম হওয়ার” সম্ভাবনা সমান। ফলে

এখন linearity: E[X] = Σ_i<j 2/(j − i + 1). d = j − i ধরলে, প্রতিটি i-এর জন্য ভেতরের যোগফল ≤ 2(1/2 + 1/3 + … + 1/n) < 2H_n, যেখানে H_n ≈ ln n হলো harmonic সংখ্যা (পর্ব ৬)। n টি i মিলিয়ে

সঠিক মান 2(n + 1)H_n − 4n. এটি প্রতিটি ইনপুটের জন্য — সাজানো, উল্টো, সব একই। Worst case এখনো O(n²) (প্রতিবার সবচেয়ে ছোটটি pivot হলে), কিন্তু তার জন্য প্রতিপক্ষকে আমাদের মুদ্রা আগে থেকে জানতে হবে। উত্তর সবসময় সঠিক সাজানো অ্যারে — তাই এটি Las Vegas. জায়গা: ছোট ভাগে আগে recursion করলে O(log n) stack.

Randomized quickselect — প্রত্যাশিত O(n)

k-তম ছোট উপাদান চাই (পর্ব ২২-এর quickselect). এলোমেলো pivot দিয়ে partition করুন, তারপর শুধু যে ভাগে k-তম আছে সেখানে recursion. সংক্ষিপ্ত যুক্তি: pivot যদি মাঝের অর্ধেকে পড়ে (rank n/4 থেকে 3n/4-এর মধ্যে; সম্ভাবনা 1/2), পরের ভাগের আকার ≤ 3n/4. তাই প্রতি “ভালো pivot” পেতে গড়ে ২টি চেষ্টা, প্রতিটি চেষ্টা O(বর্তমান আকার)। মোট প্রত্যাশিত কাজ ≤ 2c(n + 3n/4 + (3/4)²n + …) = 8cn = O(n) — একটি geometric series. Worst case O(n²); deterministic O(n) চাইলে পর্ব ২২-এর median of medians.

Treap আর skip list — এলোমেলো ভারসাম্য

পর্ব ৩৮-এ treap আর skip list দেখেছেন। একই ধরনের indicator-প্রমাণ সেখানে দেখিয়েছিল যে x_i হলো x_k-এর ancestor সম্ভাবনা 1/(|k − i| + 1)-এ, তাই প্রত্যাশিত গভীরতা < 2 ln n. লক্ষ করুন, quicksort-এর “তুলনা হয়” আর treap-এর “ancestor হয়” একই ঘটনা — treap আসলে quicksort-এর recursion tree. দুটোই Las Vegas: খোঁজার উত্তর সবসময় সঠিক, শুধু tree-র আকৃতি এলোমেলো।

৪ · Monte Carlo — Freivalds-এর matrix যাচাই আর amplification

সমস্যা: তিনটি n × n matrix A, B, C দেওয়া। C = A·B কি না? সরাসরি গুণ করে মেলাতে সাধারণ পদ্ধতিতে Θ(n³), আর সবচেয়ে উন্নত তাত্ত্বিক পদ্ধতিতেও প্রায় n^2.37 (পর্ব ৪৯-এর matrix multiplication exponent). ১৯৭৭ সালে লাটভিয়ার Rūsiņš Freivalds একটি অসাধারণ সরল উপায় দেন — Freivalds' algorithm:

1. একটি এলোমেলো vector r ∈ {0, 1}ⁿ বাছুন — প্রতিটি ঘর আলাদা মুদ্রা।
2. B·r হিসাব করুন (n² গুণ), তারপর A·(B·r) (আরও n²), আর C·r (আরও n²)।
3. A·(B·r) ≠ C·r হলে বলুন “C ভুল” — এটি নিশ্চিত. সমান হলে বলুন “সম্ভবত ঠিক”.

প্রতিটি চেষ্টা Θ(n²) — পুরো গুণফলের চেয়ে অনেক কম। নিচের ডায়াগ্রামে C-তে একটিমাত্র ঘর ভুল (row 1, column 2-এ 3-এর জায়গায় 4); r-এর মানগুলো এই পর্বের প্রোগ্রামের seed থেকে আসা।

কেন ভুল ধরার সম্ভাবনা অন্তত 1/2

D = A·B − C ধরুন। C ভুল হলে D ≠ 0, তাই D-এর কোনো ঘর d_ij ≠ 0. অ্যালগরিদম ভুল ধরে না শুধু যখন D·r = 0. D·r-এর i-তম ঘর = d_ij r_j + S, যেখানে S বাকি সব r-এর ওপর নির্ভর করে। বাকি সব r আগে ঠিক হয়ে গেছে ধরুন — তখন S একটি নির্দিষ্ট সংখ্যা। r_j-এর দুটি মানের (0 আর 1) জন্য d_ij r_j + S দুটি আলাদা সংখ্যা (কারণ d_ij ≠ 0), তাই সর্বোচ্চ একটি শূন্য হতে পারে। r_j স্বাধীন মুদ্রা, তাই P(D·r = 0) ≤ 1/2. ∎ (এই যুক্তির নাম principle of deferred decisions — একটি মুদ্রা সবার শেষে ছুড়ছি ধরে নেওয়া।)

উত্তর একমুখী: “অসমান” সবসময় সঠিক, ভুল হতে পারে শুধু “সম্ভবত ঠিক”. এমন অ্যালগরিদমের ভুলকে বলে one-sided error. এখন amplification: k বার স্বাধীন r দিয়ে চালান, একবারও অসমান হলে “ভুল” বলুন। ভুল C সব k বার পার পাওয়ার সম্ভাবনা ≤ (1/2)^k — স্বাধীনতার গুণফল। k = 20-এ ≤ 2⁻²⁰ ≈ 10⁻⁶, k = 40-এ প্রায় 10⁻¹²; মোট খরচ Θ(k n²).

৫ · Miller–Rabin — বড় সংখ্যাটি মৌলিক কি না

RSA-র মতো cryptography-তে কয়েকশো অঙ্কের মৌলিক সংখ্যা লাগে। √n পর্যন্ত ভাগ করে দেখা (trial division) সেখানে অসম্ভব — 10³⁰⁰-এর বর্গমূল 10¹⁵⁰. দরকার একটি দ্রুত primality test. এখানে শুধু মূল ধারণা; modular exponentiation আর সংখ্যাতত্ত্বের ভিত্তি পর্ব ৬২-এ।

Fermat-এর ধারণা. n মৌলিক হলে, n দিয়ে বিভাজ্য নয় এমন যেকোনো a-র জন্য aⁿ⁻¹ ≡ 1 (mod n) — Fermat's little theorem. তাই কোনো a-তে এটি না মিললে n নিশ্চিত যৌগিক। এটি Fermat test. সমস্যা: কিছু যৌগিক সংখ্যা n-এর সাথে সহমৌলিক প্রতিটি a-তে এই পরীক্ষা পাস করে — Carmichael number. সবচেয়ে ছোটটি 561 = 3 × 11 × 17, আর ১৯৯৪ সালে Alford, Granville ও Pomerance প্রমাণ করেন এমন সংখ্যা অসীম।

Miller–Rabin-এর বাড়তি পরীক্ষা. মৌলিক p-তে x² ≡ 1 (mod p) হলে x ≡ ±1 — এর বাইরে 1-এর কোনো “বর্গমূল” নেই। তাই n − 1 = 2^s · d লিখুন (d বিজোড়), আর ক্রমটি দেখুন: a^d, a^2d, a^4d, …, a^2sd (mod n). n মৌলিক হলে ক্রমটি হয় 1 দিয়ে শুরু, নয়তো শেষ 1-এর ঠিক আগে কোথাও n − 1 আসে। না হলে a হলো n-এর যৌগিকতার “সাক্ষী” (witness).

Miller–Rabin-এর মূল উপপাদ্য (Rabin ১৯৮০, আর স্বাধীনভাবে Monier ১৯৮০): n বিজোড় যৌগিক হলে 2 থেকে n − 2-এর মধ্যে অন্তত 3/4 অংশ a সাক্ষী। তাই এলোমেলো a-তে একটি চেষ্টায় যৌগিক n-কে “মৌলিক” বলার সম্ভাবনা ≤ 1/4, আর k টি স্বাধীন a-তে ≤ 4^−k. খরচ প্রতি চেষ্টায় O(log n) টি modular গুণ (binary exponentiation). আবার one-sided error: “যৌগিক” সবসময় সত্য, “সম্ভবত মৌলিক” কদাচিৎ মিথ্যা।

৬ · Karger-এর min cut — এলোমেলো edge চেপে চেপে

একটি সংযুক্ত, undirected গ্রাফকে দুই ভাগ করতে সবচেয়ে কম কয়টি edge কাটতে হয়? এটি global minimum cut — পর্ব ৪৭-এর s-t min cut-এর মতো নির্দিষ্ট source আর sink নেই। নেটওয়ার্কের সবচেয়ে দুর্বল জায়গা খোঁজা, ছবিকে অংশে ভাগ করা — এমন কাজে লাগে। ১৯৯৩ সালে David Karger একটি আশ্চর্য সরল অ্যালগরিদম দেন — Karger's algorithm:

Pseudocode — one run of Karger's contraction (graph may have parallel edges)

while the graph has more than 2 vertices:
    pick an edge (u, v) uniformly at random
    contract it: merge u and v into one supernode,
                 keep all their other edges (parallel edges stay),
                 delete the edges between u and v (they would be self-loops)
return the number of edges left between the 2 supernodes

Edge contraction মানে দুটি vertex-কে একটি “supernode”-এ জুড়ে দেওয়া। ফলে একই জোড়ার মধ্যে একাধিক edge (parallel edge) তৈরি হয় — গ্রাফটি একটি multigraph. নিচে দুটি ত্রিভুজ, মাঝে একটিমাত্র সেতু-edge c-d; min cut = 1.

সাফল্যের সম্ভাবনা ≥ 2/(n(n − 1)) — প্রমাণের রূপরেখা

একটি নির্দিষ্ট min cut C ধরুন, তার আকার c. প্রতিটি vertex-এর degree ≥ c (নইলে সেই vertex-কে একা আলাদা করাই c-র চেয়ে ছোট cut হতো), তাই মোট edge ≥ nc/2. প্রথম ধাপে C-র একটি edge চাপার সম্ভাবনা ≤ c / (nc/2) = 2/n. Contraction কোনো cut-কে ছোট করে না, তাই i টি supernode থাকা অবস্থাতেও min cut ≥ c আর edge ≥ ic/2; C-র edge চাপার সম্ভাবনা ≤ 2/i. C পুরো চালানোয় বেঁচে থাকলে উত্তর ঠিক C. তাই

গুণফলটি telescoping: লব আর হর প্রায় সব কেটে যায়। n = 6-এ 2/30 = 1/15 — ওপরের গ্রাফে আসল সম্ভাবনা এর চেয়ে বেশি, কারণ সীমাটি সবচেয়ে খারাপ গ্রাফের জন্য। একটি সুন্দর উপজাত: যেকোনো গ্রাফে min cut-এর সংখ্যা ≤ n(n − 1)/2, কারণ প্রতিটি min cut-এর সম্ভাবনা ≥ 2/(n(n − 1)), আর এই সম্ভাবনাগুলোর যোগফল ≤ 1.

পুনরাবৃত্তি. একবারে সাফল্যের সম্ভাবনা ছোট, তাই N = (n(n − 1)/2) · ln n বার স্বাধীনভাবে চালিয়ে সবচেয়ে ছোট উত্তরটি নিন। সব চালানো ব্যর্থ হওয়ার সম্ভাবনা ≤ (1 − 2/(n(n−1)))^N ≤ e^{−ln n} = 1/n (কারণ 1 − x ≤ e^−x). একটি চালানো সরলভাবে O(n²) সময়ে, বা edge-গুলো এলোমেলো ক্রমে সাজিয়ে DSU দিয়ে Kruskal-এর মতো প্রায় O(m) সময়ে করা যায় — তাই মোট O(n⁴ log n) পর্যন্ত। ১৯৯৩ সালেই Karger ও Clifford Stein দেখান, শুরুর ধাপগুলোতে ঝুঁকি কম বলে সেগুলো একবার করে শুধু শেষের দিকে বারবার চালালে মোট O(n² log³ n) — Karger–Stein অ্যালগরিদম।

৭ · ন্যায্য নমুনা আর ন্যায্য shuffle — reservoir sampling আর Fisher–Yates

Reservoir sampling — অজানা দৈর্ঘ্যের stream থেকে k টি

একটি bKash-ধাঁচের সার্ভারে সারাদিন লেনদেন আসছে — একটার পর একটা, মোট কত আসবে আগে জানা নেই, আর সব মেমরিতে রাখা সম্ভব নয়। দিনশেষে ঠিক k টি লেনদেন চাই, যেখানে প্রতিটি লেনদেনের নমুনায় থাকার সম্ভাবনা সমান। এমন ডেটাকে বলে data stream, আর সমাধান reservoir sampling (Algorithm R): প্রথম k টি সরাসরি রাখুন। i-তম উপাদান (i > k) এলে 0 থেকে i − 1-এর মধ্যে একটি এলোমেলো j নিন; j < k হলে reservoir-এর j-তম ঘরে নতুনটি বসান, নইলে বাদ দিন। জায়গা O(k), প্রতি উপাদানে O(1) সময়, ডেটার ওপর মাত্র একবার হাঁটা।

প্রমাণ — প্রত্যেকে k/n সম্ভাবনায় থাকে. একটি উপাদান i (i > k) ধরুন। সে ঢোকে k/i সম্ভাবনায়। পরের প্রতিটি ধাপ t-তে (t = i + 1, …, n) সে বের হয় শুধু যদি নতুনটি ঢোকে (k/t) আর ঠিক তার ঘরটি বাছা হয় (1/k) — সম্ভাবনা 1/t; তাই টিকে থাকে (t − 1)/t সম্ভাবনায়। মুদ্রাগুলো স্বাধীন, তাই

প্রথম k টির জন্য শুরুতে সম্ভাবনা ১, তারপর একই telescoping গুণফল k/(k+1) ⋯ (n−1)/n = k/n. ∎ আরও শক্ত কথাও সত্য: n টির মধ্য থেকে যেকোনো k-সদস্যের সেট সমান সম্ভাবনায় আসে। Jeffrey Vitter ১৯৮৫ সালে এই পদ্ধতির বিশ্লেষণ আর দ্রুততর রূপ প্রকাশ করেন; Knuth Algorithm R-কে Alan Waterman-এর বলে উল্লেখ করেন।

Fisher–Yates shuffle — প্রতিটি permutation সমান সম্ভাবনায়

লুডুর গুটি, তাসের প্যাকেট, অনলাইন পরীক্ষার প্রশ্নের ক্রম বা লটারির তালিকা — এলোমেলো সাজানো চাই, যেখানে n! টি ক্রমের প্রত্যেকটির সম্ভাবনা ঠিক 1/n!. Fisher–Yates shuffle (Durstenfeld-এর ১৯৬৪ সালের কম্পিউটার রূপ): i = n − 1 থেকে 1 পর্যন্ত, 0 থেকে i-এর মধ্যে এলোমেলো j নিয়ে a[i] আর a[j] অদলবদল করুন। ঘর i একবার ঠিক হলে আর কখনো ছোঁয়া হয় না — যেন বাক্স থেকে একটি করে কাগজ তুলে সারির শেষ থেকে বসানো। Θ(n) সময়, O(1) বাড়তি জায়গা।

পরিচিত ভুল — naive shuffle. “প্রতিটি i-তে 0 থেকে n − 1-এর যেকোনো j নিয়ে অদলবদল” দেখতে নির্দোষ, কিন্তু পক্ষপাতদুষ্ট। এর nⁿ টি সমান-সম্ভাবনার পথ। n ≥ 3-এ n − 1 সংখ্যাটি n! কে ভাগ করে কিন্তু nⁿ-এর সাথে সহমৌলিক, তাই n! কখনো nⁿ-কে ভাগ করে না — পথগুলোকে n! টি সমান দলে ভাগ করা অসম্ভব। নিচের প্রোগ্রাম n = 3-এ ২৭টি পথের সবগুলো গুনে দেখায়, তারপর fixed seed-এর একটি PRNG দিয়ে ৬০০০০ বার shuffle করে।

#include <cstdint>
#include <cstdio>
#include <initializer_list>
#include <string>
#include <utility>
#include <vector>
using namespace std;
uint64_t state = 1977;                                  // fixed seed: same numbers on every run
int rnd(int n) {                                        // xorshift64, then a number in 0..n-1
    state ^= state << 13; state ^= state >> 7; state ^= state << 17; return int(state % n); }
template <class F> string naive(F choice) {             // swap position i with ANY position
    string a = "ABC";
    for (int i = 0; i < 3; ++i) swap(a[i], a[choice(i, 3)]);
    return a;
}
template <class F> string fisher_yates(F choice) {      // swap position i with a position 0..i only
    string a = "ABC";
    for (int i = 2; i > 0; --i) swap(a[i], a[choice(i, i + 1)]);
    return a;
}
void table(const string& title, const vector<string>& f1, const vector<string>& f2) {
    printf("%s\nperm   naive  fisher-yates\n", title.c_str());
    for (string p : {"ABC", "ACB", "BAC", "BCA", "CAB", "CBA"}) {
        int c1 = 0, c2 = 0;
        for (auto& s : f1) c1 += s == p;
        for (auto& s : f2) c2 += s == p;
        printf("%s  %6d  %12d\n", p.c_str(), c1, c2);
    }
}
int main() {
    vector<string> nv, fy;                              // 1) exact: every sequence of choices once
    for (int c = 0; c < 27; ++c) nv.push_back(naive([&](int i, int) { return i == 0 ? c % 3 : i == 1 ? c / 3 % 3 : c / 9; }));
    for (int a = 0; a < 3; ++a) for (int b = 0; b < 2; ++b) fy.push_back(fisher_yates([&](int i, int) { return i == 2 ? a : b; }));
    table("all choice sequences: naive 27, fisher-yates 6", nv, fy);
    const int T = 60000; nv.clear(); fy.clear();        // 2) simulation with the fixed-seed generator
    for (int t = 0; t < T; ++t) nv.push_back(naive([](int, int n) { return rnd(n); }));
    for (int t = 0; t < T; ++t) fy.push_back(fisher_yates([](int, int n) { return rnd(n); }));
    table("simulation: " + to_string(T) + " shuffles each (fair = " + to_string(T / 6) + ")", nv, fy);
}

all choice sequences: naive 27, fisher-yates 6
perm   naive  fisher-yates
ABC       4             1
ACB       5             1
BAC       5             1
BCA       5             1
CAB       4             1
CBA       4             1
simulation: 60000 shuffles each (fair = 10000)
perm   naive  fisher-yates
ABC    8870          9927
ACB   11110         10160
BAC   11236          9896
BCA   11030          9851
CAB    8771         10078
CBA    8983         10088

from itertools import product
MASK, state = (1 << 64) - 1, 1977           # fixed seed: same numbers on every run
def rnd(n):                                 # xorshift64, then a number in 0..n-1
    global state
    state ^= (state << 13) & MASK; state ^= state >> 7; state ^= (state << 17) & MASK
    return state % n
def naive(choice):                          # swap position i with ANY position
    a = list("ABC")
    for i in range(3): j = choice(i, 3); a[i], a[j] = a[j], a[i]
    return "".join(a)
def fisher_yates(choice):                   # swap position i with a position 0..i only
    a = list("ABC")
    for i in range(2, 0, -1): j = choice(i, i + 1); a[i], a[j] = a[j], a[i]
    return "".join(a)
PERMS = ["ABC", "ACB", "BAC", "BCA", "CAB", "CBA"]
def table(title, f1, f2):
    print(title); print("perm   naive  fisher-yates")
    for p in PERMS: print(f"{p}  {f1.count(p):6d}  {f2.count(p):12d}")
# 1) exact: every possible sequence of choices, each equally likely
nv = [naive(lambda i, n: c[i]) for c in product(range(3), repeat=3)]
fy = [fisher_yates(lambda i, n: c[2 - i]) for c in product(range(3), range(2))]
table(f"all choice sequences: naive {len(nv)}, fisher-yates {len(fy)}", nv, fy)
# 2) simulation with the fixed-seed generator
T = 60000
nv = [naive(lambda i, n: rnd(n)) for _ in range(T)]
fy = [fisher_yates(lambda i, n: rnd(n)) for _ in range(T)]
table(f"simulation: {T} shuffles each (fair = {T // 6})", nv, fy)

all choice sequences: naive 27, fisher-yates 6
perm   naive  fisher-yates
ABC       4             1
ACB       5             1
BAC       5             1
BCA       5             1
CAB       4             1
CBA       4             1
simulation: 60000 shuffles each (fair = 10000)
perm   naive  fisher-yates
ABC    8870          9927
ACB   11110         10160
BAC   11236          9896
BCA   11030          9851
CAB    8771         10078
CBA    8983         10088

পরীক্ষাগার — Freivalds আর reservoir একসাথে

দ্বিতীয় প্রোগ্রামটি অংশ ৪-এর matrix নিয়ে Freivalds চালায় (প্রথম তিনটি r ডায়াগ্রামের মতোই), তারপর হাজার হাজার বার চালিয়ে ভুল ধরার হার মাপে। শেষে stream 1..10 থেকে k = 3-এর reservoir sampling ৩০০০০ বার চালিয়ে দেখে প্রতিটি সংখ্যা কতবার নমুনায় এসেছে।

#include <array>
#include <cstdio>
#include <vector>
using namespace std;
unsigned long long state = 1977;                        // fixed seed: same numbers on every run
int rnd(int n) {                                        // xorshift64, then a number in 0..n-1
    state ^= state << 13; state ^= state >> 7; state ^= state << 17; return int(state % n); }
using V = array<int, 3>; using M = array<V, 3>;
const M A = {{{1, 2, 0}, {0, 1, 3}, {2, 0, 1}}}, B = {{{1, 0, 2}, {3, 1, 0}, {0, 2, 1}}};
const M good = {{{7, 2, 2}, {3, 7, 3}, {2, 2, 5}}};    // the true product A*B
const M bad  = {{{7, 2, 2}, {3, 7, 4}, {2, 2, 5}}};    // one wrong entry: row 1, column 2
V mul(const M& m, const V& v) { V o{}; for (int i = 0; i < 3; ++i) for (int j = 0; j < 3; ++j) o[i] += m[i][j] * v[j]; return o; }
bool freivalds(const M& C, int k, bool show = false) {  // true = "C looks like A*B"
    for (int t = 1; t <= k; ++t) {
        V r; for (int& x : r) x = rnd(2);               // random 0/1 vector
        bool ok = mul(A, mul(B, r)) == mul(C, r);       // O(n^2) work, no matrix product
        if (show) printf("  trial %d: r = %d %d %d -> %s\n", t, r[0], r[1], r[2], ok ? "equal" : "DIFFER: C != A*B");
        if (!ok) return false;
    }
    return true;
}
vector<int> reservoir(int n, int k) {                   // Algorithm R on the stream 1..n, k slots
    vector<int> res;                                    // the reservoir: at most k items
    for (int i = 1; i <= n; ++i)
        if (i <= k) res.push_back(i);
        else if (int j = rnd(i); j < k) res[j] = i;     // j in 0..i-1, so P(j < k) = k/i
    return res;
}
int main() {
    puts("Freivalds on the wrong C, k = 10:"); freivalds(bad, 10, true);
    int one = 0, ten = 0, rej = 0, cnt[11] = {};
    for (int i = 0; i < 10000; ++i) one += !freivalds(bad, 1);
    for (int i = 0; i < 10000; ++i) ten += freivalds(bad, 10);
    for (int i = 0; i < 1000; ++i) rej += !freivalds(good, 10);
    printf("k = 1 : wrong C caught in %d of 10000 runs (guarantee: at least half)\n", one);
    printf("k = 10: wrong C slipped through %d of 10000 runs (bound: 10000 / 2^10 = 9.8)\n", ten);
    printf("correct C rejected in %d of 1000 runs\n", rej);
    for (int i = 0; i < 30000; ++i) for (int x : reservoir(10, 3)) ++cnt[x];
    printf("reservoir k = 3 from 1..10, 30000 runs (fair = 9000 each):\n "); for (int x = 1; x <= 10; ++x) printf(" %d", cnt[x]); printf("\n");
}

Freivalds on the wrong C, k = 10:
  trial 1: r = 0 1 0 -> equal
  trial 2: r = 0 0 0 -> equal
  trial 3: r = 1 0 1 -> DIFFER: C != A*B
k = 1 : wrong C caught in 5043 of 10000 runs (guarantee: at least half)
k = 10: wrong C slipped through 8 of 10000 runs (bound: 10000 / 2^10 = 9.8)
correct C rejected in 0 of 1000 runs
reservoir k = 3 from 1..10, 30000 runs (fair = 9000 each):
  9098 9102 8995 8948 8942 9011 9109 8922 9024 8849

MASK, state = (1 << 64) - 1, 1977           # fixed seed: same numbers on every run
def rnd(n):                                 # xorshift64, then a number in 0..n-1
    global state
    state ^= (state << 13) & MASK; state ^= state >> 7; state ^= (state << 17) & MASK
    return state % n
A = [[1, 2, 0], [0, 1, 3], [2, 0, 1]]
B = [[1, 0, 2], [3, 1, 0], [0, 2, 1]]
good = [[7, 2, 2], [3, 7, 3], [2, 2, 5]]    # the true product A*B
bad = [[7, 2, 2], [3, 7, 4], [2, 2, 5]]     # one wrong entry: row 1, column 2
def mul(M, v): return [sum(M[i][j] * v[j] for j in range(3)) for i in range(3)]
def freivalds(C, k, show=False):            # True = "C looks like A*B"
    for t in range(1, k + 1):
        r = [rnd(2) for _ in range(3)]      # random 0/1 vector
        ok = mul(A, mul(B, r)) == mul(C, r) # O(n^2) work, no matrix product
        if show: print(f"  trial {t}: r = {r[0]} {r[1]} {r[2]} ->", "equal" if ok else "DIFFER: C != A*B")
        if not ok: return False
    return True
print("Freivalds on the wrong C, k = 10:"); freivalds(bad, 10, True)
one = sum(not freivalds(bad, 1) for _ in range(10000))
ten = sum(freivalds(bad, 10) for _ in range(10000))
print(f"k = 1 : wrong C caught in {one} of 10000 runs (guarantee: at least half)")
print(f"k = 10: wrong C slipped through {ten} of 10000 runs (bound: 10000 / 2^10 = 9.8)")
print(f"correct C rejected in {sum(not freivalds(good, 10) for _ in range(1000))} of 1000 runs")
def reservoir(stream, k):                   # Algorithm R: one pass, k slots
    res = []
    for i, x in enumerate(stream, 1):
        if i <= k: res.append(x)
        else:
            j = rnd(i)                      # 0..i-1, so P(j < k) = k/i
            if j < k: res[j] = x
    return res
cnt = [0] * 11
for _ in range(30000):
    for x in reservoir(range(1, 11), 3): cnt[x] += 1
print("reservoir k = 3 from 1..10, 30000 runs (fair = 9000 each):")
print(" ", *cnt[1:])

Freivalds on the wrong C, k = 10:
  trial 1: r = 0 1 0 -> equal
  trial 2: r = 0 0 0 -> equal
  trial 3: r = 1 0 1 -> DIFFER: C != A*B
k = 1 : wrong C caught in 5043 of 10000 runs (guarantee: at least half)
k = 10: wrong C slipped through 8 of 10000 runs (bound: 10000 / 2^10 = 9.8)
correct C rejected in 0 of 1000 runs
reservoir k = 3 from 1..10, 30000 runs (fair = 9000 each):
  9098 9102 8995 8948 8942 9011 9109 8922 9024 8849

৮ · Hashing, birthday paradox আর এলোমেলো সংখ্যার উৎস

এলোমেলো seed-এর hash — universal hashing

পর্ব ৩২-এ দেখেছেন: যেকোনো একটি নির্দিষ্ট hash function-এর জন্য খারাপ key-এর সেট আছে। সমাধান universal hashing (Carter ও Wegman, ১৯৭৯): প্রোগ্রাম শুরুতে একটি পরিবার থেকে এলোমেলোভাবে hash function বেছে নেয়, যেমন h(x) = ((a·x + b) mod p) mod m, যেখানে a, b এলোমেলো। তখন যেকোনো দুটি আলাদা key-র collision-এর সম্ভাবনা ≤ 1/m, আর indicator + linearity দিয়ে একটি খোঁজার প্রত্যাশিত খরচ O(1 + n/m). Python 3.3 থেকে প্রতিটি প্রক্রিয়ায় string-এর hash এলোমেলো seed-এ চলে; কনটেস্টে unordered_map-এর বিরুদ্ধে anti-hash test এড়াতেও একই কৌশল — এলোমেলো seed-এর custom hash.

Birthday paradox — collision আসে প্রত্যাশার আগেই

একটি ক্লাসে কতজন হলে দুজনের একই জন্মদিনের সম্ভাবনা অর্ধেক ছাড়ায়? উত্তর মাত্র ২৩ (প্রায় ৫০.৭%) — এটি birthday paradox. Indicator দিয়ে সহজ ব্যাখ্যা: n জন, m = 365 দিন। প্রতিটি জোড়ার জন্য একটি indicator, P = 1/m. জোড়া C(n, 2) টি, তাই প্রত্যাশিত মিলে-যাওয়া জোড়া = n(n − 1)/(2m). n = 23-এ 253/365 ≈ 0.69. জোড়ার সংখ্যা n²-এর মতো বাড়ে, তাই মিল আসে প্রায় √m জনেই। Hashing-এ এর মানে: m টি ঘরে মাত্র প্রায় 1.18√m টি key-তেই collision-এর সম্ভাবনা অর্ধেক। 64-bit hash-এ প্রায় ৫০০ কোটি key-তেই collision-এর সম্ভাবনা অর্ধেক ছাড়ায় — তাই কনটেস্টে string hashing-এ (পর্ব ৫৮) দুটি mod বা 64-bit hash নেওয়া হয়, আর একই কারণে cryptographic hash-এর আউটপুট লম্বা রাখা হয়।

PRNG বনাম cryptographic RNG — কোনটি কোথায়

কম্পিউটারের pseudorandom number generator একটি সূত্র: seed থেকে একটি লম্বা, এলোমেলো-দেখতে ধারা। Seed জানলে পুরো ধারা জানা। এই পর্বের xorshift বা C++-এর mt19937 (Mersenne Twister) অ্যালগরিদম আর simulation-এর জন্য চমৎকার — দ্রুত, আর নির্দিষ্ট seed-এ reproducible. কিন্তু নিরাপত্তার জন্য একেবারে অচল: Mersenne Twister-এর মাত্র ৬২৪টি পরপর 32-bit আউটপুট দেখে তার পুরো ভেতরের অবস্থা বের করা যায়, তারপর ভবিষ্যতের প্রতিটি সংখ্যা।

৯ · তুলনা, ইতিহাস আর তত্ত্ব

রূপান্তর. Las Vegas → Monte Carlo: Markov দিয়ে, প্রত্যাশিত সময়ের দ্বিগুণে থামিয়ে (ব্যর্থতা ≤ 1/2). Monte Carlo → Las Vegas: উত্তর দ্রুত যাচাই করা গেলে, সঠিক উত্তর না পাওয়া পর্যন্ত বারবার চালান। প্রতিবার সাফল্যের সম্ভাবনা p হলে প্রত্যাশিত চেষ্টা 1/p — geometric distribution.

১০ · AI-এর যুগে randomized অ্যালগরিদম — এলোমেলো কোড, নির্দিষ্ট যাচাই

Randomized কোডের ভুল সবচেয়ে চতুরভাবে লুকায়। পক্ষপাতদুষ্ট shuffle-ও “এলোমেলো দেখায়”; ভুল Miller–Rabin হাজার সংখ্যায় ঠিক উত্তর দিয়ে একটি বিশেষ যৌগিক সংখ্যায় ভুল করে। AI সহকারী এসব কোড মুহূর্তে লেখে, আর প্রায়ই আত্মবিশ্বাসী কিন্তু অসম্পূর্ণ দাবি করে — “এটি O(n log n)”, “এই base-গুলো যথেষ্ট”.

Randomized কোডে AI সাধারণত যেসব ভুল করে

Randomized কোড যাচাইয়ের চেকলিস্ট

• ছোট n-এ সব সম্ভাব্য পছন্দের ক্রম গুনে (exhaustive) দেখেছি প্রতিটি ফল সমান সম্ভাবনার কি না
• নির্দিষ্ট seed-এ হাজার বার চালিয়ে হার মেপেছি (frequency test), আর কয়েকটি ভিন্ন seed-এও
• Monte Carlo-র উত্তর একটি ধীর কিন্তু নিশ্চিত brute force-এর সাথে এলোমেলো ইনপুটে মিলিয়েছি — stress testing
• Edge case: n = 1, k = n, একই উপাদান বারবার, matrix-এর শূন্য সারি, n = 2 বা 3-এর মতো ছোট মৌলিক সংখ্যা
• প্রতিটি সীমায় লিখেছি “প্রত্যাশিত”, “worst case” না “ভুলের সম্ভাবনা ≤ …”
• নিরাপত্তার কাজে শুধু OS-এর CSPRNG; পরীক্ষায় নির্দিষ্ট seed, প্রোডাকশনে অনুমান-অযোগ্য seed