NP-hard সমস্যা বাস্তবে — চেনা আর মোকাবিলা

১ · “NP-hard” শুনে ভয় নয় — এর আসল মানে কী

পর্ব ৬৪-এ দেখেছেন P, NP, NP-complete আর reduction. সংক্ষেপে: একটি সমস্যা NP-hard মানে NP-এর প্রতিটি সমস্যাকে polynomial সময়ে এতে রূপান্তর করা যায়। তাই এর একটি polynomial সময়ের নির্ভুল অ্যালগরিদম পেলে P = NP প্রমাণ হয়ে যাবে — যা পঞ্চাশ বছরের বেশি সময় ধরে কেউ পারেননি। বাস্তবে কথাটির অর্থ: এই সমস্যার এমন কোনো অ্যালগরিদম জানা নেই, যা সব ইনপুটে, দ্রুত, আর নির্ভুল উত্তর দেয়।

বাক্যটিতে তিনটি দাবি আছে — “সব ইনপুটে”, “দ্রুত” (polynomial সময়ে) আর “নির্ভুল” (প্রমাণিত সর্বোত্তম)। পুরো পর্বের মূল কৌশল একটাই: এই তিনটির অন্তত একটি ছেড়ে দিন।

NP-hard প্রমাণ কেন কাজের

একজন প্রোগ্রামার সপ্তাহের পর সপ্তাহ দ্রুত নির্ভুল অ্যালগরিদম খুঁজে হয়রান হতে পারেন। সমস্যাটি NP-hard প্রমাণ করতে পারলে সেই খোঁজ বন্ধ করে সঠিক প্রশ্নে যাওয়া যায়: “কোন দাবিটি ছাড়ব?” প্রমাণের পদ্ধতি পর্ব ৬৪-এর polynomial-time reduction: একটি জানা NP-hard সমস্যা Y-এর যেকোনো ইনপুটকে আপনার সমস্যা X-এর ইনপুটে বদলে দেখান। দিকটি গুরুত্বপূর্ণ — Y থেকে X-এ। উল্টো দিকে (X থেকে Y-এ) reduction শুধু বলে X বড়জোর Y-এর মতো কঠিন, যা সহজ সমস্যার বেলাতেও সত্য।

২ · ছদ্মবেশে NP-hard — কাজের জায়গায় চেনার উপায়

কাজের জায়গায় কেউ বলবে না “একটি graph coloring সমাধান করো।” বলবে, “আগামী মাসের পরীক্ষার রুটিন বানাও।” নিচের টেবিলে চারটি পরিচিত কাজ আর তাদের ভেতরের ক্লাসিক সমস্যা:

চেনার লক্ষণ — আর চেহারায় একই, অথচ সহজ

কিছু লক্ষণ দেখলেই সতর্ক হন: সমস্যাটি “সব উপসেট / সব ক্রম / সব ভাগাভাগি থেকে সেরাটি” চায়; জোড়ায় জোড়ায় সংঘাত আছে (“এ আর ও একসাথে নয়”); “প্রতিটি জায়গায় ঠিক একবার” যেতে হয়; অথবা ধারণক্ষমতায় পূর্ণসংখ্যক জিনিস ভরতে হয়। কিন্তু সাবধান — ছোট একটি বদলেই কঠিন সমস্যা সহজ হয়ে যায়:

৩ · হাতিয়ার ১ — ছোট ইনপুট: exponential, তবু যথেষ্ট দ্রুত

NP-hard শুনলে অনেকে ভাবেন “exponential, মানে অচল।” কিন্তু প্রতি রাতে ১৫টি দোকানে মাল দেওয়া একটি ভ্যানের জন্য 2¹⁵ · 15² ≈ ৭৪ লক্ষ ধাপ — চোখের পলকে শেষ। প্রশ্ন “exponential কি না” নয়, প্রশ্ন “আমার n-এ কত ধাপ।” মোটামুটি হিসাবে একটি সাধারণ কম্পিউটার সেকেন্ডে ১০^৮-এর কাছাকাছি সরল কাজ করে। সেই হিসাবে:

খেয়াল করুন, ভিত্তি কমানোও বিশাল লাভ। n! থেকে 2ⁿ · n²-এ গেলে n = ২০-এ কাজ প্রায় ২.৪ × ১০^১৮ থেকে নেমে ৪ × ১০^৮; আর 2ⁿ থেকে 2^n/2-এ গেলে একই সময়ে প্রায় দ্বিগুণ বড় n সামলানো যায়। Exponential অ্যালগরিদমের জগতেও “ভালো” আর “খারাপ” আছে।

