Bagaimana C menghitung sin () dan fungsi matematika lainnya?

Saya telah meneliti melalui .NET disassemblies dan kode sumber GCC, tetapi tampaknya tidak dapat menemukan di mana saja implementasi aktual sin()dan fungsi matematika lainnya ... mereka sepertinya selalu merujuk pada hal lain.

Adakah yang bisa membantu saya menemukannya? Saya merasa sepertinya tidak mungkin SEMUA perangkat keras yang akan dijalankan oleh C mendukung fungsi trigonometri di perangkat keras, jadi pasti ada algoritma perangkat lunak di suatu tempat , kan?

Saya mengetahui beberapa cara fungsi dapat dihitung, dan telah menulis rutinitas saya sendiri untuk menghitung fungsi menggunakan seri taylor untuk bersenang-senang. Saya ingin tahu tentang bagaimana nyata, bahasa produksi melakukannya, karena semua implementasi saya selalu beberapa kali lipat lebih lambat, meskipun saya pikir algoritma saya cukup pintar (jelas mereka tidak pintar).

Dalam GNU libm, implementasi sinbergantung pada sistem. Oleh karena itu Anda dapat menemukan implementasinya, untuk setiap platform, di suatu tempat di subdirektori sysdeps yang sesuai .

Satu direktori termasuk implementasi dalam C, disumbangkan oleh IBM. Sejak Oktober 2011, ini adalah kode yang benar-benar berjalan ketika Anda memanggil sin()sistem Linux x86-64. Ini tampaknya lebih cepat daripada fsininstruksi perakitan. Kode sumber: sysdeps / ieee754 / dbl-64 / s_sin.c , cari __sin (double x).

Kode ini sangat kompleks. Tidak ada algoritma perangkat lunak yang secepat mungkin dan juga akurat pada seluruh jajaran nilai x , sehingga perpustakaan mengimplementasikan beberapa algoritma berbeda, dan tugas pertamanya adalah melihat x dan memutuskan algoritma mana yang akan digunakan.

• Ketika x sangat sangat dekat dengan 0, sin(x) == xadalah jawaban yang tepat.

• Sedikit lebih jauh, sin(x)menggunakan seri Taylor yang akrab. Namun, ini hanya akurat di dekat 0, jadi ...

• Ketika sudut lebih dari sekitar 7 °, algoritma yang berbeda digunakan, menghitung pendekatan deret Taylor untuk sin (x) dan cos (x), kemudian menggunakan nilai dari tabel yang dikomputasi untuk menyempurnakan aproksimasi.

• Kapan | x | > 2, tak satu pun dari algoritma di atas akan berfungsi, jadi kode dimulai dengan menghitung beberapa nilai lebih dekat dengan 0 yang dapat diumpankan ke sinatau cossebagai gantinya.

• Masih ada cabang lain untuk berurusan dengan x menjadi NaN atau infinity.

Kode ini menggunakan beberapa peretasan angka yang belum pernah saya lihat sebelumnya, meskipun untuk semua yang saya tahu mereka mungkin terkenal di kalangan ahli floating-point. Terkadang beberapa baris kode membutuhkan beberapa paragraf untuk menjelaskan. Misalnya, dua baris ini

double t = (x * hpinv + toint);
double xn = t - toint;

digunakan (kadang-kadang) dalam mereduksi x ke nilai mendekati 0 yang berbeda dari x dengan kelipatan 2/2, khususnya xn× π / 2. Cara ini dilakukan tanpa pembagian atau bercabang agak pintar. Tapi tidak ada komentar sama sekali!

GCC / glibc versi 32-bit yang lebih lama menggunakan fsininstruksi ini, yang secara mengejutkan tidak akurat untuk beberapa input. Ada posting blog menarik yang menggambarkan ini hanya dengan 2 baris kode .

Implementasi fdlibm sindalam murni C jauh lebih sederhana daripada glibc dan berkomentar dengan baik. Kode sumber: fdlibm / s_sin.c dan fdlibm / k_sin.c

