diff --git a/chastotnie-methods/lab1/fourier.py b/chastotnie-methods/lab1/fourier.py index e70aee1..cb75aab 100644 --- a/chastotnie-methods/lab1/fourier.py +++ b/chastotnie-methods/lab1/fourier.py @@ -1,5 +1,6 @@ from sys import argv from math import pi +from math import e function_formulas = [ # 1. Square wave @@ -29,6 +30,18 @@ function_formulas = [ if n > 0 else (4 / (pi ** 2 * n ** 2) if n % 2 == 0 else -24j / (pi ** 3 * n ** 3)) ) + }, + # 5. Complex func + { + 'an': lambda n: 0, + 'bn': lambda n: 0, + 'cn': lambda n: 0 if n == 0 else ( + complex(-pi * n / 4, -1 - pi * n / 4) * e ** (1j * pi * n / 4) + + (1j + 1) * e ** (-1j * pi * n / 4) + + (1j - 1) * e ** (-3j * pi * n / 4) + + (-1j - 1) * e ** (-5j * pi * n / 4) + + complex(1 + pi * n / 4, pi * n / 4) * e ** (-7j * pi * n / 4) + ) / (pi ** 2 * n ** 2 / 2) } ] diff --git a/chastotnie-methods/lab1/Егор_Капралов_5_1.typ b/chastotnie-methods/lab1/Егор_Капралов_5_1.typ index 3eaf577..3c283c8 100644 --- a/chastotnie-methods/lab1/Егор_Капралов_5_1.typ +++ b/chastotnie-methods/lab1/Егор_Капралов_5_1.typ @@ -87,15 +87,15 @@ $ b_n &= 1/(T/2) integral^3_1 f(t)sin(omega_n t)upright(d)t = integral^3_1 f(t)s 2/(pi n)\, & n = 2k - 1 ), k in ZZ $ -$ c_n &= 1/T integral^3_1 f(t)e^(i omega_n t) = 1/2 integral^3_1 f(t)e^(i pi n t) = 1/2 (integral^2_1 e^(i pi n t) + integral^3_2 2e^(i pi n t))upright(d)t \ -&= 1/2 (lr(-(i e^(i pi n t))/(pi n)|)^2_1 - lr((2i e^(i pi n t))/(pi n)|)^3_2) = -(i e^(2i pi n))/(2pi n) + (i e^(i pi n))/(2pi n) - (2i e^(3i pi n))/(2pi n) + (2i e^(2i pi n))/(2pi n) \ -&= (i e^(2i pi n))/(2pi n) - (i e^(i pi n))/(2pi n) = i/(2pi n) - (i(-1)^n)/(2pi n) = (i(1 - (-1)^n))/(2pi n) \ +$ c_n &= 1/T integral^3_1 f(t)e^(-i omega_n t) = 1/2 integral^3_1 f(t)e^(-i pi n t) = 1/2 (integral^2_1 e^(-i pi n t) + integral^3_2 2e^-(i pi n t))upright(d)t \ +&= 1/2 (lr(-(i e^(-i pi n t))/(pi n)|)^2_1 - lr((2i e^(-i pi n t))/(pi n)|)^3_2) = -(i e^(-2i pi n))/(2pi n) + (i e^(-i pi n))/(2pi n) - (2i e^(-3i pi n))/(2pi n) + (2i e^(-2i pi n))/(2pi n) \ +&= (i e^(-2i pi n))/(2pi n) - (i e^(-i pi n))/(2pi n) = i/(2pi n) - (i(-1)^n)/(2pi n) = (i(1 - (-1)^n))/(2pi n) \ &= cases( 0\, & n = 2k, i/(pi n)\, & n = 2k - 1 ), k in ZZ $ -$ a_0 &= c_0 = 3, a_n = 0, \ +$ a_0 &= 3, a_n = 0, c_0 = a_0\/2 = 1.5 \ b_n &= cases( 0\, & n = 2k, 2/(pi n)\, & n = 2k - 1 @@ -554,7 +554,7 @@ $ f(t) = cases( Период функции $T = 4$, $omega_n = (pi n)/2$. Рассмотрим функцию на промежутке $[-2; 2]$ и найдём коэффициенты $a_0, a_n, b_n, c_n$: -$ a_n &= 1/(T/2) integral^2_(-2) f(t) upright(d)t \ +$ a_0 &= 1/(T/2) integral^2_(-2) f(t) upright(d)t \ &= ... \ &= -2/3 $ @@ -597,13 +597,13 @@ function_formulas = [ 'cn': lambda n: -1/3 if n == 0 else ( (4 / (pi ** 2 * n ** 2) if n % 2 == 0 else 24j / (pi ** 3 * n ** 3)) if n > 0 else - (4 / pi ** 2 * n ** 2) if n %2 == 0 else -24j / (pi ** 3 * n ** 3)) + (4 / pi ** 2 * n ** 2) if n % 2 == 0 else -24j / (pi ** 3 * n ** 3)) ) } ... ] ... -```, caption: "Дополнение программы нечётной функцией") +```, caption: "Дополнение программы четвёртой функцией") #figure(``` $ python fourier.py 4 3 @@ -754,7 +754,7 @@ $ upright("Re")f(t) = cases( } else if tt >= t1 and tt < t2 { return 1 } else if tt >= t2 and tt < t3 { - return 4-t + return 4-tt } else { return -1 } @@ -767,8 +767,160 @@ $ upright("Re")f(t) = cases( Период функции $T = 8$, $omega_n = (pi n)/4$. Рассмотрим функцию на промежутке $[-1; 7]$ и найдём коэффициент $c_n$: -$ c_n &= 1/T integral^7_(-1) f(t)e^(i omega_n t) upright(d)t = \ -&= ... \ -&= (((omega_n-1)sin(omega_n) + omega_n cos(omega_n))(1+i e^(2i omega_n) -i e^(4i omega_n) -i e^(6i omega_n)))/(4omega_n^2) \ -&= ((((pi n)/4-1)sin((pi n)/4) + (pi n)/4 cos((pi n)/4))(1+i e^(i (pi n)/2) -i e^(i pi n) -i e^(i (3pi n)/2)))/((pi^2 n^2)/4) $ - +$ c_n &= 1/T integral^7_(-1) f(t)e^(-i omega_n t) upright(d)t \ +&= 1/8(integral^1_(-1) (1+i t)e^(-i omega_n t) upright(d)t + integral^3_1 (2-t+i)upright(d)t + integral^5_3 (-1+i(4-t))upright(d)t + integral^7_5 (-6+t-i)upright(d)t) \ +&= 1/8(integral^1_(-1) e^(-i omega_n t)upright(d)t + 2integral^3_1 e^(-i omega_n t)upright(d)t + i integral^3_1 e^(-i omega_n t)upright(d)t - integral^5_3 e^(-i omega_n t)upright(d)t + 4i integral^5_3 e^(-i omega_n t)upright(d)t \ +&- 6integral^7_5 e^(-i omega_n t)upright(d)t - i integral^7_5 e^(-i omega_n t)upright(d)t + i integral^1_(-1) t e^(-i omega_n t)upright(d)t - integral^3_1 t e^(-i omega_n t)upright(d)t - i integral^5_3 t e^(-i omega_n t)upright(d)t + integral^7_5 t e^(-i omega_n t)upright(d)t) \ +&= lr((i e^(-i omega_n t))/(8omega_n)|)^1_(-1) + lr((2i e^(-i omega_n t))/(8omega_n)|)^3_1 - lr((e^(-i omega_n t))/(8omega_n)|)^3_1 - lr((i e^(-i omega_n t))/(8omega_n)|)^5_3 - lr((4 e^(-i omega_n t))/(8omega_n)|)^5_3 - lr((6i e^(-i omega_n t))/(8omega_n)|)^7_5 + lr((e^(-i omega_n t))/(8omega_n)|)^7_5 \ +&+ lr((e^(-i omega_n t)(i-omega_n t))/(8omega_n^2)|)^1_(-1) - lr((e^(-i omega_n t)(1+i omega_n t))/(8omega_n^2)|)^3_1 - lr((e^(-i omega_n t)(i-omega_n t))/(8omega_n^2)|)^5_3 + lr((e^(-i omega_n t)(1+i omega_n t))/(8omega_n^2)|)^7_5 \ +&= ((-i-omega_n-i omega_n)e^(i omega_n) + (i+1)e^(-i omega_n) + (i-1)e^(-3i omega_n) + (-i-1)e^(-5i omega_n) + (omega_n + 1 + i omega_n)e^(-7i omega_n))/(8omega_n^2) \ +&= ((-i-(pi n)/4-i (pi n)/4)e^(i (pi n)/4) + (i+1)e^(-i (pi n)/4) + (i-1)e^(-3i (pi n)/4) + (-i-1)e^(-5i (pi n)/4) + ((pi n)/4 + 1 + i (pi n)/4)e^(-7i (pi n)/4))/((pi^2 n^2)/2) $ + +$ c_0 &= 1/T integral^7_(-1)f(t)upright(d)t \ +&= 1/8(integral^1_(-1) (1+i t) upright(d)t + integral^3_1 (2-t+i) upright(d)t + integral^5_3 (-1+i(4-t)) upright(d)t + integral^7_5 (-6+t-i) upright(d)t) \ +&= 1/8(lr((t+(i t^2)/2)|)^1_(-1) + lr((2t+i t-t^2/2)|)^3_1 + lr((-t+4i t-(i t^2)/2)|)^5_3 + lr((-6t+t^2/2-i t)|)^7_5) \ +&= 1/8(1+i/2+1-i/2+6+3i-9/2-2-i+1/2-5+20i-(25i)/2+3-12i+(9i)/2-42+49/2 \ +&-7i+30-25/2+5i) = 0 $ +#pagebreak() + + +Дополним программу @prog новой комплексной функцией: + +#figure(```python +... +function_formulas = [ + ... + # 5. Complex func + { + 'an': lambda n: 0, + 'bn': lambda n: 0, + 'cn': lambda n: 0 if n == 0 else ( + complex(-pi * n / 4, -1 - pi * n / 4) * e ** (1j * pi * n / 4) + + (1j + 1) * e ** (-1j * pi * n / 4) + + (1j - 1) * e ** (-3j * pi * n / 4) + + (-1j - 1) * e ** (-5j * pi * n / 4) + + complex(1 + pi * n / 4, pi * n / 4) * e ** (-7j * pi * n / 4) + ) / (pi ** 2 * n ** 2 / 2) + } + ... +] +... +```, caption: "Дополнение программы комплексной функцией") + +#figure(``` +$ python fourier.py 5 3 +an values: ['0: 0', '1: 0', '2: 0', '3: 0'] +bn values: ['0: 0', '1: 0', '2: 0', '3: 0'] +cn values: ['-3: -0.12737+3.9996e-17j', '-2: 1.1249e-17+1.1249e-17j', '-1: 1.1249e-17+2.2498e-17j', '0: 0', '1: 1.1463+1.3499e-16j', '2: -6.7493e-17-6.7493e-17j', '3: 3.9996e-17+1.9998e-17j'] +```, caption: "Ввывод программы для копмлексной функции") + +Из-за погрешностей чисел с плавающей точкой (float), программа выводит очень маленькое число ($~10^(-17)$) в некоторых случаях. В таких ситуациях мы будем воспринимать коэффициент $c_n$ равным нулю. + +$ G_3(t) = -0.12737e^((-3i pi t)/4) + 1.1463e^((i pi t)/4) $ + +С помощью программы из @prog, подсчитаем также коэффициенты для $N = 10$, после чего построим график $G_N (t)$ для значений $N = 1, 2, 3, 10$ и сравним их с $f(t)$: +#pagebreak() + + +#figure( + canvas(length: 1.75cm, { + plot.plot(size: (8, 4), + x-tick-step: 1, y-tick-step: 1, + x-min: -10, x-max: 10, + y-min: -3, y-max: 3, + x-grid: true, y-grid: true, + x-label: [$t$], y-label: [Re$G_N (t)$], + { + plot.add(domain: (-9.99, 9.99), t => { + let (t0, t1, t2, t3, t4) = (-1, 1, 3, 5, 7) + if t < t0 { t = -t + t0 + 1 } + let tt = calc.rem(t - t0, t4 - t0) + t0 + if tt >= t0 and tt < t1 { + return 1 + } else if tt >= t1 and tt < t2 { + return 2 - tt + } else if tt >= t2 and tt < t3 { + return -1 + } else { + return -6 + tt + } + }, line: "linear", samples: 1000, label: [Re$f(t)$#v(1em)]) + plot.add(domain: (-9.99, 9.99), t => { + 1.1463*calc.cos(calc.pi*t/4) + }, line: "spline", samples: 100, label: [Re$G_1(t)$#v(1em)]) + plot.add(domain: (-9.99, 9.99), t => { + -0.12737*calc.cos(3*calc.pi*t/4) + 1.1463*calc.cos(calc.pi*t/4) + }, line: "spline", samples: 100, label: [Re$G_2(t)$#v(1em)]) + plot.add(domain: (-9.99, 9.99), t => { + -0.12737*calc.cos(3*calc.pi*t/4) + 1.1463*calc.cos(calc.pi*t/4) + }, line: "spline", samples: 100, label: [Re$G_3(t)$#v(1em)], style: (stroke: purple)) + plot.add(domain: (-9.99, 9.99), t => { + 0.023394*calc.cos(7*calc.pi*t/4) - 0.12737*calc.cos(3*calc.pi*t/4) + 1.1463*calc.cos(calc.pi*t/4) - 0.045853*calc.cos(5*calc.pi*t/4) + 0.014152*calc.cos(9*calc.pi*t/4) + }, line: "spline", samples: 100, label: [Re$G_10(t)$#v(1em)]) + } + ) + }), caption: [Графики частичных сумм Фурье комплексной функции $G_N (t)$ и $f(t)$ (Вещественная плоскость)] +) + +#figure( + canvas(length: 1.75cm, { + plot.plot(size: (8, 4), + x-tick-step: 1, y-tick-step: 1, + x-min: -10, x-max: 10, + y-min: -3, y-max: 3, + x-grid: true, y-grid: true, + x-label: [$t$], y-label: [Im$G_N (t)$], + { + plot.add(domain: (-9.99, 9.99), t => { + let (t0, t1, t2, t3, t4) = (-1, 1, 3, 5, 7) + if t < t0 { t = -t + t0 - 3 } + let tt = calc.rem(t - t0, t4 - t0) + t0 + if tt >= t0 and tt < t1 { + return tt + } else if tt >= t1 and tt < t2 { + return 1 + } else if tt >= t2 and tt < t3 { + return 4-tt + } else { + return -1 + } + }, line: "linear", samples: 1000, label: [Im$f(t)$#v(1em)], style: (stroke: yellow)) + plot.add(domain: (-9.99, 9.99), t => { + 1.1463*calc.sin(calc.pi*t/4) + }, line: "spline", samples: 100, label: [Im$G_1(t)$#v(1em)]) + plot.add(domain: (-9.99, 9.99), t => { + 0.12737*calc.sin(3*calc.pi*t/4) + 1.1463*calc.sin(calc.pi*t/4) + }, line: "spline", samples: 100, label: [Im$G_2(t)$#v(1em)]) + plot.add(domain: (-9.99, 9.99), t => { + 0.12737*calc.sin(3*calc.pi*t/4) + 1.1463*calc.sin(calc.pi*t/4) +}, line: "spline", samples: 100, label: [Im$G_3(t)$#v(1em)], style: (stroke: purple)) + plot.add(domain: (-9.99, 9.99), t => { + -0.023394*calc.sin(7*calc.pi*t/4) + 0.12737*calc.sin(3*calc.pi*t/4) + 1.1463*calc.sin(calc.pi*t/4) - 0.045853*calc.sin(5*calc.pi*t/4) + 0.014152*calc.sin(9*calc.pi*t/4) + }, line: "spline", samples: 100, label: [Im$G_10(t)$#v(1em)]) + } + ) + }), caption: [Графики частичных сумм Фурье комплексной функции $G_N (t)$ и $f(t)$ (Комплексаня плоскость)] +) + +Из графиков видно, что сумма Фурье $G_N (t)$ сходится к исходной функции $f(t)$ при $N -> infinity$ - коэффициент $c_n$ подсчитан верно. +#pagebreak() + +Проверим выполнение равенства Парсеваля: + +$ 1/(T/2) integral^7_(-1) (f(t))^2 upright(d)t = sum^infinity_(n=-infinity) |c_n|^2 $ + +$ 1/(T/2) integral^7_(-1) (f(t))^2 upright(d)t &= 1/4(integral^1_(-1) (1+i t)^2 upright(d)t + integral^3_1 (2-t+i)^2 upright(d)t integral^5_3 (-1+i(4-t))^2 upright(d)t + integral^7_5 (-6+t-i)^2 upright(d)t) \ +&= 1/4(lr((t+i t^2-t^3/3)|)^1_(-1) + lr((3t+4i t-2t^2-i t^2+t^3/3)|)^3_1 \ +&+ lr((-15t-8i t+i t^2+4t^2-t^3/3)|)^5_3 + lr((t^3/3-6t^2-i t^2+35t+12i t)|)^7_5)\ +&= 1/4(1+i-1/3+1-i-1/3+9+12i-18-9i+9-3-4i+2+i-1/3\ +&-75-40i+25i+100-125/3+45+24i-9i-36+9+343/3-294-49i\ +&+245+84i-125/3+150+25i-175-60i) = 0 $ + +$ sum^infinity_(n=-infinity) |c_n|^2 &= sum^infinity_(n=-infinity) lr(|cases( + (8sqrt(2))/(pi^2 n^2)\, & n = 8k-7, + -(8sqrt(2))/(pi^2 n^2)\, & n = 8k-3 +)|)^2, k in ZZ = sum^infinity_(n=-infinity) ((8sqrt(2))/(pi^2 n^2))^2, n = 4k-3, k in ZZ \ +&= 128/pi^4 sum^infinity_(k=-infinity) 1/(4k-3)^2 = 128/pi^4 dot pi^4/96 = 4/3 $ + +Левая и правая часть не равны, поскольку функция $f(t)$ не соответствует условию Дирихле. +