خب سلام به همهی گلای تو خونه!
راستش همیشه میخواستم یه چیزی توی ویرگول بنویسم ولی خب خدایی هیچوقت فکرشو نمیکردم که در مورد FFT و کاربردش تو برنامهنویسی باشه اون چیز...
در هر صورت قراره در مورد Fast Fourier Transformation for polynomial multiplication صحبت کنیم.
سیگنال درسی بود که هیچوقت نفهمیدمش و صرفا همینجوری تو دوران کرونا پاسش کردم:)
بعد از چند روز گشتن توی اینترنت و با پایه ماشالله قوی که از سیگنال و سیستم داشتم یه پست خوب دیدم از همکارم در واقع ایشون Aiswarya Prakasan این خانوم هندیه.
بعد دیدم تو نتایج وب فارسی چیزی در موردش نیست، پس گفتم اونو ترجمه کنم با اندکی دخل و تصرف.
همون طور که در جریانید و میدونید تبدیل فوریه سریع یک الگوریتم پرکاربرد در علوم کامپیوتره. و معمولا ما کامپیوتریا دو درش میکنیم میره. همکارم پدرش درومده و مثه من چند روزی گشته و هر چی یاد گرفته رو اورده جمع کرده اینجا. اینجا از FFT (تبدیل فوریه سریع) برای حل مسئله ضرب سریع دو چند جمله ای استفاده می کنیم.فقط یه مشکلی که هست خیییییلی طولااانیه و دهن مام صاف میشه احتمالا تا آخرش پس یه لطفی کنید بخونیدش!
برای درک این پست معلومه که نیازی به دانستن تبدیل فوریه ندارید. با این حال، درک اساسی از موارد زیر لازمه:
خب حالا اصلا سوال چیه؟ سوال اینه که دو تا چندجملهای دادن بهمون و ما قراره ضربشون کنیم.
A(x) =a0+a1x+a2x^2+… anx^n و B(x) =b0+b1x+b2x^2+… bnx^n
و در نهایت یه چنین چیزی داریم: C(x) = A(x)*B(x)
این شکلی:
E.g., multiply (x^2 + 2x + 3)(2x^2 + 5) =
2x^2 + 0x + 5
1x^2 + 2x + 3
-------------
6 0 15
4 0 10
2 0 5
-----------------------------
2x^4 + 4x^3 + 11x^2 +10x+15
که فرم جنرالتر و کلیترش میشه این:
ci = a0*bi + a1*bi-1 + ... + ai*b0.
اگر A و B را به عنوان بردار در نظر بگیریم، آنگاه بردار C "کانولوشن" A و B نامیده می شود (به صورت {A⊗B} نشان داده می شود). فقط باید حواسمون باشه که چندجملهای های A و B با صفر پر بشن تا به درجه C که همون 2n هستش برسن. این راهحلی که گفتیم تابلوعه که دو تا حلقه for میخواد و از O(n^2) هستش.
یبذ یبر بیر تحلیل فوریه، تجزیه هر تابع ریاضی به یک سری امواج سینوسی و کسینوس (sinusoids) است. اگر کاکوها حال و حوصلشو دارید که میتونید یه سر به این لینک که جواب آقای Had Seddiqi تو quora و یا این سیبیل عزیز که میاد با رسم شکل نشون میده اون انتگرال پیوسته فرمول سری فوریه رو.
معمولا کرمی که وجود داره اینه که میایم با تبدیل فوریه یه تابع رو از یه چیزی مثلا از یه حوزهای میبریم به یه حوزه دیگه که کار کردن باهاش برامون آسونتره. مثلا از حوزه زمان میبریم به حوزه فرکانس یا برعکس یا هر چیز دیگهای.
قضیه کانولوشن: میگه که تبدیل فوریه یک کانولوشن (که با '⊗' نشان داده می شود) حاصلضرب نقطهای تبدیلهای فوریه است. فرض کنید f و g دو تابع با کانولوشن {f⊗ g} باشند. تبدیل فوریه با عملگر "ζ" نشان داده می شود، بنابراین ζ(f) و ζ(g) به ترتیب تبدیل فوریههای f و g هستند.
ζ(f⊗g) = ζ(f).ζ(g)
عمل کانولوشن به ضرب نقطهای تبدیل میشود. (بعدا میگیم که ضرب نقطهای چیه)
نمایشی از چند جملهایها در بالا دیدیم، نمایش ضریبیشون یا همون coefficientشون بود. یه روش دیگه نمایش چند جملهای، نمایش نقطه-مقدار یا point-value که ما اینجا با این سروکار داریم.
یک چند جمله ای درجه n-1 را میتوان به طور منحصر به فرد با مقادیر آن در حداقل n نقطه مختلف نشان داد. این نتیجه نظریه اساسی جبر است که بیان می کند:
"یک چند جملهای A(x) درجه n-1 به طور منحصربهفرد با داشتن n مقدار متمایز x آن مشخص میشود."
یک چند جملهای A(x) =a0+a1x+a2x^2+… anx^n را می توان به صورت نمایش داد:
A set of n pairs {(x0, y0),(x1, y1),...,(xn, yn)} such that
• for all i != j, xi != xj. ie, the points are unique.
• for every k, yk = A(xk);
در شکل نقطهای، ضرب 𝐶(𝑥) = 𝐴(𝑥)𝐵(𝑥) به شکل 𝐶(𝑥𝑘) = 𝐴(𝑥𝑘) .𝐵(𝑥𝑘) است برای هر نقطه (𝑥𝑘)
If 𝐴 and 𝐵 are of degree-bound 'n', then 𝐶 is of degree-bound '2n'.
Need to start with “extended” point-value forms for 𝐴 and 𝐵 consisting of
2𝑛 point-value pairs each.
𝐴: (𝑥0, 𝑦0) , (𝑥1, 𝑦1) , … ,( 𝑥2𝑛−1, 𝑦2𝑛−1)
𝐵: (𝑥0, 𝑦0') , (𝑥1, 𝑦1′) , … , (𝑥2𝑛−1, 𝑦2𝑛−1′)
𝐶: (𝑥0, 𝑦0𝑦0′) ,( 𝑥1, 𝑦1𝑦1′) , … ,( 𝑥2𝑛−1, 𝑦2𝑛−1𝑦2𝑛−1′)
برای مثال:
𝐴(𝑥) = 𝑥 ^3 − 2𝑥 + 1
𝐵 𝑥 = 𝑥 ^3 + 𝑥^ 2 + 1
𝑥𝑘 = (−3, −2, −1, 0, 1, 2, 3) پس ما 7 تا ضریب نیاز داریم!
𝐴: (−3,−17) ,( −2,−3) ,( −1,1) ,( 0, 1) ,( 1, 0) , (2, 5) , (3, 22)
𝐵: (−3,−20) ,(−2,−3) , (−1, 2) ,( 0, 1) ,( 1, 3) ,( 2, 13) ,( 3, 37)
𝐶: (−3,340) , (−2,9) , (−1,2) ,( 0, 1) ,( 1, 0) ,( 2, 65) ,( 3, 814)
یعنی اومدیم و به ازای xهای مشخص و یکسان، yهاشون رو در هم ضرب کردیم.
این ضرب نقطهای(point-wise) است. زمان مورد نیاز برای ضرب دو چندجملهای با درجه n در شکل و نمایش نقطه-مقدار از Θ(n) است. به تعیین ضرایب یک چندجملهای از روی نمایش نقطه-مقدار آن درونیابی(interpolation) نامیده میشود. این همون چیزیه که قضیه کانولوشن با توجه به تبدیل فوریه یک تابع بیان می کند.
پس رودمپمون تا اینجای کار این شکلیه:
پس کل اُردِر و پرفورمنس الگوریتم ما به بازده تبدیل بین دو نمایش بستگی داره:
اینجاست که از الگوریتم FFT استفاده می کنیم.(و ایناست که بددددده:))
Let m = 2n-1. [so degree of C is less than m]
1. Pick m points x_0, x_1, ..., x_{m-1} chosen cleverly.
2. Evaluate A at each of the points: A(x_0),..., A(x_{m-1}).
3. Same for B.
4. Now compute C(x_0),..., C(x_{m-1}), where C is A(x)*B(x)
5. Interpolate to get the coefficients of C.
در یک نگاه، به نظر میرسه که مراحل 2 و 3 الگوریتم به زمان O(n^2) نیاز دارند. با این حال، FFT به ما این امکان را میده که فرآیند نمایش ضریب چندجملهای به نمایش نقطه-مقدار رو سریعتر کنیم و اردشو به O(nlogn) برسونم😍✌️. فقط نکتش اینه که باید یه طور سامورایی اون m نقطه رو انتخاب کنیم. (برای m نقطه دلخواه کار نمیکنه بعدا میگیم چرا؟؟. حالا اون نقاط خاص چیان یا کیان؟ همون ریشه واحد توی اعداد مختلط که مودبشون میشه roots of unity که اگه اعداد مختلط و این سوسول بازیام یادتون رفته مثل من نگران نباشید چون خاله Aiswarya قراره کامل و مفصل برامون توضیح بده.)
وقتی میگیم "ریشه" یا "root" در اینجا منظورمون ریشههای مربع، ریشه های مکعب و هر ریشه دیگهای که ممکنه به اون نیاز داشته باشیم اشاره داره. ریشه nام یک عدد k چیه؟ عددیه که وقتی در خودش n ضرب شه، k به دست میاد. کلمه "واحد" یا "unity" هم در جریانید که معنیش میشه یک:// بنابراین یک ریشه واحد یا root of unity هر عددی است که وقتی در خودش چند بار ضرب شود، عدد 1 به دست مییاد. برای به دست آوردن هر ریشهی واحد، به غیر از 1، باید به اعداد مختلط بپردازیم.
ما میتونیم مجموعه اعداد مختلط رو به عنوان یک صفحه دو بعدی در نظر بگیریم. ما یک عدد مختلط با دو مولفه در مختصات را به همون روشی که در نمودارهای دکارتی نمایش میدادیم مشخص میکنیم: نقطه (1،1) به عدد 1+i اشاره دارد. اگر روی نقطه (0,0) باشیم، خب چه کاریه که یکی بریم راست بعد یکی بیایم بالا تا به نقطه (1،1) برسیم. به جاش، با راه رفتن 2√ واحد در زاویه 45 درجه نسبت به محور x، اوریب میریم تا برسیم بهش! به طور نمادین، وقتی از این فاصله شعاعی و جهت استفاده میکنیم، یک عدد مختلط را به صورت r * e^iθ مینویسیم، جایی که r فاصله و θ جهت یا همون زاویهست، معمولاً به جای درجه از رادیان استفاده میکنیم. در یک دایره کامل کلا 2π رادیانه، بنابراین عدد 1 معادل e^i2π نوشته میشه. زاویه 45 درجه π/4 رادیان است، بنابراین نقطه (1،1) که در بالا گفتیم به صورت 2e^iπ/4√ نوشته میشه.
به این روشِ تعیین یک نقطه با طول و جهت که بالا گفتیم مختصات قطبی یا polar coordinates میگن.
اگه نقطهای به صورت re^iθ نمایش داده بشه، مختصات x (قسمت واقعی یا همون real) با rcosθ و مختصات y (قسمت خیالی یا همون imaginary ترجمه imaginary واقعا سخته!) با rsinθ داده میشه. (که اینا بیس مثلثاته اگه اینجا مشکل دارید استاد جواد خیابونی خیلی خوب این قسمتو مثلث قائمالزاویه و وترو توضیح دادن)
polar coordinates : re^iθ
cartesian coordinates :(rcosθ, rsinθ)
ضرب با نماد قطبی re^iθ که دیگه جدیجدی آسون و شیرینه و گلابیه:
چی کار میکنیم؟ میایم و فواصل را با هم ضرب میکنیم و زوایا رو هم جمع میکنیم به همین خوشمزگی. بنابراین 5e^iπ/6 × 2e^iπ/3 اگه گفتید چند میشه؟؟؟؟؟ آفرین گلای تو خونه 10e^iπ/2 میشه. که فک کنم واضحه اومدیم و 5 * 2 کردیم که شده 10 و π/6 هم π/3 جمع کردیم که شده 3π/6 یا همون π/2.
بنابراین ضرب هم شامل انبساط یا انقباض است (این قسمت فاصلهست) و هم چرخش (این قسمت زاویهشه)
اگه یه عدد داشته باشیم مثلا re^iθ که وقتی در خودش تعداد معینی ضرب بشه، 1 به دست بیاد، اون فاصله rش دقیقا باید 1 باشه. میگی نه؟ پس چند باشه اون r؟ اگه r بزرگتر از 1 بود، مثلاً 2، اون موقع وقتی عدد رو هی در خودش ضرب کنیم، فاصلهش از (0,0) میره به سمت 2 بعد 4 بعد 8 و .... تا به سمت بینهایت میل کنه. اگه r کوچیکتر از 1 میبود، مثلاً 1/2، اون موقع میرفتیم به سمت این مقادیر: 1/2، 1/4، 1/8 و... پس 1 تنها فاصله شعاعیایه که کاملاً بالانس میمونه، پس طبق چیزی که اون توی ضرب گفتیم این یک که ثابته و هر چی ضربش کنیم هی خودشه... فقط میمونه چی؟ ایول اون زوایا رو باید هی با هم جمع کنیم یعن نقاطی داریم که روی محیط دایره هی اینور اونور میشن in other word به قول این باکلاسا.
خاله Aiswarya برای این که قضیه برای ما گلای تو خونه بهتر جا بیفته با رسم شکل واسمون توضیح داده که مام برای این که از دنیا عقب نمونیم با رسم شکل توضیح میدیم.
خب ما و خاله میدونیم که e^i2π/7 یه ریشه واحده. وقتی اونو به توان هفت برسونیم، به e^i2π میرسیم که برابر با یکه. در واقع دور دایره، وقتی هی e^i2π/7تا، e^i2π/7تا بچرخیم، یه همچین چیزی داریم:
که هفت ریشه هفتم واحد وجود داره که هر دیسک طلایی در تصویر بالا یکیشونه. برای هر عدد n میتونیم با جایگزین کردن n به جای 7 در e^i2π/7، یک nامین ریشه واحد داشته باشیم. تصاویر مربوط به سایر ریشههای واحد بسیار شبیه به شکل بالا هستن، فقط فرقشون توی تعداد دیسکهای طلاعه.
دقیقاً nتا ریشه مختلط nام واحد وجود دارد e^2𝜋𝑖𝑘/𝑛 برای k = 0، 1، …. ,𝑛 − 1
که توی اونا e^iu= cos (u) + 𝑖sin(u). در اینجا "u" یه زاویه رو به رادیان نشون میده.
حالا بریم که هشتتاییشو داشته باشیم:
2. لم Cancellation لغوکردن یا حذف کردن اون چیزا توانا یا یه همچین چیزی
3. شیرین بعدی و آخری لم نصف کردن یا Halving Lemma
هوووف تموم شد شیرینا:////
حالا FFT از کدوماش استفاده میکنه بلا، از اینا:
As we will see, the fast Fourier transform algorithm cleverly makes
use of the following properties about ωn:
𝜔𝑛^n =1
𝜔𝑛^n+𝑘 = 𝜔𝑛^𝑘
𝜔𝑛^n/2 = -1
𝜔𝑛^(n/2+𝑘) = -𝜔𝑛^𝑘
ما مسئله رو تجزیه کردیم به مسائل کوچیکتر، الان نوبت اینه که یک چندجملهای مثل A (x)(degree m) از درجه m که توی m نقطه انتخابی ما هستش رو در زمان O(mlogm) ارزیابی کنیم. فرض کنید m از توان 2عه.
اجازه بدید با یه مثال بریم جلو.
ما اگه m=8، اون موقع یه چندجملهای این شکلی داریم:
A(x) = a_0 + a_1 x + a_2 x² + a_3 x³ + a_4 x⁴ + a_5 x⁵ + a_6x⁶ + a_7x⁷
(نمایش وکتور و آرایهایشم این شکلیه: A = [a_0, a_1, …, a_7])
و ما میخوایم در هشت نقطه انتخاب خودمو بسنجیم که چه گلی به سرمون بگیریم که بتونیم از تکنیک divide and conquer یا همون تقسیم و حل استفاده کنیم. یه ایده جالب اینه که A را به دو قسمت تقسیم کنیم، ولی به چه دو قسمتی؟ به جای این که به دو قسمت چپ و راست تقسیمشون کنیم، اونا رو زوج و فردشون میکنیم(با تشکر از دوستامون در کنترل ترافیک) بنابراین، به عنوان بردار، این دو تا رو داریم:
A_even = [a_0, a_2, a_4, a_6]
A_odd = [a_1, a_3, a_5, a_7]
یا، چندجملهایش:
A_even(x) = a_0 + a_2 x + a_4 x² + a_6 x³
A_odd(x) = a_1 + a_3 x + a_5 x² + a_7 x³
خب حالا درجه هر کدومشون شد n/2. حالا چهطوری میتونیم A(x) رو بر حسب A_even و A_odd بنویسیم؟
A(x) = A_even(x²) + x A_odd(x²).
A(-x) = A_even(x²) — x A_odd(x²).
اگه دقت کنید اون A_odd مثل همون A_even فقط یه x اضافه داره که از اون x فاکتور میگیریم و به شکل بالا مینویسیمشون. اینم یه هنر دیگه:) اون A(-x) هم داستانشو جلوتر کامل متوجه میشیم و به عینه میبینیم.
خب حالا با این هنرمون (الکی در اصل هنر دوستان بوده که من خودمو به زور چسبوندم بهشون)، به جای محاسبه مقدار یک چندجملهای درجه m، باید 2 چند جمله ای درجه (m/2) رو محاسبه کنیم و بعدش اونارو حال، به جای محاسبه مقدار 1 چند جمله ای درجه m، اکنون باید 2 چند جمله ای درجه (m/2) را محاسبه کرده و بعدش اونارو با هم ترکیب کنیم.
اجازه بدین T(m) بگیریم زمان مورد نیاز برای ارزیابی یک چندجملهای درجه-m در مجموعه m نقاط ما. ما این کارو با تقسیم به دو چندجملهای درجه-m/2 با m/2 نقطه انجام میدیم، و بعدش یه O(m) هم برای ترکیب نتایج نیاز داریم. خب خیلی خوب شد دیگه، نه؟ شد T(m) = 2T(m/2) + O(m) که میشه همون O(mlogm) خودمون.
تازه یه چیز خیلی باحال و جالب ماجرا اینجاست که وقتی A(x) رو حساب کردیم، یه جورایی A(-x) رو هم داریم و دیگه نیازی به صرف هزینه جداگونه برای محاسبه اون نداریم!(اینو جلوتر مثه ماه میبینیم)
بنابراین، ما باید به طور خاص m نقاطی را انتخاب کنیم که دارای ویژگی زیر باشند:
نیمه دوم نقاط باید منفی(قرینه) نیمه اول باشه. (اسمشو میذاریم ویژگی یونسیوسی)
این نکته بالا خیلی مهمه نمیدونم اینجا چه طوری میشه هایلاتش کرد، ولی خب حواستون حسابی جمع باشه چون قراره جلوتر ضایع بشیم سر همین قضیه!
مثلا اینو ببینید چه غوچگله: {1,2,3,4,-1,-2,-3,-4}
ولی خب بدبختی اینه که همیشه به این غوچگلی(خوشگلی) و سادگی و خوشمزگی نیست که بگیم خب اوکی نصفشون میکنیم، بعدشم بازگشتی یا ریکرسیو به قول این باکلاسا حلش میکنیم و عشق میکنیم و خَلاص!
مشکل اینجاست که اون بخش بازگشتیمون برای یه سری نقاط خاص مثلا: {1, 4, 9, 16} اون ویژگی یونسیوسی رو برآورده یا ارضا کنه! مثلا طبق اون ویژگی باید دقیقا این شکلی باشه: {1, 4, -1, -4}
هاهاها. اما اینا مربعن...تازه یه بدبختی دیگهام داریم اصن گیریم که این مرحله رو با خرشانسی مطلق رد کردیم، ولی مرحله بعد که قراره به توان دو برسن با اون مربعهاشون چی کار کنیم؟ +هیچی. دیگه لزوما اون شرط قرینه بودنشون حفظ نمیشه. خب حالا چاره چیه؟ دقیقا اعداد مختلط.
مثلا اگه اینا مربع باشن:
{1, 2, i, 2i, -1, -2, -i, -2i}
نقاط اصلیمون کدومان؟
بعد مربع اینا میشه چی؟ میشه این: 1, 4, -1, -4 خب حالا دوباره این آخریارو اگه یه توان دو دیگه بدیم میشن چی؟ میشن این دو تا: 1, 16.
اما برای لول بعد ما دوباره به ویژگی یونسیوسی نیاز داریم:) و ما میخوایم {1, -1} رو اونجا داشته باشیم! این معنیش اینه که یه حلقه گمشده داریم... حلقهمون چیه؟ حلقهمون اینه که لول قبلی باید یه چنین چیزی: {1, i, -1, -i} باشه:) که معادلشو به صورت توانی بخوایم بنویسیم میشه این: {1, i, i^2, i^3} و اون موقع دیگه همه دردامون دوا میشه😊🥳
بنابراین، برای لول اصلی، اجازه بدید ω باشه مساوی ریشه هشتم اصلی واحد، و بعد مجموعه اصلی نقاط ما به صورت زیر خواهد بود:
1,ω ,ω^2,ω^3,ω^4 (= -1),ω^5 (= -ω ),ω^6 (= -ω^2),ω^7 (= -ω^3)
که همون طور که میبینید مثه ماه شده. خب حالا بریم جلو مربع اینا چی میشه؟ میشه: 1, i, i² (= -1), i³ (= -i)
و در لول و مربع بعدی داریم: 1, -1
و در نهایت: 1
فرم کلی ریشه mام اولیه واحد همکارم زحمت کشیده براتون نوشته که یه همچین برداریه:
ω = cos(2*pi/m) + i*sin(2*pi/m)
پس الگوریتم کلی FFT رو بخوایم مرور کنیم به زیبایی زیر میشه:
Fast Fourier Transform (FFT)
The problem of evaluating 𝐴(𝑥) at 𝜔𝑛^0 , 𝜔𝑛^1 , … , 𝜔𝑛^𝑛−1 reduces to
1. evaluating the degree-bound 𝑛/2 polynomials 𝐴even(𝑥) and 𝐴odd(𝑥) at the points (𝜔𝑛^0)^2 ,(𝜔𝑛^1)^2 , … , (𝜔𝑛^𝑛−1)^2.
2. combining the results by 𝐴(𝑥) = 𝐴even(𝑥2) + 𝑥𝐴odd(𝑥2).
خب کاکو حالا کی حال و حوصله داره این مقادیرو دو بار حساب کنه؟ اصن چرا خودمونو زحمت بدیم؟؟ ها اصن چه کاریه؟
The list (𝜔𝑛^0)^2 ,(𝜔𝑛^1)^2 , … , (𝜔𝑛^𝑛−1)^2 does not contain 𝑛 distinct values
میگه که درسته که این عظمت مقادیر منحصر به فرد نداره ولی اون n/2تا عدد مختلط که n/2امین ریشه واحد دارن که.... در واقع چندجملهایهای 𝐴even و 𝐴odd به طور بازگشتی توسط n/2تا عدد مختلط که n/2امین ریشه واحدن ارزیابی میشن و این زیرمسئلهها دقیقا مشابه مسئله اصلی هستند فقط با اندازه نصف یا همون n/2
فک کنم به زبان فارسی خیلی سخت گفتمش:/ پس بریم که با دو تا مثال مشتی برسیشون کنیم این چیزایی که گفتیمو و بریم که داشته باشیم:
A(x) = 3+2x+3x^2+4x^3
A(ω4^0) = A(1) = 3+2+3+4=12
A(ω4^1) = A(i) = 3+2i-3-4i=-2i
A(ω4^2) = A(-1) = 3–2+3–4=0
A(ω4^3) = A(-i) = 3–2i-3+4i=2i
خب این اون روش اصلی که خیالتون راحت بشه. حالا بریم طرح ترافیک فرد و زوجمون رو اجرا کنیم و ببینید چه قد مثل ماهه قبونش برم این الگوریتم! پس A(x) میگیریم Aeven(x^2) + xAodd(x^2) پس داریم:
A(x)= Aeven(x^2) + xAodd(x^2)
Aeven(x) = 3+3x
Aodd(x) = 2+4x
A(ω4^0)= A(1) = Aeven(1)+1.Aodd(1)= 3+3+2+4=12
A(ω4^1)= A(i) = Aeven(-1)+i.Aodd(-1)=3–3+2i-4i=-2i
A(ω4^2)= A(-1) = Aeven(1)-1.Aodd(1)= 3+3–2–4=0
A(ω4^3)= A(-i)= Aeven(-1)-i.Aodd(-1)=3–3–2i+4i= 2i
خب حالا تمام اونایی که گفتم بعدا میبینیم بعدا میگم همش اینجا و توی این مثال بود...البته همهایم نبود به اون صورت دو تا مورد بود. یکی این که به جای این که مستقیم خودِ مسئله رو حساب کنیم اون زوج و فردشو حل کنیم دقیقا به همون جوابا میرسیم و یکی دیگهام اون A(-x) بود که یه ذره بالاتر گفتم میگم چرا منهای جملات فرد(Aodd(x)) میکنیم؛ که اینجا فک کنم خیلی واضحه دیگه نه؟؟؟ به جملات A(ω4^0) و A(ω4^1) نگاه کنید! خب دیگه بسه! حالا به A(ω4^2) و A(ω4^3) نگاه کنید. خب اینم دیگه بسه!!! جمله A(ω4^2) دقیقا همون A(ω4^0) فقط این سری به جای این که Aeven(1)+1.Aodd(1) بشه، Aeven(1)-1.Aodd(1) شده!!! و در جریان رابطه ω4^0 و ω4^2 هم هستید دیگه اگه کارتزین بخوایم بگیم اون ω4^0 معادل (0,1) و ω4^2 معادل (0, -1) پس میشه همون A(x) و A(-x) خودمون. همین رابطه برای A(ω4^1) و A(ω4^3) هم برقراره.
خب بعد این همه قصه بریم که الگوریتم و سودکدش رو با هم ببینیم:
Algorithm
Here is the general algorithm in pseudo-C:
Let A be array of length m, w be primitive mth root of unity.
Goal: produce DFT F(A): evaluation of A at 1, w, w^2,...,w^{m-1}.
FFT(A, m, w)
{
if (m==1) return vector (a_0)
else {
A_even = (a_0, a_2, ..., a_{m-2})
A_odd = (a_1, a_3, ..., a_{m-1})
F_even = FFT(A_even, m/2, w^2)//w^2 is a m/2-th root of unity
F_odd = FFT(A_odd, m/2, w^2)
F = new vector of length m
x = 1
for (j=0; j < m/2; ++j) {
F[j] = F_even[j] + x*F_odd[j]
F[j+m/2] = F_even[j] - x*F_odd[j]
x = x * w
}
return F
}
در واقع این همه حرف زدیم شد همین یه تیکه سودکد:( ولی خب به جاش ایشالا کامل متوجه شدیم!)
نکته: برای اعمال این الگوریتم، «m» باید توان 2 باشه، همون چیزی که قبلا بهش اشاره کردیم که باید یه شرایطی خاصی داشته همین بود یکی از منظورمان که البته میتونیم صفرهای اضافی به چند جملهایمون اضافه کنیم برای اطمینان از اینکه طول آرایه توانی از دو هست. ( یه شرط دیگهمون که همون قرینگی بود که بالا کامل در موردش توضیح دادیم). این الگوریتم آرایهای از اعداد مختلط رو برامون برمیگردونه که ازشون برای ضرب نقطهای استفاده میکنیم.
این زیر درخت بردارهای ورودی برای فراخوانیهای بازگشتی FFT رو داریم. بردار اولیهمون 8 تا عضو داره:
خب وقتی ضرب نقطهای رو روی تبدیل فوریه A و B انجام دادیم و C (آرایهای از اعداد مختلط) رو به دست اُوردیم، باید اونو دوباره به شکل و نمایش اولیهاش یعنی همون ضرایب coefficient برگردونیم تا باهاش چندجملهایمون رو بر اساس ضرایبش نمایش بدیم.
ممکنه براتون واضح باشه یا شایدم نباشه:/// که ما می تونیم با همون الگوریتم FFT البته معکوسش به نمایش یا فرم ضریبی از نمایش نقطه-مقدار برگردیم. اگه نقاط فعلیو معکوس کنیم(توجه کنید که ω^-j هم ریشه واحد برای برخی از اعداد صحیح j هست) و با اعمال FFT روی مقادیر، ضرایب را برمی گردونیم (یه ضرب در 1/n هم داره.) یه راه خوب برای فکر کردن در این مورد، فرمولبندی مجدد FFT به عنوان یک مسئله جبر خطیه یعنی روز از نو روزی از نو. خب حالا چیا داریم چیا نداریم... ما داریم:
A = Mn * a
که Aمون برداری از مقادیرِ A(x0),A(x1),…,A(xk)، اون یکی a برداری از ضرایبه و Mn یه ماتریس تبدیل خطیه:
ماتریس Mn(x) رو بهش میگن Vandermonde matrix که یه سری خواص داره که به کار ما میاد. مثلا اگه x0,x1,…,xn−1 متمایز باشن، ماتریسمون معکوسپذیره ولی خب هاهاها خب که چی؟ ایدهی احمقانهایه چون حساب کردن معکوس ماتریس خودش از O(n3). اما شانس ما اینه که ماتریسهای Vandermonde رو میتونیم توی O(n2) معکوس کنیم. ولی خب اینم بدردمون نمیخوره... چون که شاخ شدیم و فقط O(nlogn) کار میکنیم(به قول یکی از بچهها فقط ونک به بالا) یه حرکت باحال اینه که اسکی بریم از کاری که قبلن انجام دادهبودیم که احتمالا همهی شما هم به ذهنتون رسیده بود و من الکی وقتتون رو گرفتم که ببیند FFT چه کمکی به ما میکنه. خب حالا که متوجه شدین کامل ما باید این ماتریس تبدیل خطی یا همون Mn(x)مون رو با nامین ریشه واحد بازنویسی کنیم که به یه همچین چیزی میرسیم:
نگران قسمت خیالی یا همون imaginary اعداد مختلط نباشید. وقتی که تبدیل فوریه معکوس انجام میشه، آرایهای از اعداد مختلط به دست میان که اون قسمت imaginaryشون تقریبا برابر با صفره. بنابراین، میتونیم ایگنورش کنیم.
توی ماتریس آخری میتونید ثابت کنید که ستونهای Mn(ω) متعامدن(یعنی MM∗=nI.) خب حالا دیگه کار راحت شد دیگه...الان فقط یه تبدیل مختصات ساده داریم(تغییر خطی متخصات.) پس تبدیل معکوس FFT تبدیل میشه به:
Mn(ω)−1=M(ω−1)/n
در یک ضریب 1/n ضرب میشه اون M(ω−1). حالا اون M(ω−1) رو چه طوری حساب کنیم؟؟؟
حالا برای محاسبه M(ω−1) کافیه خیلی راحت به جای ω، در ماتریس فوریه معکوس از ω^-1 استفاده کنیم طبق قوانین زیر و در نهایت اونو در 1/n ضرب کنیم.
if ω = a+ib => ω^-1 = a-ib
if ω = e^2Πi/k => ω^-1 = e^-2Πi/k
پس Mn(ω)−1 در نهایت میشه:
پس چی شد؟؟؟؟ اون چیزی که در واقع ما انجام دادیم، محاسبه A(x) در 1ω,ω2, …,ω{m-1} بود دیگه که میتونیم به شکل زیر نمایشش بدیم:
فرم ماتریسیشم که بالاتر گفتیم چه شکلی به دست میاد همون واندرمده:
پس الگوریتم نهایی این مصیبت به شکل زیر خواهد بود:
Let F_A = FFT(A, m, ω) // time O(n log n)
Let F_B = FFT(B, m, ω) // time O(n log n)For i=1 to m, let F_C[i] = F_A[i]*F_B[i] // time O(n)
Output C = 1/m * FFT(F_C, m, ω-1). // time O(n log n)
هوووووووووف بالاخره تموم شد:) کدشم که همینجوری ریخته تو اینترنت ولی خب اگه خواستید بگین تا لینک قرار بدم و در پایان با تشکر ویژه از خاله Aiswarya از هندوستان و عمو Had که این زیر لینک ایشونام قرار میدم:
لینک خاله Aiswarya که بیس اصلی بود
و لینک عمو Had که اینا آخرارو واسه معکوس فوریه کمک کرد.