Có thể cuộn phiên bản nhanh hơn đáng kể của sqrt

Question

Trong ứng dụng tôi đang lược tả, tôi thấy rằng trong một số trường hợp, chức năng này có thể chiếm hơn 10% tổng thời gian thực hiện.

Tôi đã nhìn thấy cuộc thảo luận qua nhiều năm triển khai sqrt nhanh hơn bằng cách sử dụng thủ thuật lén lút trôi nổi, nhưng tôi không biết liệu những thứ đó có lỗi thời trên các CPU hiện đại hay không.

Trình biên dịch MSVC++ 2008 đang được sử dụng, để tham khảo ... mặc dù tôi muốn giả định sqrt sẽ không thêm nhiều chi phí mặc dù.

Xem thêm tại đây để thảo luận tương tự về chức năng modf.

EDIT: để tham khảo, this là một phương pháp được sử dụng rộng rãi, nhưng thực sự có nhanh hơn nhiều không? Có bao nhiêu chu kỳ là SQRT vào những ngày này?

Answer 1

Vâng, nó có thể thậm chí không có thủ đoạn gian trá:

1) hy sinh chính xác cho tốc độ: các thuật toán sqrt được lặp đi lặp lại, tái thực hiện với sự lặp lại ít hơn.

2) bảng tra cứu: chỉ dành cho điểm bắt đầu của lần lặp lại hoặc kết hợp với nội suy để giúp bạn có tất cả các cách ở đó.

3) lưu vào bộ nhớ cache: bạn luôn luôn sqrting cùng một tập hạn chế các giá trị? nếu vậy, bộ nhớ đệm có thể hoạt động tốt. Tôi đã tìm thấy điều này hữu ích trong các ứng dụng đồ họa mà cùng một thứ đang được tính toán cho nhiều hình dạng có cùng kích thước, vì vậy kết quả có thể được lưu trữ một cách hữu ích.

Answer 2

Bạn rất có khả năng để đạt được những cải thiện tốc độ nhanh hơn bằng cách thay đổi thuật toán của bạn hơn bằng cách thay đổi họ triển khai: Cố gắng gọi sqrt() ít thay vì thực hiện cuộc gọi nhanh hơn. (Và nếu bạn nghĩ rằng điều này là không thể - những cải tiến cho sqrt() bạn đề cập chỉ là: các cải tiến của thuật toán được sử dụng để tính căn bậc hai.)

Vì nó được sử dụng rất thường xuyên, có khả năng việc triển khai thư viện chuẩn của bạn là sqrt() gần như là tối ưu cho trường hợp chung. Trừ khi bạn có miền bị hạn chế (ví dụ: nếu bạn cần ít chính xác hơn), nơi thuật toán có thể thực hiện một số phím tắt, rất khó có ai đó đưa ra một triển khai nhanh hơn. Lưu ý rằng, vì chức năng đó sử dụng 10% thời gian thực hiện của bạn, ngay cả khi bạn quản lý để triển khai thực hiện chỉ mất 75% thời gian std::sqrt(), điều này sẽ chỉ mang lại thời gian thực hiện của bạn xuống 2,5%. Đối với hầu hết các ứng dụng, người dùng thậm chí sẽ không nhận thấy điều này, trừ khi họ sử dụng đồng hồ để đo lường.

Answer 3

