使用 PyFMI 进行灵敏度分析 - for-loop 中的 FMU [英] Sensitivity Analysis using PyFMI - FMU in for-loop

查看：163
本文介绍了使用 PyFMI 进行灵敏度分析 - for-loop 中的 FMU的处理方法，对大家解决问题具有一定的参考价值，需要的朋友们下面随着小编来一起学习吧！

问题描述

主要目标

区域供热网络的敏感性分析.

方法

  1. 使用 AixLib 和 BuildingSystem 库的系统 Modelica 模型(在 Dymola 中)

  2. 将模型导出为 FMU 协同仿真

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

  4. 使用 PyFMI 在 Python 中的 for 循环中为所有单个样本运行模型(并并行化 for 循环，可能使用 JobLib 在多个处理器上执行模拟)

  5. SALib 执行基于方差的敏感性分析(

    更新

    我做了一些更多的测试，这就是我的发现:

    根据 FMU 是从 Dymola 还是从 JModelica 导出，行为会有所不同:

    使用从 Dymola 导出的 FMU:

    • 将 load_fmu 行从 for 循环中取出似乎可行
    • 即使 load_fmu 不在 for 循环中，有时也有崩溃
    • 在 model.set(...) 命令之前添加一个新行 model.reset()似乎工作正常
    • 有或没有模拟时的结果是不同的model.reset() -> 为什么?
    • model.instantiate() 而不是 model.reset() -> 不起作用.这任务管理器中的内存使用量上升到 350 MB 左右，然后<块引用>

      FMUException:无法实例化模型.查看日志以获取更多信息.

    log_level=4 的日志文件:

    FMIL: module = FMILIB, log level = 4: XML 指定 FMI 标准版本 2.0FMIL:模块 = FMILIB，日志级别 = 4:使用默认"平台类型加载win32"二进制文件FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiInstantiateModel 已完成FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiInstantiateSlaveFMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiInstantiateModel 已完成FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiInstantiateSlaveFMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiSetReal:x1 = -1.76101FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiSetReal:x2 = -2.53414FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiSetReal:x3 = 0.116583FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmi2SetupExperiment:startTime 设置为 0FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiEnterSlaveInitializationMode...FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiSetTime 为 0FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiEnterSlaveInitializationMode 已完成FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiExitSlaveInitializationMode...FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiExitSlaveInitializationMode 已完成FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiGetReal:x1 = -1.76101FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiGetReal:x2 = -2.53414FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiGetReal:x3 = 0.116583FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiGetReal:a = 7FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiGetReal:b = 0.05FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiGetReal:f = 1.29856FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiGetDerivativesFMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiGetDerivativesFMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiGetDerivativesFMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiGetDerivativesFMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiGetDerivativesFMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiGetDerivativesFMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiGetDerivativesFMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiGetDerivativesFMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiSetTime 为 0.002FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiDoStepFMIL:模块 = 型号，日志级别 = 4:[][FMU 状态:OK] fmiSetTime 为 0.004FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiDoStepFMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiSetTime 到 0.006FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiDoStepFMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiGetDerivativesFMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiSetTime 到 0.008FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiDoStepFMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiSetTime 为 0.01FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiDoStepFMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiSetTime 到 0.012FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiDoStepFMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiSetTime 到 0.014FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiDoStepFMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiSetTime 到 0.016FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiDoStepFMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiSetTime 到 0.018FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiDoStepFMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiSetTime 为 0.02FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiDoStep...FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiDoStepFMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiSetTime 到 0.99FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiDoStepFMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiSetTime 到 0.992FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiDoStepFMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiSetTime 到 0.994FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiDoStepFMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiSetTime 到 0.996FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiDoStepFMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiSetTime 到 0.998FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiDoStepFMIL:模块 = 型号，日志级别 = 4:[][FMU 状态:OK] fmiSetTime 为 1FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiDoStepFMIL:模块 = 型号，日志级别 = 1:[][FMU status:Fatal] 未找到许可证文件.使用环境变量DYMOLA_RUNTIME_LICENSE"来指定您的 Dymola 许可证文件.FMIL:模块 = 模型，日志级别 = 1:[][FMU 状态:致命] 实例化失败FMIL:模块 = 模型，日志级别 = 4:[][FMU 状态:OK] fmiFreeModelInstance

    使用从 JModelica 导出的 FMU:

    • 即使 load_fmu 在 for 循环内也能正常工作(但速度较慢)
    • 此体验与第 5.4.2 章中 JModelica 文档中给出的示例不符(

      从 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 支持，则修复包应该可用.

      Main goal

      Sensitivity analysis of a district heating network.

      Approach

      1. Modelica model of the system (in Dymola) using the AixLib and BuildingSystem libraries

      2. Export model as FMU co-simulation

      3. Use SALib (sensitivity analysis python library) to define the samples (parameter sweep)

      4. Use PyFMI to run the model in a for-loop in Python for all the individual samples (and parallelize the for loop maybe using JobLib to perfome the simulation on multiple processors)

      5. SALib to perform a variance-based sensitivity analyses (http://salib.readthedocs.io/en/latest/basics.html#an-example)

      First step

      Simple modelica model of the Ishigami function (not time dependent). This function is often used to test sensitivity analysis methods (https://www.sfu.ca/~ssurjano/ishigami.html).

      The python code (including loading the FMU with PyFMI and the parameter sweep) works fine.

      The problem

      After a certain amount of simulation we get an error. The error output looks not always the same. Sometimes we get

      FMUException: Error loading the binary. Could not load the DLL: Eine DLL-Initialisierungsroutine ist fehlgeschlagen.

      Translation: A DLL-Initilisation routine is failed.

      And sometimes we get:

      FMUException: Error loading the binary. Could not load the DLL: Für diesen Befehl ist nicht genügend Speicher verfügbar.

      Translation: There is not enough memory available for this command.

      The error occurs after around 650 simulation runs. This is not dependent of if the simulations are performed in smaller loop-blocks which are re-run one after another or if one single for loop runs through all the simulations. By restarting the python console/process, new simulations can be run again.

      Working environment:

      Windows 10, Python 2.7, PyFMI installed using pip (not JModelica), Python coding on Jupyther notebook (on Mozilla Firefox)

      We have only basic knowledge of python and PyFMI and are really struggling with this error!

      Attachment

      Below you can find

      • Modelica model used to export co-simulation FMU from Dymola (using CVode)

      • Python code as py file

      • Output scatter plot of the python code.

      I also made a post on the JModelica Forum, where you can download the files directly (FMU, Jupyter notebook, 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;

      Python code

      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)

      Output plot from python script

      Update

      I made some more tests, and this is what I found:

      Depending on if the FMU is exported from Dymola or from JModelica the behavior is different:

      Using an FMU exported from Dymola:

      • Taking the load_fmu line out of the for-loop seems to work
      • Even with the load_fmu not in the for-loop there are sometimes crashes
      • Adding a new line model.reset() before the model.set(...) command seems to work fine
      • The results are different when simulated with or without model.reset() -> Why??
      • model.instantiate() instead of model.reset() -> doesn't work. The memory usage in the task manager goes up to around 350 MB and then

        FMUException: Failed to instantiate the model. See the log for possibly more information.

      The log file with 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

      Using an FMU exported from JModelica:

      • Works fine even if the load_fmu is within the for-loop (but slower)
      • This experience does not correspond the example given within the JModelica documentation in chapter 5.4.2 (http://www.jmodelica.org/api-docs/usersguide/2.1/ch05s04.html#d0e1854) where the load_fmu command is given within the for-loop
      • The command model.reset() or model.instatiate() is required within the for-loop (contrary to Dymola FMU) -> Why??

      My questions:

      What is the right why to perform a loop, which simulates a FMU model many times with differen parameters?

      What is the difference between using model.reset(), model.instatiate() or none of them?

      Attachment

      Here's a plot showing the difference between a for-loop with model.reset() and without it.

      An FMU exported from JModelica (doesn't need any license) can be downloaded here: http://www.jmodelica.org/27925#comment-6668

      解决方案

      The right way for a Dymola FMU (and probably the same for FMUs from other vendors) would be call fmi/fmi2Instantiate outside the for loop. These functions will allocate memory and perform license check if the FMU is exported without binary export license. By calling fmiResetSlave/fmi2Reset you can reset the FMU to the instantiated state without new memory allocation.

      • fmiInstantiateSlave/fmi2Instantiate

        creates an FMU instance that can be used for simulation, multiple calls will create multiple instances, each which would require new memory allocation and proper deletion.

      • fmiReset

        resets your instance to the state after instantiate and before fmiInitializeSlave/fmi2Intialize is called. This is faster requires no new dynamic memory allocation and should be used in your case.

      In addition the license check of Dymola FMU's exported without binary export may leak memory in older Dymola versions when calling fmiFreeSalveInstance/fmi2FreeInstance. This is in most cases not a problem as you normally terminate your program when you terminate your FMU. By instantiating your FMU inside the for loop this becomes serious and your memory will finally end. Fix pack should be available if you contact Dymola support.

      这篇关于使用 PyFMI 进行灵敏度分析 - for-loop 中的 FMU的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！

查看全文
相关文章
Python最新文章
热门教程
热门工具
登录 关闭
扫码关注1秒登录
发送“验证码”获取 | 15天全站免登陆