Fungsi seperti sinus dan cosinus diimplementasikan dalam mikrokode di dalam mikroprosesor. Chip Intel, misalnya, memiliki instruksi perakitan untuk ini. Kompiler AC akan menghasilkan kode yang memanggil instruksi perakitan ini. (Sebaliknya, kompiler Java tidak akan. Java mengevaluasi fungsi trigonometri dalam perangkat lunak daripada perangkat keras, dan karenanya berjalan jauh lebih lambat.)

Chip tidak menggunakan seri Taylor untuk menghitung fungsi trigonometri, setidaknya tidak seluruhnya. Pertama-tama mereka menggunakan CORDIC , tetapi mereka juga dapat menggunakan seri Taylor pendek untuk memoles hasil CORDIC atau untuk kasus khusus seperti menghitung sinus dengan akurasi relatif tinggi untuk sudut yang sangat kecil. Untuk penjelasan lebih lanjut, lihat jawaban StackOverflow ini .

OK kiddies, waktu untuk pro .... Ini adalah salah satu keluhan terbesar saya dengan insinyur perangkat lunak yang tidak berpengalaman. Mereka datang dalam menghitung fungsi transendental dari awal (menggunakan seri Taylor) seolah-olah tidak ada yang pernah melakukan perhitungan ini sebelumnya dalam hidup mereka. Tidak benar. Ini adalah masalah yang didefinisikan dengan baik dan telah didekati ribuan kali oleh insinyur perangkat lunak dan perangkat keras yang sangat pintar dan memiliki solusi yang terdefinisi dengan baik. Pada dasarnya, sebagian besar fungsi transendental menggunakan Polinomial Chebyshev untuk menghitungnya. Adapun polinomial yang digunakan tergantung pada keadaan. Pertama, Alkitab tentang hal ini adalah sebuah buku yang disebut "Perkiraan Komputer" oleh Hart dan Cheney. Dalam buku itu, Anda dapat memutuskan apakah Anda memiliki penambah perangkat keras, pengganda, pembagi, dll, dan memutuskan operasi mana yang tercepat. mis. Jika Anda memiliki pembagi yang sangat cepat, cara tercepat untuk menghitung sinus mungkin P1 (x) / P2 (x) di mana P1, P2 adalah polinomial Chebyshev. Tanpa pembagi cepat, mungkin hanya P (x), di mana P memiliki lebih banyak istilah daripada P1 atau P2 .... jadi itu akan lebih lambat. Jadi, langkah pertama adalah menentukan perangkat keras Anda dan apa yang dapat dilakukannya. Kemudian Anda memilih kombinasi yang sesuai dari polinomial Chebyshev (biasanya dari bentuk cos (ax) = aP (x) untuk cosinus misalnya, sekali lagi di mana P adalah polinomial Chebyshev). Maka Anda memutuskan presisi desimal apa yang Anda inginkan. misalnya jika Anda ingin presisi 7 digit, Anda melihat bahwa dalam tabel yang sesuai dalam buku yang saya sebutkan, dan itu akan memberi Anda (untuk presisi = 7,33) angka N = 4 dan angka polinomial 3502. N adalah urutan polinomial (jadi p4.x ^ 4 + p3.x ^ 3 + p2.x ^ 2 + p1.x + p0), karena N = 4. Kemudian Anda mencari nilai aktual dari p4, p3, p2, p1, nilai p0 di bagian belakang buku di bawah 3502 (mereka akan berada di floating point). Kemudian Anda mengimplementasikan algoritma Anda dalam perangkat lunak dalam bentuk: (((p4.x + p3) .x + p2) .x + p1) .x + p0 .... dan ini adalah bagaimana Anda akan menghitung cosinus hingga 7 desimal tempat di perangkat keras itu.

Perhatikan bahwa sebagian besar implementasi perangkat keras dari operasi transendental dalam FPU biasanya melibatkan beberapa mikrokode dan operasi seperti ini (tergantung pada perangkat kerasnya). Polinomial Chebyshev digunakan untuk sebagian besar transendental tetapi tidak semua. misal akar kuadrat lebih cepat menggunakan iterasi ganda dari metode Newton raphson menggunakan tabel pencarian terlebih dahulu. Sekali lagi, buku "Perkiraan Komputer" itu akan memberi tahu Anda hal itu.

