Floating point — IEEE 754, normalized, denormal আর বিশেষ মান

১ · কেন floating point — বিশাল পরিসর আর বৈজ্ঞানিক নোটেশন

রসায়ন বইয়ে Avogadro সংখ্যা লেখা থাকে ৬.০২ × ১০^২৩, আর একটি ইলেকট্রনের ভর প্রায় ৯.১১ × ১০^−৩১ কেজি। এই দুটি সংখ্যা পুরো অঙ্কে লিখতে প্রথমটিতে ২৪ অঙ্ক আর দ্বিতীয়টিতে দশমিকের পরে ৩০টি শূন্য লাগত। বিজ্ঞানীরা তাই বৈজ্ঞানিক নোটেশন ব্যবহার করেন: কয়েকটি অর্থপূর্ণ অঙ্ক (৬.০২) আর একটি ঘাত (১০^২৩)। অর্থপূর্ণ অঙ্কগুলো বলে “কত নিখুঁত”, আর ঘাত বলে “কত বড়”।

পর্ব ০৫-এ দেখা integer এই কাজ পারে না। ৩২ বিটের unsigned integer-এর পরিসর ০ থেকে প্রায় ৪.৩ × ১০^৯, আর ভগ্নাংশ একেবারেই নেই। অথচ পদার্থবিজ্ঞানের simulation, ছবির রং, শব্দের তরঙ্গ, GPS-এর অবস্থান আর AI মডেলের ওজন — সবখানে দরকার ভগ্নাংশ আর বিশাল পরিসর।

প্রথম উপায়: fixed point

সবচেয়ে সরল উপায় হলো দশমিক বিন্দুকে একটি নির্দিষ্ট জায়গায় “আঁকা” ধরে নেওয়া। একে বলে fixed point. দোকানের হিসাবখাতায় যেমন সব দাম পয়সায় লেখা হয় আর মনে মনে শেষ দুই অঙ্কের আগে বিন্দু বসানো হয় — ১০৫০৭৫ মানে ৳১০৫০.৭৫। Binary-তে একই কাজ করে Q format: Q8.8 মানে ১৬ বিটের signed integer, যার নিচের ৮ বিট ভগ্নাংশ। এর পরিসর −১২৮ থেকে প্রায় ১২৭.৯৯৬, আর প্রতিটি ধাপ ঠিক ১/২৫৬ = ০.০০৩৯০৬২৫। যোগ-বিয়োগ সাধারণ integer-এর মতোই, তাই খুব সস্তা circuit-এ চলে — DSP chip আর ছোট microcontroller-এ fixed point আজও জনপ্রিয়। Q format-এ গুণ করলে বিন্দু কোথায় সরে আর কীভাবে আবার shift করতে হয়, তা পর্ব ২১-এর গুণের আলোচনায় পাবেন।

কিন্তু fixed point-এর সমস্যা হলো ধাপ সব জায়গায় সমান। Q8.8-এ ১০০-এর কাছেও ধাপ ১/২৫৬, আর ০.০০১-এর কাছেও ধাপ ১/২৫৬ — অর্থাৎ ০.০০১ রাখতে গেলে প্রায় পুরোটাই হারিয়ে যায় (সবচেয়ে কাছের মান ০ বা ০.০০৩৯)। ছোট সংখ্যার জন্য আপেক্ষিক ভুল বিশাল, আর বড় সংখ্যা (৩০০) আঁটেই না।

দ্বিতীয় উপায়: বিন্দুকে ভাসতে দেওয়া

Floating point-এর বুদ্ধি হলো, বিটগুলোকে দুই ভাগে ভাগ করা: একটি অংশ রাখে অর্থপূর্ণ অঙ্ক (significand, পুরোনো বইয়ে mantissa), আরেকটি রাখে ঘাত (exponent). Exponent বদলালে বিন্দু ডানে-বাঁয়ে সরে — “ভাসে”। ফলে ছোট সংখ্যার কাছে ধাপ ছোট, বড় সংখ্যার কাছে ধাপ বড়, কিন্তু আপেক্ষিক নির্ভুলতা প্রায় সব জায়গায় সমান।

২ · Binary বৈজ্ঞানিক নোটেশন, hidden bit আর IEEE 754-এর জন্ম

দশমিকে ১৩.৬২৫ = ১.৩৬২৫ × ১০^১. Binary-তেও একই কাজ হয়, শুধু ভিত্তি ২। পর্ব ০৪-এর বারবার ভাগ আর বারবার গুণের নিয়মে ১৩.৬২৫ = 1101.101₂. বিন্দুকে তিন ঘর বাঁয়ে সরালে:

বিন্দুর বাঁয়ে ঠিক একটি অশূন্য অঙ্ক রেখে লেখার এই রূপকে বলে normalized রূপ। দশমিকে বাঁয়ের অঙ্কটি ১ থেকে ৯-এর যেকোনোটি হতে পারে, কিন্তু binary-তে অশূন্য অঙ্ক একটিই — 1. তাই প্রতিটি normalized binary সংখ্যা 1 দিয়ে শুরু হয়। যে বিট সবসময় 1, তা মেমরিতে লিখে রাখার কোনো মানে নেই। IEEE 754 তাই সেটি বাদ দেয় — একে বলে hidden bit বা implicit leading 1. ফলে ২৩টি বিট জায়গা দিয়ে ২৪ বিটের নির্ভুলতা পাওয়া যায়: একটি বিট বিনা খরচে।

সাধারণ রূপটি তাহলে:

এখানে s চিহ্ন-বিট (0 ধনাত্মক, 1 ঋণাত্মক), f হলো বিন্দুর পরের বিটগুলো (fraction), আর e আসল সূচক। লক্ষ করুন, চিহ্ন আলাদা একটি বিট — এটি পর্ব ০৫-এর sign-magnitude চিন্তা, two's complement নয়। তাই +x আর −x-এর মধ্যে পার্থক্য কেবল একটি বিট।