আরেকটি “ছোট” আছে — সংখ্যার মান ছোট। 0/1 knapsack NP-hard, কিন্তু পর্ব ৫৩-এর DP O(nW) সময়ের, যা ধারণক্ষমতা W ছোট হলে দারুণ চলে। এমন সময়কে বলে pseudo-polynomial. ব্যাংকের টাকার হিসাবে W হয়তো কয়েক কোটি — তখন আর কাজে আসে না; কিন্তু ১০০ কেজির ভ্যানে কেজির হিসাবে W = ১০০, খুব সহজ।

৪ · হাতিয়ার ২ — গঠন কাজে লাগান: P-তে থাকা বিশেষ রূপ

NP-hard প্রমাণ হয় সাধারণ ইনপুটের জন্য। কিন্তু বাস্তবের ইনপুট প্রায়ই সাধারণ নয় — তার একটি বিশেষ আকার থাকে, আর সেই আকারে সমস্যাটি polynomial হয়ে যায়। কয়েকটি গুরুত্বপূর্ণ উদাহরণ:

• Bipartite গ্রাফ। পরীক্ষার হলের সারি-কলামের গ্রিডে (শুধু সামনে-পেছনে আর ডানে-বাঁয়ে পাশাপাশি ধরলে) গ্রাফটি bipartite: দাবার ছকের সাদা ঘরে এক বিভাগ, কালো ঘরে আরেক — দুই রঙেই শেষ। কোণাকুণি পাশাপাশিও গুনলে চারটি বিভাগ দিয়ে ২ × ২ নকশা পুনরাবৃত্তি করলেই চলে। পরিদর্শক আর হলের জোড় বানানো (bipartite maximum matching, পর্ব ৪৭) P-তে; আর König-এর উপপাদ্য বলে bipartite গ্রাফে সবচেয়ে ছোট vertex cover-এর আকার = maximum matching-এর আকার — তাই সেখানে vertex cover-ও সহজ।
• 2-SAT. প্রতিটি শর্তে বড়জোর দুটি চলক থাকলে (“রহিম আর করিম একই দলে নয়”) implication graph-এর SCC দিয়ে রৈখিক সময়ে সমাধান (পর্ব ৪৬). তিনটি হলেই NP-complete.
• Interval graph. প্রতিটি vertex যদি সময়ের একটি ব্যবধান হয় (ক্লাস ৯টা–১০:৩০), আর মিলে যাওয়া ব্যবধানে edge — তাহলে এটি interval graph. এখানে coloring মানে পর্ব ৫০-এর interval partitioning: শুরুর সময় ধরে সাজিয়ে greedy, আর রঙের সংখ্যা = একসাথে সর্বোচ্চ কয়টি ব্যবধান মেলে — যা একটি lower bound-ও, তাই সর্বোত্তম। শুধু সময়ের সংঘাত থাকলে কক্ষ বণ্টন তাই সহজ; বিষয়ভিত্তিক সংঘাত (একই শিক্ষার্থী) যোগ হলেই সাধারণ coloring.
• Tree. সাধারণ গ্রাফে minimum vertex cover আর maximum independent set NP-hard, কিন্তু গ্রাফটি tree হলে পর্ব ৫৫-এর tree DP দিয়ে O(n) সময়ে। নিচের স্টেপারে দেখুন।

Tree-তে vertex cover-এর DP: প্রতিটি node v-এর জন্য দুটি সংখ্যা — in[v] = v-কে cover-এ নিলে v-এর subtree-তে সবচেয়ে ছোট cover; out[v] = v-কে না নিলে। v নেওয়া হলে সন্তানদের স্বাধীনতা আছে: in[v] = 1 + Σ min(in[c], out[c]). v না নেওয়া হলে প্রতিটি edge v–c ঢাকতে প্রতিটি সন্তান c-কে নিতেই হবে: out[v] = Σ in[c]. Leaf-এ in = 1, out = 0.

        1
      /   \
     2     3
    / \     \
   4   5     6

কেন ঠিক: subtree-র উচ্চতা ধরে induction. v-এর subtree-র যেকোনো cover-এ হয় v আছে, নয়তো নেই। থাকলে প্রতিটি সন্তানের subtree আলাদাভাবে ঢাকতে হয়, আর প্রতিটির সেরা খরচ min(in[c], out[c]); না থাকলে edge v–c-এর জন্য c থাকতেই হবে, সেরা খরচ in[c]. এই দুই ক্ষেত্রের বাইরে কিছু নেই। জটিলতা: প্রতিটি node একবার, প্রতিটি edge একবার — সময় O(n), জায়গা O(n).

