python - 在For回路中使用PyFMI-FMU进行灵敏度分析

主要目标

区域供热网络的灵敏度分析。

方法

使用AixLib和BuildingSystem库

的系统的Modelica模型(在Dymola中)

将模型导出为FMU联合仿真

使用SALib(灵敏度分析python库)来定义样本(参数扫描)

使用PyFMI在Python中的for循环中为所有单个样本运行模型(并并行化for循环，也许使用JobLib在多个处理器上进行仿真)

SALib进行基于方差的敏感性分析(http://salib.readthedocs.io/en/latest/basics.html#an-example)

第一步

Ishigami函数的简单modelica模型(与时间无关)。此功能通常用于测试敏感性分析方法(https://www.sfu.ca/~ssurjano/ishigami.html)。

python代码(包括使用PyFMI加载FMU和参数扫描)可以正常工作。

问题

经过一定数量的仿真后，我们得到一个错误。错误输出看起来并不总是相同。有时候我们得到



转换:DLL初始化例程失败。

有时我们得到:



翻译:没有足够的内存可用于此命令。

在大约650次模拟运行之后，发生错误。
这不取决于仿真是在较小的循环块中执行的，然后逐个重新运行，还是单个for循环遍历所有仿真。通过重新启动python控制台/进程，可以再次运行新的模拟。

工作环境:

Windows 10，Python 2.7，使用pip(不是JModelica)安装的PyFMI，在Jupyther笔记本上的Python编码(在Mozilla Firefox上)

我们只有python和PyFMI的基本知识，并且确实为这个错误而苦苦挣扎!

附件

在下面你可以找到

用于从Dymola导出协同仿真FMU的Modelica模型(使用CVode)

Python代码作为py文件

输出python代码的散点图。

我还在JModelica论坛上发了一个帖子，您可以在其中直接下载文件(FMU，Jupyter笔记本等):
http://www.jmodelica.org/27925

Modelica模型

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;

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)

来自python脚本的输出图


更新

我进行了更多测试，这是我发现的结果:

取决于是从Dymola还是从JModelica导出FMU，其行为是不同的:

使用从Dymola导出的FMU:

从for循环中取出load_fmu行似乎很有效

即使load_fmu不在for循环中，有时
崩溃

在model.reset()命令之前添加新行model.set(...)似乎可以正常工作

在有或没有模拟的情况下，结果都不同model.reset() ->为什么？

model.instantiate()而不是model.reset()->不起作用。这
任务管理器中的内存使用量上升到大约350 MB，然后

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] fmiGetReal: 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

使用从JModelica导出的FMU:

即使load_fmu在for循环内(但速度较慢)

也能正常工作

这种经验与5.4.2章(http://www.jmodelica.org/api-docs/usersguide/2.1/ch05s04.html#d0e1854)中的JModelica文档中给出的示例不同，在

中给出了load_fmu命令

在内部需要model.reset()或model.instatiate()命令
for循环(与Dymola FMU相反)->为什么？

我的问题:

为什么执行一个循环，用不同的参数多次模拟FMU模型，这是正确的选择？

使用model.reset()，model.instatiate()或不使用它们之间有什么区别？

附件

这是显示带model.reset()和不带ojit_code的for循环之间的区别的图表。


从JModelica导出的FMU(不需要任何许可证)可以在这里下载:http://www.jmodelica.org/27925#comment-6668

最佳答案

Dymola FMU(可能与其他供应商的FMU相同)的正确方法是在for循环外调用fmi/fmi2Instantiate。如果在没有二进制导出许可证的情况下导出FMU，这些功能将分配内存并执行许可证检查。通过调用fmiResetSlave/fmi2Reset，可以将FMU重置为实例化状态，而无需分配新的内存。

fmiInstantiateSlave/fmi2Instantiate
创建一个可用于仿真的FMU实例，多个调用将创建多个实例，每个实例都需要新的内存分配和适当的删除。

fmiReset
在实例化之后和调用fmiInitializeSlave/fmi2Intialize之前将实例重置为状态。这样比较快，不需要新的动态内存分配，因此应在您的情况下使用。

此外，调用fmiFreeSalveInstance/fmi2FreeInstance时，不进行二进制导出的Dymola FMU导出的许可证检查可能会泄漏旧版Dymola中的内存。在大多数情况下，这不是问题，因为在终止FMU时通常会终止程序。通过在for循环内实例化FMU，这将变得很严重，并且您的内存将最终结束。
如果您与Dymola支持人员联系，则应该可以使用修订包。