আজ x86-64-এর SSE/AVX, ARM-এর FP/NEON আর RISC-V-র F ও D extension — সবাই IEEE 754 মানে। RISC-V-তে float-এর জন্য আলাদা ৩২টি register (f0–f31) আর flw/fsw, fadd.s-এর মতো instruction আছে; এদের কয়েকটি (যেমন fadd.s, fmadd.s) পর্ব ২৪-এ আসবে। এই পর্বে আমাদের প্রশ্ন একটিই: বিটগুলো ঠিক কীভাবে সাজানো।

৩ · Single আর double precision — বিটের নকশা আর হাতে encode

IEEE 754-এর দুটি প্রধান format হলো single precision (৩২ বিট, C-র float) আর double precision (৬৪ বিট, C-র double, JavaScript-এর Number, Python-এর float). দুটির নকশা একই — তিনটি field পাশাপাশি:

Exponent কেন biased — two's complement নয় কেন?

আসল সূচক e ঋণাত্মকও হয় (০.১-এর জন্য e = −৪)। IEEE 754 e-কে সরাসরি না রেখে রাখে E = e + bias — পর্ব ০৫-এর excess-K পদ্ধতি, যেখানে single-এর bias ১২৭। কেন? কারণ তাতে E একটি সাধারণ unsigned সংখ্যা হয়ে যায়, আর ছোট সূচক মানে ছোট E. ফলে দুটি ধনাত্মক float-এর বিট-প্যাটার্নকে সাধারণ unsigned integer হিসেবে তুলনা করলেই বড়-ছোট বোঝা যায় — exponent ওপরের বিটে, fraction নিচে, তাই integer comparator নিজেই আগে exponent, পরে fraction মেলায়। এর ব্যবহার অংশ ৭-এ দেখবেন। ৮ বিটে E = ০ থেকে ২৫৫; কিন্তু ০ আর ২৫৫ বিশেষ মানের জন্য সংরক্ষিত (অংশ ৫), তাই normal সংখ্যায় E = ১…২৫৪, অর্থাৎ e = −১২৬…+১২৭।

হাতে encode — −১৩.৬২৫ থেকে বিট

বিশ্ববিদ্যালয় আর চাকরির পরীক্ষায় সবচেয়ে বেশি আসা প্রশ্ন: “−১৩.৬২৫-কে IEEE 754 single precision-এ লিখুন, hex-এ উত্তর দিন।” প্রতিটি ধাপে “পরে” চাপুন:

x = -13.625  is negative
sign bit s = 1
from now on work with |x| = 13.625

integer 13 = 8 + 4 + 1   -> 1101
fraction .625 x 2 = 1.25  -> 1
         .25  x 2 = 0.5   -> 0
         .5   x 2 = 1.0   -> 1   (stops)
13.625 = 1101.101 (base 2)

1101.101 = 1.101101 x 2^3
(point moved 3 places left)
hidden bit = 1,  f = 101101...

E = e + 127 = 3 + 127 = 130
130 = 128 + 2 = 1000 0010

drop the hidden 1:  101101
pad to 23 bits:     10110100000000000000000
(no bits left over -> exact, no rounding)

s │ E        │ f
1 │ 10000010 │ 10110100000000000000000
groups of 4:
1100 0001 0101 1010 0000 0000 0000 0000

1100 0001 0101 1010 0000 ...
 C    1    5    A    0  0 0 0
-13.625f = 0xC15A0000
little-endian memory: 00 00 5A C1

৪ · Hex থেকে সংখ্যা — decode করা, আর একটি encoder-decoder

উল্টো প্রশ্নটিও সমান জনপ্রিয়: “মেমরিতে 0x40490FDB আছে; এটি float হলে মান কত?” Debugger-এ, network packet-এ বা ফাইলের hex dump-এ float এভাবেই দেখা যায়। নিয়ম একটাই — তিনটি field আলাদা করুন, তারপর সূত্র বসান।

  4    0    4    9    0    F    D    B
0100 0000 0100 1001 0000 1111 1101 1011

s = 0
E = 1000 0000
f = 100 1001 0000 1111 1101 1011

E = 1000 0000 = 128
E is not 0 and not 255  -> normal
e = E - 127 = 1

1.f = 1.10010010000111111011011 (base 2)
    = 1 + 0x490FDB / 2^23
    = 1 + 4788187 / 8388608
    = 1.5707963705...

value = (-1)^0 x 1.5707963705... x 2^1
      = 3.1415927410125732...

pi        = 3.14159265358979...
stored    = 3.14159274101257...
error     = 8.74e-8
ULP at 2  = 2^-22 = 2.38e-7  (half = 1.19e-7)
-> 0x40490FDB is the float nearest to pi

এবার হাতের দুটি হিসাব যন্ত্রে মিলিয়ে নিন। নিচের JavaScript simulator কোনো library ছাড়াই একটি দশমিক সংখ্যাকে ধাপে ধাপে single precision-এ encode করে — normalize, bias, ২৩ বিট কেটে নেওয়া আর বাড়তি বিট দেখে rounding — আর ফলটি ব্রাউজারের নিজস্ব Float32Array/DataView-এর সাথে মেলায়। দ্বিতীয় ফাংশনটি যেকোনো ৩২ বা ৬৪ বিটের hex pattern decode করে। Run চাপুন:

// IEEE 754 by hand: encode a decimal step by step, decode hex patterns
const dv = new DataView(new ArrayBuffer(8)), b2 = (v, w) => v.toString(2).padStart(w, "0");
function encode32(x) {                       // normal numbers only
  const s = x < 0 ? 1 : 0; let a = Math.abs(x), e = 0;
  while (a >= 2) { a /= 2; e++; }            // normalise: 1 <= a < 2
  while (a < 1) { a *= 2; e--; }
  const exact = (a - 1) * 2 ** 23; let f = Math.floor(exact);
  const rest = exact - f;                    // the bits that do not fit
  if (rest > 0.5 || (rest === 0.5 && f % 2 === 1)) f++;   // nearest, ties to even
  if (f === 2 ** 23) { f = 0; e++; }         // rounding carried into exponent
  const E = e + 127, bits = ((s << 31) | (E << 23) | f) >>> 0;
  dv.setFloat32(0, x);
  console.log(`${x}: sign ${s} | ${a.toString(2).slice(0, 18)}... x 2^${e} | E = ${e}+127 = ${E} = ${b2(E, 8)}`);
  console.log(`   fraction ${b2(f, 23)} (rounding rest ${rest.toFixed(2)})`);
  console.log(`   = 0x${bits.toString(16).toUpperCase()}  Float32Array says 0x${dv.getUint32(0).toString(16).toUpperCase()}`);
}
function decode(hex, bitsW) {                // bitsW = 32 or 64
  const [EW, FW, bias] = bitsW === 32 ? [8, 23, 127] : [11, 52, 1023];
  const b = BigInt(hex), s = Number(b >> BigInt(bitsW - 1));
  const E = Number((b >> BigInt(FW)) & ((1n << BigInt(EW)) - 1n));
  const f = b & ((1n << BigInt(FW)) - 1n), max = 2 ** EW - 1;
  const v = bitsW === 32 ? (dv.setUint32(0, Number(b)), dv.getFloat32(0))
                         : (dv.setBigUint64(0, b), dv.getFloat64(0));
  const kind = E === max ? (f ? "NaN" : "infinity") : E === 0 ? (f ? "subnormal" : "zero")
             : `normal, 1.f x 2^(${E}-${bias}) = 2^${E - bias}`;
  console.log(`${hex.padEnd(18)} s=${s} E=${E} f=0x${f.toString(16)} ${kind} -> ${Object.is(v, -0) ? "-0" : v}`);
}
[-13.625, 0.1, 16777217, 1e30].forEach(encode32);
["0x40490FDB", "0xC15A0000", "0x00000001", "0x80000000", "0x7F800000", "0x7FC00000"]
  .forEach(h => decode(h, 32));
["0x3FB999999999999A", "0x400921FB54442D18"].forEach(h => decode(h, 64));

-13.625: sign 1 | 1.101101... x 2^3 | E = 3+127 = 130 = 10000010
   fraction 10110100000000000000000 (rounding rest 0.00)
   = 0xC15A0000  Float32Array says 0xC15A0000
0.1: sign 0 | 1.1001100110011001... x 2^-4 | E = -4+127 = 123 = 01111011
   fraction 10011001100110011001101 (rounding rest 0.80)
   = 0x3DCCCCCD  Float32Array says 0x3DCCCCCD
16777217: sign 0 | 1.0000000000000000... x 2^24 | E = 24+127 = 151 = 10010111
   fraction 00000000000000000000000 (rounding rest 0.50)
   = 0x4B800000  Float32Array says 0x4B800000
1e+30: sign 0 | 1.1001001111100101... x 2^99 | E = 99+127 = 226 = 11100010
   fraction 10010011111001011001010 (rounding rest 0.80)
   = 0x7149F2CA  Float32Array says 0x7149F2CA
0x40490FDB         s=0 E=128 f=0x490fdb normal, 1.f x 2^(128-127) = 2^1 -> 3.1415927410125732
0xC15A0000         s=1 E=130 f=0x5a0000 normal, 1.f x 2^(130-127) = 2^3 -> -13.625
0x00000001         s=0 E=0 f=0x1 subnormal -> 1.401298464324817e-45
0x80000000         s=1 E=0 f=0x0 zero -> -0
0x7F800000         s=0 E=255 f=0x0 infinity -> Infinity
0x7FC00000         s=0 E=255 f=0x400000 NaN -> NaN
0x3FB999999999999A s=0 E=1019 f=0x999999999999a normal, 1.f x 2^(1019-1023) = 2^-4 -> 0.1
0x400921FB54442D18 s=0 E=1024 f=0x921fb54442d18 normal, 1.f x 2^(1024-1023) = 2^1 -> 3.141592653589793

৫ · বিশেষ মান — ±0, ±∞, NaN আর subnormal

Exponent field-এর দুটি মান — সব 0 আর সব 1 — IEEE 754 সাধারণ সংখ্যার জন্য রাখেনি। এই দুটি দিয়ে সে এমন জিনিস প্রকাশ করে, যা কোনো integer format পারে না। পুরো টেবিলটি মুখস্থ রাখার মতো:

Double-এ একই টেবিল, শুধু E-এর সব-1 মান ২০৪৭, আর subnormal-এর সূচক −১০২২।

±0 — শূন্যেরও চিহ্ন আছে

চিহ্ন আলাদা বিট বলে 0x00000000 (+0) আর 0x80000000 (−0) দুটিই বৈধ। তুলনায় -0.0 == 0.0 সত্য, কিন্তু আচরণে পার্থক্য আছে: 1/+0 = +∞ আর 1/−0 = −∞. কেন রাখা হলো? একটি খুব ছোট ঋণাত্মক সংখ্যা underflow হয়ে শূন্য হলে −0 অন্তত বলে দেয় “এটি নিচের দিক থেকে এসেছে”। জটিল সংখ্যার ফাংশনে (যেমন বর্গমূল বা log-এর “branch cut”-এর দুই পাশে) এই চিহ্ন সঠিক উত্তর বেছে নিতে সাহায্য করে — Kahan এ নিয়ে ১৯৮৭ সালে একটি বিখ্যাত প্রবন্ধও লিখেছিলেন। দৈনন্দিন ভাষায়: তাপমাত্রা −০.০০০১° থেকে গোল হয়ে “−০°” দেখালে বোঝা যায় বরফ গলেনি।

±∞ — overflow-এর ভদ্র উত্তর