এই ধারণার বড় রূপ treewidth: গ্রাফটি গাছের কতটা কাছাকাছি, তার একটি মাপ। Treewidth w হলে বহু NP-hard গ্রাফ সমস্যা প্রায় O(c^w · n) সময়ে সমাধান হয় — গাছের মতো ছোট “ব্যাগ” ধরে DP. রাস্তার নেটওয়ার্ক, সার্কিট আর কম্পাইলারের control-flow গ্রাফে treewidth প্রায়ই ছোট, তাই এই কৌশল বাস্তবেও চলে।

৫ · Fixed-parameter tractability — vertex cover খোঁজের গাছে

ঢাকার একটি এলাকায় প্রতিটি রাস্তায় নজর রাখতে মোড়ে মোড়ে CCTV বসবে। একটি মোড়ের ক্যামেরা সেই মোড়ে মেলা সব রাস্তা দেখে। বাজেটে বড়জোর k টি ক্যামেরা। প্রশ্ন: “k টি মোড় বেছে কি প্রতিটি রাস্তা ঢাকা যায়?” এটাই vertex cover-এর decision রূপ — NP-complete. এলাকায় হয়তো n = ১০০০ মোড়, কিন্তু k = ১০. সব উপসেট 2¹⁰⁰⁰, অসম্ভব। কিন্তু খোঁজের exponential অংশটি যদি শুধু k-তে আটকে রাখা যায়?

এটাই fixed-parameter tractability (FPT): সময় f(k) · n^c, যেখানে k একটি ছোট প্যারামিটার আর c ধ্রুবক। Vertex cover-এর জন্য একটি সুন্দর পর্যবেক্ষণই যথেষ্ট — যেকোনো edge (u, v) ঢাকতে u বা v-এর অন্তত একটি cover-এ থাকতেই হবে। তাই দুই শাখা: u নিন, নয়তো v নিন। প্রতিটি শাখায় বাজেট k এক কমে। এই কৌশলের নাম bounded search tree.

Pseudocode — does G have a vertex cover of size ≤ k?
VC(G, k):
    if G has no edge: return YES
    if k = 0: return NO
    pick any edge (u, v)
    return VC(G − u, k − 1) or VC(G − v, k − 1)   // G − x: delete x and its edges

নিচের গ্রাফে ৫টি মোড় আর ৫টি রাস্তা: 1–2, 2–3, 3–1, 3–4, 4–5. প্রথমে k = ২ দিয়ে, তারপর k = ৩ দিয়ে খোঁজ চলবে। প্রতিবার বাকি রাস্তাগুলোর প্রথমটি বাছা হয়।

সঠিকতা: যদি k আকারের একটি cover C থাকে, তাহলে বাছা edge (u, v)-এর জন্য u ∈ C বা v ∈ C; ধরুন u ∈ C. তাহলে C − {u} গ্রাফ G − u-এর একটি cover, আকার k − 1 — অর্থাৎ ওই শাখা হ্যাঁ বলবে (k-এর ওপর induction). উল্টো দিকে, কোনো শাখা হ্যাঁ বললে সেই শাখার cover-এ বাছা vertex যোগ করলেই G-এর cover. জটিলতা: ≤ 2^k+1 node, প্রতিটিতে edge তালিকা ছেঁকে নিতে O(m). আর শুরুতেই একটি সহজ পরীক্ষা: k টি vertex বড়জোর k(n − 1) টি edge ঢাকতে পারে, তাই m > k(n − 1) হলে সরাসরি “না”. ফলে m = O(kn), আর মোট সময় O(2^k · k · n) — n-এ রৈখিক। জায়গা: রিকার্শনের গভীরতা k, প্রতিটি স্তরে একটি edge তালিকা — O(k · m).

Kernel — আগে ইনপুট ছোট করুন

আরও একটি চমৎকার কৌশল kernelization. Buss-এর নিয়ম: যে মোড়ে k-এর বেশি রাস্তা মিলেছে, তাকে cover-এ নিতেই হবে — না নিলে তার সব (k-এর বেশি) প্রতিবেশীকে নিতে হতো, বাজেট ছাড়িয়ে যেত। এমন সব vertex নিয়ে নেওয়ার পর প্রতিটি vertex-এর degree ≤ k, তাই বাকি k টি vertex বড়জোর k² টি edge ঢাকতে পারে; k²-এর বেশি edge বাকি থাকলে উত্তর “না”. যা বাকি থাকে তা একটি kernel — আকার শুধু k-এর ওপর নির্ভর করে, n-এর ওপর নয়। Kernel-এর ওপর খোঁজের গাছ চালালে মোট সময় O(n + m + 2^k · k²).

কোডে — FPT খোঁজের গাছ বনাম brute force

