Дельта-оператор при проектировании цифровых систем на малоразрядных вычислителях

Каримов Тимур Искандарович ассистент, кафедра САПР, Санкт-Петербургский государственный электротехнический университет «ЛЭТИ» им. В.И. Ульянова (Ленина) (СПбГЭТУ «ЛЭТИ») 197022, Россия, г. Санкт-Петербург, ул. Профессора Попова, 5, ауд. 1152 Karimov Timur Iskandarovich Assistant, Department of CAD, Ulyanov (Lenin) St. Petersburg State Electrotechnical University "LETI"  197022, Russia, g. Saint Petersburg, ul. Professora Popova, 5, aud. 1152 carimus@gmail.com Каримов Артур Искандарович ассистент, кафедра САПР, Санкт-Петербургский государственный электротехнический университет «ЛЭТИ» им. В.И. Ульянова (Ленина) (СПбГЭТУ «ЛЭТИ») 197022, Россия, г. Санкт-Петербург, ул. Профессора Попова, 5, ауд. 1152 Karimov Artur Iskandarovich Assistant, Department of CAD, Ulyanov (Lenin) St. Petersburg State Electrotechnical University "LETI" 197022, Russia, g. Saint Petersburg, ul. Professora Popova, 5, aud. 1152 art.krmv@gmail.com Бутусов Денис Николаевич кандидат технических наук доцент, кафедра САПР, Санкт-Петербургский государственный электротехнический университет «ЛЭТИ» им. В.И. Ульянова (Ленина) (СПбГЭТУ «ЛЭТИ») 197022, Россия, г. Санкт-Петербург, ул. Профессора Попова, 5, ауд. 1152 Butusov Denis Nikolaevich PhD in Technical Science Associate Professor, Department of CAD, Ulyanov (Lenin) St. Petersburg State Electrotechnical University "LETI" 197022, Russia, g. Saint Petersburg, ul. Professora Popova, 5, aud. 1152 butusovdn@mail.ru Сольницев Ремир Иосифович доктор технических наук профессор, кафедра САПР, Санкт-Петербургский государственный электротехнический университет «ЛЭТИ» им. В.И. Ульянова (Ленина) (СПбГЭТУ «ЛЭТИ») 197022, Россия, г. Санкт-Петербург, ул. Профессора Попова, 5, ауд. 1152 Sol'nitsev Remir Iosifovich Doctor of Technical Science Professor, Department of CAD, Ulyanov (Lenin) St. Petersburg State Electrotechnical University "LETI" 197022, Russia, g. Saint Petersburg, ul. Professora Popova, 5, aud. 1152 remira70@mail.ru Булахов Артем Викторович аспирант, кафедра Систем автоматизированного проектирования, Санкт-Петербургский государственный электротехнический университет «ЛЭТИ» им. В.И. Ульянова (Ленина) (СПбГЭТУ «ЛЭТИ») 197022, Россия, г. Санкт-Петербург, ул. Профессора Попова, 5, ауд. 1151 Bulakhov Artem Viktorovich Postgraduate Student, Department of Computer-Aided Design Systems,  Ulyanov (Lenin) St. Petersburg State Electrotechnical University "LETI"  197022, Russia, g. Saint Petersburg, ul. Professora Popova, 5, aud. 1151 bulahov.a@visteh.ru

Сигаева Мария Сергеевна ассистент, кафедра ИнЯз, Санкт-Петербургский государственный электротехнический университет «ЛЭТИ» им. В.И. Ульянова (Ленина) (СПбГЭТУ «ЛЭТИ») 197022, Россия, г. Санкт-Петербург, ул. Профессора Попова, 5, ауд. 1152 Sigaeva Maria Sergeevna Assistant, InYaZ Department, Ulyanov (Lenin) St. Petersburg State Electrotechnical University "LETI" 197022, Russia, g. Saint Petersburg, ul. Professora Popova, 5, aud. 1152 msigaeva@gmail.com

Работа выполнена при финансовой поддержке Российского фонда фундаментальных исследований, грант № 17-07-00862 А. "Теория и средства проектирования цифровых генераторов хаотических сигналов"

Состояние предметной области