সবচেয়ে বড় finite float প্রায় ৩.৪ × ১০^৩৮ (0x7F7FFFFF = (২ − ২^−২৩) × ২^১২৭); double-এ প্রায় ১.৮ × ১০^৩০৮. এর চেয়ে বড় ফল হলে প্রোগ্রাম থামে না, ফল হয় ∞ — আর ∞ নিয়ে হিসাব চলতে থাকে নিয়মমতো: ∞ + ১ = ∞, ১/∞ = ০, আর ধনাত্মক সংখ্যা/০ = +∞. পর্ব ০৫-এ দেখা integer overflow-এর মতো চুপচাপ ঘুরে ঋণাত্মক হয়ে যায় না — এটি floating point-এর একটি বড় নিরাপত্তা-গুণ।

NaN — “এটি কোনো সংখ্যা নয়”

০/০, ∞ − ∞, ০ × ∞ বা ঋণাত্মক সংখ্যার বর্গমূল (বাস্তব সংখ্যায়) — এদের কোনো অর্থপূর্ণ উত্তর নেই। IEEE 754 এখানে দেয় NaN (Not a Number): E সব 1, f অশূন্য। NaN-এর তিনটি বৈশিষ্ট্য পরীক্ষায় আর ইন্টারভিউতে প্রায়ই আসে:

• সংক্রামক: NaN-এর সাথে যেকোনো হিসাবের ফল NaN — একটি ভুল পুরো হিসাবের শেষ পর্যন্ত বয়ে যায়, হারিয়ে যায় না।
• নিজের সমানও নয়: NaN == NaN মিথ্যা, NaN != NaN সত্য, আর <, > সবই মিথ্যা (“unordered”). তাই x != x হলো NaN চেনার পুরোনো কৌশল; আধুনিক কোডে isnan(x) ব্যবহার করুন।
• দুই রকম: quiet NaN নীরবে ছড়ায়; signalling NaN ব্যবহার হলেই একটি “invalid operation” exception-flag তোলে — যেমন অ-initialized মেমরি চিহ্নিত করতে। ২০০৮-এর standard সুপারিশ করে fraction-এর সবচেয়ে বড় বিট 1 মানে quiet (তাই 0x7FC00000 quiet, 0x7F800001 signalling); কিছু পুরোনো architecture (যেমন পুরোনো MIPS) উল্টো নিয়ম মানত বলে জানা যায়।

NaN-এর বাকি ২২টি fraction বিট (double-এ ৫১টি) “payload” — চাইলে এখানে তথ্য লুকানো যায়। কিছু JavaScript engine (যেমন Firefox-এর SpiderMonkey) এই কৌশলে (“NaN-boxing”) pointer আর ছোট integer-কে double-এর NaN-এর ভেতরে রাখে।

Subnormal — শূন্যের কাছে ধীরে ধীরে নামা

সবচেয়ে ছোট normal float ২^−১২৬ ≈ ১.১৮ × ১০^−৩৮. এর নিচে কী? যদি সরাসরি ০-তে লাফ দিত, তাহলে ০ আর ২^−১২৬-এর মাঝে একটি বিশাল ফাঁক থাকত — অথচ ২^−১২৬-এর ঠিক ওপরে পাশাপাশি float-এর দূরত্ব মাত্র ২^−১৪৯. IEEE 754-এর সমাধান subnormal (পুরোনো নাম denormal) সংখ্যা: E = ০ হলে hidden bit 1 নয় 0, আর সূচক স্থির −১২৬:

এর ফলে শূন্যের কাছের মানগুলো সমান দূরত্বে (২^−১৪৯) সাজানো — fixed point-এর মতো — আর নির্ভুলতা ধীরে ধীরে কমে। একে বলে gradual underflow. নিচের stepper-এ দেখুন এটি কোন সমস্যার সমাধান করে:

0x00800000 = 0 00000001 000...0
E = 1  -> 1.0 x 2^(1-127) = 2^-126
       = 1.1754944e-38

x = 0x00C00000 = 1.5 x 2^-126 = 1.7632415e-38
y = 0x00800000 = 1.0 x 2^-126 = 1.1754944e-38
x != y

x - y = 0.5 x 2^-126 = 2^-127 = 5.877472e-39
smaller than the smallest normal (2^-126)
-> no normal float can hold it

E = 0 means 0.f x 2^-126
2^-127 = 0.1 (base 2) x 2^-126
f = 100 0000 0000 0000 0000 0000
x - y = 0x00400000   (exact!)

flush-to-zero:  x - y -> 0
x != y  but  x - y == 0
if (x != y) z = 1 / (x - y);   -> z = inf !

largest subnormal  0x007FFFFF = 1.1754942e-38
smallest subnormal 0x00000001 = 2^-149 = 1.4e-45
spacing everywhere below 2^-126: 2^-149
(2^-149) / 2 -> tie -> rounds to even -> 0

Double-এ একই চিত্র: সবচেয়ে ছোট normal ২^{−১০২২} ≈ ২.২২ × ১০^−৩০৮, সবচেয়ে ছোট subnormal ২^{−১০৭৪} ≈ ৪.৯৪ × ১০^−৩২৪. মূল্যটাও জানা দরকার: অনেক x86 CPU-তে subnormal ফল তৈরি হলে hardware একটি ধীর “microcode assist” পথে যায়, যাতে শতাধিক cycle পর্যন্ত লাগতে পারে। তাই গতি-প্রধান কাজে (গেম, অডিও, GPU) প্রায়ই flush-to-zero চালু রাখা হয় — x86-এর SSE-তে MXCSR register-এর FTZ আর DAZ বিট দিয়ে। এটি standard-এর বাইরে একটি সচেতন আপস।

এবার C দিয়ে সবকিছু একসাথে দেখি। প্রোগ্রামটি প্রতিটি float-এর বিট memcpy দিয়ে একটি uint32_t-এ আনে (type punning-এর নিরাপদ উপায়, পর্ব ০৬-এ 1.0f-এর byte দেখার সময় যেমন করেছিলেন), তিনটি field আলাদা করে, শ্রেণি বলে আর শুধু field থেকে মানটি আবার বানায়:

#include <stdio.h>
#include <stdint.h>
#include <string.h>
#include <float.h>
#include <math.h>                       /* only for the INFINITY, NAN macros */
static double pow2(int k) {             /* exact 2^k without libm */
    double r = 1.0;
    for (; k != 0; k += (k > 0) ? -1 : 1) r = (k > 0) ? r * 2.0 : r / 2.0;
    return r;
}
static void show(const char *name, float x) {
    uint32_t b;
    memcpy(&b, &x, sizeof b);           /* safe type punning */
    unsigned s = b >> 31, e = (b >> 23) & 0xFFu, f = b & 0x7FFFFFu;
    const char *kind = e == 255 ? (f ? "NaN" : "infinity")
                     : e == 0   ? (f ? "subnormal" : "zero") : "normal";
    double v = e == 255 ? (f ? (double)NAN : (double)INFINITY)
             : e == 0   ? f * pow2(-149)                 /* 0.f x 2^-126  */
             : (f | 0x800000u) * pow2((int)e - 150);     /* 1.f x 2^(e-127) */
    printf("%-7s 0x%08X  s=%u e=%3u f=0x%06X %-9s -> %s%.9g\n",
           name, (unsigned)b, s, e, f, kind, s ? "-" : "", v);
}
int main(void) {
    show("1.0", 1.0f);      show("-13.625", -13.625f); show("0.1f", 0.1f);
    show("1e-40f", 1e-40f); show("INF", INFINITY);     show("NAN", NAN); show("-0.0f", -0.0f);
    printf("FLT_EPSILON  %.9g = 2^-23\nFLT_MIN      %.9g = 2^-126\n",
           (double)FLT_EPSILON, (double)FLT_MIN);
#ifdef FLT_TRUE_MIN
    printf("FLT_TRUE_MIN %.9g = 2^-149\n", (double)FLT_TRUE_MIN);
#else
    printf("FLT_TRUE_MIN %.9g = 2^-149 (computed)\n", pow2(-149));
#endif
    printf("FLT_MAX      %.9g\n", (double)FLT_MAX);
    return 0;
}

1.0     0x3F800000  s=0 e=127 f=0x000000 normal    -> 1
-13.625 0xC15A0000  s=1 e=130 f=0x5A0000 normal    -> -13.625
0.1f    0x3DCCCCCD  s=0 e=123 f=0x4CCCCD normal    -> 0.100000001
1e-40f  0x000116C2  s=0 e=  0 f=0x0116C2 subnormal -> 9.9999461e-41
INF     0x7F800000  s=0 e=255 f=0x000000 infinity  -> inf
NAN     0x7FC00000  s=0 e=255 f=0x400000 NaN       -> nan
-0.0f   0x80000000  s=1 e=  0 f=0x000000 zero      -> -0
FLT_EPSILON  1.1920929e-07 = 2^-23
FLT_MIN      1.17549435e-38 = 2^-126
FLT_TRUE_MIN 1.40129846e-45 = 2^-149
FLT_MAX      3.40282347e+38

৬ · নির্ভুলতা — significant digit, machine epsilon আর ULP

পরিসর বড় হলেও প্রতিটি float-এ অর্থপূর্ণ বিট সীমিত: single-এ ২৪টি, double-এ ৫৩টি। এর অর্থ কী, তা তিনটি ধারণায় ধরা যায়।

ULP — পাশাপাশি দুটি float-এর দূরত্ব

একই exponent-এর সব float — অর্থাৎ ২^e থেকে ২^e+১-এর মাঝের একটি “binade” — সমান দূরত্বে বসানো, কারণ শুধু fraction-এর ২৩ বিট বদলায়। Fraction-এর সবচেয়ে ছোট বিটের মান ২^−২৩, তাই সেই binade-এ ধাপ:

এই ধাপের নাম ULP (unit in the last place). ১-এর কাছে single-এ ULP = ২^−২৩ ≈ ১.১৯ × ১০^−৭; ১০২৪-এর কাছে ২^−১৩ ≈ ১.২২ × ১০^−৪; ১০ লাখের কাছে ০.০৬২৫; আর ২^২৪ = ১৬৭৭৭২১৬-এর কাছে ঠিক ২. অর্থাৎ ১৬৭৭৭২১৬-এর পরের float ১৬৭৭৭২১৮ — মাঝের ১৬৭৭৭২১৭ float-এ রাখাই যায় না। এটিই সবচেয়ে ছোট ধনাত্মক integer, যা float-এ হুবহু থাকে না (double-এ সেটি ২^৫৩ + ১)। তাই float-এ কোনো counter ১ করে বাড়াতে থাকলে ১৬৭৭৭২১৬-এ পৌঁছে সেটি আটকে যায় — f + 1.0f == f.

Machine epsilon আর আপেক্ষিক ভুলের সীমা — একটি ছোট প্রমাণ

Machine epsilon ε হলো ১ আর তার পরের float-এর দূরত্ব: single-এ ε = ২^−২৩ ≈ ১.১৯ × ১০^−৭ (C-র FLT_EPSILON), double-এ ২^−৫২ ≈ ২.২২ × ১০^−১৬ (DBL_EPSILON). এখন দাবি: যেকোনো বাস্তব সংখ্যা x (পরিসরের ভেতরে, normal অঞ্চলে) সবচেয়ে কাছের float-এ গোল করলে আপেক্ষিক ভুল সর্বোচ্চ ε/২।

প্রমাণ: ধরুন ২^e ≤ |x| < ২^e+১. ওই binade-এ পাশাপাশি float-এর দূরত্ব ২^e−২৩ = ε · ২^e. সবচেয়ে কাছের float বেছে নিলে ভুল সর্বোচ্চ অর্ধেক ধাপ, অর্থাৎ |fl(x) − x| ≤ ε · ২^e / ২. দুই দিকে |x| দিয়ে ভাগ করি, আর |x| ≥ ২^e ব্যবহার করি: |fl(x) − x| / |x| ≤ (ε · ২^e / ২) / ২^e = ε/২. ∎ এই ε/২ = ২^−২৪ ≈ ৫.৯৬ × ১০^−৮-কে বলে unit roundoff u. সাবধান: কিছু বই “machine epsilon” বলতে এই u-কেই বোঝায় — উত্তর লেখার আগে কোন সংজ্ঞা, তা উল্লেখ করুন।

