Node-SPICE: Моделирование переходных процессов в электрической сети06.06.2016 12:03

Всем привет! Сегодня я хочу рассказать об одном своем проекте, который создавался как один из инструментов получения данных для диссертации, и так как на данный момент он свою основную задачу выполнил, я хочу пустить его в GPLv3-плавание — быть может, он будет полезен кому-то еще. Однако перед тем, как отдать швартовы, я решил воспользоваться профилировщиком Intel Vtune Implifier, чтобы убедиться в том, что мой пакет имитационного моделирования древовидной сети электроснабжения оптимально расходует вычислительные ресурсы компьютера.

Под катом подробности про себя, про проект и про оптимизацию производительности (которую за полчаса удалось повысить более, чем в два раза)

Введение

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

Представьте себе обычный асинхронный электродвигатель. Тысячи их. Выглядят так:

В различных авторитетных и не очень источниках можно найти статистическую информацию о том, что до 70% вырабатываемой электроэнергии потребляется именно асинхронными электродвигателями. Не думаю что реальная цифра далека от этого значения.

Так вот, замечали когда-нибудь, что если дома запускается старый холодильник, свет моргает? Этот эффект — фликер — возникает из-за того, что при пуске электродвигатель потребляет ток в 5–7 раз больше номинального. В самые первые моменты пуска намагничивание статора отсутствует, индуктивное сопротивление минимально и сеть фактически нагружается чисто довольно малым активным сопротивлением обмотки статора. Потом, когда двигатель начинает набирать обороты, статор намагничивается, реактивное сопротивление обмотки статора увеличивается и ток уменьшается.

А теперь представьте себе электрическую сеть предприятия:

Рис. 1 — Магистральные схемы питания электроприемников:, а — с распределенными нагрузками; б — с сосредоточенными нагрузками; в — блок трансформатор — магистраль; 1 — распределительный щит подстанции; 2 — распределительный силовой пункт; 3 — электроприемник; 4 — магистраль; 5 — шинная сборка.

Это такая древовидная разветвленная электрическая сеть с множеством электроприемников. В обобщенном виде ее можно нарисовать вот так:

Рис. 2 — Обобщенная структурная схема питания электроприемников.

В схеме на рис. 2 узел Ve является точкой подключения источника питания сети (промышленная сеть переменного тока, судовой генератор, инвертор ветрогенератора и т.п.), в результате чего напряжение узла становится равным Ue. К источнику посредством активно-индуктивной питающей линии с сопротивлением Ze=Re+jLe, подключен распределительный узел V0 с напряжением U0, которое определяется как:

где I_{e-0}— ток потребляемый от узла Ve, который равен сумме токов, потребляемых нижестоящими нагрузками:

где N— число нагрузок, запитанных от данного узла. Для схемы на рис. 2 от узла Ve питаются все имеющиеся в системе узлы-нагрузки — V3 — V6. К узлу V1 подключены узлы-нагрузки V3, V4;, а к узлу V2 узлы-нагрузки V5, V6 соответственно.

Зачем создавался Node-SPICE

Если какая-то из нагрузок изменяется, изменяется ток во всей цепи до корня, следовательно, изменяется напряжение в корне, а за ним и во всех остальных узлах. И если нам нужно стабилизировать напряжение в нескольких точках цепи, то возникает задача сделать это оптимально, ибо два стабилизатора будут оказывать влияние друг на друга. Чтобы проследить это влияние на множестве вариантов необходимо произвести имитационное моделирование сети.

Схему на рис. 2 Вполне себе можно нарисовать в пакете Matlab Simulink. Но есть одна загвоздка — если схема большая, или этих схем много, то рисовать каждую схему, запускать моделирование, снимать и сохранять результаты моделирования, графики переходных процессов, чертовски муторно, и решил я, что создать свой собственный моделлер будет быстрее (фигушки) и интереснее (а вот тут я был прав).

Для того, чтобы разработка была еще более интересной и полезной, я, суровый Сишник-железячник, в качестве языка разработки решил разобраться уже наконец с C++.

Установка

Исходники представляют собой проект Visual Studio 2013 и выложены на GitHub.
Для сборки приложения необходимо скачать библиотеку линейной алгебры Eigen и указать путь к папке с библиотекой с помощью системной переменной среды $(EIGEN_DIR). Visual Studio должна будет подхватить путь к этой папке и без особых шорохов скомпилировать приложение.

