Предисловие к первой части


Моделирование паровых турбин — повседневная задача сотен людей в нашей стране. Вместо слова модель принято говорить расходная характеристика. Расходные характеристики паровых турбин используют при решении таких задач, как вычисление удельного расхода условного топлива на электроэнергию и тепло, производимые ТЭЦ; оптимизация работы ТЭЦ; планирование и ведение режимов ТЭЦ.


Мною разработана новая расходная характеристика паровой турбины — линеаризованная расходная характеристика паровой турбины. Разработанная расходная характеристика удобна и эффективна в решении указанных задач. Однако на текущий момент она описана лишь в двух научных работах:


  1. Оптимизация работы ТЭЦ в условиях оптового рынка электроэнергии и мощности России;
  2. Вычислительные методы определения удельных расходов условного топлива ТЭЦ на отпущенную электрическую и тепловую энергию в режиме комбинированной выработки.

И сейчас в своем блоге мне бы хотелось:




1. Исходные данные


Исходными данными для построения линеаризованной расходной характеристики могут быть


  1. фактические значения мощностей Q0, N, Qп, Qт измеренные в процессе функционирования паровой турбины,
  2. номограммы qт брутто из нормативно-технической документации.

Конечно, фактические мгновенные значения Q0, N, Qп, Qт являются идеальными исходными данными. Сбор таких данных трудоемок.

В тех случаях, когда фактические значения Q0, N, Qп, Qт недоступны, можно обработать номограммы qт брутто. Они, в свою очередь, были получены на основании измерений. Подробнее об испытаниях турбин читайте в Горнштейн В.М. и др. Методы оптимизации режимов энергосистем.


2. Алгоритм построения линеаризованной расходной характеристики


Алгоритм построения состоит из трех шагов.


  1. Перевод номограмм или результатов измерений в табличный вид.
  2. Линеаризация расходной характеристики паровой турбины.
  3. Определение границ регулировочного диапазона работы паровой турбины.

При работе с номограммами qт брутто первый шаг осуществляется быстро. Такую работу называют оцифровкой (digitizing). Оцифровка 9 номограмм для текущего примера заняла у меня около 40 минут.


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


Расходная характеристика для рассматриваемой турбины строилась для следующих фиксированных значений параметров режима:


  • одноступенчатый режим работы,
  • давление пара среднего давления = 13 кгс/см2,
  • давление пара низкого давления = 1 кгс/см2.


Содержание архива Chuchueva-PT-80-linear-characteristic-curve.rar.


1) Номограммы удельного расхода qт брутто на выработку электроэнергии (отмеченные красные точки оцифрованы — перенесены в таблицу):


  • PT80_qt_Qm_eq_0_digit.png,
  • PT80_qt_Qm_eq_100_digit.png,
  • PT80_qt_Qm_eq_120_digit.png,
  • PT80_qt_Qm_eq_140_digit.png,
  • PT80_qt_Qm_eq_150_digit.png,
  • PT80_qt_Qm_eq_20_digit.png,
  • PT80_qt_Qm_eq_40_digit.png,
  • PT80_qt_Qm_eq_60_digit.png,
  • PT80_qt_Qm_eq_80_digit.png.

2) Результат оцифровки (каждому файлу csv соответствует файл png):


  • PT-80_Qm_eq_0.csv,
  • PT-80_Qm_eq_100.csv,
  • PT-80_Qm_eq_120.csv,
  • PT-80_Qm_eq_140.csv,
  • PT-80_Qm_eq_150.csv,
  • PT-80_Qm_eq_20.csv,
  • PT-80_Qm_eq_40.csv,
  • PT-80_Qm_eq_60.csv,
  • PT-80_Qm_eq_80.csv.

3) Скрипт MATLAB с расчетами и построением графиков:


  • PT_80_linear_characteristic_curve.m

4) Результат оцифровки номограмм и результат построения линеаризованной расходной характеристики в табличном виде:


  • PT_80_linear_characteristic_curve.xlsx.

Шаг 1. Перевод номограмм или результатов измерений в табличный вид


1. Обработка исходных данных


Исходными данными для нашего примера являются номограммы qт брутто.


Для перевода в цифровой вид множества номограмм нужен специальный инструмент. Я многократно использовала web-приложение для этих целей. Приложение просто, удобно, однако не имеет достаточной гибкости для автоматизации процесса. Часть работы приходится делать вручную.


На данном шаге важно оцифровать крайние точки номограмм, которые задают границы регулировочного диапазона работы паровой турбины.