দশমিক অঙ্কে কতটা? ২৪টি বিট মানে ২^২৪ রকম অবস্থা; সমান সংখ্যক দশমিক অঙ্ক d হলে ১০^d = ২^২৪, তাই d = ২৪ × log_১০২ ≈ ৭.২২. Double-এ ৫৩ × log_১০২ ≈ ১৫.৯৫. ব্যবহারিক নিয়ম দুটি: যেকোনো ৬ অঙ্কের দশমিক সংখ্যা float-এ গিয়ে ফিরে এলে অক্ষত থাকে (FLT_DIG = ৬), কিন্তু একটি float-কে নির্ভুলভাবে আবার পড়তে চাইলে ৯ অঙ্ক ছাপতে হয় (FLT_DECIMAL_DIG = ৯)। Double-এ সংখ্যা দুটি ১৫ আর ১৭।

নিচের JavaScript simulator দেখায় মেমরিতে আসলে কী আছে: BigInt দিয়ে একটি double-এর হুবহু দশমিক মান বের করে (প্রতিটি double = m × ২^k, আর ২⁻ⁿ = ৫ⁿ / ১০ⁿ, তাই দশমিক রূপ সবসময় সসীম), তারপর float-এর ULP মাপে আর ১৬৭৭৭২১৭-এর পরিণতি দেখায়:

// What is REALLY stored? Exact decimal value of a double, via BigInt
function exact(x) {                          // x > 0, finite
  const dv = new DataView(new ArrayBuffer(8)); dv.setFloat64(0, x);
  const b = dv.getBigUint64(0), E = Number((b >> 52n) & 0x7FFn);
  let m = b & ((1n << 52n) - 1n), k = -1074; // value = m x 2^k
  if (E) { m |= 1n << 52n; k = E - 1075; }
  if (k >= 0) return (m << BigInt(k)).toString();
  const d = (m * 5n ** BigInt(-k)).toString().padStart(1 - k, "0"); // m*2^k = m*5^-k / 10^-k
  return (d.slice(0, d.length + k) + "." + d.slice(d.length + k)).replace(/0+$/, "");
}
for (const [name, v] of [["0.1", 0.1], ["0.2", 0.2], ["0.3", 0.3], ["0.1+0.2", 0.1 + 0.2],
                         ["float(0.1)", Math.fround(0.1)]])
  console.log(name.padEnd(10), exact(v));
console.log("0.1 + 0.2 === 0.3 ?", 0.1 + 0.2 === 0.3, "| toPrecision(17):", (0.1 + 0.2).toPrecision(17));
const f32 = new Float32Array(1), u32 = new Uint32Array(f32.buffer);
function ulp32(x) { f32[0] = x; u32[0] += 1; return f32[0] - x; }   // gap to the next float
console.log("\nx            float ULP               double ULP");
for (const x of [1, 1024, 16777216, 1e10])
  console.log(String(x).padEnd(12), String(ulp32(x)).padEnd(23), 2 ** (Math.floor(Math.log2(x)) - 52));
console.log("\nfloat(16777217) =", Math.fround(16777217), "| float(16777219) =", Math.fround(16777219));
let s = 0; for (let i = 0; i < 10; i++) s += 0.1;
console.log("ten times 0.1 =", s, "| ten times 10 poisha =", 10 * 10, "poisha = 1.00 taka");

0.1        0.1000000000000000055511151231257827021181583404541015625
0.2        0.200000000000000011102230246251565404236316680908203125
0.3        0.299999999999999988897769753748434595763683319091796875
0.1+0.2    0.3000000000000000444089209850062616169452667236328125
float(0.1) 0.100000001490116119384765625
0.1 + 0.2 === 0.3 ? false | toPrecision(17): 0.30000000000000004

x            float ULP               double ULP
1            1.1920928955078125e-7   2.220446049250313e-16
1024         0.0001220703125         2.2737367544323206e-13
16777216     2                       3.725290298461914e-9
10000000000  1024                    0.0000019073486328125

float(16777217) = 16777216 | float(16777219) = 16777220
ten times 0.1 = 0.9999999999999999 | ten times 10 poisha = 100 poisha = 1.00 taka

৭ · গোল করার নিয়ম, ০.১ + ০.২-এর রহস্য, বিটে সাজানো আর টাকার হিসাব

০.১ কেন হুবহু থাকে না, আর কোন দিকে গোল হয়

পর্ব ০৪-এ দেখেছেন, ০.১ binary-তে কখনো শেষ হয় না: 0.000110011001100…₂. কারণ ১/১০-এর হরে ৫ আছে, আর binary-তে সসীম ভগ্নাংশ মানে হর শুধু ২-এর ঘাত। Normalize করলে ০.১ = 1.100110011001100…₂ × ২^−৪. Single-এ fraction-এর জায়গা ২৩ বিট:

kept 23 bits : 10011001100110011001100 | 1100 1100 ...
next bits    : 1100...  -> more than half an ULP
round up     : 10011001100110011001101  -> 0x3DCCCCCD

IEEE 754-এর ডিফল্ট নিয়ম হলো round to nearest, ties to even: সবচেয়ে কাছের প্রতিনিধিত্বযোগ্য মান নিন; ঠিক মাঝখানে পড়লে সেটি নিন, যার শেষ বিট 0. “Even”-এর পক্ষে যুক্তি — স্কুলের “৫ হলে সবসময় ওপরে” নিয়মে অসংখ্য হিসাবে ভুল একদিকে জমে (bias তৈরি হয়), কিন্তু ties to even-এ অর্ধেক tie ওপরে আর অর্ধেক নিচে যায়, ফলে গড়ে ভুল কাটাকাটি হয়। Standard আরও তিনটি “directed” mode রাখে (শূন্যের দিকে, +∞-এর দিকে, −∞-এর দিকে); এগুলো হিসাবের (operation) ভেতরে কীভাবে কাজ করে, আর catastrophic cancellation, Kahan summation, fused multiply-add ও Patriot-Ariane ধরনের দুর্ঘটনা — সব আসছে পরের পর্ব ২৪-এ। এখানে মূল কথা: প্রতিটি দশমিক literal মেমরিতে যাওয়ার মুহূর্তেই একবার গোল হয়।

