- vừa được xem lúc

Sự thật thú vị về phép khai căn: Từ "Ai chả biết" đến những điều "miệng chữ Ô mắt chữ A"

0 0 27

Người đăng: Erics Nguyen

Theo Viblo Asia

Note: Bài viết có nhiều nhận định cá nhân, sử dụng nhiều phép tưởng tượng ,và mang yếu tố chủ đạo là seeder cho bản thân ^^, mọi gạch đá xi măng đều được trân trọng.

Quay lại thời gian hồi còn đi học, chắc hẳn các bạn đều nhận ra rằng cái gì xuôi thì lúc nào làm cũng dễ hơn là ngược, mấy phép công, phép nhân đều dễ hơn phép chia và phép trừ. Như trường hợp lớp mình ngày xưa, 99% các bạn đều chả hiểu gì về phép trừ và phép chia cho đến khi bị ăn thước vào tay ? và mình lại là ngoại lệ ^^ . Với kỷ lục 8 năm liền là trò cưng của các cô dạy Toán, cho đến năm thứ 9 mọi hình tượng xây dựng trong lòng các cô đã sụp đổ vì sau khi thử mọi loại "đồ ăn" mình vẫn chưa hiểu làm sao để tính cái phép khai căn.

Đây các bác xem:

Thế thì Cô có oái oăm không chứ.

Và sau đó thì mình đã chấp nhận với bản thân là k có công thức tính nào cả, phải học thuộc lòng kết quả của các con số căn bậc hai cơ bản như : 4\sqrt{4}, 9\sqrt{9}... Và nếu phải tính những con số khác thì áp dụng tính chất nhân để làm nhỏ giá trị: 8\sqrt{8} = 2 x 2\sqrt{2} ^^ hoặc là bấm máy tính để lấy kết quả mỗi khi thèm "ăn" gì đó.

Và đó là câu chuyện của 10 năm về trước, giờ đây sau 10 năm đã đến lúc mình báo thù, chắc chắn phải có huyền cơ trong phép tính căn bậc hai chứ, vì con máy tính cầm tay 15 nút của mẹ mình vẫn tính được cơ mà. Và đó là động lực thôi thúc giúp mình "báo thù", mình sẽ tìm ra được cách tính căn bậc hai.

Khi giờ đây mình đã là một thợ code chân chính, với sự giúp đỡ của dàn máy tính toàn RGB tâm huyết thì mình tin là mình sẽ làm được. Và đúng là trời không phụ lòng người "tài" chỉ một vài phút ngồi tĩnh tâm mình đã ngộ ra một bộ "võ công" có thể giúp mình "trả thù" và xưng bá:

Phương pháp tính số ở giữa

Thật vậy nếu X = Y\sqrt{Y} thì X sẽ nằm trong khoảng từ 0 đến Y, Thử con số ở giữa Y2\frac{Y}{2}, tính bình phương nó lên rồi so sánh với Y, nếu nó lớn hơn thì sẽ thử con số tiếp theo trong khoảng 0 đến Y2\frac{Y}{2}. lặp lại các lần thử và chọn đến khi chọn được một con số đủ làm đối thủ khâm phục.

public static void GetTheMiddleMethod(double number, double acceptError)
{ double start = 0; double end = number; double theMiddle = GetTheMiddle(start, end); double times = 0; while (Math.Abs(theMiddle*theMiddle - number) > acceptError) { if (theMiddle*theMiddle > number) { end = theMiddle; } else { start = theMiddle; } theMiddle = GetTheMiddle(start, end); times++; } Console.WriteLine(theMiddle); Console.WriteLine(times);
} double number = 10;
double acceptError = 0.00001;
GetTheMiddleMethod(number, acceptError);

Và đây là kết quả

result = 3.1622791290283203/ time = 19 ?