Bạn cần số sqrt chính xác đến mức nào? Bạn có thể nhận được xấp xỉ hợp lý rất nhanh chóng: xem chức năng tuyệt vời của Quake3 inverse square root cho cảm hứng (lưu ý rằng mã được GPL'ed, vì vậy bạn có thể không muốn tích hợp nó trực tiếp).

Answer 4

Có một bảng so sánh tuyệt vời ở đây: http://assemblyrequired.crashworks.org/timing-square-root/

câu chuyện dài ngắn, ssqrts SSE2 là khoảng 2x nhanh hơn so với FPU fsqrt, và một xấp xỉ + lặp là khoảng 4x nhanh hơn (8x tổng thể).

Ngoài ra, nếu bạn đang cố gắng thực hiện một sqrt chính xác đơn, hãy đảm bảo đó thực sự là những gì bạn đang nhận được. Tôi đã nghe nói về ít nhất một trình biên dịch có thể chuyển đổi đối số float thành một double, gọi sqrt có độ chính xác kép, sau đó chuyển trở lại float.

Answer 5

Không biết nếu bạn sửa lỗi này, nhưng tôi đã đọc về nó trước đây, và có vẻ như điều nhanh nhất cần làm là thay thế hàm sqrt bằng phiên bản lắp ráp nội tuyến;

bạn có thể xem mô tả về tải thay thế here.

Điều tốt nhất là đoạn này của ma thuật:

double inline __declspec (naked) __fastcall sqrt(double n) 
{ 
    _asm fld qword ptr [esp+4] 
    _asm fsqrt 
    _asm ret 8 
} 

Đó là về 4.7x nhanh hơn so với tiêu chuẩn sqrt gọi với độ chính xác tương tự.

Answer 6

Đây là cách nhanh chóng với bảng tra cứu chỉ có 8KB. Sai lầm là ~ 0,5% kết quả. Bạn có thể dễ dàng mở rộng bảng, do đó làm giảm sai lầm. Chạy nhanh hơn khoảng 5 lần so với sqrt thông thường()

// LUT for fast sqrt of floats. Table will be consist of 2 parts, half for sqrt(X) and half for sqrt(2X). 
const int nBitsForSQRTprecision = 11;      // Use only 11 most sagnificant bits from the 23 of float. We can use 15 bits instead. It will produce less error but take more place in a memory. 
const int nUnusedBits = 23 - nBitsForSQRTprecision;  // Amount of bits we will disregard 
const int tableSize  = (1 << (nBitsForSQRTprecision+1)); // 2^nBits*2 because we have 2 halves of the table. 
static short sqrtTab[tableSize]; 
static unsigned char is_sqrttab_initialized = FALSE;  // Once initialized will be true 

// Table of precalculated sqrt() for future fast calculation. Approximates the exact with an error of about 0.5% 
// Note: To access the bits of a float in C quickly we must misuse pointers. 
// More info in: http://en.wikipedia.org/wiki/Single_precision 
void build_fsqrt_table(void){ 
    unsigned short i; 
    float f; 
    UINT32 *fi = (UINT32*)&f; 

    if (is_sqrttab_initialized) 
     return; 

    const int halfTableSize = (tableSize>>1); 
    for (i=0; i < halfTableSize; i++){ 
     *fi = 0; 
     *fi = (i << nUnusedBits) | (127 << 23); // Build a float with the bit pattern i as mantissa, and an exponent of 0, stored as 127 

     // Take the square root then strip the first 'nBitsForSQRTprecision' bits of the mantissa into the table 
     f = sqrtf(f); 
     sqrtTab[i] = (short)((*fi & 0x7fffff) >> nUnusedBits); 

     // Repeat the process, this time with an exponent of 1, stored as 128 
     *fi = 0; 
     *fi = (i << nUnusedBits) | (128 << 23); 
     f = sqrtf(f); 
     sqrtTab[i+halfTableSize] = (short)((*fi & 0x7fffff) >> nUnusedBits); 
    } 
    is_sqrttab_initialized = TRUE; 
} 

// Calculation of a square root. Divide the exponent of float by 2 and sqrt() its mantissa using the precalculated table. 
float fast_float_sqrt(float n){ 
    if (n <= 0.f) 
     return 0.f;       // On 0 or negative return 0. 
    UINT32 *num = (UINT32*)&n; 
    short e;         // Exponent 
    e = (*num >> 23) - 127;     // In 'float' the exponent is stored with 127 added. 
    *num &= 0x7fffff;       // leave only the mantissa 

    // If the exponent is odd so we have to look it up in the second half of the lookup table, so we set the high bit. 
    const int halfTableSize = (tableSize>>1); 
    const int secondHalphTableIdBit = halfTableSize << nUnusedBits; 
    if (e & 0x01) 
     *num |= secondHalphTableIdBit; 
    e >>= 1;         // Divide the exponent by two (note that in C the shift operators are sign preserving for signed operands 

    // Do the table lookup, based on the quaternary mantissa, then reconstruct the result back into a float 
    *num = ((sqrtTab[*num >> nUnusedBits]) << nUnusedBits) | ((e + 127) << 23); 
    return n; 
}