Работа состояла в том, чтобы в каждом файле png при помощи приложения отметить точки расходной характеристики, скачать полученный csv и собрать все данные в одной таблице. Результат оцифровки можно найти в файле PT-80-linear-characteristic-curve.xlsx, лист «PT-80», таблица «Исходные данные».


2. Приведение единиц измерения к единицам мощности


Далее на листе «PT-80» рассчитываем значения Q0 по формуле

$$display$$\begin{equation} Q_0 = \frac {q_T \cdot N} {1000} + Q_П + Q_Т \qquad (1) \end{equation}$$display$$


и приводим все исходные величины к МВт. Расчеты реализованы средствами MS Excel.


Полученная таблица «Исходные данные (ед. мощности)» является результатом первого шага алгоритма.


Шаг 2. Линеаризация расходной характеристики паровой турбины


1. Проверка работы MATLAB


На данном шаге требуется установить и открыть MATLAB версии не ниже 7.3 (это старая версия, текущая 8.0). В MATLAB открыть файл PT_80_linear_characteristic_curve.m, запустить его и убедиться в работоспособности. Все работает корректно, если по итогам запуска скрипта в командной строке вы увидели следующее сообщение:


Значения считаны из файла PT_80_linear_characteristic_curve.xlsx за 1 сек
Коэффициенты: a(N) = 2.317, a(Qп) = 0.621, a(Qт) = 0.255, a0 = 33.874
Средняя ошибка = 0.006, (0.57%)
Число граничных точек регулировочного диапазона = 37

Если у вас возникли ошибки, то разберитесь самостоятельно, как их исправить.


2. Вычисления


Все вычисления реализованы в файле PT_80_linear_characteristic_curve.m. Рассмотрим его по частям.


1) Укажем название исходного файла, лист, диапазон ячеек, содержащий полученную на предыдущем шаге таблицу «Исходные данные (ед. мощности)».


XLSFileName = 'PT_80_linear_characteristic_curve.xlsx';
XLSSheetName = 'PT-80';
XLSRange = 'F3:I334';

2) Считаем исходные данные в MATLAB.


sourceData = xlsread(XLSFileName, XLSSheetName, XLSRange);
N = sourceData(:,1);
Qm = sourceData(:,2);
Ql = sourceData(:,3);
Q0 = sourceData(:,4);
fprintf('Значения считаны из файла %s за %1.0f сек\n', XLSFileName, toc);

Используем переменную Qm для расхода пара среднего давления Qп, индекс m от middle — средний; аналогично используем переменную Ql для расхода пара низкого давления Qn, индекс l от low — низкий.


3) Определим коэффициенты ?i.


Вспомним общую формулу расходной характеристики

$$display$$\begin{equation} Q_0 = f(N, Q_П, Q_Т) \qquad (2) \end{equation}$$display$$



и укажем независимые (x_digit) и зависимые (y_digit) переменные.


x_digit = [N Qm Ql ones(size(N,1),1)];        % электроэнергия N, промышленный пар Qп, теплофикационный пар Qт, единичный вектор
y_digit = Q0;                                 % расход острого пара Q0

Если вам непонятно, зачем в матрице x_digit единичный вектор (последний столбец), то читайте материалы по линейной регрессии. На тему регрессионного анализа рекомендую книгу Draper N., Smith H. Applied regression analysis. New York: Wiley, In press, 1981. 693 p. (есть на русском языке).


Уравнение линеаризованной расходной характеристики паровой турбины


$$display$$\begin{equation} Q_0 = \alpha_N \cdot N + \alpha_П \cdot Q_П + \alpha_Т \cdot Q_Т + \alpha_0 \qquad (3) \end{equation}$$display$$



является моделью множественной линейной регрессии. Коэффициенты ?i определим при помощи «большого блага цивилизации» — метода наименьших квадратов. Отдельно отмечу, что метод наименьших квадратов разработан Гауссом в 1795 году.


В MATLAB это делается одной строчкой.


A = regress(y_digit, x_digit);
fprintf('Коэффициенты: a(N) = %4.3f, a(Qп) = %4.3f, a(Qт) = %4.3f, a0 = %4.3f\n',...
    A);

Переменная A содержит искомые коэффициенты (см. сообщение в командной строке MATLAB).


Таким образом, полученная линеаризованная расходная характеристика паровой турбины ПТ-80 имеет вид


$$display$$\begin{equation} Q_0 = 2.317 \cdot N + 0.621 \cdot Q_П + 0.255 \cdot Q_Т + 33.874 \qquad (4) \end{equation}$$display$$