১২টি মোড় আর ১৬টি রাস্তার একটি নির্দিষ্ট গ্রাফে প্রোগ্রামটি k = ০, ১, ২, … ধরে খোঁজের গাছ চালায়, প্রতিটি k-এ কতটি node লাগল ছাপে, আর শেষে সব 2¹² = ৪০৯৬টি উপসেট দেখে brute force দিয়ে উত্তর মেলায়।

#include <algorithm>
#include <bit>
#include <cstdio>
#include <utility>
#include <vector>
using namespace std;
using Edges = vector<pair<int, int>>;
const int N = 12;                                   // 12 junctions, 16 roads
const Edges E = {{0, 1}, {0, 2}, {1, 3}, {2, 3}, {3, 4}, {4, 5}, {4, 6}, {5, 7},
                 {6, 7}, {7, 8}, {8, 9}, {2, 9}, {1, 5}, {9, 10}, {10, 11}, {8, 11}};
long long nodes = 0;
vector<int> cover;
bool vc(const Edges& es, int k) {                   // is there a vertex cover of size <= k?
    nodes++;
    if (es.empty()) return true;                    // every edge covered
    if (k == 0) return false;                       // edges left, budget used up
    for (int x : {es[0].first, es[0].second}) {     // edge (u, v): u or v MUST be in the cover
        Edges rest;
        for (auto [a, b] : es) if (a != x && b != x) rest.push_back({a, b});
        cover.push_back(x);
        if (vc(rest, k - 1)) return true;
        cover.pop_back();
    }
    return false;
}
int main() {
    long long total = 0;
    for (int k = 0; ; k++) {                        // try k = 0, 1, 2, ... until YES
        nodes = 0;
        bool ok = vc(E, k);
        total += nodes;
        printf("FPT k=%d: %s, nodes=%lld (limit %d)\n", k, ok ? "YES" : "NO ", nodes, (1 << (k + 1)) - 1);
        if (ok) break;
    }
    printf("cover:");
    for (int x : cover) printf(" %d", x);
    printf("\nFPT total nodes = %lld\n", total);
    int best = N;                                   // brute force: all 2^N subsets
    for (unsigned m = 0; m < (1u << N); m++)
        if (all_of(E.begin(), E.end(), [m](auto e) { return (m >> e.first & 1) || (m >> e.second & 1); }))
            best = min(best, popcount(m));
    printf("brute force: %d subsets checked, minimum cover size = %d\n", 1 << N, best);
}

FPT k=0: NO , nodes=1 (limit 1)
FPT k=1: NO , nodes=3 (limit 3)
FPT k=2: NO , nodes=7 (limit 7)
FPT k=3: NO , nodes=15 (limit 15)
FPT k=4: NO , nodes=31 (limit 31)
FPT k=5: NO , nodes=63 (limit 63)
FPT k=6: YES, nodes=124 (limit 127)
cover: 1 2 4 7 8 10
FPT total nodes = 244
brute force: 4096 subsets checked, minimum cover size = 6

N = 12                                            # 12 junctions, 16 roads
E = [(0, 1), (0, 2), (1, 3), (2, 3), (3, 4), (4, 5), (4, 6), (5, 7),
     (6, 7), (7, 8), (8, 9), (2, 9), (1, 5), (9, 10), (10, 11), (8, 11)]
nodes, cover = 0, []

def vc(es, k):                    # is there a vertex cover of size <= k?
    global nodes
    nodes += 1
    if not es: return True        # every edge covered
    if k == 0: return False       # edges left, budget used up
    for x in es[0]:               # edge (u, v): u or v MUST be in the cover
        cover.append(x)
        if vc([(a, b) for a, b in es if a != x and b != x], k - 1): return True
        cover.pop()
    return False

total, k = 0, 0
while True:                       # try k = 0, 1, 2, ... until YES
    nodes = 0
    ok = vc(E, k)
    total += nodes
    print(f"FPT k={k}: {'YES' if ok else 'NO '}, nodes={nodes} (limit {2 ** (k + 1) - 1})")
    if ok: break
    k += 1
print("cover:", *cover)
print(f"FPT total nodes = {total}")
best = min(bin(m).count("1") for m in range(1 << N)       # brute force: all 2^N subsets
           if all(m >> a & 1 or m >> b & 1 for a, b in E))
print(f"brute force: {1 << N} subsets checked, minimum cover size = {best}")

FPT k=0: NO , nodes=1 (limit 1)
FPT k=1: NO , nodes=3 (limit 3)
FPT k=2: NO , nodes=7 (limit 7)
FPT k=3: NO , nodes=15 (limit 15)
FPT k=4: NO , nodes=31 (limit 31)
FPT k=5: NO , nodes=63 (limit 63)
FPT k=6: YES, nodes=124 (limit 127)
cover: 1 2 4 7 8 10
FPT total nodes = 244
brute force: 4096 subsets checked, minimum cover size = 6

