Аннотация и ключевые слова
Аннотация (русский):
Рассмотрена возможность применения метода безградиентной n-параметрической минимизации Powell'64 в задачах мониторинга рабочего процесса судовых дизелей. Приведен пример поиска глобального минимума тестовой функции Розенброка. С помощью метода Powell'64 минимизируются функционалы МНК в задачах синхронизации и моделирования кривых сжатия/расширения в рабочем цилиндре. Показаны случаи расчета синхронизации данных для малооборотных двухтактных и среднеоборотных четырехтактных судовых дизелей. Задача синхронизации решается на базе уравнения P ' = 0, составленного для участка от начала сжатия до начала сгорания в цилиндре. Показан выбор граничных условий моделирования. Преимуществом использования метода Powell'64 является его высокая эффективность для квадратичных функционалов. В отличие от градиентных методов, метод Powell'64 не требует расчета производных и является универсальным для минимизации сложных нелинейных функционалов общего вида. Оригинальный авторский алгоритм синхронизации данных путем анализа индикаторных диаграмм, в котором используется метод Powell'64, применен в последних версиях систем мониторинга судовых дизелей D4.0H.

Ключевые слова:
безградиентная минимизация, функция Розенброка, рабочий процесс дизеля, политропа сжатия, функционал
Текст
Введение и постановка задачи Мониторинг рабочего процесса судовых дизелей сводится к задаче непрерывного циклического контроля давленияв рабочем цилиндре P в реальном времени (каждый цикл) или в режиме псевдореального времени (через несколько циклов). Контроль предполагает запись и отображение диаграмм давления газов в цилиндрах, а также расчет основных параметров рабочего цикла, заканчивающийся определением среднего индикаторного давления и индикаторной мощности цилиндра [1]. Значения давления P представлены на оси ординат. Ось абсцисс может быть представлена углами поворота коленчатого вала (ºПКВ) - P(φ) (рис. 1, б) либо соответствующими этим углам объемами надпоршневого пространства цилиндра - P(V) (рис. 1, a). Рис. 1. Диаграммы рабочего цикла дизеля в координатах P(V) - (a) и P(φ) - (б) [1] Кроме того, для анализа тепловой напряженности цилиндра и решения некоторых диагностических задач рабочий цикл может быть представлен в виде T(S)-диаграмм [2], рассмотрение которых выходит за рамки настоящей публикации. Представление рабочего цикла в P(V)-координатах было широко распространено на флоте до 2000-х гг. Его достоинством является то, что площадь замкнутой P(V)-диаграммы однозначно характеризует цикловую работу цилиндра [3]. Если рассматривать одновременно индикаторные P(V)-диаграммы всех цилиндров, то по их форме и площади хорошо визуально различимо, какой цилиндр нагружен больше, а какой меньше. В этом случае для определения перегруженного или недогруженного цилиндра не нужно даже считать среднее индикаторное давление и индикаторную мощность. Визуальный анализ распределения мощностей по цилиндрам был достаточным на практике в большинстве эксплуатационных ситуаций. Однако это является хоть и важным, но и практически единственным преимуществом P(V)-способа представления индикаторных диаграмм. Для качественного сгорания топлива в цилиндре необходима синхронизированная работа топливной системы высокого давления (ТА) и механизма газораспределения (МГР). Фазы закрытия/открытия клапанов и фазы впрыска топлива должны быть строго регламентированы. Наглядный контроль работы ТА и МГР возможен с помощью развернутых P(φ)-диаграмм (рис. 1, б), которые в последние годы все более широко применяются на флоте [4]. Смысл в том, что на фоне диаграммы рабочего цикла, представленной в ПКВ от нижней мертвой точки (НМТ) через верхнюю мертвую точку (ВМТ) и снова в НМТ, удобно показывать фазы закрытия/открытия клапанов и фазы впрыска топлива, которые получаются разнесенными по разным углам диаграммы, что удобно для визуального анализа (рис. 1, б). На индикаторную P(φ)-диаграмму можно наложить отмасштабированную диаграмму давления впрыска топлива или вибродиаграмму впрыска, а также вибродиаграмму работы МГР [1], и это улучшит визуальное восприятие индикаторной диаграммы и повысит ее диагностическую информативность. Для представления индикаторных диаграмм в P(V)-, P(φ)- или T(S)-видах необходимо иметь как можно более точное соответствие между давлениями в цилиндре и углами ПКВ. Для этого используются датчики на маховике двигателя (pick-up sensors [4]), которые измеряют фазу ПКВ синхронно с датчиками давления в цилиндре. Далее, при известной геометрии цилиндра и кривошипно-шатунного механизма (КШМ), производится пересчет ПКВ в объемы надпоршневого пространства [3]. Проблема заключается в том, что аппаратное определение фаз ПКВ с помощью pick-up sensors, настроенных в статике, будет иметь погрешности во время работы двигателя под нагрузкой. Основные составляющие этой погрешности следующие: - выборка микролюфтов КШМ при знакопеременной нагрузке; - крутильные колебания и скручивание коленчатого вала; - волновые эффекты и задержка сигнала давления при прохождении им пути от условного центра камеры сгорания до мембраны датчика в случае измерения давления через канал индикаторного крана. Последняя составляющая может быть весьма существенной (до 8º ПКВ на среднеоборотных СОД и до 12º ПКВ на высокооборотных дизелях) [5]. А в связи с тем, что все переносные системы мониторинга используют именно этот способ измерения давления и это является обычной практикой индицирования судовых дизелей, можно смело предположить, что использование pick-up sensors для переносных систем вообще нецелесообразно. Все виды расчетного определения ВМТ дают заведомо меньшую погрешность, нежели погрешность от задержки сигнала [6]. Состояние проблемы Анализ последних исследований и публикаций показывает в первую очередь, что задача определения положения ВМТ является важной составляющей мониторинга рабочего процесса, поскольку от точности ее решения существенно зависит расчетное значение среднего индикаторного давления и индикаторной мощности цилиндра. В общем виде зависимость погрешности расчета индикаторной мощности от погрешности ВМТ показана в работе M. Tazerout - на каждый градус погрешности приходится до 9 % расчетного значения мощности [6]. Термодинамические методы расчета и корректировки положения ВМТ предложны в работах [5, 6]. Они основаны, в частности, на расчете энтропии при известных значениях давления в цилиндре. Комментируя их, можно отметить, что наличие цифровых и аналоговых шумов при измерении давления в цилиндре может значительно снизить точность сложных расчетов. Численное определение первой, а тем более второй производной не имеет смысла без специальной обработки исходных данных (например, LowPass FFT Filter [7]). Интересное свойство индикаторной диаграммы, точнее ее первой производной, показано в работе [8]. С учетом допущений, принятых на практике, показано, что в точке максимальной скорости роста давления на участке сжатия объем цилиндра равен где Таким образом, объем цилиндра над поршнем в указанной точке может быть рассчитан с помощью известных геометрических размеров цилиндра и значений Р(φ) и P'(φ), полученных из индикаторной диаграммы. Далее может быть определен объем камеры сжатия и положение ВМТ, как предложено в [9, 10]. Красивый теоретический вывод данного положения, к сожалению, также сложно использовать на практике, т. к. точность расчета объема камеры сжатия и положения ВМТ сильно зависит от шумов записанной диаграммы давления. Применение LowPass фильтра [7] необходимо обосновывать для каждого конкретного случая, учитывая параметры аналого-цифрового преобразования, дискретность и аналоговые шумы при измерении давления в рабочих цилиндрах с помощью датчиков разных фирм. Целью данной публикации является анализ решения задачи определения ВМТ с помощью моделирования кривой скорости сжатия. При анализе используется часть данных индикаторной диаграммы до начала сгорания, по которым определяются коэффициенты модели P'(φ). При этом минимизируются функционалы, составленные согласно требованиям метода наименьших квадратов (МНК): где l1, l2 - границы моделирования, определяемые от начала сжатия до начала сгорания; P'i - значения первых производных, полученные с помощью методов численного дифференцирования индикаторной диаграммы; σi - весовые функции. Минимизацию указанного функционала предлагается произвести с помощью метода безградиентной минимизации Powell'64 [11, 12]. Этот метод позволяет определить минимум нелинейной функции n-переменных путем проведения последовательных поисков вдоль системы сопряженных направлений [11]. Метод Powell'64 не использует производных для осуществления поиска, что весьма удобно в практических расчетах. Кроме того, он является эффективным не только для квадратичных функций, но и для нелинейных n-параметрических функций общего вида. Пример использования метода Powell'64 на тестовой функции Покажем вначале работу нелинейного метода n-параметрической безградиентной минимизации по сопряженным направлениям Powell'64 на примере тестовой функции Розенброка [11]: Известная нелинейная двухпараметрическая функция Розенброка часто используется для тестирования методов минимизации в связи с ее специфической формой. Ее легко представить визуально (рис. 2) - она представляет собой изогнутый овраг с крутыми склонами и длинным слабо изменяющимся криволинейным дном. Ее минимум (Z = 0) находится в точке [x1 = 1, x2 = 1]. Рис. 2. Поиск минимума функции Розенброка [11] с помощью метода Powell'64 (интерфейс ПО) Специфика поиска минимума функции Розенброка заключается в том, что если задать начальную точку поиска [x1 = -2, x2 = 2] (Z = 3609), которая относительно точки минимума находится с другой стороны криволинейного оврага, то метод минимизации должен обойти все криволинейное дно функции и выйти на точку глобального минимума. От того, насколько быстро будет найден минимум функции Розенброка из неудобной начальной точки поиска [x1 = -2, x2 = 2], зависит эффективность того или иного метода минимизации. Учитывается общее количество итераций метода и количество расчетов самой функции Розенброка. Использование сопряженных [11, 12] направлений поиска является эффективным алгоритмом для квадратичных n-параметрических функций и, соответственно, для решения задач МНК. В специальной литературе по анализу методов минимизации [11] указано, что ряд градиентных методов, в частности метод наискорейшего спуска, оказывается неэффективным для некоторых квадратичных функций. Например, в данном случае, при начальной точке поиска [x1 = -2, x2 = 2], градиентные методы могут производить большое количество вычислений функции Розенброка, определяя локальные минимумы между боковыми склонами оврага и при этом медленно двигаясь в сторону глобального минимума. Применяемый в методе Powell'64 алгоритм поиска по сопряженным направлениям лишен этого недостатка. Из каждого локального минимума на следующей итерации, при необходимости, происходит поворот осей в сопряженных направлениях [12], и движение в сторону глобального минимума происходит большими шагами по дну оврага функции (рис. 3). Пошаговый алгоритм метода Powell'64 описан в [11]. В данной публикации используется модификация метода Powell'64, которая предусматривает использование метода «золотого сечения» (GM) [11] при поиске локальных минимумов в направлении осей. В оригинальном методе используется поиск минимума в осевых направлениях с помощью квадратичной функции [12]. Оказалось, что эффективность метода «золотого сечения» выше, чем у метода последовательного деления пополам, и существенно превосходит эффективность поиска непоследовательными методами [11]. Поиск GM основан на разбиении отрезка прямой на две части, Φ1 и Φ2, при этом отношение длины всего отрезка к большей части равно отношению большей части к меньшей. Отрезки пропорциональны числам Фибоначчи: Рис. 3. Этапы поиска функции Розенброка [11] в районе «точка минимума» По мнению авторов, использование метода «золотого сечения» более эффективно для нелинейных функций общего вида, для которых расчет с помощью оригинального алгоритма [12] часто приводит к переполнению регистров памяти компьютера, что требует необходимого специального контроля. В данном случае количество общих итераций при поиске минимума функции Розенброка с поворотом осей в сопряженных направлениях составило 31 (см. рис. 2). На рис. 2 видно, что метод Powell'64 уже на второй итерации быстро определяет положение дна оврага функции и затем осуществляет большие эффективные шаги в сторону глобального минимума по дну оврага (рис. 3). Сравнительный анализ методов n-параметрической нелинейной минимизации, приведенный в [11], показывает, что метод Powell'64 является наиболее эффективным методом поиска минимумов n-параметрических нелинейных функций общего вида без использования производных. В общем случае градиентные методы дают большую скорость поиска, но необходимость задания частных производных по параметрам поиска усложняет оптимизационную задачу и может служить причиной трудновыявляемых ошибок. Часто минимизируемый функционал представляет собой сложную нелинейную функцию общего вида либо является результатом работы алгоритма со сложной условно изменяющейся структурой. В этом случае для поиска минимума такого функционала задаются граничные условия по независимым переменным, в рамках которых должен быть найден минимум. Использование первых, а тем более вторых производных по независимым переменным часто оказывается неудобным для исследователей, т. к. почти всегда исследователи модифицируют свой алгоритм, и численное приближение производных усложняет процесс исследований, являясь источником трудновыявляемых ошибок. Использование же алгоритма Powell'64 не связано с производными и обладает эффективностью, не уступающей градиентным методам [11, 12]. Решение задачи определения ВМТ с помощью уравнения Pʹ(φ) = 0 Используется тот факт, что при отсутствии сгорания в цилиндре скорость изменения давления в ВМТ равна нулю [5], за вычетом термодинамического смещения, связанного с передачей тепла в стенки цилиндра (рис. 4) [3, 10]: . Рис. 4. Моделирование на участке сжатия: Pc' - давление начала сгорания; Pc - давление в конце сжатия; Pz - максимальное давление сгорания топлива; Pexp - давление в точке 36º за ВМТ С учетом термодинамического фазового смещения в ВМТ (при φ = 0) Тогда выражение для скорости изменения давления на участке сжатия можно записать как (1) где Pa - давление в начале сжатия; - объем цилиндра: ; - объем, описываемый полным ходом поршня; - объем камеры сжатия; - отношение радиуса кривошипа к длине шатуна; Va - объем цилиндра в начале сжатия; n1 - показатель политропы сжатия. После несложных преобразований правая часть уравнения (1) преобразуется в нелинейное трансцендентное уравнение относительно угла ПКВ: (2) Действительный корень уравнения (2) определяется с помощью численного метода (например, BISECT [11]) в диапазоне между абсциссами максимумов P'(φ) на участках сжатия и сгорания ([P'max1, P''max2] (рис. 5). Этот корень представляет собой расчетное значение ВМТ. Рис. 5. Графики кривой сжатия/расширения СОД ЧН25/34 В процессе расчета алгоритм Powell'64 используется дважды: сначала для предварительной оценки положения ВМТ с помощью синусоидальной модели, затем окончательный расчет путем моделирования скорости изменения давления dp / dφ. Метод Powell'64 используется для минимизации функционала, записанного согласно требованиям МНК: , где m - количество точек на участке до . На рис. 5 представлен график скорости изменения давления, полученный с помощью методов численного дифференцирования кривой . На графике второй производной (рис. 5) показаны контрольные точки, определяющие положение ВМТ (точка m), а также граничные точки для построения линейной и синусной моделей. В отдельных случаях поиск коэффициентов модели (1) удобнее производить путем поиска минимума второй производной . Выражение для расчета второй производной запишем в виде Обозначим Запишем выражения для производных Квадрат объема , поскольку объем положительный. Ограничение учитывается , где , Должно быть выполнено условие Поскольку объем в цилиндре всегда положительный, неравенство можно сократить: или . Запишем выражение в квадратных скобках для производной (С)': С учетом производной в квадратных скобках запишем выражение для (С)' и произведем необходимые действия: Выражение для расчета второй производной от давления запишем в виде Для случая анализа произвольного ряда значений давления, где представлены несколько подряд записанных в память рабочих циклов, предварительную оценку значений ВМТ (точка m) каждого цикла удобней производить путем поиска локальных минимумов (рис. 5). Кроме того, при расчете предварительной синхронизации в синусоидальной и линейной моделях используются координаты точек B и C, где , а . Эти координаты рассчитываются как положительные корни нелинейного уравнения , где область определения корней где Δφ ≤ 20° ПКВ. Для определения левой границы в предварительной синхронизации по синусоидальной модели используется левая координата максимума второй производной. В этом случае производится поиск локального максимума на области определения значений углов φ: где Δφ ≤ 30° ПКВ. В конечном итоге уточнялись численные значения коэффициентов , рассчитанные по записанным точкам участка кривой скорости изменения давления до начала сгорания, и уточнялось истинное значение ВМТ как угол, при котором скорость изменения давления равна нулю. Пример решения задачи определения ВМТ с помощью метода Powell'64 Определение положения ВМТ с помощью моделирования кривой скорости сжатия P'(φ) и последующего решения уравнения P'(φ) = 0 является, по мнению авторов, наиболее эффективным на практике методом. Для решения задачи используется тот факт, что при отсутствии сгорания в цилиндре скорость изменения давления в ВМТ равна нулю. Цифровые и аналоговые шумы в исходных данных усложняют моделирование участка сжатия на индикаторной диаграмме, но эта задача решается с помощью, например, Butterworth LowPass FFT Filter [7]. На рис. 6. представлено решение задачи определения ВМТ вышеуказанным алгоритмом для судового среднеоборотного дизеля Yanmar 6EY18ALW. Рис. 6. Индикаторная диаграмма P (φ) и модели P', P'' среднеоборотного дизеля Yanmar 6EY18ALW Начало сгорания топлива в данном случае - до ВМТ, и этот случай более сложен для расчета, чем случай с расчетом малооборотного дизеля с поздним впрыском топлива. Интервал моделирования кривой сжатия: от начала сжатия (закрытия клапанов МГР) до начала впрыска топлива в цилиндр. Фазу начала сгорания Pc' было бы удобно рассматривать в качестве правой границы, т. к. для этого достаточно наличия только самой диаграммы P(φ) [13]. Однако если рассматривать правую границу как фазу начала сгорания, то в участок моделирования попадает часть данных после впрыска топлива, которые уже не представляют собой политропу сжатия вследствие испарения топлива, отбора на это части энергии и частичного снижения давления в цилиндре [14]. При использовании методов вибродиагностики топливной аппаратуры высокого давления можно определить фазу начала впрыска топлива в цилиндр [5]. Во всех остальных случаях необходимо учитывать задержку самовоспламенения топлива и сдвигать на нее правую границу моделирования кривой сжатия относительно точки начала сгорания (Pc', рис. 6). Заключение Алгоритм программной синхронизации данных рабочего процесса был разработан и впервые применен в переносных системах DEPAS 2.34 [5]. Последняя модификация алгоритма применяется в системах DEPAS D4.0HT и получила название PLS-алгоритм (PhaseLess Synchronization algorithm). Во всех случаях, где производился поиск локальных минимумов и максимумов нелинейных функций от двух и более переменных, использовался метод минимизации по сопряженным направлениям, не использующий производные - метод Powell'64. Метод оказался устойчивым даже для случая 7 % «зашумленности» данных белым шумом (оценка уровня шумов для данных, полученных системой D4.0H «полезный сигнал/шум» оценивается в 3 %) [5, 13]. Несмотря на то, что, по сути, метод Powell'64 - это метод поиска с условным окончанием итераций, удалось сформулировать такие начальные условия, при которых общее время решения задачи синхронизации оказалось приемлемо малым (единицы миллисекунд) для большинства современных контроллеров. Накопленный опыт работы в области мониторинга рабочего процесса судовых дизелей позволяет утверждать следующее: использование PLS-алгоритма при практическом индицировании СДВС более предпочтительно по следующим соображениям: 1) автоматический учет погрешностей определения ВМТ (установка фазового датчика и маркировка маховика выполняется на остановленном двигателе. Во время работы дизеля ВМТ смещается из-за скручивания коленчатого вала, пропорционально нагрузке увеличивающегося на дальних от маховика цилиндрах; из-за крутильных колебаний; зазоров в КШМ и других факторов, которые невозможно учесть в «статике»): - PLS-алгоритм автоматически учитывает влияние скручивания коленчатого вала на нагруженном двигателе; - PLS-алгоритм автоматически учитывает влияние несоответствия между истинным положением ВМТ и маркировки на маховике, возникающего вследствие возможной неточной маркировки маховика, влияния зазоров в деталях КШМ и других эксплуатационных факторов; - PLS-алгоритм автоматически учитывает влияние конечной скорости прохождения волны давления в канале индикаторного крана (от камеры сгорания до мембраны датчика давления); 2) возможность производить индицирование без предварительной подготовки двигателя (во время работы с системами мониторинга рабочего процесса, использующими аппаратную синхронизацию, наиболее трудоемкая и длительная часть настройки системы - установка фазового датчика и маркировка маховика): - в случае применения PLS-алгоритма отсутствуют временные и финансовые затраты на установку фазовых датчиков; - применение PLS-алгоритма позволяет производить индицирование дизеля непосредственно в процессе его эксплуатации без принудительного изменения режима его работы и без специальной подготовки, которая необходима при аппаратной синхронизации данных; - PLS-алгоритм позволяет создать универсальную переносную малогабаритную систему диагностического контроля рабочего процесса транспортных дизелей [1, 5, 13].
Список литературы

1. Varbanets R., Karianskiy A. Analyse of marine diesel engine performance // Journal of Polish CIMAC. Energetic Aspects. Gdansk: Faculty of Ocean Engineering and Ship Technology Gdansk University of Technology. 2012. Vol. 7. No. 1. P. 269-275.

2. Staś M. An Universally Applicable Thermodynamic Method for T.D.C. Determination // SAE Technical Paper 2000-01-0561, 2000. URL: http://papers.sae.org/2000-01-0561/. DOI: 10.4271/2000-01-0561.

3. Heywood J. B. Internal Combustion Engine Fundamentals: First Edition. McGraw-Hill Education, 1988. 930 p.

4. CCM - COMBUSTION MONITORING SYSTEM. URL: http://www.imes.de/cms.html.

5. Варбанец Р. А. Диагностический контроль рабочего процесса судовых дизелей в эксплуатации: дис. … д-ра техн. наук. Одесса: Одес. национ. морск. акад., 2010. 314 с.

6. Tazerout M., Le Corre O., Rousseau S. TDC Determination in IC Engines Based on the Thermodynamic Analysis of the Temperature-Entropy Diagram // SAE Technical Paper 1999-01-1489, 1999. URL: http://papers.sae.org/1999-01-1489/. DOI: 10.4271/1999-01-1489.

7. Butterworth filter. URL: https://en.wikipedia.org/wiki/Butterworth_filter.

8. Ерыганов А. В., Варбанец Р. А. Определение степени сжатия дизеля по результатам индицирования рабочего процесса // Вестн. Астрахан. гос. техн. ун-та. Сер.: Морская техника и технология. 2017. № 1. С. 44-49.

9. Семенов В. С. Современные проблемы теории судовых дизелей. М.: В/О Мортехинформреклама, 1991. 112 с.

10. Ваншейдт В. А. Судовые двигатели внутреннего сгорания. Л.: Судостроение, 1977. 392 с.

11. Химмельблау Д. Прикладное нелинейное программирование. М.: Мир, 1975. 525 с.

12. Powell M. J. D. An efficient method for finding the minimum of a function of several variables without calculating derivatives // Computer J. 1964. No. 7. P. 155.

13. Варбанец Р. А., Головань А. И., Кучеренко Ю. Н. Мониторинг частотных параметров судового дизеля с турбонаддувом // Вестн. Астрахан. гос. техн. ун-та. Сер.: Морская техника и технология. 2013. № 1. С. 103-110.

14. Белоусов Е. В. Топливные системы современных судовых дизелей: учебн. пособ. Херсон: ХГМА, 2014. 267 с.