4) Оценим ошибку линеаризации полученной расходной характеристики.


y_model = x_digit * A;
err = abs(y_model - y_digit) ./ y_digit;
fprintf('Средняя ошибка = %1.3f, (%4.2f%%)\n\n', mean(err), mean(err)*100);

Ошибка линеаризации равна 0,57% (см. сообщение в командной строке MATLAB).


Для оценки удобства использования линеаризованной расходной характеристики паровой турбины решим задачу вычисления расхода пара высокого давления Q0 при известных значениях нагрузки N, Qп, Qт.


Пусть N = 82.3 МВт, Qп = 55.5 МВт, Qт = 62.4 МВт, тогда


$$display$$\begin{equation} Q_0 = 2.317 \cdot 82,3 + 0.621 \cdot 55,5 + 0.255 \cdot 62,4 + 33.874 = 274,9 \qquad (5) \end{equation}$$display$$


Напомню, что средняя ошибка вычислений составляет 0,57%.


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


  1. Вычислите величину Q0 с указанной точностью с использованием номограмм и ваших глаз.
  2. Автоматизируйте процесс расчета Q0 с использованием номограмм.

Очевидно, что в первой задаче определение значений qт брутто на глаз чревато грубыми ошибками.


Вторая задача громоздка для автоматизации. Поскольку значения qт брутто нелинейны, то для такой автоматизации число оцифрованных точек в десятки раз больше, чем в текущем примере. Одной оцифровки недостаточно, также необходимо реализовать алгоритм интерполяции (нахождения значений между точками) нелинейных значений брутто.


Шаг 3. Определение границ регулировочного диапазона работы паровой турбины


1. Вычисления


Для вычисления регулировочного диапазона воспользуемся другим «благом цивилизации» — методом выпуклой оболочки, convex hull.


В MATLAB это делается следующим образом.


indexCH = convhull(N, Qm, Ql, 'simplify', true);
index = unique(indexCH);
regRange = [N(index) Qm(index) Ql(index)];
regRangeQ0 = [regRange ones(size(regRange,1),1)] * A;
fprintf('Число граничных точек регулировочного диапазона = %d\n\n', size(index,1));

Метод convhull() определяет граничные точки регулировочного диапазона, заданного значениями переменных N, Qm, Ql. Переменная indexCH содержит вершины треугольников, построенных при помощи триангуляции Делоне. Переменная regRange содержит граничные точки регулировочного диапазона; переменная regRangeQ0 — значения расхода пара высокого давления для граничных точек регулировочного диапазона.


Результат вычислений можно найти в файле PT_80_linear_characteristic_curve.xlsx, лист «PT-80-result», таблица «Границы регулировочного диапазона».


Линеаризованная расходная характеристика построена. Она представляет собой формулу (4) и 37 точек, задающих границы (оболочку) регулировочного диапазона в соответствующей таблице.


2. Проверка


При автоматизации процессов расчета Q0 необходимо проверять, находится ли некоторая точка со значениями N, Qп, Qт внутри регулировочного диапазона или за его пределами (режим технически не реализуем). В MATLAB это можно делать следующим образом.


Задаем значения N, Qп, Qт, которые мы хотим проверить.


n = 75;
qm = 120;
ql = 50;

Проверяем.


in1 = inpolygon(n, qm, regRange(:,1),regRange(:,2));
in2 = inpolygon(qm, ql, regRange(:,2),regRange(:,3));
in = in1 && in2;
if in
    fprintf('Точка N = %3.2f МВт, Qп = %3.2f МВт, Qт = %3.2f МВт находится внутри регулировочного диапазона\n', n, qm, ql);
else
    fprintf('Точка N = %3.2f МВт, Qп = %3.2f МВт, Qт = %3.2f МВт находится снаружи регулировочного диапазона (технически недостижима)\n', n, qm, ql);
end

Проверка осуществляется в два шага:


  • переменная in1 показывает, попали ли значения N, Qп внутрь проекции оболочки на оси N, Qп;
  • аналогично переменная in2 показывает, попали ли значения Qп, Qт внутрь проекции оболочки на оси Qп, Qт.

Если обе переменные равны 1 (true), то искомая точка находится внутри оболочки, задающей регулировочный диапазон работы паровой турбины.


Иллюстрация полученной линеаризованной расходной характеристики паровой турбины


Наиболее «щедрые блага цивилизации» нам достались в части иллюстрации результатов расчетов.