Дискретные модели динамических систем традиционно описываются с использованием z‑преобразования, иначе называемого оператором задержки. Синтез моделей в z‑области основан на хорошо изученных методах дискретизации непрерывных систем и в настоящее время полностью формализован. Однако при реализации дискретных систем на ограниченной машинной сетке из-за округления коэффициентов разностных уравнений погрешность вычислений может стать недопустимо высокой. Для повышения точности приходится выбирать большую длину машинного слова, что негативно сказывается на энергоэффективности, быстродействии и стоимости конечного решения ^[8].

В работах ^{[2, 12]} в качестве альтернативы оператору задержки предложен δ‑оператор. Собственные числа дискретной модели, полученной с помощью δ-преобразования, при малых периодах квантования не стремятся к единице или нулю, в отличие от модели, полученной z-преобразованием ^[2]. Применение δ-оператора позволяет во многих случаях добиться существенного улучшения параметров дискретных систем. Например, в работе ^[10] описано применение δ-оператора при решении задачи нахождения оптимального скользящего режима управления. Авторы отмечают, что полученные ими модели значительно меньше подвержены автоколебаниям по сравнению с системами, синтезированными традиционным способом. В работе ^[6] исследователи сравнивают реализацию ПИД‑регуляторов на основе δ- и z-оператора. Показано, что коэффициенты разностных уравнений в случае применения δ-оператора значительно менее чувствительны к ограничению длины разрядной сетки. Работа ^[11] демонстрирует возможности δ-оператора при реализации LQR-регулятора. При длине мантиссы машинного слова в 5 бит авторам удалось добиться точности, соответствующей реализации z-модели в 32-битных числах типа single с мантиссой длиной 23 бита.

Несмотря на то, что в последнее время наблюдается возрастающая публикационная активность по тематике δ-оператора, практические аспекты его применения освещаются крайне скудно. Остается открытым вопрос, в каком случае к системе должен применяться δ-оператор, а когда можно использовать классический подход на основе z-оператора. Эксперименты с произвольными дискретными системами дают противоречивые результаты: на одних системах разница в точности представления модели существенна, в то время как на других практически незаметна. Ранние работы предсказывали постепенный переход от оператора сдвига к δ-оператору ^[7], чего не наблюдается в действительности. Во многих источниках ^{[13, 14]} авторы делают акцент лишь на математическом определении δ-оператора и рассмотрении частных случаев его использования.

Таким образом, можно заключить, что разработка принципов выбора между z-оператором и δ-оператором является актуальной задачей. В данной работе авторы предлагают методику и критерии выбора дискретного оператора, рассматривая практические аспекты построения дискретных систем на платформах малой разрядности с плавающей запятой.

1. Критерий предпочтительности δ-оператора

Наиболее известно определение δ-оператора, изложенное в работах ^{[1, 3, 4]}. Пусть— определение z-оператора через оператор Лапласа s = σ + jω, а — определение δ-оператора через z-оператор. Здесь ∆ — некоторый параметр, с математической точки зрения способный принимать любые значения. В ранних публикациях по теме δ-оператора авторы предлагали установить ∆ = T_s, где T_s — период дискретизации ^[12]. Позднее в ряде работ было показано, что более целесообразен выбор ∆ = 2^-n, n ϵ N ^[9].

На практике переход от передаточной функции системы, заданной в z-области, к системе в δ-области состоит в алгебраическом преобразовании коэффициентов передаточной функции с последующей реализацией системы таким образом, что вместо оператора задержки применяется δ-оператор. Рассмотрим выражения, используемые для преобразования коэффициентов системы:

, (1)

. (2)

Здесь a_i, α_i — коэффициенты знаменателя z‑системы и δ‑системы соответственно, b_i, β_i— коэффициенты числителя, n— порядок системы, k— порядок операторной переменной при коэффициенте, а выражения вида — биномиальные коэффициенты, рассчитываемые как:

.

Для звена второго порядка, заданного в канонической форме II (Direct Form II) и реализованного с помощью δ‑оператора, алгоритм вычисления нового состояния будет следующим:

Сравним исследуемые операторы на примере проектирования цифрового БИХ‑фильтра, реализуемого в арифметике с типом данных single. Синтезируем дискретный эллиптический фильтр нижних частот 2‑го порядка со следующими параметрами: уровень колебаний в полосе пропускания 10^-4дБ, уровень подавления сигнала в полосе задерживания 40 дБ, приведенные частоты среза фильтра λ_s = 0,03 (случай 1) и λ_s = 0,002 (случай 2). Опустив для краткости математическое описание фильтра, построим графики погрешности АЧХ моделей фильтра, реализованных с помощью z и δ‑операторов, относительно аналитически рассчитанной АЧХ.

а) б)

Рисунок 1 — Сравнение АЧХ различных реализаций тестового фильтра

В первом случае (рис. 1, а) δ‑система демонстрирует более низкую точность во всем диапазоне частот. Особенно критично это может повлиять на сигнал в полосе пропускания. Во втором случае (рис. 1, б) соотношение погрешностей свидетельствует о существенном улучшении характеристик модели фильтра, полученной с помощью δ‑оператора.

Из данного примера можно сделать вывод, что критерий целесообразности применения δ‑оператора зависит от нормализованной частоты среза фильтра. Общая методика применения δ‑оператора, предложенная авторами ранее в работе ^[4], была такова:

1. Создание z‑модели с требуемыми свойствами.

2. Проверка, является ли рабочая частота системы много меньшей, чем частота Найквиста (λ_раб<<1). Если условие выполняется, есть основания для применения δ‑оператора.

3. Оценка величины S = || E_z || / || E_δ ||, где || · || — некоторая норма, E_z, E_δ — погрешности АЧХ. Если S ≥ S*, где S* — заранее известная величина, то выбор делается в пользу δ‑оператора.

В работе ^[4] пункт 3 методики был переформулирован следующим образом: дискретная z‑система, комплексные собственные числа матрицы A пространства состояния которой лежат в области предпочтительности, с вероятностью 95% покажет лучшую точность при реализации ее с помощью δ‑оператора. Область предпочтительности δ‑оператора P_δ приближенно имеет форму полуэллипса на z‑плоскости с центром в точке 0 и полуосями (0,09; 0,16) — см. рис. 2. Если система обладает вещественными собственными числами, то ее область предпочтительности — отрезок (–0,09; 0) на действительной оси.

Рисунок 2 — Область предпочтительности δ‑оператора (белая окружность) на z‑плоскости. Черным цветом обозначена единичная окружность

Применение подобного критерия выбора дискретного оператора на практике оказалось не всегда удобным, особенно при реализации систем, заданных в форме передаточных функций. Работа ^[4] также оставляет открытым вопрос о влиянии на целесообразность применения δ‑оператора точности и формата представления чисел.

В последующих разделах статьи мы исследуем влияние погрешности округления коэффициентов в арифметике с плавающей запятой в случае применения z- и δ‑операторов, а также конкретизируем условие λ_раб << 1 пункта 2 методики.

2. Погрешность представления коэффициентов

Пусть E_r — погрешность машинного представления вещественного числа r. Она может быть вычислена по формуле

, (3)

где ε — машинный эпсилон числа r, а запись [r] обозначает машинное представление числа.

Ввод формулы (3) обоснован тем, что погрешность является смещением функции mod(r, ε): она возникает при отбрасывании выходящих за машинную сетку разрядов r, то есть как остаток от деления его на ε — минимально представимое число на данной сетке, соответствующее значению младшего бита мантиссы. При этом максимальная погрешность `E_r= r - [r]` равна

.

Рисунок 3 – Пилообразный характер погрешности машинного представления коэффициентов

В ходе исследования было сделано предположение, что зависимости между предпочтительностью вида дискретного оператора и разрядностью машинных чисел не наблюдается. Докажем это утверждение.

Лемма о верхних оценках погрешности АЧХ для z- и δ-моделей

Обозначим , а ε — машинный эпсилон числа 0,5. Пусть , причем погрешности числителя и знаменателя , как для z‑модели, так и для δ‑модели.

Тогда справедливы следующие формулы верхних оценок погрешностей АЧХ z- и δ-моделей:

, (4)

. (5)

Здесь и далее E — оценка погрешности, A = |H(λ)|—точная АЧХ системы.

Доказательство леммы основано на предельных свойствах

,

и анализе поведения погрешностей типовых звеньев 1-го и 2-го порядков.