Jika Anda berencana untuk mengimplementasikan fungsi-fungsi ini, saya akan merekomendasikan kepada siapa pun bahwa mereka mendapatkan salinan buku itu. Ini benar-benar Alkitab untuk jenis algoritma ini. Perhatikan bahwa ada banyak cara alternatif untuk menghitung nilai-nilai ini seperti cordics, dll, tetapi ini cenderung yang terbaik untuk algoritma tertentu di mana Anda hanya membutuhkan presisi rendah. Untuk menjamin ketepatan setiap waktu, polinomial chebyshev adalah jalan yang harus ditempuh. Seperti saya katakan, masalah didefinisikan dengan baik. Telah dipecahkan selama 50 tahun sekarang ..... dan begitulah caranya.

Sekarang, yang dikatakan, ada teknik dimana polinomial Chebyshev dapat digunakan untuk mendapatkan hasil presisi tunggal dengan polinomial derajat rendah (seperti contoh untuk cosinus di atas). Kemudian, ada teknik lain untuk menginterpolasi antara nilai-nilai untuk meningkatkan akurasi tanpa harus pergi ke polinomial yang jauh lebih besar, seperti "Gal's Accurate Tables Method". Teknik yang terakhir ini adalah apa yang merujuk posting ke literatur ACM. Tetapi pada akhirnya, Polinomial Chebyshev adalah apa yang digunakan untuk mendapatkan 90% dari perjalanan ke sana.

Nikmati.

Untuk sinkhusus, menggunakan ekspansi Taylor akan memberi Anda:

sin (x): = x - x ^ 3/3! + x ^ 5/5! - x ^ 7/7! + ... (1)

Anda akan terus menambahkan istilah sampai selisih di antara keduanya lebih rendah dari tingkat toleransi yang diterima atau hanya untuk jumlah langkah yang terbatas (lebih cepat, tetapi kurang tepat). Contohnya akan seperti:

float sin(float x)
{
  float res=0, pow=x, fact=1;
  for(int i=0; i<5; ++i)
  {
    res+=pow/fact;
    pow*=-1*x*x;
    fact*=(2*(i+1))*(2*(i+1)+1);
  }

  return res;
}

Catatan: (1) berfungsi karena dosa aproximation (x) = x untuk sudut kecil. Untuk sudut yang lebih besar, Anda perlu menghitung lebih banyak istilah untuk mendapatkan hasil yang dapat diterima. Anda dapat menggunakan argumen sementara dan melanjutkan untuk akurasi tertentu:

double sin (double x){
    int i = 1;
    double cur = x;
    double acc = 1;
    double fact= 1;
    double pow = x;
    while (fabs(acc) > .00000001 &&   i < 100){
        fact *= ((2*i)*(2*i+1));
        pow *= -1 * x*x; 
        acc =  pow / fact;
        cur += acc;
        i++;
    }
    return cur;

}

Ya, ada algoritma perangkat lunak untuk menghitung sinjuga. Pada dasarnya, menghitung barang-barang semacam ini dengan komputer digital biasanya dilakukan dengan menggunakan metode numerik seperti mendekati deret Taylor yang mewakili fungsi tersebut.

Metode numerik dapat memperkirakan fungsi ke jumlah akurasi sembarang dan karena jumlah akurasi yang Anda miliki dalam angka mengambang terbatas, mereka sesuai dengan tugas-tugas ini dengan cukup baik.

Gunakan seri Taylor dan cobalah untuk menemukan hubungan antara syarat-syarat seri sehingga Anda tidak menghitung hal-hal berulang kali

Ini adalah contoh untuk cosinus:

double cosinus(double x, double prec)
{
    double t, s ;
    int p;
    p = 0;
    s = 1.0;
    t = 1.0;
    while(fabs(t/s) > prec)
    {
        p++;
        t = (-t * x * x) / ((2 * p - 1) * (2 * p));
        s += t;
    }
    return s;
}

menggunakan ini kita bisa mendapatkan istilah baru dari jumlah menggunakan yang sudah digunakan (kita menghindari faktorial dan x ^2p )

Ini adalah pertanyaan yang kompleks. CPU mirip Intel dari keluarga x86 memiliki implementasi sin()fungsi perangkat keras , tetapi merupakan bagian dari FPU x87 dan tidak digunakan lagi dalam mode 64-bit (di mana register SSE2 digunakan sebagai gantinya). Dalam mode itu, implementasi perangkat lunak digunakan.

Ada beberapa implementasi seperti itu di luar sana. Satu di fdlibm dan digunakan di Jawa. Sejauh yang saya tahu, implementasi glibc berisi bagian dari fdlibm, dan bagian lain yang dikontribusikan oleh IBM.

Implementasi perangkat lunak dari fungsi transendental seperti sin()biasanya menggunakan perkiraan oleh polinomial, sering diperoleh dari seri Taylor.

Polinomial Chebyshev, sebagaimana disebutkan dalam jawaban lain, adalah polinomial di mana perbedaan terbesar antara fungsi dan polinomial sekecil mungkin. Itu awal yang bagus.

Dalam beberapa kasus, kesalahan maksimum bukanlah apa yang Anda minati, tetapi kesalahan relatif maksimum. Misalnya untuk fungsi sinus, kesalahan dekat x = 0 harus jauh lebih kecil daripada untuk nilai yang lebih besar; Anda ingin kesalahan relatif kecil . Jadi Anda akan menghitung polinomial Chebyshev untuk sin x / x, dan kalikan polinomial dengan x.

Selanjutnya Anda harus mencari cara untuk mengevaluasi polinomial. Anda ingin mengevaluasinya sedemikian rupa sehingga nilai-nilai perantara kecil dan oleh karena itu kesalahan pembulatan kecil. Kalau tidak, kesalahan pembulatan mungkin menjadi jauh lebih besar daripada kesalahan dalam polinomial. Dan dengan fungsi seperti fungsi sinus, jika Anda ceroboh maka mungkin saja hasil yang Anda hitung untuk dosa x lebih besar daripada hasil untuk dosa y bahkan ketika x <y. Jadi pilihan hati-hati dari urutan perhitungan dan perhitungan batas atas untuk kesalahan pembulatan diperlukan.

Misalnya, sin x = x - x ^ 3/6 + x ^ 5/120 - x ^ 7/5040 ... Jika Anda menghitung secara sinis x = x * (1 - x ^ 2/6 + x ^ 4 / 120 - x ^ 6/5040 ...), maka fungsi dalam tanda kurung menurun, dan itu akan terjadi bahwa jika y adalah angka berikutnya yang lebih besar untuk x, maka kadang-kadang sin y akan lebih kecil dari sin x. Sebaliknya, hitung sin x = x - x ^ 3 * (1/6 - x ^ 2/120 + x ^ 4/5040 ...) di mana ini tidak dapat terjadi.

Saat menghitung polinomial Chebyshev, Anda biasanya perlu membulatkan koefisien menjadi presisi ganda, misalnya. Tetapi sementara polinomial Chebyshev optimal, polinomial Chebyshev dengan koefisien dibulatkan menjadi presisi ganda bukanlah polinomial optimal dengan koefisien presisi ganda!