Предварительно нужно сказать, что пространство, в котором мы строим графики, т. е. пространство с осями x – N, y – Qт, z – Q0, w – Qп, называем режимным пространством (см. Оптимизация работы ТЭЦ в условиях оптового рынка электроэнергии и мощности России

). Каждая точка этого пространства определяет некоторый режим работы паровой турбины. Режим может быть


  • технически реализуемым, если точка находится внутри оболочки, задающей регулировочный диапазон,
  • технически не реализуемым, если точка находится за пределами этой оболочки.

Если говорить о конденсационном режиме работы паровой турбины (Qп = 0, Qт = 0), то линеаризованная расходная характеристика представляет собой отрезок прямой. Если говорить о турбине Т-типа, то линеаризованная расходная характеристика представляет собой плоский многоугольник в трехмерном режимном пространстве с осями x – N, y – Qт, z – Q0, который легко визуализировать. Для турбины ПТ-типа визуализация наиболее сложная, поскольку линеаризованная расходная характеристика такой турбины представляет плоский многоугольник в четырехмерном пространстве (пояснения и примеры см. в Оптимизация работы ТЭЦ в условиях оптового рынка электроэнергии и мощности России, раздел Линеаризация расходной характеристики турбины).


1. Иллюстрация полученной линеаризованной расходной характеристики паровой турбины


Построим значения таблицы «Исходные данные (ед. мощности)» в режимном пространстве.




Рис. 3. Исходные точки расходной характеристики в режимном пространстве с осями x – N, y – Qт, z – Q0


Поскольку построить зависимость в четырехмерном пространстве мы не можем, до такого блага цивилизации еще не дошли, оперируем значениями Qп следующим образом: исключаем их (рис. 3), зафиксируем (рис. 4) (см. код построения графиков в MATLAB).


Зафиксируем значение Qп = 40 МВт и построим исходные точки и линеаризованную расходную характеристику.




Рис. 4. Исходные точки расходной характеристики (синие точки), линеаризованная расходная характеристика (зеленый плоский многоугольник)


Вернемся к полученной нами формуле линеаризованной расходной характеристики (4). Если зафиксировать Qп = 40 МВт МВт, то формула будет иметь вид


$$display$$\begin{equation} Q_0 = 2.317 \cdot N + 0.255 \cdot Q_Т + 58.714 \qquad (6) \end{equation}$$display$$


Данная модель задает плоский многоугольник в трехмерном пространстве с осями x – N, y – Qт, z – Q0 по аналогии с турбиной Т-типа (его мы и видим на рис. 4).


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


$$display$$\begin{equation} Q_0(N) = Q_э = Q_0 - Q_Т - Q_П \qquad (7) \end{equation}$$display$$


Вычли из расхода пара высокого давления Q0 расходы паров Qт, Qп и отнесли полученную разницу Q0(N) = Qэ на выработку электроэнергии. Полученную величину Q0(N) = Qэ поделили на N и перевели в ккал/кВт·ч, получив удельный расход qт брутто. Данный расчет не соответствует законам термодинамики.


Дорогие читатели, может, именно вы знаете неведомую причину? Поделитесь ею!


2. Иллюстрация регулировочного диапазона паровой турбины


Посмотрим оболочку регулировочного диапазона в режимном пространстве. Исходные точки для его построения представлены на рис. 5. Это те же самые точки, которые мы видим на рис. 3, однако теперь исключен параметр Q0.




Рис. 5. Исходные точки расходной характеристики в режимном пространстве с осями x – N, y – Qп, z – Qт


Множество точек на рис. 5 является выпуклым. Применив функцию convexhull(), мы определили точки, которые задают внешнюю оболочку этого множества.


Триангуляция Делоне (набор связанных треугольников) позволяет нам построить оболочку регулировочного диапазона. Вершины треугольников являются граничными значениями регулировочного диапазона рассматриваемой нами паровой турбины ПТ-80.




Рис. 6. Оболочка регулировочного диапазона, представленная множеством треугольников


Когда мы делали проверку некоторой точки на предмет попадания внутрь регулировочного диапазона, то мы проверяли, лежит ли эта точка внутри или снаружи полученной оболочки.


Все представленные выше графики построены средствами MATLAB (см. PT_80_linear_characteristic_curve.m).


Перспективные задачи, связанные с анализом работы паровой турбины при помощи линеаризованной расходной характеристики


