В настоящее время я пытаюсь смоделировать оптимизацию топологии сетевого потока в OpenMDAO. Одной из проектных переменных является массовый поток из источника в сеть. Я использую сохранение массы на каждом узле сети, чтобы определить потоки в сети. Сохранение массы моделируется с неявным компонентом.
Визуализация возможной структуры графа
На добавленном изображении показан примерный график, где источник выделен оранжевым цветом, цели — синим, а соединения — серым. Всего имеется 70 узлов, соединенных 72 ребрами. Массовый поток источника и целей служит входными данными для неявного компонента, а массовые потоки 72 ребер должны быть определены с 70 остатками.
Ниже приведен минимальный пример неявного компонента:
import numpy as np
import openmdao.api as om
class MassBalance(om.ImplicitComponent):
def setup(self):
# Declare input variables
self.add_input('m_p', val=1)
self.add_input('m_c', val=np.ones(2))
# Declare output variables
self.add_output('m_i', val=np.ones(5))
# Define Residuals
self.add_residual('m_nodal_balance', shape=(8, 1), desc = "Residual for nodal mass balance")
def setup_partials(self):
self.declare_partials(of='m_nodal_balance', wrt=['m_i', 'm_p', 'm_c'])
def apply_nonlinear(self, inputs, outputs, residuals):
m_p = inputs['m_p']
m_c = inputs['m_c']
m_i = outputs['m_i']
residuals['m_nodal_balance'] = a_i.dot(m_i) + a_p.dot(m_p) + a_c.dot(m_c)
def linearize(self, inputs, outputs, partials):
# Declare Options
s = self.options['s']
# Define partials -> mass conservation in every node
partials['m_nodal_balance', 'm_i'] = a_i
partials['m_nodal_balance', 'm_p'] = a_p
partials['m_nodal_balance', 'm_c'] = a_c
a_i = np.array([
[1, 1, 0, 0, 0],
[-1, 0, -1, -1, 0],
[0, 0, 1, 0, -1],
[0, -1, 0, 1, 1]
])
a_p = np.array([[-1], [0], [0], [0]])
a_c = np.array([[0, 0], [1, 0], [0, 1], [0, 0]])
p = om.Problem()
p.model.add_subsystem('MassBalance', MassBalance(), promotes=['*'])
p.setup()
что дает следующую ошибку: Количество остатков (8) не соответствует количеству выходов (5). Если какие-либо остатки добавляются с помощью add_residuals, их общий размер должен соответствовать общему размеру выходных данных.
Следуя инструкциям, я сначала попытался смоделировать сохранение массы с помощью явного компонента и BalanceComp, но затем получил следующую ошибку: «Если вам нужен BalanceComp в OpenMDAO с различным количеством переменных состояния и остатков, вам придется создать собственный компонент».
Если я моделирую сохранение массы как отдельный неявный компонент, мне неясно, как мы можем определить выходные данные, поскольку вектор невязки должен иметь ту же длину, что и выходной вектор (даже с помощью add_residuals).
Я мог бы увеличить количество остатков с помощью 0 (или количества выходов, в зависимости от конфигурации сети), пока количество остатков не будет совпадать с количеством выходных данных, но это излишне усугубит проблему. Есть ли возможность моделирования неявного компонента, где размер остатков не соответствует количеству переменных состояния в OpenMDAO, и есть ли в руководстве примеры, которые я пропустил?
Я думаю, что основной источник путаницы носит математический, а не программный характер. Длины переменных, которые вы определили в своем примере, не соответствуют тому, какими должны быть математические расчеты.
Чтобы иметь четко определенную систему, которую можно решить с помощью линейного или нелинейного решателя, у вас должно быть такое же количество переменных состояния, как и остатков. Это, вероятно, звучит очевидно, когда формулируется абстрактно, но иногда легко запутаться, когда вы записываете уравнения. В вашем примере есть 5 выходов (я думаю, узловые массовые потоки) и 3 входа (разделенные между двумя переменными). Затем вы определяете 8 остатков вместо 5.
Однако, когда я проверяю размеры вашего фактического уравнения невязки, я получаю результат размером 4. Таким образом, ни одно из чисел в сумме не дает того, что должно.
Вместо этого я сделаю несколько предположений относительно ваших намерений и попытаюсь рассуждать исходя из этого.
Вы говорите: «Массовые потоки 72 ребер должны быть определены с 70 остатками», что нарушает правило равного количества известных и неизвестных. Нам нужно каким-то образом свести задачу только к 70 неизвестным массовым потокам. Я думаю, что 72 получается в результате добавления двух дополнительных узлов: одного в качестве источника и одного в качестве приемника. Однако источник и сток не являются неизвестными. У вас там есть известные ценности. Таким образом, эти два узла моделируются неявными функциями (т.е. с ними не связаны никакие остатки), а явными.
Эта проблема очень похожа на анализ напряжения узла в электрических цепях, который имеет существующий пример OpenMDAO. В электрическом анализе каждый узел имеет напряжение, для которого рассчитывается (за исключением узлов источника/приемника, представляющих батарею и землю). Здесь вам понадобится какая-то похожая функция потенциального типа на узлах, которая, вероятно, будет давлением.
Таким образом, если имеется 70 неизвестных узлов, будет 70 переменных состояния давления. Затем нам нужно структурировать остаток. По сути, мы хотим, чтобы чистый массовый поток через узел был равен 0 в устойчивом состоянии. Итак, 0 = mass_in + mass_out
(где я предположил, что либо mass_in
, либо mass_out
— отрицательное значение).
Таким образом, концептуально, каждое ребро должно иметь свой собственный массовый расход. Это необходимо будет вычислить на основе относительных давлений на ребре (т. е. переменных состояния в каждом узле).
Похоже, ваша матрица a_i
предназначена для сопоставления ребер, и форма этой матрицы будет зависеть от того, сколько у вас ребер. Таким образом, для n_nodes и m_edges у вас будет что-то размером (m_edges, n_nodes), которое вы можете умножить на вектор состояния, чтобы получить массовый поток через каждое ребро. Тогда это будет вектор длины m_edge, который затем можно будет умножить на другую матрицу размера (n_nodes, m_deges), которая будет обрабатывать суммирование потока массы ребра для каждого узла.
Таким образом, вы получите остаточный вектор равной длины вашему вектору состояния. Единственная сложность заключается в том, как обрабатывать граничные условия узлов источника и приемника. Здесь существуют разные подходы. Некоторым людям нравится сжимать их в основные матрицы. Другие предпочитают добавлять строки/столбцы с дополнительными «фальшивыми» переменными состояния, которые устанавливаются на входах, но которые удаляются перед передачей в вектор остатков (эти две формы называются слабыми и сильными граничными условиями).
Поэтому, прежде чем мы сможем структурировать компонент для этого, нам нужно получить лучшее определение уравнений, которые вы хотите решить. Если вы можете выразить четко определенную (т. е. количество состояний соответствует количеству остатков) систему уравнений, то я смогу лучше помочь вам сопоставить ее с компонентом. Но я считаю, что сначала нужно начать с математики.
Вы никогда не сможете решить недоопределенную или переопределенную задачу с помощью решателя Ньютона-Рапсона, по крайней мере, не делая специальных действий, чтобы превратить их в четко определенные задачи. Математика решателя Ньютона-Рапсона требует, чтобы количество неизвестных (состояний) было такое же, как и остатков. Поэтому вам нужно переформулировать свою проблему до тех пор, пока это не станет тем случаем, если вы хотите использовать этот метод решения.
Извините за мой немного запутанный вопрос. A_i описывает связь между двумя узлами графа (каждый столбец представляет собой ребро). A_p и a_c указывают расположение источника и стоков. В сочетании эти уравнения образуют 4 остатка для 5 неизвестных (к сожалению, в моем примере форма была неправильной). Разве я не мог решить эту неопределенную систему с помощью решателя Ньютона-Рафсона? Зная, что этот метод может не привести к однозначному решению. С другой стороны, может ли OpenMDAO работать с переопределенными системами? Потенциальной функцией на узлах является давление. Спасибо за это предложение.