Float-কে integer হিসেবে সাজানো — biased exponent-এর পুরস্কার

অংশ ৩-এ বলেছিলাম, exponent biased বলেই ধনাত্মক float-এর বিট-প্যাটার্ন unsigned integer হিসেবে একই ক্রমে সাজে: ১.৫ = 0x3FC00000 < ২.০ = 0x40000000 < ১০^৩০ = 0x7149F2CA. ঋণাত্মক সংখ্যায় sign-magnitude বলে ক্রম উল্টো (বড় magnitude মানে বেশি ঋণাত্মক)। একটি ছোট কৌশলে পুরো ক্রম ঠিক করা যায় — চিহ্ন-বিট 1 হলে সব বিট উল্টে দিন, নইলে শুধু চিহ্ন-বিট উল্টান:

key = (b & 0x80000000) ? ~b : (b | 0x80000000);   /* unsigned compare of keys = float order */

এর ফলে database index, GPU-র radix sort আর কিছু hardware comparator floating-point তুলনা করে সাধারণ integer circuit দিয়ে। একই কারণে “পরের float” বের করা মানে বিট-প্যাটার্নে ১ যোগ — C-র nextafterf ভেতরে প্রায় এটিই করে, আর simulator-এর ulp32() ঠিক এভাবেই ULP মেপেছে। (NaN আর −0/+0 আলাদা করে সামলাতে হয়।)

টাকা কখনো float-এ নয় — integer আর decimal floating point

৳১০৫০.৭৫ binary float-এ হুবহু থাকে না, ৳০.১০ও না। একটি bKash-ধাঁচের অ্যাপে লাখো লেনদেনে প্রতিবার সামান্য ভুল জমলে হিসাব মেলে না, আর হিসাবের আইনি নিয়ম (কোন অঙ্কে কোন দিকে গোল হবে) binary float মানতে পারে না। পেশাদার সমাধান দুটি:

1. ক্ষুদ্রতম এককে integer: ৳১০৫০.৭৫ = ১০৫০৭৫ পয়সা, একটি ৬৪-বিটের integer-এ। যোগ-বিয়োগ হুবহু; শুধু সুদ বা শতাংশের হিসাবে গোল করার নিয়মটি স্পষ্টভাবে লিখতে হয়।
2. Decimal floating point: IEEE 754-2008 ভিত্তি-১০-এর format যোগ করেছে — decimal32, decimal64 (১৬ দশমিক অঙ্ক) আর decimal128 (৩৪ অঙ্ক) — যেখানে ০.১ হুবহু থাকে। IBM-এর কিছু POWER আর mainframe processor-এ এর hardware সমর্থন আছে বলে জানা যায়; অন্যরা software-এ চালায়: Python-এর decimal, Java-র BigDecimal, আর C23-এ ঐচ্ছিক _Decimal64.

C আর Python-এর দিক থেকে এই অভ্যাসগুলো ধাপে ধাপে দেখানো আছে C কোর্সের data representation পর্বে আর Python কোর্সের সংখ্যা ও Decimal পর্বে।

৮ · ছোট format — FP16, bfloat16, FP8 আর AI hardware

একই নকশা (চিহ্ন + biased exponent + hidden bit-সহ fraction) যেকোনো মাপে বানানো যায়। Exponent-এর বিট ঠিক করে পরিসর, fraction-এর বিট ঠিক করে নির্ভুলতা — মোট বিট স্থির রাখলে একটি বাড়ালে অন্যটি কমে। AI-এর যুগে এই দর-কষাকষি হঠাৎ খুব গুরুত্বপূর্ণ হয়ে উঠেছে:

• FP16 (IEEE 754-2008-এর binary16): ছবি, গ্রাফিক্স আর মোবাইল GPU-তে বহুদিনের পরিচিত। নির্ভুলতা প্রায় ৩.৩ দশমিক অঙ্ক, কিন্তু পরিসর ছোট — ৬৫৫০৪-এর ওপরে ∞. AI training-এ ছোট gradient-গুলো প্রায়ই এর নিচে পড়ে শূন্য হয়ে যায়, তাই FP16-এ training করতে “loss scaling” নামের কৌশল লাগে।
• bfloat16 (“brain float”, Google-এর গবেষকদের প্রস্তাব; training-এর জন্য বানানো দ্বিতীয় প্রজন্মের TPU থেকে ব্যাপক ব্যবহার — প্রথম TPU ছিল ৮ বিট integer-এর): float-এর ঠিক ওপরের ১৬ বিট। ৮ বিট exponent বলে পরিসর float-এর সমান, কিন্তু মাত্র ৭ বিট fraction — নির্ভুলতা প্রায় ২.৪ দশমিক অঙ্ক। Neural network সামান্য নির্ভুলতা হারানো সহ্য করে, কিন্তু পরিসরের বাইরে যাওয়া সহ্য করে না — তাই training-এ এটি জনপ্রিয়। আজ অনেক GPU, ARM আর x86 CPU-তে bfloat16 instruction আছে।
• FP8: দুটি প্রধান রূপ — E4M3 (বেশি নির্ভুলতা, সাধারণত ওজন আর activation-এর জন্য) আর E5M2 (বেশি পরিসর, gradient-এর জন্য)। শিল্পের কয়েকটি কোম্পানি মিলে ২০২৩-এর দিকে OCP (Open Compute Project)-এর মাধ্যমে একটি spec প্রকাশ করে; সেখানে E4M3-এ ∞ নেই আর NaN শুধু একটি pattern-এ (exponent আর fraction-এর সব বিট 1, চিহ্ন-বিট যা-ই হোক), যাতে বড় মান ৪৪৮ পর্যন্ত পৌঁছানো যায়, আর E5M2 IEEE-র নিয়মের কাছাকাছি। বিভিন্ন hardware-এর সংস্করণে খুঁটিনাটি পার্থক্য আছে, আর IEEE-ও (P3109 নামে) ছোট format standard করার কাজ করছে — তাই নির্দিষ্ট chip-এর জন্য তার ডকুমেন্টেশন দেখুন।