৬ · হাতিয়ার ৩ — নিশ্চয়তাসহ approximation

ইনপুট বড়, গঠনও সাধারণ — তখন “নির্ভুল” দাবিটি একটু ছাড়ুন, কিন্তু নিয়ন্ত্রিতভাবে। একটি approximation algorithm polynomial সময়ে চলে আর প্রমাণ দেয়: উত্তর কখনো সর্বোত্তমের নির্দিষ্ট গুণের বেশি খারাপ নয়। পর্ব ৬৬ পুরোটাই এর জন্য; এখানে শুধু স্বাদ নিন।

• Vertex cover, ২ গুণ: যেকোনো না-ঢাকা edge নিন, তার দুই প্রান্তই cover-এ দিন, তাদের ছোঁয়া সব edge মুছুন, আবার করুন। বাছা edge-গুলো কোনো vertex ভাগ করে না (একটি matching), আর যেকোনো cover-কে প্রতিটি থেকে অন্তত একটি প্রান্ত নিতে হয় — তাই সর্বোত্তম ≥ বাছা edge-এর সংখ্যা, আর আমাদের উত্তর ঠিক তার দ্বিগুণ।
• Metric TSP, ১.৫ গুণ: দূরত্ব যদি triangle inequality মানে (সরাসরি যাওয়া কখনো ঘুরে যাওয়ার চেয়ে লম্বা নয়), Christofides-এর ১৯৭৬ সালের অ্যালগরিদম MST আর matching মিলিয়ে সর্বোত্তমের বড়জোর ১.৫ গুণ দেয়। কিন্তু দূরত্ব যেমন খুশি হলে কোনো ধ্রুবক গুণের নিশ্চয়তাই সম্ভব নয়, যদি না P = NP.
• Bin packing: বড় থেকে ছোট সাজিয়ে প্রতিটি জিনিস প্রথম যে ট্রাকে ধরে সেখানে দিন (First Fit Decreasing) — ট্রাকের সংখ্যা মোটামুটি সর্বোত্তমের ১১/৯ গুণের বেশি নয়, সাথে একটি ছোট ধ্রুবক।

৭ · হাতিয়ার ৪ — heuristic আর local search: 2-opt

বাস্তবের বেশিরভাগ বড় রুটিং আর সময়সূচির সফটওয়্যার শেষ পর্যন্ত heuristic চালায়: দ্রুত, প্রায়ই চমৎকার, কিন্তু কোনো প্রমাণিত নিশ্চয়তা নেই। কাঠামোটি সাধারণত দুই ধাপের:

1. Greedy construction: একটি সম্পূর্ণ সমাধান দ্রুত বানানো — যেমন TSP-এর nearest neighbour (সবসময় সবচেয়ে কাছের অদেখা দোকানে যাওয়া), বা coloring-এ সবচেয়ে বেশি সংঘাতের বিষয়কে আগে রং দেওয়া।
2. Improvement: local search — সমাধানটিতে একটি ছোট পরিবর্তন (move) করে দেখা, ভালো হলে রাখা, যতক্ষণ উন্নতি হয়।

TSP-এর সবচেয়ে বিখ্যাত move হলো 2-opt: tour থেকে দুটি edge a–b আর c–e সরিয়ে a–c আর b–e জুড়ুন। এতে b থেকে c পর্যন্ত অংশটি উল্টো দিকে চলে। লাভ = d(a,b) + d(c,e) − d(a,c) − d(b,e); ধনাত্মক হলেই move নিন। সমতলে দুটি edge পরস্পর কাটলে (crossing) এই move সবসময় লাভজনক — triangle inequality থেকে আসে, কারণ কাটা দুই কর্ণের যোগফল সামনের দুই বাহুর যোগফলের চেয়ে বড়। তাই 2-opt শেষে tour-এ কোনো crossing থাকে না। নিচে ৮টি দোকানের উদাহরণ (স্থানাঙ্ক কিলোমিটারে, দোকান ০ থেকে শুরু):

জটিলতা: একটি পূর্ণ পাসে O(n²) জোড়া (i, j) পরীক্ষা, প্রতিটির লাভ O(1); move নিলে অংশ উল্টাতে O(n). কয়টি পাস লাগবে তার ভালো সীমা নেই — কিছু কৃত্রিম ইনপুটে 2-opt exponential সংখ্যক move নিতে পারে, যদিও বাস্তবে n-এর ছোট গুণিতক পাসেই থামে। বড় n-এ প্রতিটি শহরের শুধু কাছের ১০–২০টি প্রতিবেশী (neighbour list) দেখে কাজ অনেক কমানো হয়। জায়গা O(n).

