CC BY
98
РАДИОЭЛЕКТРОНИКА
УДК 621.513.6
В. И. Д ж и г а н
ПАРАЛЛЕЛЬНЫЕ ВЫЧИСЛЕНИЯ В RLS-АЛГОРИТМАХ АДАПТИВНОЙ ФИЛЬТРАЦИИ
Рассмотрен способ представления RLS-алгоритмов адаптивной фильтрации со скользящим окном и регуляризацией корреляционной матрицы адаптивного фильтра, позволяющий реализовывать такие алгоритмы и их частные случаи с помощью параллельных вычислений. Данные алгоритмы математически эквивалентны известным одноименным последовательным и параллельным алгоритмам адаптивной фильтрации. Результаты моделирования демонстрируют работоспособность рассмотренных алгоритмов.
Цифровая обработка сигналов (ЦОС) — одна из областей современной радиотехники — имеет в своем составе ряд научно-технических направлений. Одно из таких направлений — это адаптивная фильтрация [1]. Адаптивные фильтры используются для решения радиотехнических задач широкого класса [2-4]. Качество решения этих задач определяется применяемыми в адаптивных фильтрах алгоритмами. В большинстве таких фильтров в основном используются простейшие (с точки зрения числа арифметических операций) градиентные алгоритмы по критерию наименьшего среднеквадратического отклонения (Least Mean Squares, LMS) или нормализованные LMS (Normalized LMS, NLMS) алгоритмы. Такие алгоритмы привлекательны своей низкой вычислительной сложностью (числом арифметических операций, необходимых для выполнения одной итерации алгоритма в течение одного периода дискретизации обрабатываемых сигналов), но имеют некоторые недостатки, среди которых медленная сходимость, проявляющаяся при обработке нестационарных сигналов, и зависимость длительности переходного процесса и остаточной ошибки на выходе адаптивного фильтра от шага сходимости. Такие недостатки отсутствуют в рекурсивных алгоритмах по критерию наименьших квадратов (Recursive Least Squares, RLS) [5]. Однако вычислительная сложность RLS-алгоритмов значительно превышает вычислительную сложность простых алгоритмов адаптивной фильтрации, что на протяжении многих лет служило весомым препятствием на пути применения сложных алгоритмов в приложениях.
Возможность реализации того или иного алгоритма ЦОС определяется производительностью используемых при этом цифровых сигнальных процессоров (ЦСП). Производительность современных ЦСП уже позволяет эффективно реализовывать вычислительно сложные алгоритмы. Кроме того, современные большие интегральные схемы (БИС) могут содержать несколько ЦСП в одном кристалле [6], что позволяет создавать компактные устройства ЦОС, функционирующие на основе использования параллельных вычислений. Алгоритмы адаптивной фильтрации также могут строиться с использованием параллельных вычислительных процедур. Некоторые из таких алгоритмов рассмотрены в работах [7-14]. Эти алгоритмы имеют достаточно сложное математическое описание.
В настоящей статье рассмотрен простой способ описания RLS-и линейно-ограниченных (Linearly Constrained, LC) RLS-алгоритмов адаптивной фильтрации, позволяющий реализовывать такие алгоритмы с помощью параллельных вычислений, для многоканальных фильтров с неодинаковым числом комплексных весовых коэффициентов в каналах. Возможность использования неодинакового числа весовых коэффициентов в каналах позволяет уменьшать вычислительную нагрузку на ЦСП, если в соответствии с природой решаемой задачи одинаковое число весовых коэффициентов в каналах адаптивного фильтра не требуется. Оценка корреляционной матрицы входных сигналов рассматриваемого фильтра осуществляется на скользящем окне (Sliding Window, SW) с экспоненциальным взвешиванием. Кроме того, при обращении этой матрицы (математической операции, присутствующей в явном или неявном виде во всех RLS-алгоритмах) используется динамическая регуляризация [15] в целях обеспечения устойчивости алгоритма. Многоканальные адаптивные алгоритмы без линейных ограничений, нерегуляризированные алгоритмы со скользящим окном, регуляризированные и нерегуляризированные алгоритмы с бесконечным окном (Prewindowed, PW), одноканальные версии таких алгоритмов или алгоритмы с действительными весовыми коэффициентами — все это частные случаи данного общего решения.
Получение вычислительных процедур большинства RLS-алгорит-мов основано на использовании леммы об обращении матриц [1]. С ее помощью осуществляется рекурсивное обращение корреляционной матрицы адаптивного фильтра. Обеспечение следящих свойств таких фильтров при обработке нестационарных сигналов часто осуществляется с помощью экспоненциального взвешивания сигналов и (или) применения скользящего окна. В силу ограниченного числа отсчетов, используемых при оценке корреляционной матрицы адаптивного фильтра на скользящем окне, SW RLS-алгоритмы в ряде случаев мо-
гут становиться нестабильными. Как уже упоминалось, эффективным способом обеспечения стабильной работы адаптивных фильтров является динамическая регуляризация корреляционной матрицы. Регуляризация может применяться во всех RLS-алгоритмах: одноканальных и многоканальных, с бесконечным или скользящим окном, включая быстрые (вычислительно эффективные) RLS-алгоритмы. Скользящее окно увеличивает вычислительную сложность RLS-алгоритмов примерно в 2 раза по сравнению с алгоритмами с бесконечным окном. Динамическая регуляризация также увеличивает примерно в 2 раза вычислительную сложность алгоритмов по сравнению с алгоритмами без регуляризации.
В случае использования нескольких ЦСП, вычислительную нагрузку при реализации сложных алгоритмов адаптивной фильтрации можно уменьшить с помощью параллельных вычислений. Наиболее привлекательным является распараллеливание вычислений, связанных с существующими потоками независимых данных. При этом отсутствуют дополнительные вычислительные затраты для разделения данных на параллельные потоки. В SW RLS- и в регуляризированных PW RLS-алгоритмах имеется два (F = 2), а в регуляризированных SW RLS-алгоритмах — четыре (F = 4) независимых потоков входных данных, обусловленных модификацией корреляционной матрицы адаптивного фильтра. В известных последовательных RLS-алгоритмах в результате последовательного применения леммы об обращении матрицы эти потоки обрабатываются последовательно. Рассматриваемые далее параллельные RLS-алгоритмы базируются на возможности независимой обработки таких потоков.
Основные приемы получения параллельных RLS-алгоритмов. Для получения RLS-алгоритмов, допускающих параллельные вычисления, могут быть использованы приемы, приведенные в работе [7]. В этой работе был рассмотрен способ построения одноканальных SW RLS- и быстрых SW RLS-алгоритмов адаптивной фильтрации с действительными весовыми коэффициентами, которые могут быть реализованы с помощью двух параллельных процессоров. Подобные многоканальные SW RLS- и LC SW RLS-алгоритмы с комплексными весовыми коэффициентами были получены в работе [9], а многоканальные регуляризированные PW RLS-алгоритмы — в работе [10]. Многоканальные регуляризированные SW RLS-алгоритмы для параллельной реализации с помощью четырех процессоров были рассмотрены в работах [8, 11-14]. Все эти алгоритмы имеют достаточно сложное математическое описание, являющееся следствием использования сложного выражения леммы для обращения матрицы, модифицируемой на каждой итерации алгоритма адаптивной фильтрации за счет двух или
четырех потоков обрабатываемых данных. Такая лемма получается путем последовательного двукратного или четырехкратного применения леммы об обращении матрицы [1], модифицируемой за счет одного потока обрабатываемых данных.
Далее рассматривается простой способ описания вычислительных процедур параллельных ЯЬБ-алгоритмов адаптивной фильтрации. При таком описании вычислительные процедуры регуляризирован-ных алгоритмов с бесконечным окном, нерегуляризированных или регуляризированных алгоритмов со скользящим окном имеют один и тот же вид, близкий к виду нерегуляризированных КЬ8-алгоритмов с бесконечным окном. Особенностью представляемого описания параллельных ЯЬБ-алгоритмов является замена ряда скалярных и векторных переменных, свойственных последовательным алгоритмам, на векторные и матричные переменные, размерность которых пропорциональна числу ^ независимых потоков входных данных адаптивного фильтра. Вычисления элементов векторов и столбцов таких матриц не зависят друг от друга, так как определяются только одним из раннее отмеченных независимых потоков входных данных адаптивного фильтра, а значит, могут быть выполнены с помощью ^ параллельных процессоров.
Блок-схема М-канального адаптивного фильтра, постановка задачи адаптивной фильтрации по критерию наименьших квадратов, а также основные обозначения и определения переменных, используемых в настоящей статье, могут быть найдены в работе [8]. Такой фильтр может иметь неодинаковое число комплексных весовых коэффициентов в каждом из каналов. Наиболее общий случай адаптивной фильтрации с помощью такого фильтра — это ЬС-фильтрация, при которой накладываются дополнительные ограничения вида С^Ь^ (к) = fJ, где С^ и fJ — матрица и вектор линейных ограничений [16-18], И-(к) — вектор весовых коэффициентов адаптивного фильтра.
Решением данной задачи является вектор весовых коэффициентов адаптивного фильтра [18]
Ь- (к) = К^(к)г- (к)+
+ К-^С^р^К-^С-^-1^ - С^К^(к)г-(к)], (1)
где К- (к) и г- (к) — корреляционная матрица обрабатываемых сигналов и вектор взаимной корреляции между компонентами вектора
входных сигналов адаптивного фильтра х- (к) и требуемым сигналом
м
^(к); к — номер итерации алгоритма; N = ^^ Nm — общее число
т=1
весовых коэффициентов адаптивного фильтра; Nm — число весовых коэффициентов в m-м канале. Здесь и далее один нижний индекс N, J или F используется для обозначения размерности векторов с числом элементов N, J, F и квадратных матриц с числом элементов N х N, J х J и F х F, а два индекса NJ и NF — для обозначения размерности прямоугольных (нетранспонированных) матриц с числом элементов N х J и N х F.
Первое слагаемое в уравнении (1) является решением задачи RLS адаптивной фильтрации при отсутствии ограничений. Второе слагаемое обусловлено линейными ограничениями. В соответствии с этим, RLS-алгоритмы, с помощью которых осуществляется рекурсивное вычисление вектора весовых коэффициентов адаптивного фильтра (1), также состоят из двух вычислительных процедур. При этом, в зависимости от числа слагаемых в уравнениях оценки матрицы RN (k) и вектора rN (k) [8], а также интервалов накопления обрабатываемых данных, участвующих в формировании этих оценок, уравнение (1) может быть использовано для получения PW RLS-алгоритмов, регу-ляризированных PW RLS-алгоритмов, SW RLS-алгоритмов или регу-ляризированных SW RLS-алгоритмов. Последовательное применение леммы [1] для обращения матрицы RN (k) позволяет получить соответствующие последовательные RLS-алгоритмы [5], а использование приемов [7] — параллельные алгоритмы, ориентированные на реализацию с помощью двух или четырех ЦСП [7-14].
Лемма об обращении матриц [1] имеет вид
R-1 = B-1 - B-1CA-1DB-1,
где A = DB-1C + 1 и R = B + CD, а C и D — векторы (это позволяет использовать данную лемму для рекурсивного обращения матрицы Rn(k)).
Общий вид леммы
R-1 = B-1 - B-1CA-1DB-1, (2)
где A = DB-1C + S, а C и D — матрицы, приведен в работе [19]. Уравнение (2) лежит в основе рассматриваемых далее RLS-алгоритмов. Для использования этого уравнения необходимо сформировать матрицы C = [y, x, z, v] = Xnf(k) = [m0'5Xn(k - L), Xn(k),/'5?PN(k - L), £pN(k)]; D = CH и S = diag(-1, 1, -1, 1). Число столбцов F матрицы XNF и их вид определяют возможные разновидности RLS-алгоритмов.
Уравнение (2) позволяет получить параллельный RLS-алгоритм адаптивной фильтрации.
Параллельный RLS-алгоритм
Шаг 0:
Инициализация : xN (0) = XN (0 - L + 1) = 0n,
pN (0) = 0n,•••, Pn (0 - L + 1) = 0n ,
d(0) = 0,...,d(0 - L + 1) = 0, Xnf(0) = Onf,
RN1 (0) = i-2ÄN, hN (0) = 0n,
An = diag(1,A...,ANl-1,...,1,A,...,ANM-1)
Шаг 1
Шаг 2 Шаг 3 Шаг 4
R-1 (k - 1)Xnf(k)
For k = 1, 2,...,K
Gnf(k) = AsF + XNf(k)RN(k - 1)Xnf(k) RN1 (k) = A-1[R-1(k - 1) - Gnf(k)NF(k)RN1(k - 1)] aF (k) = dF(k) - hH(k - 1)Xnf(k) hN(k) = hN(k - 1) + Gnf(k)aH(k) End for k
Данный алгоритм отличается от PW ЯЬ8-алгоритма тем, что в параллельном алгоритме Х-Р(к), С-Р(к) — это матрицы, а не векторы; ^(к) и ар(к) являются векторами, а не скалярными величинами. Кроме того, знаменатель в уравнении (шаг 2) является не скалярной величиной, а квадратной матрицей с числом элементов ^ х ^. Эта матрица обеспечивает математическую эквивалентность данного параллельного алгоритма и соответствующего последовательного регу-ляризированного SW ЯЬБ-алгоритма.
Эквивалентность означает одинаковое функционирование на всех итерациях алгоритмов при одинаковых входных сигналах и одинаковых параметрах адаптивного фильтра. Здесь р-(к) — вектор сигналов динамической регуляризации; — параметр динамической регуляризации и параметр £2 — начальной регуляризации корреляционной матрицы [8], а Ь — длина скользящего окна, т.е. число выборок обрабатываемых сигналов, участвующих в оценке матрицы К- (к) и вектора г-(к).
В представленном алгоритме вектор ^(к) определяется как ^(к) = ¿(к — Ь),^(к), 0, 0], а сигнал ошибки на выходе адаптивного фильтра — как а-1Х(к) = ¿(к) — Ь^(к — 1)х-(к) = ар^к), где ар^к) — второй элемент вектора ар (к). Рассматриваемые в работе векторы, за исключением специально упоминаемых, являются вектор-столбцами. Векторы ^ (к) и ар(к) — это векторы-строки. Матрица С-Р(к) представляет собой матрицу коэффициентов Калмана С-р (к) = К-1(к)Х-р (к).
Вычислительная сложность такого алгоритма пропорциональна О^2^) арифметических операций, необходимых для выполнения одной итерации алгоритма.
В силу небольшого числа элементов матрицы в знаменателе уравнения (шаг 2) и аналогичных матриц в других алгоритмах, рассматриваемых далее, вычислительная сложность обращения этих матриц не сказывается на общей вычислительной сложности RLS-алгоритмов при N ^ 4. Например, при F = 4 использование стандартной процедуры обращения матриц с учетом структуры матрицы позволяет выполнять это обращение с вычислительной сложностью: 88 умножений, 47 сложений и одно деление за одну итерацию алгоритма адаптивной фильтрации. Использование модифицированного приема [20] позволяет выполнять обращение со следующей вычислительной сложностью: 80 умножений, 64 сложения и 4 деления, а использование леммы об обращении клеточных матриц [19, 21], позволяет выполнить такое обращение, используя 54 умножения, 32 сложения и 4 деления за одну итерацию алгоритма адаптивной фильтрации.
Параллельные быстрые RLS-алгоритмы. Рассмотренный способ может быть использован и для получения параллельных быстрых многоканальных RLS-алгоритмов адаптивной фильтрации с вычислительной сложностью O(NF) арифметических операций. Вычислительные процедуры таких алгоритмов получаются путем использования теории линейного предсказания [22] и представления потоков обрабатываемых сигналов в виде матрицы. Полученный таким образом так называемый быстрый алгоритм Калмана (Fast Kalman, FK) приведен ниже.
Параллельный FK-алгоритм
Шаг 0
Инициализация: xn(0) = On,..., XN(0 - L + 1) = On, pN(0 - L + 1) = On,d(0) = 0,...,d(0 - L + 1) = 0, Xnf(0) = Onf, hN (0) = On,E#m)(0) = ¿2,
/ (m).
hN
GN?(1) = 0
(0) = On, hNm)(0) = On ,m =1: M,
NF
Шаг 1 Шаг 2 Шаг 3 Шаг 4 Шаг 5
Шаг 6
For k = 1, 2,...,K
For m = M,M - 1,
h
/ (m)
F
b(m), F / (m) I
N f(m)
F
(k) = xFm)(k) - hN(m)H (k - 1)xNF (k)
(k) = xFm)(k - Nm) - h^(k - 1)XNTij(k)
b(m)H /
-(m-1),
f(m) lN
— rv(m)/
N
(m) I
(k) = h^;(k - 1) + GNF(k)«F(m)H(k)
(k) = [xr (k) - hN(m)H (k)xNF (k)]SF
-(m)/
E
/(m)m _
N
(k) = AEN(m)(k - 1) + e£m;(k)a
f (m)/
f (m)H
G
(m)
(N+1)
(k) =
h
/ (m)
N
(k)
,(/(m)
F
(k)
F
(k)
E
/ (m)
N
(k)
+
OF
G(m)
gnf
(k)
1
Шаг 7 : Шаг 8 :
Шаг 9 :
Шаг 10 Шаг 11 Шаг 12
S(m) T(m)T/G(m) SN+1 TN+1 {G(N +1)F
Л =
G
(m-1)m -
NF
(k) =
Q NF(k)
Q NF(k)
Lq Fm)(k)
h
b(m)(k - 1)q Fm) (k)
x
If - aF(m)H (k)q/m)(k)
h
b(m)
N
X
(k) = hNm)(k - 1) + GNmF-1)(k)aF(m)H(k)
End for m
«f(k) = dF(k) - hH (k - 1)Xnf(k) hN(k) = hN(k - 1) + GNF(k)aH(k)
IN (
.(M),
GNF(k +1) = GNF(k)
End for k
В данном алгоритме IF — единичная матрица, а векторы xF")(k) = = [/0°'5Xm(k - L),Xm(k),/^Pm(k - L),^m(k)], xFm) (k - Nm) = = [^°'5Xm(k-Nm -L),Xm(k-Nm) ,^°'5^Pm(k-Nm -L)),£pm(k-Nm)],
а
/(m),
(k), eF(m)(k), (k) и eF(k) являются вектор-строками. Матри-
(m)
F
цы xNm)(k) = [/-5xNm)(k - L), xNm)(k),^°'5?pNm)(k - L),£p<m)'
формируются из векторов:
xi°)(k) = Xn (k)>
x!1)(k) =
x
(k - 1),xN2(k),...,
XNm (k), •••, X
Nm (k)
T
xNn)(k)= [xNx(k - 1), xN2(k - 1),...,
T
x
T
Nm (k - 1), xNm+i (k)xNm (k)
xNM)(k) = [xNi(k - 1),xN2(k - 1),..
xNm (k - 1),...,xNm (k - 1)' T
xN°)(k - L) = Xn (k - L), x!1)(k - L) =
x
T
xNm (k - L)'.'
Ni(k - L - 1), xN2(k - L),
т
T
' xNm
(k - L)
x5m}(k - L) = [xN(k - L - 1), xN2(k - L - 1),...,
x
(k - L - 1), xNm+1 (k - L),...,xNm (k - L)
т
XM}(k - L) = xNi(k - L - 1), x^2(k - L - 1),
x
(k - L - 1),
т
,xNm
(k - L - 1)
т
Pw(k) = Pn (k)
Pw(k) = [pNi (k - 1)> PN2 (k)>--->Pwm (k)>---> pnm (k)]T,
pT (k) = [pNi (k - 1), PN2 (k - 1),..., pNro (k - 1) PNm+1 (k),---,pNM (k)]T,
pNT} (k) = [pNi (k - 1), p^2 (k - 1),...,
pNm (k - 1),...,pNm (k - 1)]T;
p50}(k - L) = р^(k - L),
PN (k - L) = [pNi (k - L - 1), pN2 (k - L)
T
pNm (k - L),...,pNm (k - L)]
T
T
pft0(k - L) = [pNi(k - L - 1), p^2(k - L - 1),..., pNm(k - L - 1), pNm+i(k - L),..., pNm(k - L)]T,
pN") (k - L) = [pNi (k - L - 1), p^2 (k - L - 1),...,
pNm(k - L - 1),...,pNm(k - L - 1)]
T
T
В быстром КЬ8-алгоритме матрицы В^! и являются
перестановочными. Умножение векторов на эти матрицы не требует арифметических операций. Некоторые детали аналогичного способа
-i(m)
m
m
построения быстрых многоканальных RLS-алгоритмов адаптивной фильтрации с неодинаковым числом весовых коэффициентов в каналах адаптивного фильтра с использованием перестановочных матриц рассмотрены в работе [23]. Перестановочные матрицы для случая M = 3 приведены в работе [8].
Вычислительные процедуры параллельных версий других быстрых RLS-алгоритмов приведены ниже: Fast Transversal Filter (FTF), Fast a Posteriori Error Sequential Technique (FAEST) и стабилизированного FAEST-алгоритма [5, 22, 24].
Параллельный FTF-алгоритм
Шаг 0
Инициализация: xN(0) = On,•••, XN (0 — L + 1) = On,
Pn (0) = OnPn (0 — L + 1) = On,
d(o) = 0,...,d(0 — L + 1) = 0, Xnf (0) = Onf ,
hN (0) = On X(m)(0) = ¿2,ENm)(0) = ¿2A-Nm, h
f (m)/
N (0) = 0n, hNm) (0)= ON,m =1: M,
-\(m)
t nf (i) = onf , ) (i) = s
F
Шаг 1 : Шаг 2 :
Шаг 3 : Шаг 4 :
Шаг 5 :
Шаг 6 Шаг 7 Шаг 8 Шаг 9
Шаг 10
Шаг 11 Шаг 12 Шаг 13 Шаг 14
For k = 1, 2,...,K
For m = M, M - 1,
a
f (m)
(k) = xFm)(k) - h^ (k - 1)XNF (k)
1
f (m)H,
F
f(m) F
(m)
-(m) .
(k) = a
N (m)
F(m)(k^r (k)
f(m)
nFm)(k) =
a
F
(k)
EN(m )(k -1)
(m)
-(N+1)F
(k) =
-hN(m)(k - 1)
G (m) T(m)Tr T(m) (k)\ = SN+1TN+1 {T(N+1)F (k)} =
nFm)(k)
(m)
O
(m)
_T NF (k).
Q NF (k)
q Fm)(k)
h
f(m)
n (k) = hN(m)(k - 1) + A
1
a
b(m)/n — ^(m)
F
(k) = q Г (k)E$Jm)(k - 1)
(m) NF
(k)eF(m)H (k)
ф (m-1b ,4 ~ (mb ,4
TNF (k) = QNF(k)
f(m)
'N (k) = AEN
hNNm)(k - 1)qГ (k)
(m)
EN(m)(k) = AEN(m)(k - 1) + e£m;(k)a
f (m)H
ф Fm)(k) = $Fm)
(k)
a
f (m)H
f (k)
F
F
(k)eF(m)(k)
фFm-1)(k) = [If - А-1ФFm)(k)qFm)H(k)a
f (m)
N
-(m)H,
EN"(k)
b(m),
F
nb(m),
(k) = aF(m)(k^Fm_i ;(k)
(m-1)
b(m)
F
(,)]-1Ф Fm)(k)
N
b(m),
nb(m) N b(m)
N
E^^(k) = AE^j(k - 1) + eF(m)(k)aFmj"(k)
b(m)H,
hrj(k) = hrj(k - 1) + A-1TN"F-1)(k)e^" (k) End for k
Параллельный ПГ-алгоритм
Шаг 15 Шаг 16 Шаг 17 Шаг 18
«f(k) = dF(k) - hH (k - 1)Xnf(k) eF (k) = «f (к)Ф^0)(к) hw(k) = hw(k - 1) + A-1TNF(k)eH(k)
(0)
-NF 1
End for k
T NF (k +1) = T NF (k), Ф^ +1) = Ф^ (k)
(0)
Параллельный FAEST-алгоритм отличается от параллельного FTF-
алгоритма тем
что в первом рекурсивно вычисляются матрицы
[Ф^(к)]-1, а во втором — матрицы Ф^")(к). Аналогичные матрицы в последовательных быстрых И^-алгоритмах отсутствуют. Они становятся скалярными величинами, известными в теории адаптивной фильтрации как отношения правдоподобия.
Параллельный FAEST-алгоритм
Шаг 0 :
Инициализация: xN (0) = 0n,•••, XN (0 - L + 1) = 0
Pn (0) = 0n , •••, Pn (0 - L + 1) = 0n, d(o) = 0,---,d(0 - L + 1) = 0, Xnf(0) = Onf,
N,
hN (0) = 0n ,eNN(m) (0) = ¿2,er(0) = ¿2 a
hN(m)(0) = 0n, hNNm)(0) = 0n,m =1: m,
~ (m)
Nm
(1) = 0nf, ФFM)(1) = S
NF
F
Шаг 1 : Шаг 2 :
Шаг 3 : Шаг 4 :
Шаг 5 :
Шаг 6 Шаг 7 Шаг 8 Шаг 9 Шаг 10 Шаг 11 Шаг 12
For k = 1,
For m = M,M - 1,^,1
aF(m)(k) = xFm)(k) - hN(m)H(k - l)xNF(k)
J (m)
F
nm)(k) =
(k) = а
J(m) F
а
J(m) F
(k^F
(m)
(k)
(k)
EN(m)(k -1)
T(m+1)F(k) —
-hN(m)(k - l)
S(m) T(m)TrtT(m) (kU _ SN+1TN+1 lT(N+1)F(k)} =
nFm)(k)
(m)
0
T F (m)
NF
(k)
J(m)
Q NF (k)
Lq Fm)(k)
N
(k) = hN(m)(k - l) + A-1TNF(k)e
h
a?m)(k) = q (m^Ef^ - l)
f (m)H
F
(k)
(m-1)
(m)
-NF (k) = Q NF
[Ф Fm)(k)]-1 = [ФГ (k)]-1 + A
(k) + h
6(m)(k - l)qFm)(k)
N
(m)
F
f (m)H
а
(k)nFm)(k)
EN(m)(k) = AEN(m)(k - l) + eF(m)(k)aF(m)H(k) фFm-1)(k)]-1 = [ФFm)(k)]-1 - A-1aF(m)H(k)qFm)(k)
b(m) I
F
(k) = a^k^"1 '(k)
(m-1),
Шаг 13 : Шаг 14 :
Шаг 15 Шаг 16 Шаг 17 Шаг 18
E5m)(k) = AEJm)(k - 1) + eF(m)(k)a
h
JN b(m),
N
(k) = hNNm)(k - 1) + A-1TNF *y(k)eFm^ (k)
-1n(m-1),
b(m)H
F
(k)
b(m)H,
End for k
«F (k) = dF (k) - hH(k - 1)Xnf (k)
eF (k) = «f (k)фF°) (k)
hN(k) = hN(k - 1) + A-1TN"F-1)(k)eH(k)
IN (
ц(М),
(°)
TT NF (k + 1) = TT NF (k), ФFM )(k + 1) = Ф^ (k)
End for k
(°)
Параллельная версия стабилизированного FAEST-алгоритма [24]
Инициализация: xN (0) = 0n,•••, XN (0 - L + 1) = 0
Pn (0) = 0n ,•••, Pn (0 - L + 1) = 0n, d(0) = 0,...,d(0 - L +1) = 0, Xnf (0) = Onf, hN (0) = 0n,
N;
E
/(m)
N
(0) = ¿2 ,E^m)(0) = ¿2A
/(m) lN
,(M)
l —N„
hir (0) = 0n ,
N
hN(m)(0) = 0n,m = 1 : M,
T NF (1) = 0nf , ФFM )(1) = S
F
For k = 1, 2,...,K For m = M,M - 1,...,1
f(m)(,,) _(m)^) uf(m)H(k - 1)xNF(k)
«F '(k) = xF )(k) - hN
/(m)
F
(k) = а
/(m) F
nFm)(k) =
а
/(m)
(k^F
(m)
(k)
F
(k)
EN(m)(k - 1)
(m)
-(N + 1)F
(k) =
-hN(m)(k - 1)
S(m) T(m)TrTn(m) (k)i = SN+1TN+1 {T(N+1)F (k)} —
h
/(m)
N
(k) = hN(m) (k - 1) + A
nFm)(k)
(m)
0
T F (m)
-NF
(k)
Q NF (k)
Fm)(k)
-1 т (mv ^
NF (k)eF
T (k)eF (k)
EN(m)(k) = AEN(m)(k - 1) + eF(m)(k)«F(m)H(k)
1 /(m)H
VN W — VV ^ CF
[Ф Fm)(k)]-1 = [фFm)(k)]-1 + A-«F
„k(m)/7 ч _ „(m)/i дГ \ ib(m)H/
(k)nF
(m)
F
«Г(k) = хГ(k - Nm) - hNmjH(k - 1)XNTij(k)
(k) (m-1)
q(m)(k) =
а
b(m) F
(k)
Ej(m)(k - 1)
а
F(m)(k) = q Fm)(k)E55m)(k -1)
Шаг 12 : aF(1)(m)(k) = K!aF(m)(k) + (1 - Ki)aF(m)(k)
Шаг 13 aF(2)(m)(k) = K2aF(m)(k) + (1 - K2)äF(m)(k)
Шаг 14 : aF(5)(m)(k) = K5«F(m)(k) + (1 - K>)ÄF(m)(k)
Шаг 15 : tFm)(k) = K4qFm)(k) + (1 - K)qFm)(k)
Шаг 16 : TNT^k) = QNF(k) + hNm)(k - 1)tFm)(k)
Шаг 17 : [Ф Fm-1)(k)]-1 = A-iTNm-1)H (k)xNm-i)(k)+Sf
Шаг 18 : [ФFm-1)(k)]-1 = [ФFm)(k)]-1 - A-1 aF(5)(m)H(k)qFm)(k)
Шаг 19 : ^Fm-1)(k)]-1 = кз[Ф Fm-1)(k)]-1 + (1 - кз)[Ф Fm-1)(k)]-1
Шаг 20 : eF(1)(m)(k)= aF(1)(m)(k)^& Fm-1) (k)
Шаг 21 : eF(2)(m)(k)= aF(2)(m)(k)<& Fm-1) (k)
Шаг 22 : EN(m)(k) = AE^k - 1) + eF(2)(m)(k)aF(2)(m)H(k)
Шаг 23 : hN(m)(k) = hNm)(k - 1) + A-1TNmF-1)(k)eF(1)(m)H(k) End for m
Шаг 24 : aF (k) = dF (k) - hH(k - 1)Xnf (k)
Шаг 25 : eF (k) = «f (k)фF0)(k)
Шаг 26 : hN(k) = hN(k - 1) + A-1TNF(k)eH(k)
Шаг 27 : T NM;(k +1) = T NF (k), ФFM )(k +1) = ФF0)(k) End for k
Рассмотренные вычислительные процедуры могут быть использованы в качестве алгоритмов безусловной адаптивной фильтрации (без ограничений) или в качестве процедур вычисления матрицы коэффициентов Калмана С^(к) = С^(к) = А-1Т^р(к)Ф^(к) в рассматриваемых далее линейно-ограниченных (ЬС) ЯЬБ-алгоритмах.
Параллельные LC RLS-алгоритмы. В параллельных ЬС ЯЬБ-алгоритмах используется лемма (2) не только для вычисления матрицы И^г1(к), но и для вычисления матриц Г^(к) = И^^С^, Ф-1 (к) = [С^Г^(к)]-1 и д^(к) = Г^7(к)Ф-1(к), обусловленных линейными ограничениями. Эти вычисления базируются на приемах, аналогичных рассмотренным в работе [18]. Таким образом, можно последовательно получить три разновидности параллельных ЬС ЯЬБ-алгоритмов, условно обозначаемых как алгоритмы № 1, № 2 и № 3.
Шаг 0 :
Инициализация: х- (0) = 0-,..., х-(0 — Ь + 1) = 0
р- (0) = 0-,...,р- (0 — Ь +1) = 0-,
ад = 0,...д0 — ь + 1) = 0, Х-р (0) = о-р, К-1(0) = ¿"2Л-,
ГNJ(0) = К^(0)С^, Ф-1(0) = [cNJГNJ(0)]-1,
Ь- (0) = Г^(0)Ф-1(0^ J, Л- = ^(1,А, ...,А-1-1 ,...,1,Л,...,Л-М-1)
>N,
Шаг 1 : Шаг 2 : Шаг 3 :
Шаг 4 : Шаг 5 : Шаг 6 : Шаг 7 : Шаг 8 :
For k = 1, 2,...,K Вычисление GNF(k)
Vjf (k) = CNjGnf (k)
NHf(k) = XNf(k)rNj(k - 1)
H
Ф-1^ - 1)Vjf(k)
О (k) = — H
A } If - NHf(k^-1(k - 1)Vjf(k) Ф-1^) = A^-^k - 1) + О/(k)NHF(^Ф-1^ - 1)]
rNj(k) = A-1[rN/(k - 1) - Gnf(k)NHF(k)]
«f(k) = dF(k) - hH (k - 1)Xnf (k)
hN(k) = hN(k - 1) + [Gnf(k) - AGn/xW^j(k)]«H(k) End for k
Модификацией алгоритма № 1 является алгоритм № 2. Отличие между этими алгоритмами заключается в выполнении шага 8, который для алгоритма № 2 определяется как Ь- (к) = Р-(к)[Ь-(к — 1) + + С-р(к)аН(к)] + а-(к), где Р-(к) = I- — Г-Д^Ф-^С-, и а- (к) = Г^(к)Ф-1(к) (см. шаги 8 и 9 алгоритма № 2). Здесь Р-(к) — динамически изменяемая проекционная матрица, благодаря которой условия ограничений С^Ь- (к) = Г ^ удовлетворяются на каждой итерации алгоритма № 2 в отличие от алгоритма № 1, где такие ограничения задаются лишь при инициализации, а в процессе функционирования могут перестать удовлетворяться в силу накопления ошибок округления. Из-за наличия матрицы Р-(к) вычислительная сложность алгоритма № 2 пропорциональна О^2), что не зависит от способа вычисления матрицы С-р (к) (быстрого или не быстрого).
Шаг 0 :
Инициализация: х- (0) = 0-,..., х- (0 — Ь + 1) = 0
р- (0) = 0-,...,р- (0 — Ь + 1) = 0-,
ад = 0,...,^(0 — Ь + 1) = 0,
Х-р (0) = о-р, И-1(0) = ¿-2Л-,
Г^(0) = и-1 (0)CNJ, Ф-1(0) = [CNJГNJ(0)]
>N,
1
hN (0) = Г^(0)Ф-1(0)/,
An = diag(1,A, ...,ANl-1,..., 1,A, ...,Anm-1)
Шаг 1 Шаг 2 Шаг 3
Шаг 4 :
Шаг 5 Шаг 6 Шаг 7 Шаг 8 Шаг 9 Шаг 10
For k = 1, 2,...,K Вычисление GNF (k)
Vjf (k) = CNjGnf (k)
NHf(k) = XNf(k)rNj(k - 1)
О/(k) =_^ - 1)VJF(k)_
JV } If - NHf(k^-1(k - 1)Vjf(k) Ф-1^) = A^k - 1) + О/(k)NHF(^(k - 1)] rNj(k) = A-1[rN/(k - 1) - Gnf(k)N^(k)] аf(k) = dF(k) - hH (k - 1)Xnf(k)
H
''NJ
Pn(k) = In - Г^^Ф-^С
aN (k) = rNj(k^-1(k)f j
hN(k) = Pn(k)[hN(k - 1) + Gnf(k)«H(k)] + aN(k) End for k
Упрощение ЬС RLS-aлгоритмa № 2 достигается в алгоритме № 3, в котором, вместо рекурсивного вычисления матриц Г^(к) и Ф-1(к), вычисляется матрица QNJ(k).
Параллельный LC RLS-алгоритм № 3
Шаг 0 : Инициализация: х-(0) = 0-,..., х-(0 — Ь + 1) = 0-,
р- (0) = 0-,...,р- (0 — Ь +1) = 0-, ^(0) = 0,...,^(0 — Ь +1) = 0,
Х-р (0) = о-р, и-1(0) = Г2Л-,
ГNJ(0) = Я-^С-/, QNJ(0) = ГNJ(0)[CNJГNJ(0)]-1,
Ь-(0) = QNJ(0)fJ,
_ Л- = Аа^1,А, ...,А-1 -1,..., 1,А, ...,А-М-1)
For к = 1,2,...,К Шаг 1 : Вычисление С-р(к) Шаг 2: VJр(к) = С^С-р(к) Шаг 3: (к) = Х-р(k)QNJ(k — 1)
Шаг 4 :
Шаг 5 Шаг 6 Шаг 7 Шаг 8
QNj(k) = [QWj(k - 1) - GWF(k)NfF(k)x
x
VjF (k)N* (k)
IJ +
H I-i- NHf(H)Vjf(k).
QNJ(k) = QNJ(k) + CNJ(CNJCNJ)- [IJ - CJVjQWj(k)
«F(k) = dF(k) - hH (k - 1)Xnf(k) h'w(k) = hw(k - 1) + Gwf(k)aH(k)
hw(k) = h'w(k) + QWj(k)[fj - CNjhN(k)] End for k
Рассмотренные LC RLS-алгоритмы также являются параллельными, поскольку в таких алгоритмах вычисление всех матриц с числом столбцов, равным F, может быть выполнено с помощью F параллельных процессоров, так как эти вычисления не зависят друг от друга и определяются независимыми элементами матриц XNF(k).
Моделирование. Поскольку все рассмотренные параллельные алгоритмы являются математически эквивалентными одноименным последовательным алгоритмам, то физическая зависимость между независимыми параллельными вычислениями на самом деле существует. Она обеспечивается с помощью квадратных матриц с размером F x F, фигурирующих в выражениях параллельного RLS-алгоритма (шаг 1), параллельного FK-алгоритма (шаг 8); параллельных FTF-алгоритма и FAEST-алгоритма (шаг 11); параллельного FAEST-алгоритма и шагов 4 параллельных LC RLS-алгоритмов № 1, № 2 и № 3 (шаг 18).
В силу эквивалентности параллельных и последовательных алгоритмов, результаты моделирования, подтверждающие работоспособность рассмотренных в статье алгоритмов, полностью совпадают с результатами, представленными в работе [5]. Эти результаты также приводятся далее. На рис. 1, а и б сравниваются многоканальный вариант LC PW RLS-алгоритма и многоканальный нерегуляризирован-ный LC SW RLS-алгоритм при решении задачи идентификации трех-канального линейного фильтра и речевых (нестационарных) входных сигналов.
Из рис. 1, а видно, что ограничение 0дБ амплитудно-частотной характеристики (АЧХ) адаптивного фильтра, задаваемое на частотах 1 и 2 кГц, обеспечивается обоими алгоритмами. Вертикальные стрелки указывают на частоты ограничений, а горизонтальная штриховая линия обозначает уровень ограничения АЧХ. Однако в алгоритмах со скользящим окном (благодаря следящим свойствам этих алгоритмов) такой параметр, как Echo Return Loss Enhancement (ERLE), достигает более высокого значения, чем в алгоритмах с бесконечным окном (см.
1Д/)1,ДБ
ERLE, дБ
10
0
-10
-20
fV
-30 -40
О 1000 2000 3000 /, Гц О 1 2 3 4 Ь<10
Рис. 1. Результаты моделирования (АЧХ — а и ERLE — б) линейно-ограниченного RLS-алгоритма без регуляризации с бесконечным (1) и со скользящим (2) окном
рис. 1). Параметр ERLE является одним из показателей качества функционирования алгоритмов адаптивной фильтрации в задачах идентификации и представляет собой отношение энергии эхо-сигнала d(k) к энергии подавленного эхо-сигнала aN>x(k), измеряемых на каждой итерации алгоритма в пределах некоторого скользящего окна B как
Из рис. 1, б видно, что в рассмотренных примерах значения параметра ERLE в случае RLS-алгоритма с бесконечным окном примерно на 20 дБ ниже, чем в случае RLS-алгоритма со скользящим окном. Это объясняется тем, что при обработке нестационарных сигналов во втором случае корреляционная матрица сигналов вычисляется более точно. В рассматриваемых примерах длина скользящего окна L в адаптивном алгоритме равнялась выборкам, что при частоте дискретизации сигналов, равной 8 кГц, примерно соответствует 30 мс (интервалу стационарности речевых сигналов). Значения ERLE также рассчитывались на скользящем окне B = 256.
Улучшение разработанных алгоритмов достигается путем использования динамической регуляризации при обращении корреляционной матрицы (рис. 2, а и б). Из рис. 2, а видно, что заданные ограничения АЧХ также выполняются и в случае регуляризированного алгоритма, а из рис. 2, б следует, что при использовании регуляризации достигаемое значение ERLE в среднем не меньше, чем при отсутствии регуляризации.
Заключение. Таким образом, в настоящей работе рассмотрен простой способ описания RLS- и LC RLS-алгоритмов в виде однотипных параллельных вычислительных процедур независимо от возможных вариантов модификации корреляционной матрицы адаптивного
а
б
фильтра. Данные алгоритмы эквивалентны одноименным последовательным алгоритмам и одноименным параллельным алгоритмам. Различие между параллельными и последовательными алгоритмами заключается в вычислительных процедурах, а также в возможности или невозможности реализации этих процедур с помощью параллельных вычислений. Полное число арифметических операций в одноименных последовательных и параллельных алгоритмах примерно одинаково. Однако при использовании процессоров для реализации параллельных алгоритмов вычислительная нагрузка на один
ЦСП уменьшается в F раз. гуляризацией 1 и без регуляризации (2) Рассмотренные алгоритмы и
их частные случаи могут быть использованы во всех традиционных приложениях адаптивной обработки сигналов. Эти алгоритмы были разработаны в качестве элементов библиотеки [25] БИС серии "МУЛЬТИКОР" [6].
СПИСОК ЛИТЕРАТУРЫ
1. S a y e d A. H. Fundamentals of adaptive filtering. - NJ, Hoboken: John Wiley and Sons, Inc., 2003. - 1125 p.
2. Benesty J., Huang Y. (Eds). Adaptive signal processing: applications to real-world problems. - Berlin, Heidelberg, New York, Springer-Verlag, 2003. -356 p.
3. Первачев С. В., Петров А. И. Адаптивная фильтрация сообщений. -М.: Радио и связь, 1991. - 160 с.
4. ПистолькорсА. А, Литвинов О. С. Введение в теорию адаптивных антенн. - М.: Наука, 1991. - 200 с.
5. Джиган В. И. Многоканальные RLS- и быстрые RLS-алгоритмы адаптивной фильтрации // Успехи современной радиоэлектроники. - 2004. - № 11. - C. 4877.
6. Петричкович Я. Я., Солохина Т. В. Цифровые сигнальные контроллеры "МУЛЬТИКОР" - новые отечественные серии систем на кристалле // Докл. 6-й Международной конф. "Цифровая обработка сигналов и ее применения (DSPA-2004)" (Москва, 31 марта-2 апреля 2004 г.). - М., 2004. - Т. 1. -С. 8-15.
1#(/)|,дб 10
о
-10
-20
-30 -40
1 л
w г П
1000
2000
3000 /, Гц
Рис. 2. Результаты моделирования (АЧХ — а и ERLE — б) линейно-ограниченного RLS-алгоритма со скользящим окном с ре-
7. Papaodysseus C. A robust, parallelizable, O(m), a posteriori recursive least squares algorithm for efficient adaptive filtering // IEEE Trans. Signal Processing. -1999. - Vol. 47. - № 9. - P. 2552-2558.
8. Джиган В. И. Многоканальный быстрыйRLS-алгоритм адаптивной фильтрации для параллельной реализации с помощью четырех процессоров // Вестник МГТУ им. Н.Э. Баумана. Серия "Приборостроение". - 2005. - № 1. - С. 83-98.
9. Д ж и г а н В. И. Параллельный регуляризированный быстрый RLS-алгоритм многоканальной адаптивной фильтрации со скользящим окном и линейными ограничениями // Тр. 10-й Международной конф. "Радиолокация, навигация, связь (RLNC-2004)" (Воронеж, 13-15 апреля 2004 г.). - Воронеж, 2004. - Т. 1. -С. 132-142.
10. Д ж и г а н В. И. Параллельные регуляризированные RLS-алгоритмы многоканальной адаптивной фильтрации // Цифровая обработка сигналов. - 2004. -№ 2. - C. 7-13.
11. Д ж и г а н В. И. RLS-алгоритм адаптивной фильтрации для параллельной реализации с помощью четырех процессоров // Цифровая обработка сигналов.
- 2004. - № 3. - C. 2-8.
12. D j i g a n V. I. Multichannel RLS adaptive filtering algorithm for parallel implementation by means of four processors // Proceedings of the 4-th International Scientific and Practical Conference "Internet-Science-Education-2004 (ISE-2004)" (Baku-Vinnytsia-Tyrnovo, September 29-October 12, 2004). - Vinnytsia, 2004. -Vol. 2. - P. 687-691.
13. D j i g a n V. I. Parallelizable multichannel SW fast RLS algorithm for implementation by means of four processors // Proceedings of the Second International Conference on Information Systems and Technology (IST-2004) (November 8-12, 2004, Minsk, Belarus). - Minsk, 2004. - Vol. 2. - P. 100-105.
14. Джиган В. И. Параллельные вычисления в адаптивной обработке сигналов // Материалы 13-й Международной науч.-техн. конф. "Проблемы передачи и обработки информации в сетях и системах телекоммуникаций". - Рязань, 2004.
- C. 40-42.
15. Gay S. L. Dynamically regularized fast RLS with application to echo cancellation // Proceedings of the International Conference on Acoustic Speech and Signal Processing (ICASSP-96), (May 7-9, 1996, Atlanta, USA). - Atlanta, 1996. - P. 957960.
16. F r o s t O. L. An algorithm for linearly constrained adaptive array processing // Proceedings of the IEEE. - 1972. - Vol. 60. - August. - P. 926-935.
17. T a k a o K., F u j i t a M., N i s h i T. An adaptive antenna array under directional constraint // IEEE Trans. Antennas and Propagation. - 1976. - Vol. 24. -№ 5. - P. 662-669.
18. R e s e n d e L. S., Romano J. M. T., B e 11 a n g e r M. G. A fast least-squares algorithm for linearly constrained adaptive filtering // IEEE Trans. Signal Processing.
- 1996. - Vol. 44. -№ 5. - P. 1168-1174.
19. Giordano A. A., Hsu F. M. Least square estimation with application to digital signal processing. - Canada, Toronto: John Wiley and Sons, Inc., 1985. - 412 p.
20. M a z u c h T., K o z a n e k J. New recurrent algorithm for a matrix inversion // Journal of Computational and Applied Mathematics. - 2001. - Vol. 136, № 1-2. -P. 219-226.
21. Л а н ш о ц К. Практические методы прикладного анализа. Пер. с англ. М.З.Кайнера / Под ред. А.М. Лопшица. - М.: Государственное издательство физико-математической литературы, 1961. - 524 с.
22. Z e l n i k e r G., Taylor F. J. Advanced digital signal processing: theory and applications. - New York: Marcel Dekker, Inc., 1994. - 666 p.
23. G l e n t i s G. A., K a l o u p t s i d i s N. Fast adaptive algorithms for multichannel filtering and system identification // IEEE Trans. Signal Processing. -1992. - Vol.40. -№№ 10. - P.2433-2458.
24. S l o c k D. T. M., K a i l a t h T. Numerically stable fast transversal filters for recursive least squares adaptive filtering // IEEE Trans. Signal Processing. - 1991. -Vol.39. -№ 1. - P.92-114.
25. Д ж и г а н В. И. Двести пятьдесят алгоритмов адаптивной фильтрации // Материалы 13-й Международной науч.-техн. конф. "Проблемы передачи и обработки информации в сетях и системах телекоммуникаций". - Рязань, 2004. - с. 42-44.
Статья поступила в редакцию 8.02.2005
Виктор Иванович Джиган родился в 1958 г., окончил в 1984 г. Московский институт электронной техники (МИЭТ). Канд. техн. наук, член РНТОРЭС им. А.С. Попова и института IEEE (США), ведущий научный сотрудник ГУП НПЦ "ЭЛВИС" (Зеленоград, г. Москва). Автор около 120 научных работ в области адаптивной обработки сигналов, цифровой обработки речи и цифровой передачи данных.
V.I. Djigan (b. 1958) graduated from the Moscow Institute of Electronic Engineering in 1984. Ph. D. (Eng.), Russian Popov's Radio Society member, IEEE senior member (M'96, SM'2004), principal researcher of "ELVEES" R&D Center (Zelenograd, Moscow region). Author of about 120 publications in the field of adaptive signal processing, digital speech processing and digital communication.
В издательстве МГТУ им. Н.Э. Баумана в 2006 г. вышла в свет книга
Суржиков С.Т.
Физическая механика газовых разрядов. - М.: Изд-во МГТУ им. Н.Э. Баумана, 2006. - 640 с.: 384 ил. (Компьютерные модели физической механики).
Рассмотрены методы компьютерного моделирования электроразрядных процессов и динамики частично ионизованных газов, которые используются в задачах физической механики, физики газовых разрядов и аэрофизики. Основное внимание уделено решению двумерных задач физической механики тлеющих разрядов в аэрокосмических приложениях.
Для научных сотрудников и инженеров, работающих в области физической газовой динамики, физики низкотемпературной плазмы и газовых разрядов, а также для студентов и аспирантов физико-технических специальностей университетов.
По вопросам приобретения обращаться по тел. 433-82-98; e-mail: surg@ipmnet.ru