কেন AI hardware এগুলো চায়? (১) বিট অর্ধেক মানে একই memory bandwidth-এ দ্বিগুণ সংখ্যা — আর AI-এর গতি প্রায়ই মেমরিতেই আটকে থাকে (পর্ব ৪৪-এর memory wall). (২) গুণক circuit-এর আকার মোটামুটি significand-এর বিটের বর্গের অনুপাতে বাড়ে, তাই ২৪ বিটের গুণকের জায়গায় অনেকগুলো ৮ বা ৪ বিটের গুণক বসানো যায়। এজন্যই tensor unit সাধারণত ছোট format-এ গুণ করে আর যোগফল জমায় FP32-এ — একে বলে mixed precision. GPU আর tensor unit-এর ভেতরের গঠন পর্ব ৫৫-এ আসবে।

bfloat16 বানানো সহজ: float-এর নিচের ১৬ বিট ফেলে দিলেই হয়। কিন্তু “ফেলে দেওয়া” (truncation) আর “গোল করা”-র পার্থক্য কত? নিচের প্রোগ্রামে দুটিই দেখুন:

#include <stdio.h>
#include <stdint.h>
#include <string.h>

static uint32_t bits(float x) { uint32_t b; memcpy(&b, &x, 4); return b; }
static float from(uint32_t b) { float x; memcpy(&x, &b, 4); return x; }

/* bfloat16 = the top 16 bits of a float: 1 sign, 8 exponent, 7 fraction */
static uint16_t bf16_trunc(float x) { return (uint16_t)(bits(x) >> 16); }
static uint16_t bf16_round(float x) {           /* round to nearest, ties to even */
    uint32_t b = bits(x);                       /* (NaN/overflow checks omitted) */
    b += 0x7FFFu + ((b >> 16) & 1u);
    return (uint16_t)(b >> 16);
}
int main(void) {
    const float v[] = { 3.14159265f, 0.1f, 1.0f / 3.0f, 123.456f, 65504.0f, 3.0e38f };
    printf("%-12s %-10s | %-6s %-13s | %-6s %-13s %s\n", "float", "bits",
           "trunc", "value", "round", "value", "rel.err");
    for (int i = 0; i < 6; i++) {
        uint16_t t = bf16_trunc(v[i]), r = bf16_round(v[i]);
        float ft = from((uint32_t)t << 16), fr = from((uint32_t)r << 16);
        double err = (fr - v[i]) / (double)v[i];
        printf("%-12.8g 0x%08X | 0x%04X %-13.8g | 0x%04X %-13.8g %+.2e\n",
               (double)v[i], (unsigned)bits(v[i]), (unsigned)t, (double)ft,
               (unsigned)r, (double)fr, err);
    }
    printf("bf16 epsilon = 2^-7 = %g, float epsilon = 2^-23 = %g\n",
           (double)(from(0x3F810000u) - 1.0f), (double)(from(0x3F800001u) - 1.0f));
    return 0;
}

float        bits       | trunc  value         | round  value         rel.err
3.1415927    0x40490FDB | 0x4049 3.140625      | 0x4049 3.140625      -3.08e-04
0.1          0x3DCCCCCD | 0x3DCC 0.099609375   | 0x3DCD 0.10009766    +9.77e-04
0.33333334   0x3EAAAAAB | 0x3EAA 0.33203125    | 0x3EAB 0.33398438    +1.95e-03
123.456      0x42F6E979 | 0x42F6 123           | 0x42F7 123.5         +3.56e-04
65504        0x477FE000 | 0x477F 65280         | 0x4780 65536         +4.89e-04
3e+38        0x7F61B1E6 | 0x7F61 2.990763e+38  | 0x7F62 3.0040553e+38 +1.35e-03
bf16 epsilon = 2^-7 = 0.0078125, float epsilon = 2^-23 = 1.19209e-07

৯ · AI-এর যুগে — float-এর বিট নিজে যাচাই করুন

IEEE 754-এর প্রশ্নে AI-এর উত্তর প্রায়ই “প্রায় ঠিক” — আর এই বিষয়ে প্রায় ঠিক মানে ভুল। একটি বিট এক ঘর সরে গেলে মান দ্বিগুণ বা অর্ধেক হয়ে যায়। সুখবর হলো, প্রতিটি দাবি এক লাইনের কোডে যাচাই করা যায়।

(উত্তর মিলিয়ে নিন: −০.৪০৬২৫ = −1.101₂ × ২^−২, E = ১২৫, তাই 0xBED00000.)

AI-এর উত্তর হাতে পেলে — যাচাইয়ের চেকলিস্ট

• প্রতিটি bit pattern নিজে চালানো কোডে (DataView, memcpy, Python-এর struct.pack) মিলিয়েছি
• bias (১২৭/১০২৩) আর hidden bit ঠিক জায়গায় আছে কি না দেখেছি, decode করে মূল মান ফেরত পেয়েছি
• “সবচেয়ে ছোট”, “epsilon”, “কত অঙ্ক”-এর মতো দাবিতে কোন সংজ্ঞা ব্যবহার হয়েছে, তা স্পষ্ট করে নিয়েছি
• FP16/bfloat16/FP8-এর তথ্য নির্দিষ্ট hardware-এর ডকুমেন্টেশনে মিলিয়েছি — vendor-ভেদে পার্থক্য থাকতে পারে
• টাকা-পয়সার কোডে AI float/double ব্যবহার করলে সেটি বদলে integer পয়সা বা decimal type ব্যবহার করেছি