Журнал вычислительной математики и математической физики, 2020, T. 60, № 12, стр. 2015-2027
Задача минимизации суммы разностей взвешенных сверток
А. В. Кельманов 1, 2, Л. В. Михайлова 1, *, П. С. Рузанкин 1, 2, **, С. А. Хамидуллин 1, ***
1 Ин-т матем. им. С.Л. Соболева
630090 Новосибирск, пр-т Акад. Коптюга, 4, Россия
2 Новосибирский гос. ун-т
630090 Новосибирск, ул. Пирогова, 2, Россия
* E-mail: mikh@math.nsc.ru
** E-mail: ruzankin@math.nsc.ru
*** E-mail: kham@math.nsc.ru
Поступила в редакцию 30.07.2019
После доработки 30.07.2019
Принята к публикации 04.08.2020
Аннотация
В работе рассматривается неизученная экстремальная задача суммирования элементов числовых последовательностей $Y$ длины $N$ и $U$ длины $q \leqslant N$. В задаче требуется минимизировать сумму разностей взвешенных сверток последовательностей переменной длины (не менее $q$). В каждой разности первая невзвешенная свертка – автосвертка растянутой на переменную длину последовательности $U$ (путем кратных повторов ее элементов), вторая – взвешенная свертка этой растянутой последовательности с подпоследовательностью из $Y$. Анализируется вариант задачи с оптимизируемым числом суммируемых разностей. Показано, что задача эквивалентна одной из проблем аппроксимации последовательности $Y$ элементом $X$ из экспоненциального по мощности множества последовательностей. Это множество объединяет все последовательности длины $N$, которые в качестве подпоследовательностей включают переменное число $M$ допустимых квазипериодических (флуктуационных) повторов последовательности $U$. Каждый квазипериодический повтор порождается допустимыми преобразованиями последовательности $U$. Этими преобразованиями являются: 1) сдвиг $U$ на переменную величину, которая между соседними повторами не превышает ${{T}_{{max}}} \leqslant N$, 2) переменное растягивающее отображение $U$ в последовательность переменной длины, которое определяется в виде повторов элементов из $U$, кратность этих повторов – переменная величина. Критерием аппроксимации является минимум суммы квадратов расстояний между элементами последовательностей. Доказано, что рассматриваемая экстремальная задача и вместе с ней задача аппроксимации разрешимы за полиномиальное время. А именно, показано, что существует точный алгоритм, который находит решение задачи за время $\mathcal{O}(T_{{max}}^{3}N)$. Если ${{T}_{{max}}}$ – фиксированный параметр задачи, то время работы алгоритма линейно. Примерами численного моделирования проиллюстрирована применимость алгоритма к решению модельных прикладных задач помехоустойчивой обработки ECG-подобных и PPG-подобных квазипериодических сигналов (electrocardiogram-like and photoplethysmogram-like signals). Библ. 13. Фиг. 4.
1. ВВЕДЕНИЕ
В работе рассматривается неизученная экстремальная задача суммирования элементов числовых последовательностей. Цель работы – анализ вычислительной сложности задачи и обоснование алгоритма, обеспечивающего решение задачи с доказуемыми оценками качества (точности и трудоемкости). Исследование мотивировано новизной математической задачи и отсутствием каких-либо эффективных (полиномиальных) вычислительных алгоритмов с априорно гарантированными оценками точности для ее решения. Практическая мотивировка обусловлена важностью задачи, в частности, для помехоустойчивого дистанционного мониторинга природных объектов, состояние которых квазипериодически флуктуирует во времени (см. разд. 3, 7).
Статья имеет следующую структуру. В разд. 2 дана формулировка задачи, приведены ее истоки и некоторые интерпретации. В следующем разделе рассмотрены частные случаи задачи и наиболее близкие по постановке задачи. В разд. 4 анализируется размер (мощность) множества допустимых решений задачи. Вспомогательные утверждения, необходимые для алгоритмического решения задачи, доказаны в разд. 5. Основной результат работы – точный полиномиальный алгоритм – представлен и обоснован в разд. 6. В разд. 7 на примерах численного моделирования продемонстрирована применимость алгоритма к решению задач помехоустойчивой обработки ECG- и PPG-подобных сигналов.
2. ФОРМУЛИРОВКА ЗАДАЧИ, ЕЕ ИСТОКИ И ТРАКТОВКИ
Рассматриваемая экстремальная задача имеет следующую формулировку.
Задача 1. Дано: числовые последовательности $Y = ({{y}_{1}},\; \ldots ,\;{{y}_{N}})$, $U = ({{u}_{1}},\; \ldots ,\;{{u}_{q}})$, и натуральные числа ${{T}_{{max}}}$, $\ell $. Найти: набор $\mathcal{M} = {\text{\{ }}{{n}_{1}},\; \ldots {{n}_{m}},\; \ldots {\text{\} }}$ номеров последовательности $Y$, набор $\mathcal{P} = {\text{\{ }}{{p}^{{(1)}}},\; \ldots ,\;{{p}^{{(m)}}},\; \ldots {\text{\} }}$ натуральных чисел, набор $\mathcal{J} = {\text{\{ }}{{J}^{{(1)}}},\; \ldots ,\;{{J}^{{(m)}}},\; \ldots {\text{\} }}$ сжимающих отображений, в котором ${{J}^{{(m)}}}:{\text{\{ }}1,\; \ldots ,\;{{p}^{{(m)}}}{\text{\} }} \to {\text{\{ }}1,\; \ldots ,\;q{\text{\} }}$, а также размерность $M$ этих наборов, которые минимизируют целевую функцию
(2.1)
$F(\mathcal{M},\mathcal{P},\mathcal{J}) = \sum\limits_{m = 1}^M {\sum\limits_{i = 1}^{{{p}^{{(m)}}}} {\left\{ {u_{{{{J}^{{(m)}}}(i)}}^{2} - 2{{y}_{{{{n}_{m}} + i - 1}}}{{u}_{{{{J}^{{(m)}}}(i)}}}} \right\}} } ,$(2.2)
$q \leqslant {{p}^{{(m)}}} \leqslant \ell \leqslant {{T}_{{max}}} \leqslant N,\quad m = 1,\; \ldots ,\;M,$(2.3)
${{p}^{{(m - 1)}}} \leqslant {{n}_{m}} - {{n}_{{m - 1}}} \leqslant {{T}_{{max}}},\quad m = 2,\; \ldots ,\;M,$(2.4)
$\begin{array}{*{20}{c}} {{{J}^{{(m)}}}(1) = 1,\quad {{J}^{{(m)}}}({{p}^{{(m)}}}) = q,} \\ {0 \leqslant {{J}^{{(m)}}}(i) - {{J}^{{(m)}}}(i - 1) \leqslant 1,\quad i = 2,\; \ldots ,\;{{p}^{{(m)}}},} \\ {m = 1,\; \ldots ,\;M,} \end{array}$Приведем несколько трактовок задачи 1 и ее истоки. Далее будем считать, что номерам элементов входных последовательностей $Y$ и $U$ соответствуют дискретно-временные отсчеты непрерывных сигналов. На фиг. 1 в виде графиков изображены пример последовательности $U$ (слева) и пример последовательности $Y$ (справа). Раскрашенные участки последовательности $U$ нам понадобятся для одной из интерпретаций задачи.
Очевидно, что задачу 1 можно трактовать как задачу об оптимальном (в смысле минимума (2.1)) суммировании элементов двух последовательностей. Эта трактовка видна непосредственно из правой части формулы (2.1), в которой оптимизируемыми переменными являются индексы только двух последовательностей $U$ и $Y$.
Целевую функцию (2.1) можно записать в эквивалентном виде:
Наконец, задачу 1 можно интерпретировать как поиск оптимального разбиения целочисленного отрезка $[1,N]$ на $M + 1$ интервалов, границы которых определяются неравенствами $1 \leqslant {{n}_{1}} < \; \ldots \; < {{n}_{m}} < \; \ldots \; < {{n}_{M}} \leqslant N$.
Источником задачи 1 является задача аппроксимации последовательности $Y$ последовательностью $X \in \mathcal{X}$, где $\mathcal{X}$ – множество допустимых аппроксимирующих последовательностей, по критерию минимума суммы квадратов расстояний между элементами $Y$ и $X$. Предполагается, что каждый элемент множества $\mathcal{X}$ порождается допустимыми квазипериодическими (флуктуационными) повторами $m = 1,\; \ldots ,\;M$, последовательности $U$. При этом каждый квазипериодический повтор порождается допустимыми преобразованиями последовательности $U$. Этими преобразованиями являются: 1) сдвиг $U$ на переменную величину, которая между соседними повторами не превышает ${{T}_{{max}}} \leqslant N$, 2) переменное растягивающее отображение $U$ в последовательность переменной длины, которое определяется в виде повторов элементов из $U$, причем кратность каждого из этих повторов – переменная величина.
Рассмотрим следующую модель порождения последовательностью $U$ последовательности $X$ и покажем, что эта модель и квадратичный критерий аппроксимации индуцируют задачу 1. Иными словами, покажем, что задачу 1 можно интерпретировать как задачу приближения по упомянутому критерию.
Обозначим через $k_{t}^{{(m)}}$, $t = 1,\; \ldots ,\;q$, $m = 1,\; \ldots ,\;M$, кратность повтора $t$-го элемента из $U$ в ее $m$-м растянутом повторе
(2.5)
$(\underbrace {{{u}_{1}},\; \ldots ,\;{{u}_{1}}}_{k_{1}^{{(m)}}},\;\underbrace {{{u}_{2}}, \ldots ,\;{{u}_{2}}}_{k_{2}^{{(m)}}},\; \ldots ,\;\underbrace {{{u}_{q}},\; \ldots ,\;{{u}_{q}}}_{k_{q}^{{(m)}}})$(2.7)
$h_{i}^{{(m)}} = \left\{ \begin{gathered} {{u}_{1}},\quad {\text{если}}\quad i = 1,\; \ldots ,\;k_{1}^{{(m)}}, \hfill \\ {{u}_{2}},\quad {\text{если}}\quad i = k_{1}^{{(m)}} + 1,\; \ldots ,\;k_{1}^{{(m)}} + k_{2}^{{(m)}}, \hfill \\ \ldots \hfill \\ {{u}_{q}},\quad {\text{если}}\quad i = k_{1}^{{(m)}} + \; \ldots \; + k_{{q - 1}}^{{(m)}} + 1,\; \ldots ,\;k_{1}^{{(m)}} + \; \ldots \; + k_{q}^{{(m)}}, \hfill \\ 0,\quad {\text{если}}\quad i < 1,\quad i > {{p}^{{(m)}}}, \hfill \\ \end{gathered} \right.$(2.8)
${{x}_{n}} = \sum\limits_{m = 1}^M {h_{{n - {{n}_{m}} + 1}}^{{(m)}}} ,\quad n = 1,\; \ldots ,\;N.$Легко видеть, что при каждом $m$ образ множества $\{ 1,\; \ldots ,\;{{p}^{{(m)}}}\} $ при отображении ${{J}^{{(m)}}}$ можно представить схематично в виде
(2.9)
$k_{t}^{{(m)}} = \left| {\left\{ {i\,|\,{{J}^{{(m)}}}(i) = t,\;i \in \{ 1,\; \ldots ,\;{{p}^{{(m)}}}\} } \right\}} \right|,\quad t = 1,\; \ldots ,\;q,$(2.10)
$\begin{gathered} h_{i}^{{(m)}} = \left\{ \begin{gathered} {{u}_{{{{J}^{{(m)}}}(i)}}},\quad {\text{если}}\quad i = 1,\; \ldots ,\;{{p}^{{(m)}}}, \hfill \\ 0,\quad {\text{если}}\quad i < 1,\quad i > {{p}^{{(m)}}}, \hfill \\ \end{gathered} \right. \\ m = 1,\; \ldots ,\;M, \\ \end{gathered} $В соответствии с (2.7) и (2.8) участок последовательности $X$ с номерами элементов от $n = {{n}_{{m - 1}}}$ до $n = {{n}_{m}}$, соответствующий $(m - 1)$-му растянутому повтору, можно представить схематично в следующем виде (в нижней строке – номера элементов последовательности $X$):
Пример последовательности $X$ приведен на фиг. 2 в виде графика. На этой фигуре представлены 5 растянутых повторов, которые порождены последовательностью $U$, изображенной на фиг. 1. Кратности повторений элементов из $U$ и интервалы между двумя растянутыми повторами генерировались датчиком случайных чисел. Чередующимся раскрашенным участкам последовательности $U$ на фиг. 1 соответствуют чередующиеся растянутые (по отношению к участкам из $U$) раскрашенные участки в последовательности $X$ на фиг. 2. Ступенчатость сигнала на этой фигуре обусловлена кратными повторениями элементов из $U$.
Фиг. 2.
Пример последовательности $X$, порожденной последовательностью $U$ и ее допустимыми преобразованиями.
Поясним смысл ограничений (2.2) и (2.3) в формулировке задачи 1. В этих формулах ${{T}_{{max}}}$ обозначает верхнюю границу на интервал между двумя соседними растянутыми повторами, $\ell $ обозначает верхнюю границу на длину ${{p}^{{(m)}}}$ растяжения внутри каждого из этих интервалов, $q$ устанавливает нижнюю границу на длину допустимых растяжений. Нижняя граница ${{p}^{{(m - 1)}}}$ в формуле (2.3) означает, что каждый интервал между двумя последовательными растянутыми повторами не может быть меньше, чем длина растянутой последовательности из этого интервала. В целом формулы (2.2) и (2.3) устанавливают соотношение верхних и нижних границ для допустимых порождающих преобразований в виде сдвигов и растяжений последовательности $U$.
Для суммы квадратов расстояний между элементами $Y$ и $X$ с учетом (2.8) и (2.10) путем несложных преобразований получим равенство
Первая сумма в правой части этого равенства не зависит от каких-либо переменных задачи 1, а вторая – совпадает с (2.1). Поэтому минимизация суммы квадратов расстояний (аппроксимационная задача) эквивалентна сформулированной задаче 1.3. БЛИЗКИЕ ПО ПОСТАНОВКЕ ЗАДАЧИ
Приложениями задачи 1 являются проблемы обработки, анализа и интерпретации дискретно-временных сигналов, типичных для мониторинга состояний наземных, подземных, подводных, аэрокосмических, медицинских и других объектов. Эта задача, очевидно, является подходящей математической моделью для решения прикладных проблем из отмеченных областей в тех случаях, когда сигналы имеют квазипериодическую структуру, в виде флуктуирующих повторов участков сигнала (например, импульсов). Поэтому задачу 1 можно трактовать в терминах каждого из отмеченных приложений.
Задача 1 является обобщением ранее исследованной задачи дискретной оптимизации. Частный случай задачи 1, в котором ${{p}^{{(m)}}} = q$, ${{J}^{{(m)}}}(i) = i$, при каждом $m = 1,\; \ldots ,\;M$, исследовался в [1]–[3]. В [2] и [3] обоснован точный алгоритм с временем работы $\mathcal{O}({{T}_{{max}}}N)$.
В [1] исследовался вариант указанного частного случая, в котором $M$ – часть входа (не является оптимизируемой переменной). В цитируемой работе построен алгоритм, гарантирующий отыскание оптимального решения за время $\mathcal{O}({{T}_{{max}}}MN)$. Этот вариант частного случая задачи 1 послужил подходящей моделью, в частности, для решения прикладных проблем геофизического и аэрокосмического мониторинга подземных и космических объектов. Алгоритм с несущественными модификациями программно-инженерного плана был дважды переизложен в [4] и [5], а затем использован в центре космических исследований NASA. Он, в частности, позволил обнаружить не одну сотню новых экзопланет по сильно зашумленным астрофизическим данным [5], [6]. Авторы цитируемых работ отмечают, что известные подходы к решению такого рода прикладных мониторинговых проблем оказались существенно менее результативными в плане помехоустойчивости. Другим результативным применением алгоритма является решение обратных геофизических задач, связанных с проблемами выявления сейсмических источников, а также прогнозированием техногенных катастроф [7].
Что касается приложений задачи 1, то для нее в качестве иллюстрации в разд. 7 приведены лишь два примера к обработке медицинских сигналов, имеющих вид квазипериодической последовательности импульсов. Как известно, в биомедицинских импульсных последовательностях интервалы между повторами импульсов и форма самих импульсов в повторах сопровождаются возможными флуктуациями, обусловленными естественными для природы квазипериодическими процессами.
4. РАЗМЕР МНОЖЕСТВА ДОПУСТИМЫХ РЕШЕНИЙ
Оценим размер (мощность) множества $\mathcal{X}$, т.е. число допустимых решений задачи 1.
Заметим, что значение $M$, очевидно, ограничено сверху, т.е. $M \leqslant {{M}_{{max}}} \leqslant N$, где ${{M}_{{max}}} = \left\lfloor {\tfrac{N}{q}} \right\rfloor $ – максимально возможное число повторов последовательности длины не менее $q$ в последовательности длины $N \geqslant q$. Из комбинаторных соображений имеем оценку
(4.11)
$\left| \mathcal{X} \right| \leqslant ({{T}_{{max}}} - q + 1){{q}^{{{{T}_{{max}}} - q}}}{{M}_{{max}}}(N - q + 1){{({{T}_{{max}}} - q + 1)}^{{{{M}_{{max}}} - 1}}}.$Из анализа (4.11) ясно, что при фиксированном $q$ имеется
5. ОСНОВЫ АЛГОРИТМА
Ниже мы формулируем две вспомогательные экстремальные задачи и доказываем утверждения, гарантирующие отыскание их точных алгоритмических решений. Эти задачи и утверждения необходимы для обоснования алгоритма решения задачи 1.
Нам потребуется следующая вспомогательная
Задача 2. Дано: таблица $\left\{ {{{w}_{{i,j}}},\;i = 1,\; \ldots ,\;p,\;j = 1,\; \ldots ,\;q} \right\}$ вещественных чисел. Найти: индексы суммирования, доставляющие минимум функции
при ограничениях(5.12)
$\begin{array}{*{20}{c}} {J(1) = 1,\quad J(p) = q,} \\ {0 \leqslant J(i) - J(i - 1) \leqslant 1,\quad i = 2,\; \ldots ,\;p,} \end{array}$Справедлива следующая
Лемма 1. Минимум целевой функции $W$ в задаче 2 находится по формуле
а значение ${{W}_{{p,q}}}$ находится в результате вычислений по рекуррентной формуле(5.14)
${{W}_{{s,t}}} = min\left\{ {{{W}_{{s - 1,t}}},{{W}_{{s - 1,t - 1}}}} \right\} + {{w}_{{s,t}}},\quad s = 1,\; \ldots ,\;p,\quad t = 1,\; \ldots ,\;q,$(5.15)
${{W}_{{s,t}}} = \left\{ \begin{gathered} 0,\quad s = 0,\quad t = 0, \hfill \\ + \infty ,\quad s = 0,\quad t = 1,\; \ldots ,\;q, \hfill \\ + \infty ,\quad s = 1,\; \ldots ,\;p,\quad t = 0. \hfill \\ \end{gathered} \right.$(5.16)
$\begin{array}{*{20}{c}} {J{\text{*}}(p) = q,} \\ {J{\text{*}}(i - 1) = \left\{ \begin{gathered} J{\text{*}}(i),\quad \;\;\;\;{\text{если}}\quad {{W}_{{i - 1,J*(i)}}} \leqslant {{W}_{{i - 1,J*(i) - 1}}}, \hfill \\ J{\text{*}}(i) - 1,\quad {\text{если}}\quad {{W}_{{i - 1,J*(i)}}} > {{W}_{{i - 1,J*(i) - 1}}}, \hfill \\ \end{gathered} \right.} \\ {i = p,\;p - 1,\; \ldots ,\;2.} \end{array}$Доказательство. Легко видеть, что ограничения (5.12) определяют множество
(5.17)
${{\mathcal{J}}_{q}}(p) = {\text{\{ }}J:{\text{\{ }}1,\; \ldots ,\;p{\text{\} }} \to {\text{\{ }}1,\; \ldots ,\;q{\text{\} }}\,|\,J(1) = 1,\;J(p) = q,\;0 \leqslant J(i) - J(i - 1) \leqslant 1,\;i = 2,\; \ldots ,\;p{\text{\} }}$Обозначим задачу 2 через $\left\langle {p,q} \right\rangle $ и погрузим эту задачу в семейство подзадач $\left\langle {s,t} \right\rangle $ поиска минимальных значений
Из принципа оптимальности Беллмана и (5.17) следует, что значение ${{W}_{{s,t}}}$ может быть выбрано только из двух возможных вариантов:
(5.18)
${\text{\{ }}{{W}_{{s - 1,t}}} + {{w}_{{s,t}}},\;{{W}_{{s - 1,t - 1}}} + {{w}_{{s,t}}}{\text{\} }},$(5.19)
$\begin{array}{*{20}{c}} {{{W}_{{s - 1,t}}} = \mathop {min}\limits_{J \in {{\mathcal{J}}_{t}}(s - 1)} \sum\limits_{i = 1}^{s - 1} {{{w}_{{i,J(i)}}}} ,} \\ {{{W}_{{s - 1,t - 1}}} = \mathop {min}\limits_{J \in {{\mathcal{J}}_{{t - 1}}}(s - 1)} \sum\limits_{i = 1}^{s - 1} {{{w}_{{i,J(i)}}}} .} \end{array}$Из (5.18) и (5.19) следует справедливость формул (5.14) и (5.13). Справедливость начального и граничных условий в формуле (5.15) следует из (5.17). Правило (5.16) является очевидным следствием формул (5.13), (5.14) и (5.15). Лемма 1 доказана.
Кроме того, нам потребуется следующая вспомогательная
Задача 3. Дано: совокупность $\{ {{g}_{p}}(n),\;p = q,\; \ldots ,\;\ell ,\;n = 1,\; \ldots ,\;N - p + 1\} $ числовых последовательностей и натуральное число ${{T}_{{max}}} \geqslant \ell $. Найти: набор ${\text{\{ }}({{n}_{1}},{{p}_{1}}),\; \ldots ,\;({{n}_{M}},{{p}_{M}}){\text{\} }}$ пар натуральных чисел, в котором ${{n}_{m}} \in {\text{\{ }}1,\; \ldots ,\;N - {{p}_{m}} + 1{\text{\} }}$, ${{p}_{m}} \in {\text{\{ }}q,\; \ldots ,\;\ell {\text{\} }}$, $m = 1,\; \ldots ,\;M$, а также размерность $M$ этого набора, доставляющие минимум целевой функции
(5.20)
$G(({{n}_{1}},{{p}_{1}}),\; \ldots ,\;({{n}_{M}},{{p}_{M}})) = \sum\limits_{m = 1}^M {{{g}_{{{{p}_{m}}}}}} ({{n}_{m}}),$Справедлива следующая
Лемма 2. Минимум целевой функции $G$ в задаче 3 находится по формуле
(5.21)
$G{\kern 1pt} * = \mathop {min}\limits_{p \in {\text{\{ }}q, \ldots ,\ell {\text{\} }}} \mathop {min}\limits_{n \in {\text{\{ }}1, \ldots ,N - p + 1{\text{\} }}} {{G}_{p}}(n),$(5.22)
${{G}_{p}}(n) = \left\{ \begin{gathered} {{g}_{p}}(n),\quad p = q,\; \ldots ,\;\ell ,\quad n = 1,\; \ldots ,\;q, \hfill \\ min\left\{ {0,\mathop {min}\limits_{i \in {\text{\{ }}q, \ldots ,min{\text{\{ }}n - 1,\ell {\text{\} \} }}} \mathop {min}\limits_{j \in {{\gamma }_{i}}(n)} {{G}_{i}}(j)} \right\} + {{g}_{p}}(n),\quad p = q,\; \ldots ,\;\ell ,\quad n = q + 1,\; \ldots ,\;N - p + 1, \hfill \\ \end{gathered} \right.$(5.23)
${{\gamma }_{i}}(n) = \left\{ {k\,|\,max{\text{\{ }}n - {{T}_{{max}}},1{\text{\} }} \leqslant k \leqslant n - i} \right\},\quad n = i + 1,\; \ldots ,\;N,\quad i = q,\; \ldots ,\;\ell .$Доказательство. Опираясь на ограничения (2.2) и (2.3), определяем множество
(5.24)
$\begin{gathered} \Phi = \left\{ {{\text{\{ }}({{n}_{1}},{{p}_{1}}),\; \ldots ,\;({{n}_{M}},{{p}_{M}}){\text{\} }}} \right.\,|\,{{p}_{m}} \in {\text{\{ }}q,\; \ldots ,\;\ell {\text{\} }},\;{{n}_{m}} \in {\text{\{ }}1,\; \ldots ,\;N - {{p}_{M}} + 1{\text{\} }},\;m = 1, \ldots ,\;M; \\ {{p}_{{m - 1}}} \leqslant {{n}_{m}} - {{n}_{{m - 1}}} \leqslant {{T}_{{max}}},\;\left. {m = 2,\; \ldots ,\;M} \right\} \\ \end{gathered} $Для упрощения записи доказательства положим
В соответствии с (5.20), (5.24) и (5.25), для оптимального значения целевой функции $G$ имеемОпределим подмножество
(5.27)
${{\Phi }_{p}}(n) = \left\{ {\varphi \in \Phi \,|\,({{n}_{M}},{{p}_{M}}) = (n,p)} \right\},\quad p = q,\; \ldots ,\;\ell ,\quad n = 1,\; \ldots ,\;N - p + 1,$Легко видеть, что
(5.28)
$\Phi = \bigcup\limits_{p = q}^\ell {\bigcup\limits_{n = 1}^{N - p + 1} {{{\Phi }_{p}}} } (n).$Заметим, что для искомой размерности $M$ наборов возможны два варианта: $M = 1$ и $M \geqslant 1$. Поэтому при каждом $p = q,\; \ldots ,\;\ell $ множество (5.27) в зависимости от $M$ может быть представлено в виде
(5.29)
${{\Phi }_{p}}(n) = \left\{ \begin{gathered} {\text{\{ }}\varphi \in {{\Phi }_{p}}(n)\,|\,M = 1{\text{\} }},\quad n = 1,\; \ldots ,\;p - 1, \hfill \\ {\text{\{ }}\varphi \in {{\Phi }_{p}}(n)\,|\,M = 1{\text{\} }} \cup \{ {\kern 1pt} \varphi \in {{\Phi }_{p}}(n)\,|\,M > 1{\text{\} }},\quad n = p,\; \ldots ,\;N - p + 1. \hfill \\ \end{gathered} \right.$Для индексов $p = q,\; \ldots ,\;\ell $, $n = 1,\; \ldots ,\;N - p + 1$ справедливо очевидное равенство
(5.30)
${\text{\{ }}\varphi \in {{\Phi }_{p}}(n)\,|\,M = 1{\text{\} }} = {\text{\{ }}(n,p){\text{\} }}.$Рассмотрим множество ${\text{\{ }}\varphi \in {{\Phi }_{p}}(n)\,|\,M > 1{\text{\} }}$, $p = q,\; \ldots ,\;\ell $, $n = p,\; \ldots ,\;N - p + 1$. Из ограничений (2.2), (2.3) и определения (5.23) следует, что если значения ${{n}_{m}}$ и ${{p}_{m}}$ фиксированы, то при каждом $m = 2,\; \ldots ,\;M$ для значений компонент $(m - 1)$-й пары справедливо
(5.31)
$\begin{gathered} {\text{\{ }}\varphi \in {{\Phi }_{p}}(n)\,|\,M > 1{\text{\} }} = {\text{\{ }}{\kern 1pt} {\text{\{ }}({{n}_{1}},{{p}_{1}}),\; \ldots ,\;({{n}_{M}},{{p}_{M}}){\text{\} }}\,|\,{{n}_{M}} = n,\;{{p}_{M}} = p, \\ {{p}_{{m - 1}}} \in {\text{\{ }}q,\; \ldots ,\;min{\text{\{ }}\ell ,{{n}_{m}} - 1{\text{\} }}{\kern 1pt} {\text{\} }},\;{{n}_{{m - 1}}} \in {{\gamma }_{{{{p}_{{m - 1}}}}}}({{n}_{m}}),\;m = 2,\; \ldots ,\;M{\text{\} }} = \\ = \;{\text{\{ }}{\kern 1pt} {\text{\{ }}({{n}_{1}},{{p}_{1}}),\; \ldots ,\;({{n}_{{M - 1}}},{{p}_{{M - 1}}}),(n,p){\text{\} }}\,|\,{\text{\{ }}({{n}_{1}},{{p}_{1}}), \ldots ,\;({{n}_{{M - 1}}},{{p}_{{M - 1}}}){\text{\} }} \in {{\Phi }_{i}}(j), \\ i \in {\text{\{ }}q,\; \ldots ,\;min{\text{\{ }}\ell ,n - 1{\text{\} }}{\kern 1pt} {\text{\} }},\;j \in {{\gamma }_{i}}(n){\text{\} }}. \\ \end{gathered} $Докажем формулу (5.22), используя полученные равенства (5.30) и (5.31). Рассмотрим семейство задач $\left\langle {p,n} \right\rangle $ вычисления оптимальных значений
(5.32)
${{G}_{p}}(n) = \mathop {min}\limits_{\varphi \in {{\Phi }_{p}}(n)} G(\varphi ),\quad p = q,\; \ldots ,\;\ell ,\quad n = 1,\; \ldots ,\;N - p + 1.$Два варианта вычислений, представленных в (5.22), следуют из двух вариантов в формуле (5.29). Для первого варианта формулы (5.22), когда $p = q,\; \ldots ,\;\ell $, $n = 1,\; \ldots ,\;p - 1$, из (5.32), (5.29) и (5.30) имеем
(5.33)
${{G}_{p}}(n) = \mathop {min}\limits_{\varphi \in {{\Phi }_{p}}(n),M = 1} G(\varphi ) = \mathop {min}\limits_{\varphi \in {\text{\{ }}(n,p){\text{\} }}} G(\varphi ) = {{g}_{p}}(n).$Для второго варианта формулы (5.22), когда $p = q,\; \ldots ,\;\ell $, $n = p,\; \ldots ,\;N - p + 1$, из (5.32) и (5.29) имеем
(5.34)
${{G}_{p}}(n) = min\left\{ {\mathop {min}\limits_{\varphi \in {{\Phi }_{p}}(n),M = 1} G(\varphi ),\mathop {min}\limits_{\varphi \in {{\Phi }_{p}}(n),M > 1} G(\varphi )} \right\}.$Далее, в формуле (5.34) для первого минимума справедлива цепочка равенств (5.33). Второй минимум в формуле (5.34) в соответствии с (5.31) и принципом оптимальности Беллмана выбирается из совокупности
(5.35)
$\mathop {min}\limits_{\varphi \in {{\Phi }_{p}}(n),M > 1} {{G}_{p}}(n) = \mathop {min}\limits_{i \in {\text{\{ }}q, \ldots ,min{\text{\{ }}n - 1,\ell {\text{\} \} }}} \mathop {min}\limits_{j \in {{\gamma }_{i}}(n)} {{G}_{i}}(j) + {{g}_{p}}(n).$Для отыскания компонент оптимальных наборов определим две функции:
(5.36)
$\pi (n) = \left\{ \begin{gathered} 0,\quad {\text{если}}\quad n = 1, \ldots ,q, \hfill \\ 0,\quad {\text{если}}\quad \mathop {min}\limits_{i \in \{ q, \ldots ,min\{ n - 1,\ell \} {\kern 1pt} \} } \mathop {min}\limits_{j \in {{\gamma }_{i}}(n)} {{G}_{i}}(j) \geqslant 0,\quad n = q + 1,\; \ldots ,\;N - q + 1, \hfill \\ arg\mathop {min}\limits_{i \in \{ q, \ldots ,min\{ n - 1,\ell \} {\kern 1pt} \} } \left\{ {\mathop {min}\limits_{j \in {{\gamma }_{i}}(n)} {{G}_{i}}(j)} \right\},\quad {\text{если}}\quad \hfill \\ \mathop {min}\limits_{i \in \{ q, \ldots ,min\{ n - 1,\ell \} {\kern 1pt} \} } \mathop {min}\limits_{j \in {{\gamma }_{i}}(n)} {{G}_{i}}(j) < 0,\quad n = q + 1,\; \ldots ,\;N - q + 1, \hfill \\ \end{gathered} \right.$(5.37)
$\begin{gathered} I(n) = \left\{ \begin{gathered} - 1,\quad {\text{если}}\quad \pi (n) = 0, \hfill \\ arg\mathop {min}\limits_{j \in {{\gamma }_{{\pi (n)}}}(n)} {{G}_{{\pi (n)}}}(j),\quad {\text{если}}\quad \pi (n) > 0, \hfill \\ \end{gathered} \right. \\ n = 1,\; \ldots ,\;N - q + 1. \\ \end{gathered} $Следствие 1. Пусть выполнены условия леммы 2. Допустим, кроме того, что
(5.38)
${{\pi }_{1}} = arg\mathop {min}\limits_{p \in {\text{\{ }}q, \ldots ,\ell {\text{\} }}} \left\{ {\mathop {min}\limits_{n \in {\text{\{ }}1, \ldots ,N - p + 1{\text{\} }}} {{G}_{p}}(n)} \right\},$(5.39)
${{\nu }_{1}} = arg\mathop {min}\limits_{n \in {\text{\{ }}1, \ldots ,N - {{\pi }_{1}} + 1{\text{\} }}} {{G}_{{{{\pi }_{1}}}}}(n),$(5.40)
${{\pi }_{m}} = \pi ({{\nu }_{{m - 1}}}),\quad {{\nu }_{m}} = I({{\nu }_{{m - 1}}}),\quad m = 2,\; \ldots ,\;M,$6. АЛГОРИТМ
Опираясь на вспомогательные задачи и утверждения из разд. 5, формулируем алгоритм решения задачи 1.
Алгоритм $\mathcal{A}$
Вход: $Y$, $U$, ${{T}_{{max}}}$ и $\ell $.
Прямой ход алгоритма.
Шаг 1 (построение семейства решений задач 2). Для каждого $n = 1,\; \ldots ,\;N - q + 1$ выполним:
Для каждого $p = q,\; \ldots ,\;min{\text{\{ }}\ell ,N - n + 1{\text{\} }}$ выполним (Решение задачи 2):
(1) Вычислим
(6.41)
${{w}_{{s,t}}} = u_{t}^{2} - 2{{y}_{{n + s - 1}}}{{u}_{t}},\quad s = 1,\; \ldots ,\;p,\quad t = 1,\; \ldots ,\;q.$(2) Найдем значение $W{\text{*}}$ по формулам (5.13), (5.14), (5.15) и последовательность $J{\text{*}}(1),\; \ldots ,\;J{\text{*}}(p)$ по формулам (5.16).
(3) Положим $W{\text{*}}(n,p) = W{\text{*}}$; $J{\text{*}}(n,p) = {\text{\{ }}J{\text{*}}(i),\;i = 1,\; \ldots ,\;p{\text{\} }}$.
Конец построения семейства решений задач 2.
Шаг 2 (определение входа задачи 3). Для каждого $n = 1,\; \ldots ,\;N - p + 1$ положим ${{g}_{p}}(n) = W{\text{*}}(n,p),$ $p = q,\; \ldots ,\;\ell $.
Шаг 3 (решение задачи 3). Вычислим совокупность значений ${{G}_{p}}(n)$, $p = q,\; \ldots ,\;\ell $, $n = 1,\; \ldots ,\;N - p + 1$, по формулам (5.22), (5.23), а так же значение $G{\text{*}}$ по формуле (5.21). Положим ${{F}_{A}} = G{\text{*}}$.
Обратный ход алгоритма.
Шаг 4. Вычислим $\pi (n)$ и $I(n)$, $n = 1,\; \ldots ,\;N - q + 1$, по формулам (5.36) и (5.37).
Шаг 5. Найдем компоненты наборов $({{\pi }_{1}},\; \ldots ,\;{{\pi }_{M}})$ и $({{\nu }_{1}},\; \ldots ,\;{{\nu }_{M}})$ и их размерность $M$ по формулам (5.38), (5.39) и (5.40).
Шаг 6. Положим ${{M}_{A}} = M$, ${{\mathcal{M}}_{A}} = {\text{\{ }}{{\nu }_{M}},\; \ldots ,\;{{\nu }_{1}}{\text{\} }}$, ${{\mathcal{P}}_{A}} = {\text{\{ }}{{\pi }_{M}},\; \ldots ,\;{{\pi }_{1}}{\text{\} }}$, ${{\mathcal{J}}_{A}} = {\text{\{ }}J{\text{*}}({{\nu }_{M}},{{\pi }_{M}}),\; \ldots ,\;J{\text{*}}({{\nu }_{1}},{{\pi }_{1}}){\text{\} }}{\kern 1pt} .$
Выход: ${{M}_{A}}$, ${{\mathcal{M}}_{A}}$, ${{\mathcal{P}}_{A}}$, ${{\mathcal{J}}_{A}}$ и ${{F}_{A}}$.
Справедлива следующая
Теорема 1. Алгоритм $\mathcal{A}$ находит точное решение задачи 1 за время $\mathcal{O}(T_{{max}}^{3}N)$.
Доказательство. Для оптимального значения $F{\text{*}}$ целевой функции задачи 1 справедлива следующая цепочка равенств:
(6.42)
$\begin{gathered} F{\text{*}}\;\mathop = \limits^1 \;\mathop {min}\limits_{\mathcal{M},\mathcal{P},\mathcal{J}} \sum\limits_{m = 1}^M {\sum\limits_{i = 1}^{{{p}_{m}}} {\left\{ {u_{{{{J}^{{(m)}}}(i)}}^{2} - 2{{y}_{{{{n}_{m}} + i - 1}}}{{u}_{{{{J}^{{(m)}}}(i)}}}} \right\}} } \;\mathop = \limits^2 \;\mathop {min}\limits_{\mathcal{M},\mathcal{P}} \mathop {min}\limits_{{{J}^{{(1)}}}, \ldots ,{{J}^{{(M)}}}} \sum\limits_{m = 1}^M {\sum\limits_{i = 1}^{{{p}_{m}}} {\left\{ {u_{{{{J}^{{(m)}}}(i)}}^{2} - 2{{y}_{{{{n}_{m}} + i - 1}}}{{u}_{{{{J}^{{(m)}}}(i)}}}} \right\}} } = \\ \mathop = \limits^3 \,\mathop {min}\limits_{\mathcal{M},\mathcal{P}} \sum\limits_{m = 1}^M {\left\{ {\mathop {min}\limits_{{{J}^{{(m)}}} \in {{\mathcal{J}}_{q}}({{p}_{m}})} \sum\limits_{i = 1}^{{{p}_{m}}} {\left\{ {u_{{{{J}^{{(m)}}}(i)}}^{2} - 2{{y}_{{{{n}_{m}} + i - 1}}}{{u}_{{{{J}^{{(m)}}}(i)}}}} \right\}} } \right\}} = \\ \mathop = \limits^4 \,\mathop {min}\limits_{\mathcal{M},\mathcal{P}} \sum\limits_{m = 1}^M {\left\{ {\mathop {min}\limits_{J \in {{\mathcal{J}}_{q}}({{p}_{m}})} \sum\limits_{i = 1}^{{{p}_{m}}} {\left\{ {u_{{J(i)}}^{2} - 2{{y}_{{{{n}_{m}} + i - 1}}}{{u}_{{J(i)}}}} \right\}} } \right\}} \;\mathop = \limits^5 \;\mathop {min}\limits_{\mathcal{M},\mathcal{P}} \sum\limits_{m = 1}^M {W{\text{*}}} ({{n}_{m}},{{p}_{m}})\;\mathop = \limits^6 \;{{F}_{A}}. \\ \end{gathered} $В цепочке равенств (6.42) равенство 1 – определение оптимального решения задачи 1. Равенство 2 очевидно. Справедливость равенства 3 следует из того, что каждый элемент внутренней суммы в левой части равенства 2 зависит только от ${{J}^{{(m)}}}$ при каждом $m = 1,\; \ldots ,\;M$. Равенство 4 очевидно. Равенство 5 следует из того, что выражение во внешних скобках в левой части равенства 5 есть оптимальное решение (в соответствии с леммой 1) задачи 2, полученное при выполнении шага 1 алгоритма. Равенство 6 справедливо, так как левая часть этого равенства, согласно лемме 2, есть оптимальное значение целевой функции задачи 3 с входными данными, заданными на шаге 2 алгоритма. Таким образом, оптимальность решения, найденного алгоритмом, следует из (6.42).
Оценим трудоемкость алгоритма, используя его пошаговую запись. На шаге 1 при каждом $n = 1,\; \ldots ,\;N - p + 1$ задача 2 решается, очевидно, за время $\mathcal{O}(q{{T}_{{max}}})$, так как $q \leqslant p \leqslant \ell \leqslant {{T}_{{max}}}$. Последовательность решений этих задач находится за время $\mathcal{O}(T_{{max}}^{2}N)$, так как при каждом $n = 1,\; \ldots ,\;N - p$ имеет место вложенность таблиц:
Шаг 2 выполняется за время $\mathcal{O}({{T}_{{max}}}N)$, так как $\ell \leqslant {{T}_{{max}}}$. На шаге 3 при каждом $n = 1,\; \ldots ,\;N - p + 1$ и каждом $p = q,\; \ldots ,\;\ell $ вычисление ${{G}_{p}}(n)$ по формуле (5.22) выполняется за время $\mathcal{O}(T_{{max}}^{2})$, так как $q \leqslant p \leqslant \ell \leqslant {{T}_{{max}}}$ и $\mathop {max}\limits_{i,n} \left| {{{\gamma }_{i}}(n)} \right| \leqslant {{T}_{{max}}}$ в соответствии с (5.23). Следовательно, итоговая трудоемкость шага 3 оценивается величиной $\mathcal{O}(T_{{max}}^{3}N)$.
На шаге 4 вычисление по формулам (5.36) и (5.37), как и на шаге 1, выполняется за время $\mathcal{O}(T_{{max}}^{2}N)$. Шаги 5 и 6 выполняются за время $\mathcal{O}(N)$.
Суммируя временные затраты на всех шагах алгоритма, получаем итоговую трудоемкость алгоритма $\mathcal{A}$. Теорема 1 доказана.
Замечание 1. Если ${{T}_{{max}}}$ – фиксированный параметр, то время работы алгоритма равно $\mathcal{O}(N)$.
В следующем вспомогательном демонстрационном разделе мы приводим примеры применения полученного результата к решению двух модельных задач прикладного плана.
7. ПРИМЕРЫ ЧИСЛЕННОГО МОДЕЛИРОВАНИЯ
Очевидно, что с математической точки зрения не имеет значения, импульсу какой формы соответствует входная последовательность $U$. Тем не менее ниже для иллюстрации возможного применения алгоритма в биомедицинских приложениях приведены два примера численного моделирования процесса обработки временных рядов (последовательностей), которые можно трактовать как квазипериодические последовательности флуктуирующих ECG-подобных (фиг. 3) и PPG-подобных импульсов (фиг. 4), искаженных аддитивной помехой.
На каждой фигуре в верхней части слева изображена последовательность $U$, соответствующая некоторому ECG- или PPG-импульсу. Описание формы этих импульсов, включая их характерные участки и значимые точки, широко представлено в литературе и интернете (см., например, [8]–[13]). На фиг. 3 этими участками являются P-Q-R-S-T-U, а на фиг. 4 – систолическая точка, дикротическая точка, диастолическая точка. Характерные участки или интервалы между значимыми точками импульса выделены раскраской.
Справа от этого импульса изображена модельная программно сгенерированная квазипериодическая последовательность, соответствующая некоторой импульсной последовательности, порожденной квазипериодическими флуктуирующими повторами импульса $U$. В модельной последовательности каждый квазипериодический повтор последовательности $U$ сопровождается ее флуктуирующими растяжениями, а интервалы между повторами не постоянны. Квазипериодичность этих повторов видна по растянутым характерным участкам импульса, которые выделены цветом.
В средней части каждой фигуры (под модельной последовательностью) изображена последовательность $Y$. Эта последовательность генерировалась как результат сложения модельной последовательности с последовательностью независимых одинаково распределенных гауссовских случайных величин с нулевым математическим ожиданием. Высокий уровень помехи (дисперсии) моделировался специально для демонстрации способности алгоритма как инструмента для обработки сильно зашумленных данных (как известно, высокий уровень помехи нетипичен для биомедицинских приложений). Входными последовательностями для алгоритма являются $U$ и $Y$.
Последовательность ${{X}_{A}}$, приведенная в нижней части фигуры, является оценкой для модельной (ненаблюдаемой) последовательности в верхней части фигуры. Элементы ${{X}_{A}}$ легко находятся по входной последовательности $U$ и трем выходным наборам алгоритма по формулам (2.9), (2.7), (2.8) разд. 2 при $M = {{M}_{A}}$, $\mathcal{M} = {{\mathcal{M}}_{A}}$, $\mathcal{P} = {{\mathcal{P}}_{A}}$, $\mathcal{J} = {{\mathcal{J}}_{A}}$. По этим же формулам можно найти оценки для характерных растянутых участков и точек. Данные на фиг. 3 получены при $q = 128$, ${{T}_{{max}}} = 370$, $N = 2000$, а на фиг. 4 при $q = 120$, ${{T}_{{max}}} = 280$, $N = 1200$.
Сравнивая визуально верхние и нижние графики (временные ряды), можно оценить качество получаемых решений как вполне приемлемое для практики помехоустойчивой обработки данных в виде квазипериодических последовательностей флуктуирующих импульсов. В настоящей работе мы ограничиваемся лишь приведенными примерами, так как прикладные задачи находятся вне рамок настоящей работы. Здесь нас интересовали лишь новая дискретная экстремальная задача и алгоритм ее решения как математический инструмент для робастного решения сходных прикладных проблем. Вместе с этим отметим, что из любопытства алгоритм был опробован на целом ряде модельных импульсов различной формы и показал такую же высокую результативность.
ЗАКЛЮЧЕНИЕ
В работе установлена полиномиальная разрешимость не изученной ранее экстремальной задачи. Построен точный алгоритм, который имеет линейное время выполнения при фиксации одного из входов задачи как параметра.
В численных экспериментах мы продемонстрировали, что предложенный алгоритм послужит подходящим инструментом для решения некоторых прикладных проблем помехоустойчивой обработки квазипериодических импульсных последовательностей, в частности проблем обработки, анализа и интерпретации ECG и PPG сигналов.
Практический интерес представляет применение предложенного алгоритма к решению задачи распознавания, в которой имеется не одна последовательность $U$, а алфавит последовательностей, содержащий последовательности, соответствующие эталонным импульсам, например, каких-либо заболеваний.
С математической точки зрения значительный интерес представляет исследование дискретной экстремальной задачи, в которой по входной последовательности $Y$, имеющей квазипериодическую структуру, требуется найти последовательность $U$ фиксированной длины. Не меньший математический интерес представляет модифицированная задача 1, в которой число суммируемых разностей сверток является частью входа задачи. Исследование этих задач представляется делом ближайшей перспективы.
Список литературы
Кельманов А.В., Хамидуллин С.А. Апостериорное обнаружение заданного числа одинаковых подпоследовательностей в квазипериодической последовательности // Ж. вычисл. матем. и матем. физ. 2001. Т. 41. № 5. С. 807–820. Kel’manov A.V., Khamidullin S.A. Posterior detection of a given number of identical subsequences in a quasi-periodic sequence // Comp. Мath. and Math. Phys. 2001. V. 41. № 5. P. 762–774. (перевод)
Kel’manov A.V., Khamidullin S.A., Okol’nishnikova L.V. A Posteriori Detection of Identical Subsequences in a Quasiperiodic Sequence // Pattern Recognition and Image Analysis. 2002. V. 12. № 4. P. 438–447.
Кельманов А.В., Хамидуллин С.А., Окольнишникова Л.В. Апостериорное обнаружение одинаковых подпоследовательностей-фрагментов в квазипериодической последовательности // Сиб. ж. индустр. математики. 2002. Т. 5. № 2(10). С. 94–108.
Kel’manov A.V., Jeon B. A posteriori joint detection and discrimination of pulses in a quasiperiodic pulse train // IEEE Transactions on Signal Processing. 2004. V. 52. № 3. P. 645–656.
Carter J.A., Agol E., et al. Kepler-36: a pair of planets with neighboring orbits and dissimilar densities // Science. 2012. V. 337. № 6094. P. 556–559.
Carter J.A., Agol E. The quasiperiodic automated transit search algorithm // The Astrophysical Journal. 2013. V. 765. № 2. P. 132.
Voskoboynikova G., Khairetdinov M. Numerical modeling of posteriori algorithms for the geophysical monitoring // Communications in Computer and Information Science. 2015. V. 549. P. 190–200.
https://www.stepwards.com/page_id=24147
Zhang D., Zuo W., Wang P. Computational Pulse Signal Analysis. Springer Nature Singapore Pte Ltd, 2018. 28 p.
Rajni R., Kaur I. Electrocardiogram signal analysis – an overview // Int. J. Comput. Appl. 2013. V. 84. № 7. P. 22–25.
https://www.comm.utoronto.ca/~biometrics/PPG_Dataset/index.html
Shelley K., Shelley S. Pulse Oximeter Waveform: Photoelectric Plethysmography. In: Carol Lake, Hines R., Blitt C. (eds). Clinical Monitoring. W.B. Saunders Company. 2001. P. 420–428.
Elgendi M. On the analysis of fingertip photoplethysmogram signals // Current Cardiology Reviews. 2012. V. 8. № 1. P. 14–25.
Дополнительные материалы отсутствуют.
Инструменты
Журнал вычислительной математики и математической физики
Пн.-пт.: 10.00-19.00