Processing math: 0%
Skip to main content

Гладкая сшивка кусочно-постоянной функции.

Накануне мне понадобилось найти гладкую функцию, определенную на \textbf{R}[a,c], равную единице на одном отрезке числовой оси, скажем, [a, b-\epsilon] и нулю на другом, скажем [b+\epsilon, c]. Для того, чтобы гладко склеить эти два кусочка, я решил найти монотонно убывающий многочлен f(x), равный f(0) = 1 и f(1) = 0, а потом "вставить" его в отрезок [b-\epsilon, b+\epsilon]. Для гладкости сшивки должны выполняться условия f'(0) = f'(1) = 0. Я решил пойти далее и положить f^{(n)}(0) = f^{(n)}(1) = 0 для всех n = 1,2,\cdots,N. То есть первые N производных должны быть равны нулю в т. x = 0 и x = 1.

Понятно, что это легко сделать при наперед известном N, но как найти общую формулу? Я некоторое время промучился над этой задачей, пока не решил её следующим способом.

Пусть искомая функция f(x) = 1 + x^{N+1}L_{1}(x), где L_{1}(x) — некий многочлен. Тогда у нас выполняется условие f(0) = 1 и f^{(n)}(0) = 0 (так как первая производная будет иметь нуль порядка N).

Пусть теперь та же функция равна f(x) = (x-1)^{N+1}L_{2}(x) мы получаем выполнение условий в т. x = 1

Приравниваем эти два определения и получаем уравнение:
1 + x^{N+1}L_{1}(x) = (x-1)^{N+1}L_{2}(x)

Положим теперь для простоты N — нечетное, тогда уравнение можно переписать в виде:

x^{N+1}L'_{1}(x) + (1-x)^{N+1}L_{2}(x) = 1, где L'_{1}(x) = -L_{1}(x)

Заметим, что наибольший общий делитель у многочленов (1-x)^{N+1} и x^{N+1} равен единице, а значит уравнение выше есть ни что иное, как уравнение Безу. Есть алгоритм его решения — расширенный алгоритм Евклида, однако именно для этого случая известен ответ в явном виде (из теории вейвлетов). Перепишем уравнение следующим образом:

x^{N+1}P(1-x) + (1-x)^{N+1}P(x) = 1
Одним из решений (с наименьшей степенью) этого уравнения будет многочлен P(x) = \sum_{n=0}^{N}C_{N+n}^{n}x^n

Теперь найдем f(x) пользуясь первым определением. Получаем L_{1}(x) = -L'_{1}(x) = -P(1-x). Тогда f(x) = 1 - x^{N+1}\sum_{n=0}^{N}C_{N+n}^{n}(1-x)^n
Это и есть ответ. Замечу, что ничего не изменится, если N — четное, только решение уравнения Безу в явном виде придется поискать.

UPD 25 jan 2019:
Приведенное решение, оказывается, годится и для четного N. Для того, чтобы в этом убедится достаточно представить f(x) в следующем виде:
f(x) = 1 - x^{N+1}P(1-x) = (1-x)^{N+1}P(x)
Конец апдейта

Приведу пару примеров. Для N = 1 получаем f(x) = 2x^3 - 3x^2+1.
Для N = 3 получаем f(x) = 20x^7 - 70x^6 + 84x^5-35x^4+1.
Вот график этих функций:
Чего тут нет — так это доказательства монотонности (т.е. нет локальных максимумов и минимумов внутри интервала). Ну и Б-г с ним.

Аналогично получается второе семейство функций — тригонометрических многочленов — обладающих теми же свойствами.

Так же как и в первый раз даем два определения f(\omega):
f(\omega) = 1 - (\frac{1 - \cos \pi \omega}{2})^N Q_{1}(\cos \pi \omega) и
f(\omega) = (\frac{1 + \cos \pi \omega}{2})^N Q_{2}(-\cos \pi \omega)

Первое определение удовлетворяет условию слева (x = 0), второе ­— справа (x = 1).