Для вывода и сохранения графиков приложение использует пакет gnuplot с модулем cairo — gnuplot должен уметь сохранять изображения в формате PNG. Проверить это можно выполнив в консоли gnuplot команду set terminal png. Gnuplot не должен ругаться на неверный аргумент — последним грешил gnuplot, идущий в комплекте с octave. Путь к gnuplot должен быть указан в $(PATH).

Архитектура приложения

Приложение должно было состоять из независимых друг от друга модулей (рис 3), но что-то пошло не так:

Рис. 3 — Структурная схема программы

Основными модулями системы являются:

1. Вычислительный модуль Board. В данном модуле производится создание рабочих столов Workbench, в которых непосредственно осуществляется построение схем узлов нагрузки. Кроме того этот модуль отвечает за процесс моделирования в целом.
2. Модуль Clock. Отвечает за тактирование вычислений. На данный момент реализовано тактирование по принципу «Фиксированный шаг». Входит в состав модуля Board
3. Модуль Open. Отвечает за чтение конфигурационного файла и фалов данных в случае если таковые имеются. Входит в состав модуля Board
4. Модуль Save. Используется для сохранения результатов моделирования в файлах в сыром виде или в формате изображений. Входит в состав модуля Board
5. Модуль Plot. Отвечает за построение графиков результата.

Интерфейс программы консольный — типы и параметры электроприемников, а также конфигурация узла нагрузки описываются в конфигурационных файлах.

Команда запуска выглядит следующим образом:

node-spice.exe -f {путь к конфигурационному файлу}

Формат конфигурационного файла — текстовый, состоящий из строк вида:

command -a -b -c 1 -d 2 -e 3

где a, b, c, d, e — ключи параметров, часть из которых (a, b) имеют булев тип данных — активную или неактивную опцию или режим. Другая часть, например c, d, e — имеющие текстовое или числовое значение параметра.

Конфигурационный файл, в котором в трехфазному источнику напряжения через анализатор качества подключен электромотор и несимметричная нагрузка выглядит следующим образом:

//длительность моделирования 5 секунд. 
//Частота тактирования математического процессора 8192 Гц
setup -Off 2 -f 8192
//Создаем новый рабочий стол с именем wb0
load -t workbench -name wb0
//Создаем трехфазный источник напряжения с фазным амплитудным напряжением
//Uф = 310В и частотой 50Гц. 
//Внутреннее сопротивление источника R = 0,1 Ом, индуктивность L = 0,01Гн
load -t source -name ideal3f -f 50 -Ua 310 -R 0.1 -L 0.01
//Создаем анализатор качества электроэнергии
//анализатор будет считывать действующие значения тока(ключ -I), 
//напряжения(ключ -U), полной мощности (ключ -S), 
//коэффициента мощности (ключ -Phi),
//активной и реактивной мощности (ключи -P и -Q), 
//а также вычислять потребление электроэнергии (ключ -E)
//при этом вычисление действующего значения производится каждые 0,02c(ключ -tRMS)
//а номинальное напряжение (для регистрации провалов напряжения) 220В
load -t analyzer -name analyzer //-I -U -S -Phi -tRMS 0.02 -Unom 220 -P -Q -E
//создаем электродвигатель  4A80B4Y3, указав параметры его схемы замещения
//ключ -saveGraph активирует режим построения графиков скорости и момента
load -t acmotor -name 4A80B4Y3 -Rs 5.85 -Rr 3.0 -Ls 0.015 -Lr 0.023 -L0 0.350 -J 0.1 -p 2 -saveGraph
//создаем несимметричную нагрузку
//она будет подключена в t=1c(ключ -On 1) и отключена в t=2c(ключ -Off 2).
load -t rlc -name rl1 -On 1 -Off 2 -Ra 100 -Rb 100 -Rc 100 -La 0.01 -Lb 0.01 -Lc 0.01
//к выходу источника напряжения подключаем анализатор качества
link -output ideal3f -input  analyzer
//к выходу анализатора качества подключаем электродвигатель
link -output analyzer -input 4A80B4Y3
//к выходу анализатора подключаем несимметричную нагрузку
link -output analyzer -input rl1
//запускаем математический процессор на выполнение
solve
//по завершении строим графики и сохраняем их в виде изображений
graph

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

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

Вне зависимости от своего типа узлы имеют общий интерфейс, позволяющий создавать различные конфигурации оборудования. Добавление элементарного узла осуществляется с помощью команды load.

Общий вид команды load:

load -n {имя узла} -t {тип узла} [-ключ значение]

Существуют следующие общие для всех узлов ключи конфигурации:

Реализованые типы элементарных узлов Node: Рабочий стол -t workbench.
Предполагается, что каждый рабочий стол представляет собой некую схему и должна существовать возможность создавать вложенные схемы, т. е. вложенные рабочие столы. Эта возможность заложена в тестовой версии программы (но, естественно, не реализована :)). Уникальные ключи для рабочего стола отсутствуют. Так как может существовать несколько рабочих столов, после введения второго и более рабочего стола для узлов следует указывать, к какому рабочему столу они относятся. Если ключ -wb не указать, то элементарный узел будет размещен на последнем созданном рабочем столе.Трехфазный источник напряжения -t acsource
В текущей версии программного комплекса может быть только один источник напряжения, что несколько ограничивает возможности программы, но является достаточным для моей задачи.

Есть у меня мыслишки взять все и переписать, используя комплексное исчисление, любое число источников и приемников электроэнергии любой конфигурации, но я слезно умоляю себя если и садиться за это, то ПОСЛЕ защиты диссертации. Пока держусь.

На рисунке 5 показан процесс моделирования источника напряжения без нагрузки:

Рис. 5 — Графики тока и напряжения источника напряжения, работающего в режиме холостого хода

Анализатор качества -t analyzer
Анализатор качества потребления включается в любой участок системы и анализирует различные параметры потребления. Данный узел отвечает за построение графиков.

После проведения имитационного моделирования данный узел с помощью модуля Plot выводит требуемые графики и сохраняет их на диске в виде изображений.Асинхронный электродвигатель -t acmotor
Данный элементарный узел реализует математическую модель асинхронного электродвигателя.

На рисунке 6 показан процесс пуска асинхронного электродвигателя. В момент времени 1 с. к валу прикладывается момент 700 Н*м² и двигатель переходит в рабочий режим.

Рис. 6 — Графики частоты вращения вала двигателя, а также момента на валу и статического момента при пуске двигателя

Параллельная RLC — нагрузка -t rlc
Данный элементарный узел представляет собой параллельное соединение активного сопротивления, индуктивности и емкости. В зависимости от параметров позволяет производить моделирование следующих штатных и нештатных режимов воздействия на источник напряжения: одно- и двухфазная нагрузка, несимметричная нагрузка, короткое замыкание по фазе краткое и длительное по времени, короткое замыкание на землю по всем фазам, краткое и длительное во времени.

Моделируем кратковременное КЗ в сети:

load -t acmotor -Rs 0.02 -Rr 0.02 -Ls 0.0008 -Lr 0.0002 -Lm 0.00015 -J 3 -p 2 -Ms 700 -Tload 1
load -t rlc -Ra 0.2 -Rb 0.2 -Rc 0.2 -On 1.5 -Off 1.6

КЗ 0,1 с. Скорость не успевает упасть ниже критической, двигатель восстанавливает скорость после снятия КЗ.

load -t acmotor -Rs 0.02 -Rr 0.02 -Ls 0.0008 -Lr 0.0002 -Lm 0.00015 -J 3 -p 2 -Ms 700 -Tload 1
load -t rlc -Ra 0.2 -Rb 0.2 -Rc 0.2 -On 1.5 -Off 2

КЗ 0,5 с, двигатель успевает затормозиться и после включения момент двигателя становится меньше момента на валу и происходит аварийный останов двигателя

load -t acmotor -Rs 0.02 -Rr 0.02 -Ls 0.0008 -Lr 0.0002 -Lm 0.00015 -J 3 -p 2 -Ms 700 -Tload 1
load -t rlc -Ra 0.2 -On 1.5

Замыкание в Фазе А. Скорость практически не проседает, из-за особенностей работы асинхронного электродвигателя ему достаточно двух фаз. Вращающееся магнитное поле в зазоре принимает овальную форму и вал начинает вибрировать с частотой питающей сети.

Оптимизация кода

Вообще, как оказалось в результате, сам основной процесс моделирования написан достаточно аккуратно и по результатам моделирования каких-либо архитектурных изменений сделано не было. Но дьявол кроется в деталях.
Открываем Intel Vtune Amplifier, создаем новый проект:

Указываем путь к нашей программе и ключи запуска. Неплохо будет воспользоваться кнопками Binary/Symbol Search и Source Search и указать пути к исходному коду и бинарникам с Debud-символами — потом будет удобнее перемещаться по проекту и исходному коду.

Используем следующий конфиг:

//create new solve system:
setup -Off 10 -f 3200 //128 ticks per period
load -t workbench -name wb0
load -t acsource -name ideal3f -f 50 -Ud 220 -R 0.1 //-L 0.001
load -t motor -name motor5 -On 0.5 -Off 4 -Rs 2 -Rr 0.8 -Ls 0.00991 -Lr 0.00991 -Lm 0.008419 -J 0.5 -p 2 -Ms 50 -Tload 2 -saveGraph//15kW
load -t analyzer -name analyzer1 -tRMS 0.02 -U -I -P -E -Collect
link -output ideal3f -input analyzer1
link -output analyzer1 -input motor5
solve
graph 

Все приведенные конфиг-файлы есть в папке /doc проекта.

Начнем с самого простого basic hotspot с интервалом 1ms

И запускаем.

Top Hotspots:

Святые нейтроны… Я конечно знал что iostream работает медленно, но чтобы настолько… Это, кстати, с отключением синхронизации с

stdio ios_base::sync_with_stdio(false);

20 секунд процессорного времени из общих 35 секунд. Больше 50% времени. Это не лезет ни в какие ворота.
Больше о том, насколько медленны потоки можно прочитать здесь. Имеет смысл переписать все на бронепаровозный fprintf (). Еще меня заинтересовало что в таблице функция cout фигурирует дважды. И точно — прослойка для gnuplot создает временные файлы, а потом удаляет их. Добавим ключик -raw к node для сохранения сырых файлов графиков. Есть ключи — сохранил, нет, не сохранил.
Запускаем профилировщик. Ха!

Top Hotspots:

В лидерах по-прежнему файловый вывод, но потребляющий уже меньше 5% процессорного времени. Серьезный успех! Смотрим Bottom-Up three

Второе и третье место занимают указатели и итераторы:

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

Данный код писался как проверка концепции скользящего режима измерений. Как видно из кода, каждый новый шаг солвера сопряжен со сдвигом небольшого (64–128 символов), но все же массива. Имеет смысл воспользоваться кольцевым буфером для решения данной задачи. Тогда операция добавления нового элемента будет иметь стоимость О (1) вместо О (N).

«Зачем это надо?» скажете вы, мол, анализатор качества один в системе, моторов лучше добавь в конфиг. И окажетесь наполовину правы — моторы мы обязательно добавим, вот только анализаторов в системе может быть ровно столько сколько в системе узлов — это фишка моей диссертации такая.
Глянем заодно что там такого с GetVoltage и GetCurrent нехорошего:

Хм, как насчет воспользоваться ссылками?

Перезапускаем профилирование:

Top Hotspots:

Bottom-Up three показывает, что первый в списке опять таки наш fprintf и pango, вылезающий из-под gnuplot — в них лезть уже не будем (хотя стоило бы).

А что действительно радует, так это то что NewStep, от которого пара шагов до Solve вырвался в лидеры. Запустим моделирование на 40 секунд и посмотрим как изменится картина:

Эффект масштабируется, так что здесь нам пока делать нечего.

Подведем итог

Неплохо для получаса работы?
Добавим в систему пяток моторов:

//create new solve system:
setup -Off 10 -f 3200 //64 ticks per period
load -t workbench -name wb0
load -t acsource -name ideal3f -f 50 -Ud 220 -R 0.1 //-L 0.001
load -t motor -name motor1 -On 0.5 -Off 5 -Rs 2 -Rr 0.8 -Ls 0.00991 -Lr 0.00991 -Lm 0.008419 -J 0.5 -p 2 -Ms 50 -Tload 2 -saveGraph//15kW
load -t motor -name motor2 -On 1 -Off 6 -Rs 2 -Rr 0.8 -Ls 0.00991 -Lr 0.00991 -Lm 0.008419 -J 0.5 -p 2 -Ms 50 -Tload 2 -saveGraph//15kW
load -t motor -name motor3 -On 1.5 -Off 7 -Rs 2 -Rr 0.8 -Ls 0.00991 -Lr 0.00991 -Lm 0.008419 -J 0.5 -p 2 -Ms 50 -Tload 2 -saveGraph//15kW
load -t motor -name motor4 -On 2 -Off 8 -Rs 2 -Rr 0.8 -Ls 0.00991 -Lr 0.00991 -Lm 0.008419 -J 0.5 -p 2 -Ms 50 -Tload 2 -saveGraph//15kW
load -t motor -name motor5 -On 2.5 -Off 9 -Rs 2 -Rr 0.8 -Ls 0.00991 -Lr 0.00991 -Lm 0.008419 -J 0.5 -p 2 -Ms 50 -Tload 2 -saveGraph//15kW
load -t analyzer -name analyzer1 -tRMS 0.02 -U -I -P -E -Collect
link -output ideal3f -input analyzer1
link -output analyzer1 -input motor1
link -output analyzer1 -input motor2
link -output analyzer1 -input motor3
link -output analyzer1 -input motor4
link -output analyzer1 -input motor5
solve
graph 

