Análise de sensibilidade usando PyFMI - FMU no loop for
Objetivo principal
Análise de sensibilidade de uma rede de aquecimento urbano.
Aproximação
Modelo Modelica do sistema (em Dymola) usando as bibliotecas AixLib e BuildingSystem
Modelo de exportação como co-simulação de FMU
Use SALib (biblioteca python de análise de sensibilidade) para definir as amostras (varredura de parâmetro)
Use o PyFMI para executar o modelo em um loop for no Python para todas as amostras individuais (e paralelize o loop for talvez usando o JobLib para realizar a simulação em vários processadores)
SALib para executar uma análise de sensibilidade baseada em variância (http://salib.readthedocs.io/en/latest/basics.html#an-example)
Primeiro passo
Modelo modelica simples da função Ishigami (não depende do tempo). Essa função é frequentemente usada para testar métodos de análise de sensibilidade (https://www.sfu.ca/~ssurjano/ishigami.html)
O código python (incluindo o carregamento da FMU com PyFMI e a varredura de parâmetros) funciona bem.
O problema
Após uma certa quantidade de simulação, obtemos um erro. A saída de erro nem sempre é a mesma. Às vezes temos
FMUException: Erro ao carregar o binário. Não foi possível carregar a DLL: Eine DLL-Initialisierungsroutine ist fehlgeschlagen.
Tradução: Falha na rotina de inicialização da DLL.
E às vezes temos:
FMUException: Erro ao carregar o binário. Não foi possível carregar a DLL: Für diesen Befehl ist nicht genügend Speicher verfügbar.
Tradução: Existememória insuficiente disponível para este comando.
O erro ocorredepois de cerca de 650 simulações. Isso não depende de se as simulações são executadas em blocos de loop menores que são executados um após o outro ou se um único loop for executado em todas as simulações. Ao reiniciar o console / processo python, novas simulações podem ser executadas novamente.
Ambiente de trabalho:
Windows 10, Python 2.7, PyFMI instalado usando pip (não JModelica), codificação Python no notebook Jupyther (no Mozilla Firefox)
Temos apenas conhecimento básico de python e PyFMI e estamos realmente lutando com esse erro!
Anexo
Abaixo você pode encontrar
Modelo Modelica usado para exportar FMU de co-simulação do Dymola (usando CVode)
Código Python como arquivo py
Gráfico de dispersão de saída do código python.
Também fiz um post no Fórum JModelica, onde você pode baixar os arquivos diretamente (FMU, notebook Jupyter, etc.):http://www.jmodelica.org/27925
Modelica model
model IshigamiFunction
final parameter Real a = 7;
final parameter Real b = 0.05;
parameter Real x1 = 1;
parameter Real x2 = 1;
parameter Real x3 = 1;
Real f;
equation
f = sin(x1) + a * sin(x2)^2 + b * x3^4 * sin(x1);
end IshigamiFunction;
Código Python
import numpy as np
import pylab as pl
from pyfmi import load_fmu
from SALib.sample import saltelli
from SALib.analyze import sobol
from ipywidgets import FloatProgress
from IPython.display import display
n = 100
problem = {
'num_vars': 3,
'names': ['x1', 'x2', 'x3'],
'bounds': [[-np.pi, np.pi],
[-np.pi, np.pi],
[-np.pi, np.pi]]
}
param_values = saltelli.sample(problem, n)
fmu = 'Model\IshigamiFunction\IshigamiFunction.fmu'
n_sim = param_values.shape[0]
# Progress bar
f = FloatProgress(min = 0, max = n_sim, description='Progress:')
display(f)
# Numpy array to save results
y = np.zeros([param_values.shape[0]])
x1 = np.zeros([param_values.shape[0]])
x2 = np.zeros([param_values.shape[0]])
x3 = np.zeros([param_values.shape[0]])
for i, X in enumerate(param_values):
model = load_fmu(fmu)
model.set(problem['names'], X)
res = model.simulate(final_time = 1)
y[i] = res['f'][-1]
x1[i] = res['x1'][-1]
x2[i] = res['x2'][-1]
x3[i] = res['x3'][-1]
f.value += 1
# Scatter plots
fig = pl.figure(figsize=(20, 5))
pl.clf()
pl.subplot(1,3,1)
pl.plot(x1, y, 'or')
pl.ylabel('x1')
pl.xlabel('f')
pl.subplot(1,3,2)
pl.plot(x2, y, 'ob')
pl.ylabel('x2')
pl.xlabel('f')
pl.subplot(1,3,3)
pl.plot(x3, y, 'og')
pl.ylabel('x3')
pl.xlabel('f')
pl.suptitle('Scatter plots')
pl.show()
# Sensitivity analysis
Si = sobol.analyze(problem, y, print_to_console=True)
Gráfico de saída do script python
AtualizarFiz mais alguns testes e foi isso que descobri:
Dependendo se a FMU é exportada do Dymola ou do JModelica, o comportamento é diferente:
Usando uma FMU exportada do Dymola:
Pegando oload_fmu
linha fora do loop for parece funcionarMesmo com oload_fmu
não no loop for, às vezes há falhasAdicionando uma nova linhamodel.reset()
antes de omodel.set(...)
comando parece funcionar bemOs resultados são diferentes quando simulados com ou semmodel.reset()
-> Por quê ??model.instantiate()
ao invés demodel.reset()
-> não funciona. O uso de memória no gerenciador de tarefas sobe para cerca de 350 MB e depoisFMUException: falha ao instanciar o modelo. Veja o log para possivelmente mais informações.
O arquivo de log com log_level = 4:
FMIL: module = FMILIB, log level = 4: XML specifies FMI standard version 2.0
FMIL: module = FMILIB, log level = 4: Loading 'win32' binary with 'default' platform types
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiInstantiateModel completed
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiInstantiateSlave
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiInstantiateModel completed
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiInstantiateSlave
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiSetReal: x1 = -1.76101
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiSetReal: x2 = -2.53414
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiSetReal: x3 = 0.116583
FMIL: module = Model, log level = 4: [][FMU status:OK] fmi2SetupExperiment: startTime is set to 0
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiEnterSlaveInitializationMode...
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiSetTime to 0
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiEnterSlaveInitializationMode completed
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiExitSlaveInitializationMode...
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiExitSlaveInitializationMode completed
FMIL: module = Model, log level = 4: [][FMU status:OK] fmi,GetReal: x1 = -1.76101
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiGetReal: x2 = -2.53414
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiGetReal: x3 = 0.116583
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiGetReal: a = 7
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiGetReal: b = 0.05
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiGetReal: f = 1.29856
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiGetDerivatives
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiGetDerivatives
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiGetDerivatives
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiGetDerivatives
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiGetDerivatives
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiGetDerivatives
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiGetDerivatives
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiGetDerivatives
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiSetTime to 0.002
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiDoStep
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiSetTime to 0.004
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiDoStep
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiSetTime to 0.006
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiDoStep
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiGetDerivatives
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiSetTime to 0.008
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiDoStep
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiSetTime to 0.01
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiDoStep
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiSetTime to 0.012
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiDoStep
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiSetTime to 0.014
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiDoStep
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiSetTime to 0.016
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiDoStep
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiSetTime to 0.018
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiDoStep
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiSetTime to 0.02
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiDoStep
...
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiDoStep
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiSetTime to 0.99
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiDoStep
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiSetTime to 0.992
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiDoStep
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiSetTime to 0.994
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiDoStep
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiSetTime to 0.996
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiDoStep
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiSetTime to 0.998
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiDoStep
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiSetTime to 1
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiDoStep
FMIL: module = Model, log level = 1: [][FMU status:Fatal] The license file was not found. Use the environment variable "DYMOLA_RUNTIME_LICENSE" to specify your Dymola license file.
FMIL: module = Model, log level = 1: [][FMU status:Fatal] Instantiation failed
FMIL: module = Model, log level = 4: [][FMU status:OK] fmiFreeModelInstance
Usando uma FMU exportada do JModelica:
Funciona bem, mesmo se oload_fmu
está dentro do loop for (mas mais lento)Essa experiência não corresponde ao exemplo fornecido na documentação do JModelica no capítulo 5.4.2 (http://www.jmodelica.org/api-docs/usersguide/2.1/ch05s04.html#d0e1854) onde oload_fmu
comando é fornecido dentro do loop forO comandomodel.reset()
oumodel.instatiate()
é necessário dentro do loop for (ao contrário do FMU Dymola)-> Por quê ??Minhas perguntas:
Qual é o motivo certo para executar um loop, que simula um modelo de UMF muitas vezes com parâmetros diferentes?
Qual é a diferença entre usarmodel.reset()
, model.instatiate()
ou nenhum deles?
Anexo
Aqui está um gráfico mostrando a diferença entre um loop for commodel.reset()
e sem ele.
Um FMU exportado do JModelica (não precisa de licença) pode ser baixado aqui:http://www.jmodelica.org/27925#comment-6668