Mengapa GCC tidak mengoptimalkan a * a * a * a * a * a to (a * a * a) * (a * a * a)?

Saya melakukan beberapa optimasi numerik pada aplikasi ilmiah. Satu hal yang saya perhatikan adalah bahwa GCC akan mengoptimalkan panggilan pow(a,2)dengan mengkompilasinya a*a, tetapi panggilan pow(a,6)tersebut tidak dioptimalkan dan benar-benar akan memanggil fungsi perpustakaan pow, yang sangat memperlambat kinerja. (Sebaliknya, Intel C ++ Compiler , dapat dieksekusi icc, akan menghilangkan panggilan perpustakaan pow(a,6).)

Yang saya ingin tahu adalah ketika saya diganti pow(a,6)dengan a*a*a*a*a*amenggunakan GCC 4.5.1 dan opsi " -O3 -lm -funroll-loops -msse4", ia menggunakan 5 mulsdinstruksi:

movapd  %xmm14, %xmm13
mulsd   %xmm14, %xmm13
mulsd   %xmm14, %xmm13
mulsd   %xmm14, %xmm13
mulsd   %xmm14, %xmm13
mulsd   %xmm14, %xmm13

sedangkan jika saya menulis (a*a*a)*(a*a*a), itu akan menghasilkan

movapd  %xmm14, %xmm13
mulsd   %xmm14, %xmm13
mulsd   %xmm14, %xmm13
mulsd   %xmm13, %xmm13

yang mengurangi jumlah instruksi kalikan ke 3. iccmemiliki perilaku serupa.

Mengapa kompiler tidak mengenali trik pengoptimalan ini?

Karena Floating Point Math bukan asosiatif . Cara Anda mengelompokkan operan dalam perkalian floating point memiliki efek pada akurasi numerik jawabannya.

Akibatnya, sebagian besar kompiler sangat konservatif dalam menyusun ulang perhitungan titik apung kecuali mereka dapat yakin bahwa jawabannya akan tetap sama, atau kecuali jika Anda memberi tahu mereka Anda tidak peduli dengan akurasi numerik. Sebagai contoh: pada -fassociative-mathpilihan dari gcc yang memungkinkan gcc untuk operasi floating point reassociate, atau bahkan -ffast-mathpilihan yang memungkinkan bahkan pengorbanan lebih agresif akurasi terhadap kecepatan.

Lambdageek dengan benar menunjukkan bahwa karena associativity tidak berlaku untuk angka floating-point, "optimasi"a*a*a*a*a*auntuk(a*a*a)*(a*a*a)dapat mengubah nilai. Inilah mengapa itu dilarang oleh C99 (kecuali secara khusus diizinkan oleh pengguna, melalui flag compiler atau pragma). Secara umum, asumsi adalah bahwa programmer menulis apa yang dia lakukan karena suatu alasan, dan kompiler harus menghargai itu. Jika Anda mau(a*a*a)*(a*a*a), tulis itu.

Itu bisa jadi menyusahkan untuk menulis; mengapa kompiler tidak bisa melakukan [apa yang Anda anggap sebagai] hal yang benar ketika Anda gunakan pow(a,6)? Karena itu akan menjadi hal yang salah untuk dilakukan. Pada platform dengan perpustakaan matematika yang baik, pow(a,6)jauh lebih akurat daripada salah satu a*a*a*a*a*aatau (a*a*a)*(a*a*a). Hanya untuk memberikan beberapa data, saya menjalankan percobaan kecil pada Mac Pro saya, mengukur kesalahan terburuk dalam mengevaluasi ^ 6 untuk semua angka mengambang presisi tunggal antara [1,2):

worst relative error using    powf(a, 6.f): 5.96e-08
worst relative error using (a*a*a)*(a*a*a): 2.94e-07
worst relative error using     a*a*a*a*a*a: 2.58e-07

Menggunakan powbukannya pohon multiplikasi mengurangi kesalahan terikat oleh faktor 4 . Kompiler tidak boleh (dan umumnya tidak) membuat "optimisasi" yang meningkatkan kesalahan kecuali diizinkan oleh pengguna (misalnya melalui -ffast-math).