কোডে — nearest neighbour + 2-opt, তারপর brute force দিয়ে যাচাই

প্রোগ্রামটি ওপরের ৮টি দোকানে nearest neighbour বানায়, first-improvement 2-opt চালায় (প্রথম যে লাভজনক move পায়, সেটিই নেয়), প্রতিটি move ছাপে, আর শেষে দোকান ০ স্থির রেখে বাকি ৭টির সব ৭! = ৫০৪০টি ক্রম দেখে সর্বোত্তম বের করে।

#include <algorithm>
#include <cmath>
#include <cstdio>
#include <numeric>
#include <utility>
#include <vector>
using namespace std;
const vector<pair<int, int>> P = {{12, 3}, {16, 3}, {18, 8}, {5, 10}, {7, 8}, {9, 4}, {17, 14}, {4, 8}};
const int n = P.size();                             // 8 shops (x, y) in km
double d(int a, int b) {
    double dx = P[a].first - P[b].first, dy = P[a].second - P[b].second;
    return sqrt(dx * dx + dy * dy);
}
double length(const vector<int>& t) {
    double s = 0.0; for (int i = 0; i < n; i++) s += d(t[i], t[(i + 1) % n]);
    return s;
}
void show(const char* name, const vector<int>& t) {
    printf("%-9s", name); for (int x : t) printf("%d-", x);
    printf("0  length %.2f\n", length(t));
}
int main() {
    vector<int> t = {0};                            // nearest neighbour from shop 0
    while ((int)t.size() < n) {
        int c = -1;
        for (int j = 0; j < n; j++)
            if (find(t.begin(), t.end(), j) == t.end() && (c < 0 || d(t.back(), j) < d(t.back(), c))) c = j;
        t.push_back(c);
    }
    show("NN", t);
    for (bool improved = true; improved; ) {        // 2-opt, first improvement
        improved = false;
        for (int i = 0; i < n - 1; i++) for (int j = i + 2; j < n; j++) {
            int a = t[i], b = t[i + 1], c = t[j], e = t[(j + 1) % n];
            double gain = d(a, b) + d(c, e) - d(a, c) - d(b, e);
            if (a == e || gain <= 1e-9) continue;   // replace a-b, c-e by a-c, b-e
            reverse(t.begin() + i + 1, t.begin() + j + 1); improved = true;
            printf("  2-opt: drop %d-%d,%d-%d add %d-%d,%d-%d gain %.2f\n", a, b, c, e, a, c, b, e, gain);
        }
    }
    show("2-opt", t);
    vector<int> p(n), best; iota(p.begin(), p.end(), 0);
    double bl = INFINITY; long long cnt = 0;
    do { cnt++; if (length(p) < bl) { bl = length(p); best = p; } }   // exact: all (n-1)! tours
    while (next_permutation(p.begin() + 1, p.end()));
    show("optimum", best);
    printf("checked %lld tours; 2-opt gap = %.2f%%\n", cnt, 100 * (length(t) - bl) / bl);
}

NN       0-5-4-3-7-1-2-6-0  length 49.25
  2-opt: drop 5-4,7-1 add 5-7,4-1 gain 0.77
  2-opt: drop 4-1,6-0 add 4-6,1-0 gain 6.72
2-opt    0-5-7-3-4-6-2-1-0  length 41.76
optimum  0-1-2-6-3-7-4-5-0  length 40.99
checked 5040 tours; 2-opt gap = 1.88%

import itertools, math
P = [(12, 3), (16, 3), (18, 8), (5, 10), (7, 8), (9, 4), (17, 14), (4, 8)]   # 8 shops (x, y) in km
n = len(P)

def d(a, b):
    dx, dy = P[a][0] - P[b][0], P[a][1] - P[b][1]
    return math.sqrt(dx * dx + dy * dy)

def length(t):
    s = 0.0
    for i in range(n): s += d(t[i], t[(i + 1) % n])
    return s

def show(name, t): print(f"{name:<9}" + "-".join(map(str, t)) + f"-0  length {length(t):.2f}")

t = [0]                                           # nearest neighbour from shop 0
while len(t) < n:
    t.append(min((j for j in range(n) if j not in t), key=lambda j: d(t[-1], j)))