19 lần tính toán để ra kết quả với "sai số ảo" 0.00001 thì thật là đáng thất vọng, ta phải cải tiến nó thôi. ra rà xem nào:

  • hầu hết Y\sqrt{Y} đều từ 0 cho đến Y2\frac{Y}{2}, trừ vài trường hợp ngao ngáo ra như đám nhỏ hơn 4. ? ta thêm if else vào là nhanh hơn rồi.

...

Nhưng đó chỉ là suy nghĩ nông cạn, khi áp dụng nó cũng chỉ giảm con số từ 19 xuống 18 ?. Thế này thì còn báo cái gì nữa. Kiểu này phải tự phế thế này đi mà luyện thế mới thôi ?

Áp dụng phương pháp tính Newton-Raphson

Cũng mất 1 thơi gian khá lâu để quên đi nỗi tự nhục ^^ và sau một vài đường cơ bản search google thì chúng ta có thể sẵn sàng luyện thế thứ 2 này. Về cơ bản ý tưởng sẽ đến từ lý thuyết của môn Toán, và chúng ta lại thấy rằng toán đã cứu cánh chúng ta như thế nào trong các tình huống mà nghĩ hoài ko ra cách. Nói qua một chút về phương pháp tính Newton-Raphson, đó là cách giúp chúng ta có thể tìm được nghiệm xấp xỉ của phương trình. Ở đây mình có một hình vẽ khá dễ hiểu cho việc này:

Giải thích về phương pháp tính Newton-Raphson (ref: Wikipedia)

Giải thích:

Chúng ta sẽ cần tìm điểm x0x_{0} - cái điểm mà hàm f(x)f(x) cắt với trục OxOx khi đó f(x)f(x) =00xx = x0x_{0} cũng chính là nghiệm của phương trình f(x)f(x) =00. Ngắn gọi ý tưởng là mỗi lần kẻ đạo hàm chúng ta sẽ tìm được điểm xnx_{n} sát với x0x_{0} hơn. Đạo hàm f(x)f'(x) thì là hệ số góc của tiếp tuyến của hàm f(x)f(x) tại điểm xx.

Công thức như sau:

  • f(x)f'(x) = tan(Θtan(\Theta) = f(x)xxn\frac{f(x)}{x - x_{n}} => xnx_{n} = xx - f(x)f(x)\frac{f(x)}{f'(x)}

Quay lại với phép tính căn bậc hai, về cơ bản chúng ta cần tính

  • XX = Y\sqrt{Y} => X2X^2 = YY

=> Phương trình cần tính:

  • f(x)f(x) = x2x^2 - yy = 00