Perhatikan bahwa GCC menyediakan __builtin_powi(x,n)sebagai alternatif pow( ), yang seharusnya menghasilkan pohon multiplikasi sebaris. Gunakan itu jika Anda ingin menukar akurasi untuk kinerja, tetapi tidak ingin mengaktifkan matematika cepat.

Kasus serupa lain: sebagian besar kompiler tidak akan mengoptimalkan a + b + c + duntuk (a + b) + (c + d)(ini adalah optimasi sejak ekspresi kedua dapat pipelined lebih baik) dan mengevaluasinya seperti yang diberikan (yaitu sebagai (((a + b) + c) + d)). Ini juga karena kasus sudut:

float a = 1e35, b = 1e-5, c = -1e35, d = 1e-5;
printf("%e %e\n", a + b + c + d, (a + b) + (c + d));

Ini output 1.000000e-05 0.000000e+00

Fortran (dirancang untuk komputasi ilmiah) memiliki operator daya bawaan, dan sejauh yang saya tahu kompiler Fortran biasanya akan mengoptimalkan peningkatan daya integer dengan cara yang mirip dengan yang Anda gambarkan. C / C ++ sayangnya tidak memiliki operator listrik, hanya fungsi perpustakaan pow(). Ini tidak mencegah kompiler pintar memperlakukan powkhusus dan menghitungnya dengan cara yang lebih cepat untuk kasus khusus, tetapi tampaknya mereka melakukannya lebih jarang ...

Beberapa tahun yang lalu saya mencoba membuatnya lebih nyaman untuk menghitung kekuatan bilangan bulat secara optimal, dan muncul sebagai berikut. Ini C ++, bukan C, dan masih tergantung pada kompiler yang agak pintar tentang cara mengoptimalkan / inline sesuatu. Ngomong-ngomong, semoga bermanfaat dalam praktik:

template<unsigned N> struct power_impl;

template<unsigned N> struct power_impl {
    template<typename T>
    static T calc(const T &x) {
        if (N%2 == 0)
            return power_impl<N/2>::calc(x*x);
        else if (N%3 == 0)
            return power_impl<N/3>::calc(x*x*x);
        return power_impl<N-1>::calc(x)*x;
    }
};

template<> struct power_impl<0> {
    template<typename T>
    static T calc(const T &) { return 1; }
};

template<unsigned N, typename T>
inline T power(const T &x) {
    return power_impl<N>::calc(x);
}

_{Klarifikasi untuk yang penasaran: ini tidak menemukan cara optimal untuk menghitung kekuatan, tetapi karena menemukan solusi optimal adalah masalah NP-complete dan ini hanya layak dilakukan untuk kekuatan kecil (bukan menggunakan pow), tidak ada alasan untuk meributkan dengan detail.}

Kemudian gunakan saja sebagai power<6>(a).

Ini membuatnya mudah untuk mengetikkan kekuatan (tidak perlu menjelaskan 6 adetik dengan parens), dan memungkinkan Anda memiliki optimasi semacam ini tanpa -ffast-mathseandainya Anda memiliki sesuatu yang bergantung pada ketelitian seperti penjumlahan terkompensasi (contoh di mana urutan operasi sangat penting) .

Anda mungkin juga bisa lupa bahwa ini adalah C ++ dan hanya menggunakannya dalam program C (jika dikompilasi dengan kompiler C ++).

Semoga ini bisa bermanfaat.

EDIT:

Ini yang saya dapatkan dari kompiler saya:

Untuk a*a*a*a*a*a,

    movapd  %xmm1, %xmm0
    mulsd   %xmm1, %xmm0
    mulsd   %xmm1, %xmm0
    mulsd   %xmm1, %xmm0
    mulsd   %xmm1, %xmm0
    mulsd   %xmm1, %xmm0

Untuk (a*a*a)*(a*a*a),

    movapd  %xmm1, %xmm0
    mulsd   %xmm1, %xmm0
    mulsd   %xmm1, %xmm0
    mulsd   %xmm0, %xmm0