Misalnya untuk sin (x), di mana Anda memerlukan koefisien untuk x, x ^ 3, x ^ 5, x ^ 7 dll. Anda melakukan hal berikut: Hitung perkiraan sin terbaik x dengan polinomial (ax + bx ^ 3 + cx ^ 5 + dx ^ 7) dengan presisi lebih tinggi dari dua kali lipat, kemudian putaran ke presisi ganda, memberikan A. Perbedaan antara a dan A akan cukup besar. Sekarang hitung perkiraan terbaik (sin x - Axe) dengan polinomial (bx ^ 3 + cx ^ 5 + dx ^ 7). Anda mendapatkan koefisien yang berbeda, karena mereka beradaptasi dengan perbedaan antara a dan A. Putaran b untuk menggandakan presisi B. Kemudian perkiraan (sin x - Ax - Bx ^ 3) dengan polinomial cx ^ 5 + dx ^ 7 dan seterusnya. Anda akan mendapatkan polinomial yang hampir sebagus polinomial Chebyshev asli, tetapi jauh lebih baik daripada Chebyshev dibulatkan menjadi presisi ganda.

Selanjutnya Anda harus memperhitungkan kesalahan pembulatan dalam pemilihan polinomial. Anda menemukan polinomial dengan kesalahan minimum di polinomial mengabaikan kesalahan pembulatan, tetapi Anda ingin mengoptimalkan polinomial plus kesalahan pembulatan. Setelah Anda memiliki polinomial Chebyshev, Anda dapat menghitung batas untuk kesalahan pembulatan. Katakanlah f (x) adalah fungsi Anda, P (x) adalah polinomial, dan E (x) adalah kesalahan pembulatan. Anda tidak ingin mengoptimalkan | f (x) - P (x) |, Anda ingin mengoptimalkan | f (x) - P (x) +/- E (x) | Anda akan mendapatkan polinomial yang sedikit berbeda yang mencoba untuk menjaga kesalahan polinomial di mana kesalahan pembulatan besar, dan sedikit merelaksasi kesalahan polinomial di mana kesalahan pembulatan kecil.

Semua ini akan membuat Anda dengan mudah membulatkan kesalahan sebanyak 0,55 kali bit terakhir, di mana +, -, *, / memiliki kesalahan pembulatan paling banyak 0,50 kali bit terakhir.

Mengenai fungsi trigonometri seperti sin(), cos(), tan()belum ada menyebutkan, setelah 5 tahun, dari aspek penting dari fungsi trigonometri kualitas tinggi: Rentang pengurangan .

Langkah awal dalam fungsi-fungsi ini adalah untuk mengurangi sudut, dalam radian, hingga rentang interval 2 * π. Tetapi π tidak rasional, pengurangan sederhana seperti x = remainder(x, 2*M_PI)memperkenalkan kesalahan sebagai M_PI, atau pi mesin, adalah perkiraan π. Jadi, bagaimana caranya x = remainder(x, 2*π)?

Pustaka awal menggunakan ketelitian yang diperluas atau pemrograman yang dibuat untuk memberikan hasil yang berkualitas tetapi masih dalam kisaran terbatas double. Ketika nilai besar diminta seperti sin(pow(2,30)), hasilnya tidak berarti atau 0.0dan mungkin dengan flag kesalahan diatur ke sesuatu seperti TLOSSkehilangan total presisi atau PLOSShilangnya sebagian presisi.

Pengurangan rentang nilai besar yang baik ke interval seperti -π ke π adalah masalah yang menantang yang menyaingi tantangan fungsi trigonometri dasar, seperti sin(), itu sendiri.

Laporan yang bagus adalah reduksi argumen untuk argumen besar: Baik sampai akhir (1992). Ini mencakup masalah dengan baik: membahas kebutuhan dan bagaimana hal-hal pada berbagai platform (SPARC, PC, HP, 30+ lainnya) dan menyediakan algoritma solusi yang memberikan hasil berkualitas untuk semua double dari -DBL_MAXhingga DBL_MAX.

Jika argumen asli dalam derajat, namun mungkin bernilai besar, gunakan fmod()dulu untuk meningkatkan presisi. Barang fmod()akan memperkenalkan tidak ada kesalahan dan karenanya memberikan pengurangan rentang yang sangat baik.

// sin(degrees2radians(x))
sin(degrees2radians(fmod(x, 360.0))); // -360.0 < fmod(x,360) < +360.0

Berbagai identitas pemicu dan remquo()menawarkan peningkatan lebih besar lagi. Contoh: sind ()