show("NN", t)
improved = True
while improved:                                   # 2-opt, first improvement
    improved = False
    for i in range(n - 1):
        for j in range(i + 2, n):
            a, b, c, e = t[i], t[i + 1], t[j], t[(j + 1) % n]
            if a == e: continue
            gain = d(a, b) + d(c, e) - d(a, c) - d(b, e)
            if gain > 1e-9:                       # replace a-b, c-e by a-c, b-e
                t[i + 1:j + 1] = t[i + 1:j + 1][::-1]; improved = True
                print(f"  2-opt: drop {a}-{b},{c}-{e} add {a}-{c},{b}-{e} gain {gain:.2f}")
show("2-opt", t)
best, bl, cnt = None, math.inf, 0
for p in itertools.permutations(range(1, n)):     # exact: all (n-1)! tours
    cnt += 1; tour = [0, *p]; l = length(tour)
    if l < bl: bl, best = l, tour
show("optimum", best)
print(f"checked {cnt} tours; 2-opt gap = {100 * (length(t) - bl) / bl:.2f}%")

NN       0-5-4-3-7-1-2-6-0  length 49.25
  2-opt: drop 5-4,7-1 add 5-7,4-1 gain 0.77
  2-opt: drop 4-1,6-0 add 4-6,1-0 gain 6.72
2-opt    0-5-7-3-4-6-2-1-0  length 41.76
optimum  0-1-2-6-3-7-4-5-0  length 40.99
checked 5040 tours; 2-opt gap = 1.88%

৮ · Hill climbing থেকে simulated annealing, tabu আর genetic

ওপরের 2-opt আসলে hill climbing: শুধু ভালো move নিন, আর না থাকলে থামুন। কুয়াশার মধ্যে পাহাড়ে ওঠার মতো — সবসময় ওপরের দিকে পা ফেলুন, আর যেখানে চারদিকেই নিচে, সেখানে থামুন। সমস্যা: সেটি হয়তো একটি ছোট টিলা, পাশেই বড় পাহাড়। Local optimum থেকে বেরোনোর কৌশলগুলোকে একসাথে বলে metaheuristic.

Simulated annealing-এর বুদ্ধি ধাতু গলিয়ে ধীরে ঠান্ডা করার (annealing) প্রক্রিয়া থেকে: গরম অবস্থায় পরমাণু ছোটাছুটি করে, ধীরে ঠান্ডা হলে সুশৃঙ্খল কম-শক্তির কাঠামোয় বসে। অ্যালগরিদমে একটি “তাপমাত্রা” T থাকে। একটি এলোমেলো প্রতিবেশী নিন; খরচ বাড়ার পরিমাণ Δ. Δ ≤ 0 (উন্নতি) হলে নিন; খারাপ হলেও (Δ > 0) নিন e^−Δ/T সম্ভাবনায়। T ধীরে কমে — এই নিয়মের নাম cooling schedule.