Untuk power<6>(a),

    mulsd   %xmm0, %xmm0
    movapd  %xmm0, %xmm1
    mulsd   %xmm0, %xmm1
    mulsd   %xmm0, %xmm1

GCC tidak benar-benar mengoptimalkan a*a*a*a*a*auntuk (a*a*a)*(a*a*a)saat adalah bilangan bulat. Saya mencoba dengan perintah ini:

$ echo 'int f(int x) { return x*x*x*x*x*x; }' | gcc -o - -O2 -S -masm=intel -x c -

Ada banyak bendera gcc tetapi tidak ada yang mewah. Maksudnya: Baca dari stdin; gunakan level optimisasi O2; daftar bahasa rakitan keluaran bukan biner; daftar harus menggunakan sintaksis bahasa assembly Intel; input dalam bahasa C (biasanya bahasa disimpulkan dari ekstensi file input, tetapi tidak ada ekstensi file saat membaca dari stdin); dan menulis ke stdout.

Inilah bagian penting dari output. Saya telah menjelaskannya dengan beberapa komentar yang menunjukkan apa yang terjadi dalam bahasa assembly:

; x is in edi to begin with.  eax will be used as a temporary register.
mov  eax, edi  ; temp = x
imul eax, edi  ; temp = x * temp
imul eax, edi  ; temp = x * temp
imul eax, eax  ; temp = temp * temp

Saya menggunakan sistem GCC di Linux Mint 16 Petra, turunan Ubuntu. Ini versi gcc:

$ gcc --version
gcc (Ubuntu/Linaro 4.8.1-10ubuntu9) 4.8.1

Seperti yang telah dicatat oleh poster lain, opsi ini tidak dimungkinkan dalam floating point, karena aritmatika floating point tidak asosiatif.

Karena angka floating-point 32-bit - seperti 1,024 - bukan 1,024. Di komputer, 1.024 adalah interval: dari (1.024-e) ke (1.024 + e), di mana "e" mewakili kesalahan. Beberapa orang gagal untuk menyadari hal ini dan juga percaya bahwa * dalam a * adalah singkatan dari penggandaan angka presisi arbitrer tanpa ada kesalahan yang melekat pada angka-angka itu. Alasan mengapa beberapa orang gagal untuk menyadari ini mungkin adalah perhitungan matematika yang mereka lakukan di sekolah dasar: bekerja hanya dengan angka ideal tanpa kesalahan terpasang, dan percaya bahwa boleh saja mengabaikan "e" saat melakukan penggandaan. Mereka tidak melihat "e" tersirat dalam "float a = 1.2", "a * a * a" dan kode C serupa.

Jika sebagian besar programmer mengenali (dan dapat menjalankan) gagasan bahwa ekspresi C * a * a * a * a * a sebenarnya tidak bekerja dengan angka ideal, maka kompiler GCC kemudian akan GRATIS untuk mengoptimalkan "a * a * a * a * a * a "into say" t = (a * a); t * t * t "yang membutuhkan lebih banyak perkalian. Namun sayangnya, kompiler GCC tidak tahu apakah programmer yang menulis kode berpikir bahwa "a" adalah angka dengan atau tanpa kesalahan. Jadi, GCC hanya akan melakukan seperti apa kode sumbernya - karena itulah yang dilihat GCC dengan "mata telanjang".

... begitu Anda tahu programmer seperti apa Anda , Anda dapat menggunakan tombol "-fast-math" untuk memberi tahu GCC bahwa "Hei, GCC, saya tahu apa yang saya lakukan!". Ini akan memungkinkan GCC untuk mengubah a * a * a * a * a * a menjadi bagian teks yang berbeda - itu terlihat berbeda dari * a * a * a * a * a - tetapi masih menghitung angka dalam interval kesalahan a * a * a * a * a * a. Ini OK, karena Anda sudah tahu Anda bekerja dengan interval, bukan angka ideal.