Implementasi fungsi perpustakaan yang sebenarnya tergantung pada kompiler dan / atau penyedia perpustakaan tertentu. Apakah itu dilakukan dalam perangkat keras atau perangkat lunak, apakah itu ekspansi Taylor atau tidak, dll., Akan bervariasi.

Saya sadar itu sama sekali tidak membantu.

Mereka biasanya diimplementasikan dalam perangkat lunak dan tidak akan menggunakan perangkat keras yang sesuai (yaitu, biasanya) panggilan dalam banyak kasus. Namun, seperti yang ditunjukkan Jason, ini adalah implementasi spesifik.

Perhatikan bahwa rutin perangkat lunak ini bukan bagian dari sumber kompiler, tetapi lebih suka ditemukan di perpustakaan correspoding seperti clib, atau glibc untuk kompilator GNU. Lihat http://www.gnu.org/software/libc/manual/html_mono/libc.html#Trig-Functions

Jika Anda ingin kontrol yang lebih besar, Anda harus hati-hati mengevaluasi apa yang sebenarnya Anda butuhkan. Beberapa metode khas adalah interpolasi dari tabel pencarian, panggilan rakitan (yang sering lambat), atau skema aproksimasi lainnya seperti Newton-Raphson untuk akar kuadrat.

Jika Anda menginginkan implementasi dalam perangkat lunak, bukan perangkat keras, tempat untuk mencari jawaban pasti untuk pertanyaan ini adalah Bab 5 dari Resep Numerik . Salinan saya ada di dalam kotak, jadi saya tidak bisa memberikan detail, tetapi versi singkatnya (jika saya ingat ini benar) adalah bahwa Anda menganggap tan(theta/2)operasi primitif Anda dan menghitung yang lain dari sana. Perhitungan dilakukan dengan aproksimasi deret, tetapi ini adalah sesuatu yang menyatu jauh lebih cepat daripada deret Taylor.

Maaf saya tidak bisa mengingat lagi tanpa mendapatkan buku itu.

Tidak ada yang seperti memukul sumber dan melihat bagaimana seseorang benar-benar melakukannya di perpustakaan yang umum digunakan; mari kita lihat satu implementasi C library pada khususnya. Saya memilih uLibC.

Inilah fungsi dosa:

http://git.uclibc.org/uClibc/tree/libm/s_sin.c

yang kelihatannya menangani beberapa kasus khusus, dan kemudian melakukan pengurangan argumen untuk memetakan input ke kisaran [-pi / 4, pi / 4], (membelah argumen menjadi dua bagian, bagian besar dan ekor) sebelum menelepon

http://git.uclibc.org/uClibc/tree/libm/k_sin.c

yang kemudian beroperasi pada dua bagian tersebut. Jika tidak ada ekor, jawaban perkiraan dihasilkan menggunakan polinomial derajat 13. Jika ada ekor, Anda mendapatkan tambahan korektif kecil berdasarkan prinsip bahwasin(x+y) = sin(x) + sin'(x')y

Setiap kali fungsi tersebut dievaluasi, maka pada tingkat tertentu kemungkinan besar ada:

• Tabel nilai yang diinterpolasi (untuk aplikasi yang cepat dan tidak akurat - misalnya grafik komputer)
• Evaluasi seri yang konvergen ke nilai yang diinginkan --- mungkin bukan seri taylor, lebih mungkin sesuatu yang didasarkan pada quadrature mewah seperti Clenshaw-Curtis.

Jika tidak ada dukungan perangkat keras maka kompiler mungkin menggunakan metode yang terakhir, hanya memancarkan kode assembler (tanpa simbol debug), daripada menggunakan pustaka ac --- membuatnya sulit bagi Anda untuk melacak kode aktual di debugger Anda.

Seperti yang ditunjukkan oleh banyak orang, ini tergantung pada implementasi. Tapi sejauh yang saya mengerti pertanyaan Anda, Anda tertarik pada implementasi perangkat lunak nyata dari fungsi matematika, tetapi tidak berhasil menemukannya. Jika demikian, maka di sini Anda berada:

• Unduh kode sumber glibc dari http://ftp.gnu.org/gnu/glibc/
• Lihat file yang dosincos.cterletak di folder glibc root \ sysdeps \ ieee754 \ dbl-64
• Demikian pula Anda dapat menemukan implementasi dari sisa perpustakaan matematika, hanya mencari file dengan nama yang sesuai

Anda juga dapat melihat file-file dengan .tblekstensi, isinya tidak lebih dari tabel besar dari nilai-nilai yang dikomputasi dari berbagai fungsi dalam bentuk biner. Itulah sebabnya implementasinya sangat cepat: alih-alih menghitung semua koefisien dari seri apa pun yang mereka gunakan, mereka hanya melakukan pencarian cepat, yang jauh lebih cepat. BTW, mereka menggunakan seri Penjahit untuk menghitung sinus dan cosinus.

Saya harap ini membantu.

Saya akan mencoba menjawab untuk kasus sin()dalam program C, dikompilasi dengan kompiler C GCC pada prosesor x86 saat ini (katakanlah Intel Core 2 Duo).

Dalam bahasa C Perpustakaan C Standar mencakup fungsi matematika umum, tidak termasuk dalam bahasa itu sendiri (misalnya pow, sindan cosuntuk daya, sinus, dan cosinus masing-masing). Header yang termasuk dalam math.h .

Sekarang pada sistem GNU / Linux, fungsi perpustakaan ini disediakan oleh glibc (GNU libc atau GNU C Library). Tetapi kompiler GCC ingin Anda menautkan ke perpustakaan matematika ( libm.so) menggunakan -lmflag kompiler untuk memungkinkan penggunaan fungsi-fungsi matematika ini. Saya tidak yakin mengapa itu bukan bagian dari pustaka C standar. Ini akan menjadi versi perangkat lunak dari fungsi floating point, atau "soft-float".

Selain: Alasan memisahkan fungsi matematika adalah bersejarah, dan hanya dimaksudkan untuk mengurangi ukuran program yang dapat dieksekusi dalam sistem Unix yang sangat lama, mungkin sebelum perpustakaan bersama tersedia, sejauh yang saya tahu.

Sekarang kompiler dapat mengoptimalkan fungsi pustaka C standar sin()(disediakan oleh libm.so) untuk diganti dengan panggilan ke instruksi asli ke fungsi sin () FSINbawaan CPU / FPU Anda, yang ada sebagai instruksi FPU ( untuk x86 / x87) pada prosesor yang lebih baru seperti seri Core 2 (ini benar sejauh i486DX). Ini akan tergantung pada flag optimasi yang dilewatkan ke kompiler gcc. Jika kompiler diminta untuk menulis kode yang akan dijalankan pada prosesor i386 atau yang lebih baru, itu tidak akan membuat optimasi seperti itu. The -mcpu=486flag akan menginformasikan compiler bahwa itu aman untuk membuat optimasi tersebut.