Pseudocode — simulated annealing (minimisation)
S ← initial solution (e.g. nearest neighbour);  best ← S;  T ← T0
repeat until T < Tmin or time is up:
    S' ← random neighbour of S                 // e.g. a random 2-opt move
    Δ ← cost(S') − cost(S)
    if Δ ≤ 0 or random() < exp(−Δ / T): S ← S'
    if cost(S) < cost(best): best ← S
    T ← α · T                                  // geometric cooling, α ≈ 0.95 … 0.999
return best

নিচের স্টেপারে একটি geometric schedule: T = ১৬ থেকে শুরু, প্রতিবার অর্ধেক (α = ০.৫ — বোঝানোর জন্য খুব দ্রুত; বাস্তবে α ০.৯৫-এর ওপরে)। দুটি খারাপ move ধরুন: সামান্য খারাপ (Δ = +১ km) আর বেশ খারাপ (Δ = +৪ km).

একই পরিবারের আরও দুই সদস্য। Tabu search প্রতিবার সেরা প্রতিবেশীতে যায় — খারাপ হলেও — কিন্তু সাম্প্রতিক move-গুলো একটি “নিষিদ্ধ তালিকা”য় রাখে, যাতে সাথে সাথে ফিরে না আসে; এভাবে টিলা থেকে নেমে অন্য পাহাড়ে ওঠে। Genetic algorithm একসাথে অনেক সমাধানের “জনসংখ্যা” রাখে, ভালোগুলোর অংশ মিলিয়ে (crossover) আর একটু বদলে (mutation) নতুন প্রজন্ম বানায়। সবগুলোর সাধারণ সত্য: কোনো নিশ্চয়তা নেই — কতক্ষণ চালালে কতটা ভালো, তার প্রমাণ নেই, আর টিউনিং (T0, α, তালিকার দৈর্ঘ্য, জনসংখ্যা) ইনপুট ধরে বদলায়। এলোমেলোতা লাগে বলে একটি নির্দিষ্ট seed-এর PRNG ব্যবহার করুন, যাতে ফল পুনরায় পাওয়া যায় আর বাগ ধরা যায়।

৯ · হাতিয়ার ৫ — solver: SAT, ILP আর constraint programming

অনেক সময় সবচেয়ে ভালো অ্যালগরিদম হলো নিজে অ্যালগরিদম না লেখা। কয়েক দশক ধরে হাজারো গবেষকের পরিশ্রমে তৈরি তিন ধরনের solver আছে; আপনার কাজ শুধু সমস্যাটি তাদের ভাষায় লেখা, যাকে বলে formulation.

সমস্যা লেখা — পরীক্ষার রুটিন একটি integer program হিসেবে

Pseudocode — exam timetabling (graph colouring) as an integer program
variables:  x[e][s] ∈ {0,1}   exam e is held in slot s
            y[s]    ∈ {0,1}   slot s is used
minimise    Σ_s y[s]                                   // fewest slots
subject to  Σ_s x[e][s] = 1          for every exam e  // each exam exactly once
            x[e][s] ≤ y[s]           for every e, s    // a used slot is counted
            x[e][s] + x[f][s] ≤ 1    for every clash (e, f), every slot s

শেষ শর্তটিই coloring-এর নিয়ম: একই শিক্ষার্থী আছে এমন দুই বিষয় e আর f একই স্লটে নয়। একই সমস্যা SAT-এ লেখা যায় clause (¬x_e,s ∨ ¬x_f,s) দিয়ে, আর CP-তে “প্রতিটি বিষয়ের স্লট একটি চলক, সংঘাতী জোড়ায় মান আলাদা” লিখে। নতুন শর্ত যোগ করা খুব সহজ — “কোনো শিক্ষার্থীর একদিনে দুটির বেশি পরীক্ষা নয়” মানে আরও কয়েক সারি linear শর্ত — যা নিজে লেখা heuristic-এ পুরো কোড বদলে দিতে পারে।

Benchmark আর “যথেষ্ট ভালো”

পদ্ধতি বাছাইয়ের আগে মাপুন। তুলনার জন্য প্রচলিত benchmark আছে — TSP-এর TSPLIB, integer programming-এর MIPLIB, আর প্রতি বছরের SAT Competition. নিজের প্রতিষ্ঠানের আসল ডেটার কয়েকটি নমুনা রাখুন, আর প্রতিটি পদ্ধতির সময় আর gap একটি টেবিলে লিখুন। প্রায়ই দেখবেন: শেষ ১%-এর জন্য সময় দশ গুণ বাড়ে, অথচ ইনপুট নিজেই অনিশ্চিত — ঢাকার রাস্তায় কাল বিকেলের যাত্রার সময় ১০–২০% এদিক-ওদিক হতেই পারে। ১% gap-এর রুট তখন বাস্তবে “সর্বোত্তম”-এর চেয়ে আলাদা কিছু নয়।

১০ · AI-এর যুগে NP-hard সমস্যা — দাবি যাচাই করুন

“আমার ডেলিভারি রুটের সেরা সমাধান দাও” — AI আত্মবিশ্বাসী সুরে একটি greedy কোড দেবে আর বলবে “এটি সর্বোত্তম রুট দেয়, O(n²) সময়ে।” এই এক বাক্যে দুটি দাবি, আর দুটিই একসাথে সত্য হলে P = NP প্রমাণ হয়ে যায়। NP-hard সমস্যায় AI-এর উত্তর তাই বিশেষভাবে যাচাই করতে হয়।

এই বিষয়ে AI সাধারণত যেসব ভুল করে

NP-hard সমস্যার সমাধান যাচাইয়ের চেকলিস্ট

• পদ্ধতিটি exact, approximation না heuristic — লিখে রেখেছি, আর “optimal” শব্দটি শুধু exact-এর জন্য ব্যবহার করেছি
• n ≤ ৮-১০-এর শত শত এলোমেলো ইনপুটে brute force-এর সাথে মিলিয়েছি — stress testing (exact হলে উত্তর সমান, heuristic হলে কখনো কম নয়)
• একটি ছোট ইনপুটে প্রথম কয়েকটি ধাপ (খোঁজের গাছের শাখা, 2-opt move) হাতে trace করেছি
• Edge case: কোনো edge নেই, k = ০, সব বিন্দু এক রেখায়, সমান দূরত্ব, n = ১ বা ২
• এলোমেলোতা থাকলে নির্দিষ্ট seed দিয়েছি, আর ফলাফল কয়েকটি seed-এ মিলিয়ে দেখেছি
• বড় ইনপুটে একটি lower bound দিয়ে gap মেপেছি, আর সময়সীমা দিয়ে benchmark করেছি