Вспоминаем, что \cos 2x = \cos^2 x - \sin^2 x = 1 - 2\sin^2 x = 2\cos^2 x -1, откуда получаем \sin^2 x = \frac{1-\cos 2x}{2}, \cos^2 x = \frac{1+\cos 2x}{2}

Обозначим x = \sin^2 \frac{\pi}{2} \omega и перепишем определения f(\omega) в термах x, приравняв их друг к другу.

1 - x^N Q_{1}(1-2x) = (1-x)^N Q_{2}(2x-1) или
x^N Q_{1}(1-2x) + (1-x)^N Q_{2}(2x-1) = 1

Пусть теперь P(x) = Q_{2}(2x-1), а P(1-x) = Q_{1}(1-2x). По сути, Q_1 и Q_2 — одни и те же многочлены, обозначены они по-разному из-за того, что мы не знаем этого заранее.

Получаем уравнение:
x^N P(1-x) + (1-x)^N P(x) = 1

Решение мы знаем:
 P(x) = \sum_{n=0}^{N-1}C_{N+n-1}^{n}x^n

Из P(x) находим Q_{1}(x) и Q_{1}(\cos \pi \omega).

Пару примеров. При N = 1 получаем тривиальное f(\omega) = \frac{1+\cos \pi \omega}{2}.

При N = 2 получаем P(x) = 1+2x, P(1-x) = 1+2(1-x) = 3 - 2x, а Q_{1}(1-2x)  = 2 + (1-2x) = 2 + \cos \pi \omega.

Тогда f(\omega) = 1 - (\frac{1- \cos \pi \omega}{2})^{2} (2 + \cos \pi \omega) = \frac{1}{16}(-\cos 3 \pi \omega + 9 \cos \pi \omega + 8)

График этих функций:
Вот так были выведены два семейства полезных многочленов для склейки двух постоянных "сегментов" функции так, чтобы результирующая функция имела N первых непрерывных производных.

Comments

Popular posts from this blog

На что влияет размер (глубина) семпла?

В этом посте я разберу ещё одну характеристику представления аудио на компьютере — размер или глубину семпла. В теореме Шеннона, с которой мы ознакомились в первом посте, сказано, что семплы нужно брать с определенной частотой, чтобы восстановить исходный сигнал. Однако, в ней ничего не сказано о том, сколько памяти компьютера нужно выделить на один семпл. По сути  ней семплы — это обычные вещественные числа. В компьютере же семплы чаще всего представляются целыми числами в диапазоне от -2^{n-1} до 2^{n-1}-1 (целые числа со знаком), где n — количество бит на один семпл (bps, bits per sample). С количеством бит на семпл связана другая величина — динамический диапазон (или отношение сигнал-шум, что в данном случае то же самое). Он измеряется в децибелах определяется так: DR = 20 \lg (2^{n}) \approx 6 n Динамический диапазон показывает логарифм отношения самого сильного сигнала к самому слабому. Теперь, что будет, если мы выберем n маленьким? Если взять синусоидальную ...

Ресемплинг аудио.

В этом посте расскажу о ресемплинге аудио (в частности, о даунсемплинге), под чем я подразумеваю понижение частоты дискретизации, с которой был записан звук. Известно, что в компьютере непрерывный аудио сигнал представлен в виде семплов (по-русски, измерений), сделанных в равно отстоящие друг от друга промежутки времени. Например, если звук записан с частотой дискретизации 44.1 кГц, то за одну секунду было сделано 44100 семплов. Здесь я расскажу о том, как восстановить исходный непрерывный сигнал по взятым семплам и как выбирать и менять частоту семплирования (или, что то же самое, частоту дискретизации). Преобразование Фурье. Преобразованим Фурье F(\xi) функции f(x) называется следующая функция: F(\xi) = \int_{-\infty}^{\infty} f(x) e^{-ix\xi}dx обратное преобразование определено как f(x) = \frac{1}{2\pi}\int_{-\infty}^{\infty} F(\xi) e^{ix\xi}d\xi Вопрос о том, когда существует преобразование Фурье и обратимо ли оно (т.е. работает ли вторая формула), очень сложе...