Cứ thế ốp phương pháp tính vào thôi:

  • xn+1x_{n+1} = xnx_{n} - f(xn)f(xn)\frac{f(x_n)}{f'(x_n)} = xn2\frac{x_n}{2} + Y2xn\frac{Y}{2x_n}

Kaka, implement thôi:

public static void ByTheNewtonMethod(double number, double acceptError)
{ double root = number; int times = 0; while (Math.Abs(root*root - number) > acceptError) { root = root / 2 + number / (2 * root); times++; } Console.WriteLine(root); Console.WriteLine(times);
} double number = 10;
double acceptError = 0.00001;
ByTheNewtonMethod(number, acceptError);

kaka và đây là kết quả:

result = 3.1622776604441363/ time = 4 ???

Mê luôn mọi người ạ, Đúng là có chút kiến thức toán bỏ vào thôi nó vẫn khác, Đó đây là 1 ví dụ nhỏ về việc toán sẽ đóng vai trò gì trong lập trình đó các bạn!

Note: Về việc chứng minh phương pháp Newton Raphson:

Phương pháp này có thể chứng minh từ định lý Taylor:

  • f(x2)f(x_2) = f(x1)f(x_1) + f(x1)f'(x_1) (x2x_2 - x1x_1) + 12\frac{1}{2} f(x1)f''(x_1) (x2x1)2(x_2 - x_1)^2
  • Giả sử có x2x_2 \approx x1x_1 thì có thể rút ngọn lại là:
  • f(x2)f(x_2) = f(x1)f(x_1) + f(x1)f'(x_1) (x2x_2 - x1x_1) (1)

Giả sử ta có:

  • x3x_3 = x2x_2 + δ\deltaf(x3)f(x_3) = 0 theo (1) => f(x3)f(x_3) = f(x2)f(x_2) + f(x2)f'(x_2) δ\delta

    => f(x1)f(x_1) + f(x1)f'(x_1) (x2x1)(x_2 - x_1) + f(x2)f'(x_2) δ\delta = 0

    => (x2x1)(x_2-x_1) f(x1)f'(x_1) = - f(x1)f(x_1) - f(x2)f'(x_2) δ\delta

    => x2x_2 = x1x_1 - f(x1)f(x1)\frac{f(x_1)}{f'(x_1)} - δ\delta f(x2)f(x1)\frac{f'(x_2)}{f'(x_1)}

    OK, Nếu δ\delta (x2x_2 càng sát vs nghiệm x3x_3 của pt f(x3)f(x_3) = 0) đủ nhỏ thì ta sẽ có:

    x2x_2 = x1x_1 - f(x1)f(x1)\frac{f(x_1)}{f'(x_1)}

    Đây chính là phương pháp Newton Raphson của chúng ta.

Kaka, tầm này là hành tẩu jang hồ đi trả thù được rồi ^^ (*)

Fast inverse square root

(*) Đó là tôi của ngày hôm qua...

Đây là phát chốt của mình, một bí kíp thất truyền được mình học đc từ sư phụ gu gự, và sự phụ vi qui. đó bí kíp Quake 3 (tên ngắn gọn hơn của fast inverse square root, bị ae fb hiểu lầm là 3 qu** e chứ ?) Tác giả được biết đến là kỹ sư John Carmack.

Một bí kíp võ công có thể tính căn bậc hai trong 1 nốt nhạc, ko cần lặp chả cần thử. 2 bí kíp bên trên độ phức tạp thuật toán có thể tạm gọi là O(n) (tạm gọi vì mình thấy O(n) không tuyến tính theo n). tuy nhiên với bí kíp này chỉ là O(1) thôi ? Chúng ta xem qua nó nhé:

{ long i; float x2, y; const float threehalfs = 1.5F; x2 = number * 0.5F; y = number; i = * ( long * ) &y; // evil floating point bit level hacking i = 0x5f3759df - ( i >> 1 ); // what the fuck?  y = * ( float * ) &i; y = y * ( threehalfs - ( x2 * y * y ) ); // 1st iteration //	y = y * ( threehalfs - ( x2 * y * y ) ); // 2nd iteration, this can be removed return y;
}

Để cho anh em k phải chờ lâu mình sẽ đi thằng vào vấn đề luôn: để tu luyện bí kíp này ta cần 5 bước:

  • Theo tiêu chuẩn số phẩy động IEEE754 (ai tríp bồ i) được sử dụng để biểu diễn số float trong máy tính ta có:

    • Một số float đc biểu diễn dạng 32 bit:

      • S là bít đầu biểu thị dấu của số
      • E là 8 bít kế, biểu diễn phần mũ
      • m là 23 bít cuối, biểu diễn phần trị

      Công thức liên hệ giữa số float và các bit E, m biểu diễn như sau:

      aa = (1 + mam_a) x 2Ea1272^{E_a-127} (1)

      Có một lưu ý là E tính theo hệ cơ số 2 bình thường, còn m tính hệ cơ số 2 ngược, tức là bit đầu là 212^{-1}, bít thứ n là 2n2^{-n}, Mình lấy ví dụ là số 0.15625 nhé khi biểu diễn dạng số phẩy động thì ta có:

      E = 124, m = 222^{-2} = 0.25

      => công thức liên hệ sẽ là: 0.15625 = (1+0.25) x 21241272^{124-127}

  • Khai triển công thức (1) ta có:

    aa = (1 + mam_a) x 2Ea1272^{E_a-127}

    => log2(a)log_2(a) = log2(1+ma)log_2(1 + m_a) + log2(2Ea127)log_2(2^{E_a-127})

    => log2(a)log_2(a) = log2(1+ma)log_2(1 + m_a) + Ea127E_a-127 (2)

    • Chúng ta sẽ bắt đầu với sự kì diệu đầu tiên:

      • Chúng ta có : log2(1+m)log_2(1 + m) \approx mm với tất cả mm \in [0,1) Đúng vậy m của chúng ta tính bằng dãy 212^{-1} ... 2232^{-23} luôn thuộc khoảng [0,1) nên có thể áp dụng phương trình trên để phá log của (2).

      • log2(1+m)log_2(1 + m) \approx mm + δ\delta. với δ\delta được tính toán trung bình bằng 0.0450466.

      Từ đó ta tính được:

      • (2) <=> log2(a)log_2(a) = mam_a + Ea127E_a-127 + δ\delta (3)
  • Điều kỳ diệu của các dãy số:

    Mình sẽ để sẵn một cái ảnh ở đây nhé:

    ps: con số ví dụ của mình sai: 0.15625 chứ k phải là 0.51625

    Việc của chúng ta giờ là phải làm sao loại bỏ bớt các biến đi cho phương trình gọn gàng và dễ tính hơn. Điều gì sẽ xẩy ra khi chúng ta chuyển dãy số biển diễn số phẩy động sang dạng long? Well các bạn thấy trên hình rồi đó:

    • aIa_I = E >> 23 + m >> 23

    Tiện thể thì mình đang dùng phép dịch bit trái, công thức là: M << n = M x 2n2^n

    • Kết hợp với (3) ta có:

    log(a)log(a) = M223\frac{M}{2^{23}} + E223223\frac{E*2^{23}}{2^{23}} -127 + δ\delta = aI223\frac{a_I}{2^{23}} -127 + δ\delta (4)

  • khử log2(a)log_2(a) :

    • Giả sử ta đang cần tính x = y\sqrt{y}, sử dụng log2log_2 cho cả hai vế của phương trình ta sẽ có:

      • log2(x)log_2(x) = 12\frac{1}{2} log2(y)log_2(y)
      • Kết hợp với (3) ta có: xI223\frac{x_I}{2^{23}} -127 + δ\delta = 12\frac{1}{2} (yI223\frac{y_I}{2^{23}} -127 + δ\delta )

    Ok ở đoạn này ta có thể tính đuợc xIx_I rồi convert về xfx_f (x ở dạng số phẩy động) để ra kết quả. Tuy nhiên tại sao ta lại cần tính 1y\frac{1}{\sqrt{y}} ? Vì thuật toán có tên là fast inverse square root chứ còn sao nữa ?

    Đùa đó, thực ra nó liên quan đến xử lý 3D. Để tạo ra được khung hình đẹp, game của bạn sẽ phải xử lý các vấn đề liên quan đến xử lý ánh sáng, đổ bóng, phản chiếu ... các phép xử lý này đều là ứng dụng của các phép tính liên quan đến vector, cụ thể là chuẩn hóa vector. được sử dụng rất rất nhiều lần (hàng triệu lần 1s), Công thức chuẩn hóa vector thì sẽ như thế này:

    A\vec A = (x,y,z) khi chuẩn hóa sẽ thành:

    A||A|| = (xx1(x2+y2+z2)\frac{1}{\sqrt(x^2 + y^2 + z^2)}, yy1(x2+y2+z2)\frac{1}{\sqrt(x^2 + y^2 + z^2)}, zz1(x2+y2+z2)\frac{1}{\sqrt(x^2 + y^2 + z^2)})

    Ảnh: Ref: slideshare

    Phép chuẩn hóa vector cần dùng 1x\frac{1}{\sqrt x} chứ ko phải x\sqrt x. Và đó là lý do tại sao chúng ta sẽ khử log2(a)log_2(a) theo cách khác:

    • Giả sử ta đang cần tính x = 1(y)\frac{1}{\sqrt(y)},đặt log2log_2 cho cả 2 vế của phương trình ta sẽ có:
      • log2(x)log_2(x) = log2(1y)log_2(\frac{1}{ \sqrt{y}}) <=> log2(x)log_2(x) = log2(y12)log_2(y^\frac{-1}{2}) <=> log2(x)log_2(x) = 12\frac{-1}{2} log2(y)log_2(y)
      • Kết hợp với (3) ta có: xI223\frac{x_I}{2^{23}} -127 + δ\delta = 12\frac{-1}{2} (yI223\frac{y_I}{2^{23}} -127 + δ\delta )
      • => xIx_I = 32\frac{3}{2} 2232^{23} (127 - δ\delta) - 12\frac{1}{2} yIy_I
      • 32\frac{3}{2} 2232^{23} (127 - δ\delta) là một hằng số không phụ thuộc vào xIx_I hay yIy_I, Tác giả của thuật toán đã tính ra được giá trị của nó là: 0x5f3759df
      • Kết hợp lại ta có : xIx_I = 0x5f3759df - (yIy_I >> 1) // dịch bit phải để chia hai
  • Kết hợp phương pháp tính Newton Raphson 1 lần nữa:

    • Vì sao? Vì con số 0x5f3759df của chúng ta là giá trị ước lượng dựa theo δ\delta một giá trị tính toán trung bình thế nên chúng ta sẽ sử dụng phương pháp Newton Raphson 1 lần nữa để tìm ra nghiệm chính xác hơn là đủ:
    • Phương trình: xx = 1y\frac{1}{\sqrt y}
    • Áp dụng phương pháp tính Newton Raphson: xx' = xx - 1/x2y2y3\frac{1/x^2 - y}{-2 y^{-3}}
    • => xx' = xx (1.5F - 0.5 x2x^2 yy) (5)
  • Kết hớp từ (1) đến (5) vào code thôi:

static unsafe void Fivs(float number)
{ float yF = number; long iI = *(long*) &yF; // hack ép kiểu từ float sang long long yI = 0x5f3759df - ( iI >> 1 ); float xF = * ( float * ) &yI; // hack ép kiểu từ long sang float float x = xF * (1.5F - 0.5F * yF * xF * xF); float real = (float) (1 / Math.Sqrt(number)); Console.WriteLine(x); Console.WriteLine(real); Console.WriteLine(100*(Math.Abs(x) - real)/real);
}
Fivs(10);

Đây là kết quả: -0.31568578 - 0.31622776 - 0.17139071 Sai số của chúng ta là 0.17% khá là tuyệt vời!! ?

Giờ thì ko ngán bố con nào r nhé ???

Kết luận

Các bạn thấy thế nào, những con số nghe khô khan những lại kết hợp với nhau thật kỳ diệu phải không, hẹn gặp lại các bạn vào những bài viết sau nhé: Bonus cho các bạn về các bàn luận về comment trong code:

  • The reason why "The Evil Bit Hack" is indeed Evil, is because according to the C standard, this is actually Undefined Behavior.
  • "Premature optimisation is the fast inverse square root of all evil" - Donald Knuth, probably

Ref:

Seeder:

Đây là video mình làm về bài viết, các bạn có thể tham khảo thêm, many thanks?

Bình luận

Bài viết tương tự

- vừa được xem lúc

Lập trình viên có nên sử dụng nhiều ngôn ngữ không ?

Câu hỏi kỳ này mà Du code lởm là bạn có hay không nên code nhiều ngôn ngữ nhé hì . .

0 0 5

- vừa được xem lúc

Mình và quá trình khi code một thuật toán Deep Learning

Update 1 chút thì mình đã kiếm được việc rồi ạ, sau đâu đó 10 ngày thất nghiệp Hiện tại mình là AI developer tại Kaopiz Solution nhé mọi người =))). Đợt vừa qua mình có theo đuổi 1 thuật toán và đến h

0 0 25