Belum ada poster yang menyebutkan kontraksi ekspresi mengambang (standar ISO C, 6.5p8 dan 7.12.2). Jika FP_CONTRACTpragma diatur ke ON, kompiler diperbolehkan untuk menganggap ekspresi seperti a*a*a*a*a*aoperasi tunggal, seolah dievaluasi secara tepat dengan pembulatan tunggal. Sebagai contoh, sebuah kompiler dapat menggantinya dengan fungsi daya internal yang lebih cepat dan lebih akurat. Ini sangat menarik karena perilaku sebagian dikendalikan oleh programmer secara langsung dalam kode sumber, sementara opsi kompiler yang disediakan oleh pengguna akhir kadang-kadang dapat digunakan secara tidak benar.

Keadaan standar FP_CONTRACTpragma ditentukan oleh implementasi, sehingga seorang kompiler diizinkan untuk melakukan optimasi tersebut secara default. Jadi kode portabel yang harus benar-benar mengikuti aturan IEEE 754 harus secara eksplisit mengaturnya OFF.

Jika kompiler tidak mendukung pragma ini, ia harus konservatif dengan menghindari optimasi seperti itu, seandainya pengembang telah memilih untuk mengaturnya OFF.

GCC tidak mendukung pragma ini, tetapi dengan opsi default, ia menganggapnya ON; jadi untuk target dengan FMA perangkat keras, jika seseorang ingin mencegah transformasi a*b+cke fma (a, b, c), kita perlu memberikan opsi seperti -ffp-contract=off(untuk secara eksplisit mengatur pragma ke OFF) atau -std=c99(untuk memberitahu GCC agar sesuai dengan beberapa Versi standar C, di sini C99, jadi ikuti paragraf di atas). Di masa lalu, opsi terakhir tidak mencegah transformasi, yang berarti bahwa GCC tidak sesuai pada titik ini: https://gcc.gnu.org/bugzilla/show_bug.cgi?id=37845

Seperti yang ditunjukkan oleh Lambdageek, perkalian float bukan asosiatif dan Anda bisa mendapatkan akurasi yang lebih sedikit, tetapi juga ketika mendapatkan akurasi yang lebih baik, Anda dapat menentang optimasi, karena Anda menginginkan aplikasi deterministik. Misalnya dalam simulasi permainan klien / server, di mana setiap klien harus mensimulasikan dunia yang sama yang Anda inginkan perhitungan floating point menjadi deterministik.

Fungsi perpustakaan seperti "pow" biasanya dibuat dengan hati-hati untuk menghasilkan kesalahan seminimal mungkin (dalam kasus umum). Ini biasanya mencapai fungsi perkiraan dengan splines (menurut komentar Pascal implementasi yang paling umum tampaknya menggunakan algoritma Remez )

secara fundamental operasi berikut:

pow(x,y);

memiliki kesalahan inheren kira-kira sama besarnya dengan kesalahan dalam setiap perkalian atau pembagian tunggal .

Sedangkan operasi berikut:

float a=someValue;
float b=a*a*a*a*a*a;

memiliki kesalahan bawaan yang lebih besar dari 5 kali kesalahan satu perkalian atau pembagian (karena Anda menggabungkan 5 perkalian).

Compiler harus benar-benar berhati-hati dengan jenis optimasi yang dilakukannya:

1. jika mengoptimalkan pow(a,6)untuk a*a*a*a*a*aitu mungkin meningkatkan kinerja, tapi secara drastis mengurangi akurasi untuk angka floating point.
2. jika mengoptimalkan a*a*a*a*a*a untuk pow(a,6)itu benar-benar dapat mengurangi akurasi karena "a" adalah beberapa nilai khusus yang memungkinkan perkalian tanpa kesalahan (kekuatan 2 atau beberapa nomor bilangan bulat kecil)
3. jika mengoptimalkan pow(a,6)ke (a*a*a)*(a*a*a)atau (a*a)*(a*a)*(a*a)ada masih bisa menjadi kehilangan akurasi dibandingkan dengan powfungsi.