Из Bottom-Up three мало что уже понятно:

Вот если заглянуть в Caller counter то можно увидеть куда деваются ресурсы. На решение матричных уравнений при расчете мат. модели электродвигателя — большую часть времени работает библиотека Eigen.

В библиотеку мы не полезем, лучше заменим моторы на rl-нагрузки. Они для меня намного важнее — можно создавать всякие разные перекосы фаз, КЗ, возмущения и прочие радости.

Так как на один тик толком считать ничего не надо, увеличим частоту тактирования солвера, да и нагрузок доведем до 10 штук.

//create new solve system:
setup -Off 10 -f 6400 //128 ticks per period
load -t workbench -name wb0
load -t acsource -name ideal3f -f 50 -Ud 220 -R 0.1 //-L 0.001
load -t rlc -name rl1 -On 1 -Off 34 -Ra 100 -Rb 100 -Rc 100 -La 0.01 -Lb 0.01 -Lc 0.01
load -t rlc -name rl2 -On 2 -Off 35 -Ra 100 -Rb 100 -Rc 100 -La 0.01 -Lb 0.01 -Lc 0.01
load -t rlc -name rl3 -On 3 -Off 36 -Ra 100 -Rb 100 -Rc 100 -La 0.01 -Lb 0.01 -Lc 0.01
load -t rlc -name rl4 -On 4 -Off 37 -Ra 100 -Rb 100 -Rc 100 -La 0.01 -Lb 0.01 -Lc 0.01
load -t rlc -name rl5 -On 5 -Off 38 -Ra 100 -Rb 100 -Rc 100 -La 0.01 -Lb 0.01 -Lc 0.01
load -t rlc -name rl11 -On 11 -Off 24 -Ra 100 -Rb 100 -Rc 100 -La 0.01 -Lb 0.01 -Lc 0.01
load -t rlc -name rl21 -On 12 -Off 25 -Ra 100 -Rb 100 -Rc 100 -La 0.01 -Lb 0.01 -Lc 0.01
load -t rlc -name rl31 -On 13 -Off 26 -Ra 100 -Rb 100 -Rc 100 -La 0.01 -Lb 0.01 -Lc 0.01
load -t rlc -name rl41 -On 14 -Off 27 -Ra 100 -Rb 100 -Rc 100 -La 0.01 -Lb 0.01 -Lc 0.01
load -t rlc -name rl51 -On 15 -Off 28 -Ra 100 -Rb 100 -Rc 100 -La 0.01 -Lb 0.01 -Lc 0.01
load -t analyzer -name analyzer1 -tRMS 0.02 -U -I -P -E -Collect
link -output ideal3f -input analyzer1
link -output analyzer1 -input rl1
link -output analyzer1 -input rl2
link -output analyzer1 -input rl3
link -output analyzer1 -input rl4
link -output analyzer1 -input rl5
link -output analyzer1 -input rl11
link -output analyzer1 -input rl21
link -output analyzer1 -input rl31
link -output analyzer1 -input rl41
link -output analyzer1 -input rl51
solve
graph 

Fprintf мы не трогаем, а вот основной виновник:

Здесь мы копируем векторы double[4] друг в друга. Как видно, копирование вектора средствами самого вектора не очень оптимально. Забабахаем-ка мы цикл — для 4-х элементов особо изгаляться не стоит:

И последний раз

Выводы:

А нету их. Я решил для себя, что негоже в OpenSource выкладывать тормозные приложения и посидел немножко с удобным и мощным инструментом профилирования. В отличие от расстановки таймстампов внутри кода, Vtune, что называется, «мордой тычет» в медленный код, намекая на то, что неплохо бы тот или иной кусок переписать.

Мое приложение, на самом деле, можно бесконечно оптимизировать — ибо костыль на костыле. Можно выкинуть Eigen и переписать Acmotor используя Boost, можно на том же Boost написать вывод графиков, можно переписать кучу мест используя векторные инструкции (здесь кстати будет кстати профилировщик intel Advisor), переписать программу используя многопоточность (TBB, OpenMP, openCL) и т.д.

Кстати вот здесь можно получить бесплатную версию Intel parallel Studio для студенческих и обучающих нужд.

© Habrahabr.ru