Если вы делаете диплом или диссертацию, то могу предложить вам несколько задач, научную новизну которых вы легко сможете доказать всему миру. Кроме того, вы сделаете отличную и полезную работу.


Задача 1


Покажите, как изменится плоский многоугольник при изменении давления пара низкого давления Qт.


Задача 2


Покажите, как изменится плоский многоугольник при изменении давления в конденсаторе.


Задача 3


Проверьте, можно ли представить коэффициенты линеаризованной расходной характеристики в виде функций дополнительных параметров режима, а именно:


$$display$$\begin{equation} \alpha_N = f(p_{0},...); \\ \alpha_П = f(p_{П},...); \\ \alpha_Т = f(p_{Т},...); \\ \alpha_0 = f(p_{2},...). \end{equation}$$display$$



Здесь p0 — давление пара высокого давления, pп — давление пара среднего давления, pт — давление пара низкого давления, p2 — давление отработанного пара в конденсаторе, все единицы измерения кгс/см2.


Обоснуйте результат.


Ссылки


Чучуева И.А., Инкина Н.Е. Оптимизация работы ТЭЦ в условиях оптового рынка электроэнергии и мощности России // Наука и образование: научное издание МГТУ им. Н.Э. Баумана. 2015. № 8. С. 195-238.



Чучуева И.А. Вычислительные методы определения удельных расходов условного топлива ТЭЦ на отпущенную электрическую и тепловую энергию в режиме комбинированной выработки // Наука и образование: научное издание МГТУ им. Н.Э. Баумана. 2016. № 2. С. 135-165.



Линеаризованная расходная характеристика паровой турбины.


Поделиться с друзьями
-->

Комментарии (5)


  1. stDistarik
    23.05.2017 12:44

    Может кому-то будет интересно почитать про историю паровых машин.


  1. propell-ant
    23.05.2017 15:00
    +1

    0.57%… это много или мало…
    Вы не проводили анализ, с какой точностью нанесены кривые номограмм в исходных файлах?
    И что являлось исходными данными для построения этих номограмм?

    Просто мы тут доли процента считаем, а на кривых всё так гладко, плавно, на сырых данных обычно всё иначе.
    Кто-то проводил большое количество измерений, все данные подвергал статистической обработке, потом сводил воедино и рисовал номограммы.
    Это я к тому, что анализ точности вычисления каждой из точек оболочки мог быть еще более интересен.

    И, кстати, при всем глубочайшем уважении к Гауссу, вы могли и не использовать переход от табличной функции к аналитической, а ограничиться интерполяцией. По хорошему, рядом с формулой (4) нужно указывать границы допустимых значений аргументов, а то расчетчик может и забыть, откуда она взялась. С номограммами и интерполяцией выход за границы допустимых значений проверяется просто.

    А вообще статья мне понравилась, и материал и манера изложения, спасибо.


    1. mbureau
      23.05.2017 15:49

      Вы не проводили анализ, с какой точностью нанесены кривые номограмм в исходных файлах?

      Хороший вопрос. Это нормативно-техническая документация реальной станции, она по этим графикам работает и отчитывается. К сожалению, в России даже нет норматива на ошибку удельных в НТД. Я писала об этом уже Проблемы оценки экономичности работы турбины.

      Кто-то проводил большое количество измерений, все данные подвергал статистической обработке, потом сводил воедино и рисовал номограммы.
      Чтобы узнать, как это обычно делается, читайте Горнштейн В.М. и др. Методы оптимизации режимов энергосистем. Но как делали тут, я не могу ответить. Для меня эта НТД — исходные данные.
      По хорошему, рядом с формулой (4) нужно указывать границы допустимых значений аргументов, а то расчетчик может и забыть, откуда она взялась.

      Не очень внимательно читаете, я ниже пишу о вычислении регулировочного диапазона, он и есть граничные значения. А работать с этими горе-номограммами плохо: сложно, неточно (когда на глаз) и, что самое главное, поперек законов физики!


  1. Frimen3
    23.05.2017 15:54

    Спасибо за статью.
    Скажите, а центробежные воздушные компрессоры средствами MATLAB проектируются? Возможно ли это?
    Может быть у Вас был такой опыт?


    1. mbureau
      23.05.2017 15:59

      Пожалуйста!

      Скажите, а центробежные воздушные компрессоры средствами MATLAB проектируются?


      Не знаю. Это не мой профиль, такого опыта у меня нет. А вообще, MATLAB — инструмент очень мощный, думаю, что в нем проектировать (рассчитывать) можно, но ведь и понимать предметную область нужно.