Теорема 1. Пусть АЧХ на частоте λ равна H_δ(λ) = β/α, где α, β — комплексные значения числителя и знаменателя передаточной функции системы H_δ на данной частоте; аналогично H_z(λ) = b/a. Тогда при |α| >Kε, ε<ε*, где K, ε* — известные вещественные числа, величина

не зависит от машинного ε, то есть от разрядности типа данных.

Для доказательства рассмотрим четыре возможных случая.

1. Случай:

При |α| >Kε верно, что |a| + 1,5ε·sgn(a) ≈|a| в силу малости ε, откуда, используя формулу (4),

,

а при K >>ε* >ε верно 1 +ε·sgn(a) ≈ 1, откуда

.

Следовательно,

. (6)

Формула (6) не содержит ε, поэтому величина S₁ от него не зависит.

2. Случай .

Из-за неточности представления коэффициентов реальные значения АЧХ A_δ≠ 0,

A_z≠ 0.

Для z-системы

,

Eсли A_z<<ε, что весьма вероятно в силу условия A = 0.

Для δ-системы

.

В реальной системе |β₂| ≠ 0, но

.

Поэтому искомое отношение равно

и не зависит от ε.

3. Случай. Характерен для фильтров низких и высоких частот, когда коэффициенты числителя стремятся к целым числам. Тогда

,

,

и искомое отношение

,

и не зависит от ε.

4. Случай отвечает тривиальному условию H ≡ 0, и для него оценка S₁ не имеет смысла.

Теорема доказана.

Выводы из данной теоремы имеют важное прикладное значение.

1. Поскольку во многих случаях |α| >1,5ε (фактически это значит, что a₁ отличается от -2 больше чем на единицу младшего разряда, а a₂ — отличается от 1 больше чем на единицу младшего разряда), следовательно, отношение верхних оценок погрешностей z- и δ‑систем при малых λ не зависит от разрядности машинного представления коэффициентов системы ε. Прирост точности z‑модели за счет увеличения разрядности эквивалентен приросту точности δ‑модели. Имея данное λ или T, а также зная разрядность вычислителя, следует выбрать δ‑модель, если она даст лучшие результаты на той же разрядности.

2. Если используются высокие частоты дискретизации или низкая разрядность, то есть |α| >1,5ε, то δ‑модель демонстрирует сохранение точности, пропорциональное ε

,

в то время как z-модель показывает погрешность, сопоставимую с реальной величиной АЧХ, становясь при этом совершенно непригодной для практической реализации

.

В этом случае δ-оператор позволяет строить дискретные системы, не реализуемые стандартным способом.

3. Зависимость критерия предпочтительности от нормализованной частоты и типа фильтра

Рассмотрим зависимость значений критериев предпочтительности S и S₁ от частоты среза λ для фильтра Баттерворта нижних частот 2-го порядка. Особенностью данного фильтра является точное задание коэффициентов числителя z-модели, оба из которых стремятся к -1, поэтому используем оценку

. (7)

Формула (7) получена путем нахождения отношения выражений (4) и (5) и сокращения общих множителей. График отношения критериев S и S₁ от нормированной частоты среза представлен на рис. 4. Квадратом отмечена точка |α| =1,5ε, положение которой рассчитано для разрядности машинного слова 16 бит. При |α| >1,5ε и возрастании λ график почти линейно убывает. При |α| <1,5ε график горизонтален: S₁ практически не зависит от λ. При этом оценка S₁ коррелирует с S = E_z(λ)/E_δ(λ), но при |α| >1,5ε она явно завышена. Поэтому оценку S₁ нельзя рекомендовать для практического применения при расчете реальных погрешностей цифровых систем.

Рисунок 4 — Зависимость S и S₁ от нормированной частоты среза λ

Исследуем погрешность АЧХ для z- и δ-систем, а также погрешность машинного представления коэффициентов знаменателей их передаточных функций в окрестности λ= 0,02 (рис. 5). Рассматриваемая система — цифровой режекторный фильтр 2-го порядка с шириной полосы подавления 0,1λ_s.

Рисунок 5 — Погрешность представления машинных коэффициентов z- и δ‑реализаций режекторного фильтра

На основании анализа рис. 5 можно сделать ряд выводов:

1. На произвольном фрагменте интервала частот дискретизации длиной M∆λ гарантировано нахождение хотя бы одного максимума и минимума S. В зависимости от того, какая из погрешностей меньше — погрешность z‑системы в точке минимума, или δ‑системы в точке максимума — следует выбрать соответствующую форму представления дискретной системы.