Sekarang jika program mengeksekusi versi perangkat lunak dari fungsi sin (), itu akan dilakukan berdasarkan pada CORDIC (Koordinat Rotasi Komputer Digital) atau algoritma BKM , atau lebih mungkin perhitungan tabel atau rangkaian daya yang umum digunakan sekarang untuk menghitung fungsi transendental seperti itu. [Src: http://en.wikipedia.org/wiki/Cordic#Application]

Versi gcc terbaru (karena sekitar 2,9x) juga menawarkan versi dosa bawaan, __builtin_sin()yang akan digunakan untuk mengganti panggilan standar ke versi pustaka C, sebagai optimisasi.

Saya yakin itu sejelas lumpur, tetapi mudah-mudahan memberi Anda lebih banyak informasi daripada yang Anda harapkan, dan banyak titik melompat untuk belajar lebih banyak sendiri.

Jika Anda ingin melihat implementasi GNU sebenarnya dari fungsi-fungsi di C, periksa trunk terbaru dari glibc. Lihat GNU C Library .

Jangan gunakan seri Taylor. Polinomial Chebyshev keduanya lebih cepat dan lebih akurat, seperti yang ditunjukkan oleh beberapa orang di atas. Berikut ini adalah implementasinya (berasal dari ZX Spectrum ROM): https://albertveli.wordpress.com/2015/01/10/zx-sine/

Komputasi sinus / cosinus / tangen sebenarnya sangat mudah dilakukan melalui kode menggunakan seri Taylor. Menulis sendiri membutuhkan waktu 5 detik.

Seluruh proses dapat disimpulkan dengan persamaan ini di sini:

Berikut adalah beberapa rutinitas yang saya tulis untuk C:

double _pow(double a, double b) {
    double c = 1;
    for (int i=0; i<b; i++)
        c *= a;
    return c;
}

double _fact(double x) {
    double ret = 1;
    for (int i=1; i<=x; i++) 
        ret *= i;
    return ret;
}

double _sin(double x) {
    double y = x;
    double s = -1;
    for (int i=3; i<=100; i+=2) {
        y+=s*(_pow(x,i)/_fact(i));
        s *= -1;
    }  
    return y;
}
double _cos(double x) {
    double y = 1;
    double s = -1;
    for (int i=2; i<=100; i+=2) {
        y+=s*(_pow(x,i)/_fact(i));
        s *= -1;
    }  
    return y;
}
double _tan(double x) {
     return (_sin(x)/_cos(x));  
}

Versi kode yang ditingkatkan dari jawaban Blindy

#define EPSILON .0000000000001
// this is smallest effective threshold, at least on my OS (WSL ubuntu 18)
// possibly because factorial part turns 0 at some point
// and it happens faster then series element turns 0;
// validation was made against sin() from <math.h>
double ft_sin(double x)
{
    int k = 2;
    double r = x;
    double acc = 1;
    double den = 1;
    double num = x;

//  precision drops rapidly when x is not close to 0
//  so move x to 0 as close as possible
    while (x > PI)
        x -= PI;
    while (x < -PI)
        x += PI;
    if (x > PI / 2)
        return (ft_sin(PI - x));
    if (x < -PI / 2)
        return (ft_sin(-PI - x));
//  not using fabs for performance reasons
    while (acc > EPSILON || acc < -EPSILON)
    {
        num *= -x * x;
        den *= k * (k + 1);
        acc = num / den;
        r += acc;
        k += 2;
    }
    return (r);
}

Inti dari bagaimana hal ini terletak pada kutipan dari Analisis Numerik Terapan oleh Gerald Wheatley:

Ketika program perangkat lunak Anda meminta komputer untuk mendapatkan nilai atau , pernahkah Anda bertanya-tanya bagaimana ia bisa mendapatkan nilai jika fungsi yang paling kuat yang dapat ia hitung adalah polinomial? Itu tidak terlihat dalam tabel dan interpolasi! Sebaliknya, komputer mendekati setiap fungsi selain polinomial dari beberapa polinomial yang dirancang untuk memberikan nilai yang sangat akurat.

Beberapa poin untuk disebutkan di atas adalah bahwa beberapa algoritma melakukan infa interpolasi dari sebuah tabel, meskipun hanya untuk beberapa iterasi pertama. Juga perhatikan bagaimana ia menyebutkan bahwa komputer menggunakan polinomial yang mendekati tanpa menentukan jenis polinom yang kira-kira sama. Seperti yang orang lain utarakan, polinomial Chebyshev lebih efisien daripada polinomial Taylor dalam kasus ini.

jika Anda ingin sinmaka

 __asm__ __volatile__("fsin" : "=t"(vsin) : "0"(xrads));

jika Anda ingin cosmaka

 __asm__ __volatile__("fcos" : "=t"(vcos) : "0"(xrads));

jika Anda ingin sqrtmaka

 __asm__ __volatile__("fsqrt" : "=t"(vsqrt) : "0"(value));

jadi mengapa menggunakan kode yang tidak akurat ketika instruksi mesin akan dilakukan?