Secara umum Anda tahu bahwa untuk nilai floating point arbitrer "pow" memiliki akurasi yang lebih baik daripada fungsi apa pun yang pada akhirnya dapat Anda tulis, tetapi dalam beberapa kasus khusus multiplikasi mungkin memiliki akurasi dan kinerja yang lebih baik, tergantung pada pengembang yang memilih apa yang lebih tepat, akhirnya mengomentari kode sehingga tidak ada orang lain yang akan "mengoptimalkan" kode itu.

Satu-satunya hal yang masuk akal (pendapat pribadi, dan tampaknya pilihan di GCC tanpa optimasi atau kompiler flag tertentu) untuk mengoptimalkan harus mengganti "pow (a, 2)" dengan "a * a". Itu akan menjadi satu-satunya hal waras yang harus dilakukan vendor kompiler.

Saya tidak berharap kasus ini dioptimalkan sama sekali. Tidak dapat sering di mana ekspresi berisi subekspresi yang dapat dikelompokkan ulang untuk menghapus seluruh operasi. Saya berharap penulis kompiler menginvestasikan waktu mereka di bidang-bidang yang lebih mungkin menghasilkan peningkatan yang nyata, daripada menutupi kasus tepi yang jarang ditemui.

Saya terkejut mengetahui dari jawaban lain bahwa ungkapan ini memang dapat dioptimalkan dengan sakelar kompiler yang tepat. Entah optimasi itu sepele, atau itu adalah kasus tepi dari optimasi yang jauh lebih umum, atau penulis kompiler sangat teliti.

Tidak ada yang salah dengan memberikan petunjuk kepada kompiler seperti yang Anda lakukan di sini. Itu adalah bagian normal dan yang diharapkan dari proses optimasi mikro untuk mengatur ulang pernyataan dan ekspresi untuk melihat perbedaan apa yang akan mereka bawa.

Sementara kompiler dapat dibenarkan dalam mempertimbangkan dua ekspresi untuk memberikan hasil yang tidak konsisten (tanpa saklar yang tepat), Anda tidak perlu terikat oleh batasan itu. Perbedaannya akan sangat kecil - begitu banyak sehingga jika perbedaan itu penting bagi Anda, Anda seharusnya tidak menggunakan aritmatika floating point standar di tempat pertama.

Sudah ada beberapa jawaban yang baik untuk pertanyaan ini, tetapi demi kelengkapan saya ingin menunjukkan bahwa bagian yang berlaku dari standar C adalah 5.1.2.2.3 / 15 (yang sama dengan bagian 1.9 / 9 dalam C ++ 11 standar). Bagian ini menyatakan bahwa operator hanya dapat dikelompokkan kembali jika mereka benar-benar asosiatif atau komutatif.

gcc sebenarnya dapat melakukan optimasi ini, bahkan untuk angka floating-point. Sebagai contoh,

double foo(double a) {
  return a*a*a*a*a*a;
}

menjadi

foo(double):
    mulsd   %xmm0, %xmm0
    movapd  %xmm0, %xmm1
    mulsd   %xmm0, %xmm1
    mulsd   %xmm1, %xmm0
    ret

dengan -O -funsafe-math-optimizations. Penataan ulang ini melanggar IEEE-754, jadi, itu membutuhkan flag.

Bilangan bulat yang ditandatangani, seperti yang ditunjukkan oleh Peter Cordes dalam komentar, dapat melakukan pengoptimalan ini tanpa -funsafe-math-optimizationssejak itu berlaku tepat ketika tidak ada overflow dan jika ada overflow Anda mendapatkan perilaku yang tidak terdefinisi. Jadi kamu mengerti

foo(long):
    movq    %rdi, %rax
    imulq   %rdi, %rax
    imulq   %rdi, %rax
    imulq   %rax, %rax
    ret

dengan adil -O. Untuk bilangan bulat tak bertanda, itu bahkan lebih mudah karena mereka bekerja dengan kekuatan mod 2 dan sehingga dapat disusun ulang secara bebas bahkan dalam menghadapi overflow.