2. Ввиду негладкости S нецелесообразно использовать стандартные алгоритмы оптимизации. Задача поиска оптимального λ может быть упрощена, поскольку тактовая частота целевой платформы может принимать ограниченное количество значений, линейно связанных с частотой встроенного генератора.

3. При достаточно большой разрядности представления данных z-система может быть точнее δ-системы даже при малом λ. С уменьшением разрядности целесообразность применения δ-оператора возрастает.

Таким образом, при определенном значении частоты среза λ < λ*с высокой долей вероятности предпочтительным становится δ-оператор. Для экспериментального подтверждения этого предположения были поставлены опыты по нахождению λ* для различных типов фильтров 2-го порядка. На рис. 6 приведены соотношения погрешностей АЧХ z- и δ-реализаций фильтра Баттерворта.

Рисунок 6 — Погрешности АЧХ z- и δ-реализаций тестового фильтра Баттерворта

Пользуясь критерием, сформулированным в работе ^[4], построим обобщенный график зависимости R(λ) для фильтров различного типа (рис. 7). Анализ графиков показывает, что эллиптический ФВЧ практически всегда выгодно реализовывать с помощью δ-оператора.

Рисунок 7 — Зависимость R(λ) для различных типов цифровых фильтров

4. Обобщенный критерий применения дельта-оператора

Чтобы визуализировать область предпочтительности δ-оператора, найдем значение величины

для систем, являющихся дискретным преобразованием непрерывной системы с двумя комплексно-сопряженными полюсами s = σ ±jω и передаточной функцией вида

.

Перейдем к H(z) с помощью преобразования Тастина ^[5], и далее найдем H(δ). Нормируем коэффициенты таким образом, чтобы при старшей степени знаменателя был единичный коэффициент, и уменьшим точность представления мантиссы до 8 бит. Затем построим АЧХ всех дискретных систем, нормированная частота дискретизации которых лежит на интервале (0; 1), и сравним со значением АЧХ этих дискретных систем в случае, когда их коэффициенты заданы в формате double.

В результате мы можем построить поверхность S(σT, ωT) над комплексной плоскостью (σT;ωT), где T — период дискретизации. Поверхность ограничена плоскостью на уровне . Внутри области S>1 δ‑оператор является более предпочтительным. На рис. 8(а) эта область закрашена белым.

а) б)

Рисунок 8 — Область предпочтительности δ-оператора на непрерывной плоскости (а) и пространственная визуализация критерия S (б)

Анализируя полученные результаты, можно предложить модифицированный алгоритм выбора дискретного оператора:

1. При необходимости декомпозировать исходную непрерывную систему на звенья первого и второго порядков.

2. Проверить корни характеристического полинома системы на принадлежность области предпочтительности S. Если условие соблюдается, то делается выбор в пользу δ-оператора.

3. Если возможна вариация частоты дискретизации системы λ в допустимой окрестности λ*, то следует численно отыскать локальные минимум и максимум определенного критерия S. Тогда, если

,

то выбор делается в пользу δ‑оператора, иначе — в пользу z‑оператора. Если точно найти S невозможно и используется оценочное значение S, то выбор в пользу δ‑оператора делается по условию S > 1.

Заключение

В работе сравниваются реализации дискретных линейных систем с помощью z- и δ-операторов в арифметике с плавающей запятой. Было установлено, что точность машинного представления чисел (разрядность вычислителя) не оказывает влияния на отношение погрешностей z- и δ-моделей: ключевым параметром, оказывающим влияние на предпочтительность выбора дискретного оператора, является частота дискретизации. Сформулирован обобщенный критерий, позволяющий по виду непрерывной системы и требуемой частоте дискретизации выбрать дискретный оператор. В случае когда, согласно критерию, предпочтительным оказывается δ-оператор, выигрыш в точности по сравнению с реализацией на основе z-оператора может составлять несколько порядков (десятки дБ соотношения сигнал/машинный шум) при незначительном усложнении структуры фильтра.

Вопросы масштабирования переменных состояния дискретной системы и ее коэффициентов оставлены в данной работе без внимания, мх решение является предметом дальнейших исследований.