From cfdb620a341e0381bae0526f9d97205b6ca5a87a Mon Sep 17 00:00:00 2001 From: Adrian Pop Date: Tue, 20 Oct 2026 14:12:07 +0200 Subject: [PATCH 1/3] calibrate SSP parameters against reference data (#1674) Add OMSimulator.calibration: fits the parameters of a Calibration activity of the OpenSCALING UQ layered standard to the reference data of a CSV file by bounded nonlinear least squares (scipy.optimize.least_squares). The UQ file is the one the SSD references by ssc:MetaData of type application/x-uncertainty-quantification. * each evaluation resets the instantiated model, sets the parameters, initializes and steps through the reference time points, no result file * the model starts at its own start time; a calibration window only selects the compared samples * RMSE per observed variable as ErrorMetric in the ResultSet, calibrated values as _calibrated.ssv next to the UQ file * column mapping, weights, time window, scaling and optimizer settings, not yet in the standard, from a KeyValues element of type org.openmodelica.omsimulator.calibration in the SimulationSetting * command line: OMSimulator --calibrate [ACTIVITY] [--calibrationOutput] [--uqFile] model.ssp * pip extra OMSimulator[calibration] for numpy and scipy * test calibration/HelloWorldCalibration.py fits a and x_start of HelloWorld to exp(-t), over all samples and over t >= 0.5 Co-Authored-By: Claude Opus 5.5 --- doc/UsersGuide/source/OMSimulatorPython3.rst | 82 +++ src/OMSimulatorPython/CMakeLists.txt | 1 + src/OMSimulatorPython/calibration.py | 504 ++++++++++++++++++ src/OMSimulatorPython/cli.py | 37 ++ src/pip/pyproject.toml | 5 + testsuite/tests/CMakeLists.txt | 1 + testsuite/tests/calibration/CMakeLists.txt | 3 + .../calibration/HelloWorldCalibration.py | 110 ++++ 8 files changed, 743 insertions(+) create mode 100644 src/OMSimulatorPython/calibration.py create mode 100644 testsuite/tests/calibration/CMakeLists.txt create mode 100644 testsuite/tests/calibration/HelloWorldCalibration.py diff --git a/doc/UsersGuide/source/OMSimulatorPython3.rst b/doc/UsersGuide/source/OMSimulatorPython3.rst index f7ee53e48..751cb7f40 100644 --- a/doc/UsersGuide/source/OMSimulatorPython3.rst +++ b/doc/UsersGuide/source/OMSimulatorPython3.rst @@ -109,6 +109,16 @@ refer to :ref:`OMSimulator Flags `. * - ``--validate`` - ``.fmu`` | ``.ssp`` - Only validate the file against its schema; do not simulate. + * - ``--calibrate`` + - ``[ACTIVITY]`` + - Fit the parameters of a Calibration activity instead of simulating + (``.ssp`` only), see :ref:`oms3-calibration`. + * - ``--calibrationOutput`` + - ``FILE`` + - SSP file the calibrated model is written to. + * - ``--uqFile`` + - ``PATH`` + - UQ file inside the SSP, if the SSP does not reference it. .. index:: OMSimulator; Command Line Interface Examples @@ -134,6 +144,78 @@ Examples OMSimulator model.fmu --stripRoot OMSimulator model.ssp --stripRoot + # Calibrate the parameters of activity cal_01 against its reference data + OMSimulator model.ssp --calibrate cal_01 --calibrationOutput calibrated.ssp + +.. index:: OMSimulatorPython3; Parameter Calibration + +.. _oms3-calibration : + +Parameter Calibration +--------------------- + +The ``OMSimulator.calibration`` module fits parameters of an SSP model to +reference data. The problem is described by a ``Calibration`` activity of the +OpenSCALING UQ layered standard, in a UQ file inside the SSP that the SSD +references with an ``ssc:MetaData`` element of type +``application/x-uncertainty-quantification``: + +- ``CalibrationTargets``: the parameters to fit, with initial value and bounds +- ``ObservedVariables``: the simulated variables compared to the reference data +- ``ReferenceData``: a CSV file; the first column is the time +- ``ParameterSet``: further parameter values to set +- ``SimulationSetting``: stop time and tolerance + +The fit is a bounded nonlinear least squares problem solved with +``scipy.optimize.least_squares``, so it needs numpy and scipy +(``pip install OMSimulator[calibration]``). Each evaluation resets the model, +sets the parameters, initializes it and steps through the time points of the +reference data, so the FMUs are instantiated only once. + +The standard does not yet define how observed variables map to the reference +data, their weights, the calibration time window or optimizer options. These +are read from an optional ``KeyValues`` element of type +``org.openmodelica.omsimulator.calibration`` in the ``SimulationSetting``: + +.. list-table:: Calibration options + :header-rows: 1 + :widths: 30 70 + + * - Key + - Description + * - ``startTime``, ``stopTime`` + - Calibration time window (default: the reference data up to the + ``stopTime`` of the ``SimulationSetting``). + * - ``dataName:`` + - Reference data column of an observed variable (default: the column + named like the variable, optionally followed by ``|unit``). + * - ``weight:`` + - Residual weight of an observed variable (default: 1). + * - ``scaling`` + - ``none`` (default), ``range`` or ``std``: divide the residuals of + each observed variable by the range or standard deviation of its + reference data. + * - ``relativeTolerance`` + - Optimizer tolerance on cost, step and gradient (default: 1e-8). + * - ``maxEvaluations`` + - Maximum number of model evaluations. + * - ``diffStep`` + - Relative finite difference step for the Jacobian. + +The calibrated values are written to ``_calibrated.ssv`` next to +the UQ file and the ``ResultSet`` of the activity gets one ``ErrorMetric`` +of type ``RMSE`` per observed variable, referring to that parameter set. + +.. code-block:: python + + from OMSimulator.calibration import calibrateSSP + + result = calibrateSSP("model.ssp", "cal_01", output="calibrated.ssp") + print(result.parameters, result.rmse) + +``CalibrationProblem.fromXML`` and ``Calibration`` run the same fit on an +already instantiated model, e.g., to set further values before calibrating. + Quick start Example ------------------- The following example demonstrates how to diff --git a/src/OMSimulatorPython/CMakeLists.txt b/src/OMSimulatorPython/CMakeLists.txt index 9979efd2c..a3a87a627 100644 --- a/src/OMSimulatorPython/CMakeLists.txt +++ b/src/OMSimulatorPython/CMakeLists.txt @@ -41,6 +41,7 @@ configure_file("${CMAKE_CURRENT_SOURCE_DIR}/OMSimulator.bat.in" "${CMAKE_CURRENT install(FILES "${CMAKE_CURRENT_BINARY_DIR}/__init__.py" "${CMAKE_CURRENT_SOURCE_DIR}/__main__.py" + "${CMAKE_CURRENT_SOURCE_DIR}/calibration.py" "${CMAKE_CURRENT_SOURCE_DIR}/cli.py" "${CMAKE_CURRENT_BINARY_DIR}/capi.py" "${CMAKE_CURRENT_SOURCE_DIR}/component.py" diff --git a/src/OMSimulatorPython/calibration.py b/src/OMSimulatorPython/calibration.py new file mode 100644 index 000000000..695ef2be7 --- /dev/null +++ b/src/OMSimulatorPython/calibration.py @@ -0,0 +1,504 @@ +# This file is part of OpenModelica. +# +# Copyright (c) 1998-2026, Open Source Modelica Consortium (OSMC), +# c/o Linköpings universitet, Department of Computer and Information Science, +# SE-58183 Linköping, Sweden. +# +# All rights reserved. +# +# THIS PROGRAM IS PROVIDED UNDER THE TERMS OF AGPL VERSION 3 LICENSE OR +# THIS OSMC PUBLIC LICENSE (OSMC-PL) VERSION 1.8. +# ANY USE, REPRODUCTION OR DISTRIBUTION OF THIS PROGRAM CONSTITUTES +# RECIPIENT'S ACCEPTANCE OF THE OSMC PUBLIC LICENSE OR THE GNU AGPL +# VERSION 3, ACCORDING TO RECIPIENTS CHOICE. +# +# The OpenModelica software and the OSMC (Open Source Modelica Consortium) +# Public License (OSMC-PL) are obtained from OSMC, either from the above +# address, from the URLs: +# http://www.openmodelica.org or +# https://github.com/OpenModelica/ or +# http://www.ida.liu.se/projects/OpenModelica, +# and in the OpenModelica distribution. +# +# GNU AGPL version 3 is obtained from: +# https://www.gnu.org/licenses/licenses.html#GPL +# +# This program is distributed WITHOUT ANY WARRANTY; without +# even the implied warranty of MERCHANTABILITY or FITNESS +# FOR A PARTICULAR PURPOSE, EXCEPT AS EXPRESSLY SET FORTH +# IN THE BY RECIPIENT SELECTED SUBSIDIARY LICENSE CONDITIONS OF OSMC-PL. +# +# See the full OSMC Public License conditions for more details. + +''' +Parameter calibration of SSP models against reference data. + +The calibration problem is described by a Calibration activity of the +OpenSCALING UQ layered standard (namespace UQ_NS): the parameters to fit +(CalibrationTargets), the observed variables, the reference data (a CSV +file), the simulation setting and fixed parameter values (ParameterSet). +The standard does not (yet) define how observed variables map to the +reference data, their weights, the calibration time window or optimizer +options. These are read from an optional KeyValues element of type +KEY_VALUES_TYPE in the SimulationSetting of the activity: + + startTime, stopTime calibration time window (default: the reference data up + to the stopTime of the SimulationSetting) + scaling none (default) | range | std, per observed variable + maxEvaluations maximum number of model evaluations + relativeTolerance optimizer tolerance on cost, step and gradient (default 1e-8) + diffStep relative finite difference step for the Jacobian + weight: residual weight of an observed variable (default 1) + dataName: reference data column of an observed variable + +Without dataName, the column whose header is the variable name, optionally +followed by "|unit", is used. The first column is the time. + +The fit is a bounded nonlinear least squares problem solved with +scipy.optimize.least_squares; numpy and scipy are only needed here, so they +are imported on first use. +''' + +import csv +import math +import time +from dataclasses import dataclass, field +from pathlib import Path + +from lxml import etree as ET + +from OMSimulator.capi import Capi, Status +from OMSimulator.cref import CRef + +UQ_NS = 'http://openscaling.org/UQ1/UncertaintyQuantification' +SSV_NS = 'http://ssp-standard.org/SSP1/SystemStructureParameterValues' +SSC_NS = 'http://ssp-standard.org/SSP1/SystemStructureCommon' +UQ_MIME_TYPE = 'application/x-uncertainty-quantification' +KEY_VALUES_TYPE = 'org.openmodelica.omsimulator.calibration' +SCALINGS = ('none', 'range', 'std') + +# residual of a simulation that failed, per sample; keeps the optimizer away from there +_FAILED_RESIDUAL = 1e10 + + +def _uq(tag: str) -> str: + return f'{{{UQ_NS}}}{tag}' + + +def _float(element, name: str, default=None): + value = element.get(name) + return default if value is None or value == '' else float(value) + + +@dataclass +class CalibrationTarget: + name: str + initialValue: float + lowerBound: float = -math.inf + upperBound: float = math.inf + unit: str | None = None + + +@dataclass +class ObservedVariable: + name: str + unit: str | None = None + dataName: str | None = None + weight: float = 1.0 + + +@dataclass +class CalibrationProblem: + '''A Calibration activity of a UQ layered standard file.''' + activityId: str | None + targets: list[CalibrationTarget] + observed: list[ObservedVariable] + referenceFile: Path + parameters: dict[str, float | int | bool | str] = field(default_factory=dict) + startTime: float | None = None + stopTime: float | None = None + tolerance: float | None = None + scaling: str = 'none' + maxEvaluations: int | None = None + relativeTolerance: float = 1e-8 + diffStep: float | None = None + experimentPoint: str | None = None + + @staticmethod + def fromXML(filename, activityId: str | None = None) -> 'CalibrationProblem': + '''Read a Calibration activity; activityId may be omitted if the file has only one. + + Relative file references are resolved against the directory of filename.''' + filename = Path(filename) + base = filename.parent + root = ET.parse(str(filename)).getroot() + activities = root.findall(f'.//{_uq("Calibration")}') + if activityId is not None: + activities = [a for a in activities if a.get('id') == activityId] + if not activities: + raise ValueError(f"No Calibration activity with id '{activityId}' in {filename}") + elif len(activities) != 1: + ids = ', '.join(str(a.get('id')) for a in activities) + raise ValueError(f"{filename} has {len(activities)} Calibration activities ({ids}); select one by id") + activity = activities[0] + + targets = [] + for target in activity.findall(f'{_uq("CalibrationTargets")}/{_uq("Target")}'): + name = target.get('parameterName') + initial = _float(target, 'initialValue') + lower = _float(target, 'lowerBound', -math.inf) + upper = _float(target, 'upperBound', math.inf) + if initial is None: + if math.isinf(lower) or math.isinf(upper): + raise ValueError(f"Calibration target '{name}' needs an initialValue or both bounds") + initial = 0.5 * (lower + upper) + if not lower <= initial <= upper: + raise ValueError(f"Initial value {initial} of calibration target '{name}' is outside [{lower}, {upper}]") + targets.append(CalibrationTarget(name, initial, lower, upper, target.get('unit'))) + if not targets: + raise ValueError(f"Calibration activity '{activity.get('id')}' has no CalibrationTargets") + + observed = [ObservedVariable(v.get('name'), v.get('unit')) + for v in _elements(activity, 'ObservedVariables', 'ObservedVariable', base)] + if not observed: + raise ValueError(f"Calibration activity '{activity.get('id')}' has no ObservedVariables") + + parameters = {} + for parameter in _elements(activity.find(_uq('ParameterSet')), 'Parameters', 'Parameter', base): + value = _ssvValue(parameter) + if value is not None: + parameters[parameter.get('name')] = value + + referenceFile = None + for source in activity.findall(f'{_uq("ReferenceData")}/{_uq("DataSource")}'): + if source.get('file') and (source.get('type') in (None, 'text/csv') or source.get('file').endswith('.csv')): + referenceFile = base / source.get('file') + break + if referenceFile is None: + raise ValueError(f"Calibration activity '{activity.get('id')}' has no CSV ReferenceData/DataSource") + + problem = CalibrationProblem(activity.get('id'), targets, observed, referenceFile, parameters) + + point = activity.find(f'.//{_uq("ExperimentPoints")}//{_uq("Point")}') + if point is not None: + problem.experimentPoint = point.get('id') + + setting = activity.find(_uq('SimulationSetting')) + if setting is not None: + problem.stopTime = _float(setting, 'stopTime') + problem.tolerance = _float(setting, 'tolerance') + for keyValues in setting.findall(_uq('KeyValues')): + if keyValues.get('type') == KEY_VALUES_TYPE: + for kv in keyValues.findall(_uq('KeyValue')): + problem._setOption(kv.get('key'), kv.get('value')) + return problem + + def _setOption(self, key: str, value: str): + observed = {v.name: v for v in self.observed} + if ':' in key: + option, name = key.split(':', 1) + if name not in observed: + raise ValueError(f"Calibration option '{key}' refers to an unknown observed variable '{name}'") + if option == 'weight': + observed[name].weight = float(value) + elif option == 'dataName': + observed[name].dataName = value + else: + raise ValueError(f"Unknown calibration option '{key}'") + elif key in ('startTime', 'stopTime', 'relativeTolerance', 'diffStep'): + setattr(self, key, float(value)) + elif key == 'maxEvaluations': + self.maxEvaluations = int(value) + elif key == 'scaling': + if value not in SCALINGS: + raise ValueError(f"Calibration scaling must be one of {', '.join(SCALINGS)}, got '{value}'") + self.scaling = value + else: + raise ValueError(f"Unknown calibration option '{key}'") + + +def _elements(parent, container: str, child: str, base: Path): + '''The child elements of parent/container, read from container/@source if that is given.''' + if parent is None: + return [] + result = [] + for element in parent.findall(_uq(container)): + result += element.findall(_uq(child)) + if element.get('source'): + external = ET.parse(str(base / element.get('source'))).getroot() + result += external.findall(f'.//{_uq(child)}') + return result + + +def _ssvValue(parameter): + for type_, convert in (('Real', float), ('Integer', int), ('Boolean', lambda v: v == 'true'), ('String', str)): + element = parameter.find(f'{{{SSV_NS}}}{type_}') + if element is not None and element.get('value') is not None: + return convert(element.get('value')) + return None + + +@dataclass +class ReferenceData: + time: list[float] + values: dict[str, list[float]] + + @staticmethod + def fromCSV(filename, observed: list[ObservedVariable], startTime=None, stopTime=None) -> 'ReferenceData': + '''Read the columns of the observed variables; the first column is the time. + + Rows with a repeated time point (events) keep the last value.''' + with open(filename, newline='', encoding='utf-8') as file: + rows = [row for row in csv.reader(file) if row] + header = [h.strip() for h in rows[0]] + + columns = {} + for variable in observed: + if variable.dataName is not None: + candidates = [variable.dataName] + else: + candidates = [h for h in header if h == variable.name or h.split('|')[0] == variable.name] + if not candidates or candidates[0] not in header: + raise ValueError(f"No column for observed variable '{variable.name}' in {filename}; " + f"columns are {', '.join(header)}; set it with the dataName:{variable.name} option") + columns[variable.name] = header.index(candidates[0]) + + data = {} + for row in rows[1:]: + t = float(row[0]) + if (startTime is not None and t < startTime) or (stopTime is not None and t > stopTime): + continue + data[t] = [float(row[i]) for i in columns.values()] + times = sorted(data) + values = {name: [data[t][i] for t in times] for i, name in enumerate(columns)} + if len(times) < 2: + raise ValueError(f"{filename} has less than two samples in the calibration time window") + return ReferenceData(times, values) + + +@dataclass +class CalibrationResult: + parameters: dict[str, float] + initialParameters: dict[str, float] + success: bool + message: str + evaluations: int + failedEvaluations: int + cost: float + rmse: dict[str, float] + time: list[float] + simulated: dict[str, list[float]] + reference: dict[str, list[float]] + elapsed: float + + +class Calibration: + '''Fits the targets of a CalibrationProblem on an instantiated model.''' + + def __init__(self, model, problem: CalibrationProblem, callback=None): + '''callback(evaluation: int, parameters: dict, cost: float) is called after each model evaluation.''' + self.model = model + self.problem = problem + self.callback = callback + # the model starts at its own start time; the window only selects the samples that are compared + startTime = model.getStartTime() + if problem.startTime is not None: + startTime = max(startTime, problem.startTime) + self.reference = ReferenceData.fromCSV(problem.referenceFile, problem.observed, startTime, problem.stopTime) + self.evaluations = 0 + self.failures = 0 + self._crefs = {} + for name in [t.name for t in problem.targets] + [v.name for v in problem.observed] + list(problem.parameters): + self._crefs[name] = self._resolve(name) + + def _resolve(self, name: str) -> CRef: + '''The CRef of a dotted variable name, with or without the root system name in front.''' + root = self.model.system.name + full = name if name.startswith(root + '.') else f'{root}.{name}' + if full in self.model.mappedCrefs and '.' not in full[len(root) + 1:]: + return CRef(*full.split('.')) + # the element is the longest proper prefix that is known; the rest is the FMU variable name + parts = full.split('.') + for k in range(len(parts) - 1, 1, -1): + element = '.'.join(parts[:k]) + if element in self.model.mappedCrefs: + _, status = Capi.getVariableType(f"{self.model.mappedCrefs[element]}.{'.'.join(parts[k:])}") + if status == Status.ok: + return CRef(*parts[:k], '.'.join(parts[k:])) + break + raise KeyError(f"Unknown variable '{name}' in the SSP model") + + def simulate(self, values) -> dict[str, list[float]]: + '''Simulate with the targets set to values and return the observed variables on the reference time grid.''' + model = self.model + if model.fmuInitialized: + model.reset() + for name, value in self.problem.parameters.items(): + model.setValue(self._crefs[name], value) + for target, value in zip(self.problem.targets, values): + model.setValue(self._crefs[target.name], float(value)) + model.initialize() + result = {v.name: [] for v in self.problem.observed} + for t in self.reference.time: + if t > model.getTime(): + model.stepUntil(t) + for v in self.problem.observed: + result[v.name].append(model.getValue(self._crefs[v.name])) + return result + + def _scale(self, values: list[float]) -> float: + import numpy as np + if self.problem.scaling == 'range': + scale = np.ptp(values) + elif self.problem.scaling == 'std': + scale = np.std(values) + else: + return 1.0 + return scale if scale > 0.0 else 1.0 + + def residuals(self, values): + import numpy as np + self.evaluations += 1 + n = len(self.reference.time) + try: + simulated = self.simulate(values) + r = np.concatenate([v.weight * (np.asarray(simulated[v.name]) - self.reference.values[v.name]) / self._scales[v.name] + for v in self.problem.observed]) + if not np.all(np.isfinite(r)): + raise ArithmeticError('non-finite simulation result') + except (RuntimeError, ArithmeticError): + self.failures += 1 + r = np.full(n * len(self.problem.observed), _FAILED_RESIDUAL) + if self.callback: + self.callback(self.evaluations, self._named(values), 0.5 * float(r @ r)) + return r + + def _named(self, values) -> dict[str, float]: + return {t.name: float(v) for t, v in zip(self.problem.targets, values)} + + def run(self) -> CalibrationResult: + import numpy as np + from scipy.optimize import least_squares + + problem = self.problem + start = time.time() + model = self.model + model.setResultFile('') + model.setStopTime(self.reference.time[-1]) + if problem.tolerance is not None: + model.setTolerance(problem.tolerance) + self._scales = {v.name: self._scale(self.reference.values[v.name]) for v in problem.observed} + + x0 = np.array([t.initialValue for t in problem.targets]) + bounds = ([t.lowerBound for t in problem.targets], [t.upperBound for t in problem.targets]) + tol = problem.relativeTolerance + fit = least_squares(self.residuals, x0, bounds=bounds, method='trf', x_scale='jac', + ftol=tol, xtol=tol, gtol=tol, max_nfev=problem.maxEvaluations, + diff_step=problem.diffStep) + + simulated = self.simulate(fit.x) + rmse = {v.name: float(np.sqrt(np.mean((np.asarray(simulated[v.name]) - self.reference.values[v.name]) ** 2))) + for v in problem.observed} + return CalibrationResult(self._named(fit.x), self._named(x0), bool(fit.success), fit.message, + self.evaluations, self.failures, float(fit.cost), rmse, self.reference.time, simulated, + self.reference.values, time.time() - start) + + +def writeResults(filename, problem: CalibrationProblem, result: CalibrationResult, ssvFile=None): + '''Store the result in the ResultSet of the Calibration activity in filename. + + The RMSE of each observed variable becomes an ErrorMetric. If ssvFile is + given, the calibrated values are written there as an SSP parameter set and + the Results element refers to it; a relative ssvFile is taken relative to + the directory of filename.''' + filename = Path(filename) + tree = ET.parse(str(filename)) + activity = next(a for a in tree.getroot().iter(_uq('Calibration')) + if problem.activityId is None or a.get('id') == problem.activityId) + + resultSet = activity.find(_uq('ResultSet')) + if resultSet is None: + # appended last: it takes over the tail of the previous last child, which gets the sibling indentation + last = activity[-1] if len(activity) else None + resultSet = ET.SubElement(activity, _uq('ResultSet')) + if last is not None: + resultSet.tail, last.tail = last.tail, (activity.text or '\n') + for old in resultSet.findall(_uq('Results')): + resultSet.remove(old) + results = ET.SubElement(resultSet, _uq('Results')) + if problem.experimentPoint is not None: + results.set('point', problem.experimentPoint) + if ssvFile is not None: + ssvFile = Path(ssvFile) + path = ssvFile if ssvFile.is_absolute() else filename.parent / ssvFile + writeSSV(path, problem, result) + results.set('source', ssvFile.as_posix()) + results.set('type', 'application/x-ssp-parameter-values') + for v in problem.observed: + ET.SubElement(results, _uq('ErrorMetric'), type='RMSE', + value=repr(result.rmse[v.name]), observedVariable=v.name) + + # indent only what was added, the rest of the file keeps its formatting + ET.indent(resultSet, space=' ', level=len(list(resultSet.iterancestors()))) + tree.write(str(filename), encoding='utf-8', xml_declaration=True) + + +def writeSSV(filename, problem: CalibrationProblem, result: CalibrationResult): + '''Write the calibrated values as an SSP parameter set.''' + root = ET.Element(f'{{{SSV_NS}}}ParameterSet', nsmap={'ssc': SSC_NS, 'ssv': SSV_NS}, + version='2.0', name=problem.activityId or 'calibration') + parameters = ET.SubElement(root, f'{{{SSV_NS}}}Parameters') + for target in problem.targets: + parameter = ET.SubElement(parameters, f'{{{SSV_NS}}}Parameter', name=target.name) + # no unit: in an SSV it has to refer to a unit defined in the same file + ET.SubElement(parameter, f'{{{SSV_NS}}}Real', value=repr(result.parameters[target.name])) + Path(filename).parent.mkdir(parents=True, exist_ok=True) + ET.ElementTree(root).write(str(filename), encoding='utf-8', xml_declaration=True, pretty_print=True) + + +def findUQFiles(ssp) -> list[str]: + '''The UQ files referenced by ssc:MetaData of the active variant of an SSP: of the SSD, its systems or its components.''' + files = [] + pending = [ssp.activeVariant, ssp.activeVariant.system] + while pending: + element = pending.pop(0) + for meta in getattr(element, 'metaDataResources', []): + if meta.get('type') == UQ_MIME_TYPE and meta.get('source') and meta['source'] not in files: + files.append(meta['source']) + pending += list(getattr(element, 'elements', {}).values()) + return files + + +def calibrateSSP(ssp, activityId: str | None = None, uqFile: str | None = None, output=None, callback=None) -> CalibrationResult: + '''Calibrate an SSP and write the result into its UQ file. + + ssp is an SSP object or the path to an .ssp file. uqFile is the path of the + UQ file inside the SSP; by default it is the one referenced by the + ssc:MetaData of type UQ_MIME_TYPE. The calibrated values are stored next to + the UQ file in _calibrated.ssv. If output is given, the updated + SSP is exported there.''' + from OMSimulator.ssp import SSP + if not isinstance(ssp, SSP): + ssp = SSP(str(ssp)) + + if uqFile is not None: + path = Path(ssp.temp_dir) / uqFile + else: + candidates = findUQFiles(ssp) + if not candidates: + raise ValueError(f"The SSP references no UQ file (ssc:MetaData with type '{UQ_MIME_TYPE}')") + path = Path(ssp.temp_dir) / candidates[0] + problem = CalibrationProblem.fromXML(path, activityId) + + model = ssp.instantiate() + try: + result = Calibration(model, problem, callback).run() + finally: + if model.fmuInitialized: + model.terminate() + model.delete() + + writeResults(path, problem, result, ssvFile=f'{problem.activityId or "calibration"}_calibrated.ssv') + if output is not None: + ssp.export(str(output)) + return result diff --git a/src/OMSimulatorPython/cli.py b/src/OMSimulatorPython/cli.py index 4b6f9fcfe..e4f9de74d 100644 --- a/src/OMSimulatorPython/cli.py +++ b/src/OMSimulatorPython/cli.py @@ -50,6 +50,11 @@ simulation e.g.: OMSimulator --validate test.ssp OMSimulator --validate model.fmu + +For an SSP with a UQ layered standard file (see OMSimulator.calibration), +pass --calibrate to fit the parameters of a Calibration activity to its +reference data instead of simulating, e.g.: + OMSimulator --calibrate cal_01 --calibrationOutput calibrated.ssp test.ssp ''' import argparse @@ -76,6 +81,10 @@ def _runSSP(path: Path, args: argparse.Namespace) -> None: return ssp = SSP(str(path)) + if args.calibrate is not None: + _calibrateSSP(ssp, path, args) + return + model = ssp.instantiate() if args.startTime is not None: @@ -92,6 +101,28 @@ def _runSSP(path: Path, args: argparse.Namespace) -> None: model.terminate() model.delete() +def _calibrateSSP(ssp: SSP, path: Path, args: argparse.Namespace) -> None: + from OMSimulator.calibration import calibrateSSP + + best = [float('inf')] + def progress(evaluation: int, parameters: dict, cost: float) -> None: + if cost < best[0]: + best[0] = cost + _logInfo(f"evaluation {evaluation}: cost {cost:.6g}") + + output = args.calibrationOutput or path.with_name(f"{path.stem}.calibrated.ssp") + result = calibrateSSP(ssp, args.calibrate or None, args.uqFile, output, progress) + _logInfo( + '*** Calibration Info ***\n' + f'- status : {"converged" if result.success else "not converged"} ({result.message})\n' + f'- evaluations: {result.evaluations} ({result.failedEvaluations} failed)\n' + f'- cost : {result.cost:.6g}\n' + f'- time : {result.elapsed:.2f}s\n' + '- parameters :\n' + '\n'.join(f' {name} = {value:.6g}' for name, value in result.parameters.items()) + '\n' + '- RMSE :\n' + '\n'.join(f' {name} = {value:.6g}' for name, value in result.rmse.items()) + '\n' + f'- output : {output}' + ) + def _logInfo(msg: str) -> None: '''Print a message using the "info: " style of the OMSimulator log output.''' lines = msg.split('\n') @@ -185,6 +216,9 @@ def main(argv=None) -> int: parser = argparse.ArgumentParser(prog='OMSimulator', description=__doc__) parser.add_argument('model', type=Path, help='Path to a .ssp or .fmu file to simulate') parser.add_argument('--validate', action='store_true', help='Only validate the file against its schema; do not simulate') + parser.add_argument('--calibrate', nargs='?', const='', default=None, metavar='ACTIVITY', help='Fit the parameters of a Calibration activity of the UQ file of the SSP to its reference data instead of simulating; ACTIVITY is the activity id, needed if there are several (.ssp only)') + parser.add_argument('--calibrationOutput', type=Path, default=None, help='SSP file to write the calibrated model to (default: the input name plus .calibrated.ssp)') + parser.add_argument('--uqFile', default=None, help='Path of the UQ file inside the SSP (default: the one referenced by ssc:MetaData of type application/x-uncertainty-quantification)') parser.add_argument('--version', action='version', version=f'{Capi.getVersion()}') _addOption(parser, '--addParametersToCSV', 'false', 'Export parameters to a .csv file', action=argparse.BooleanOptionalAction) _addOption(parser, '--algLoopSolver', 'kinsol', 'Specifies the loop solver method used for algebraic loops spanning multiple components', choices=['fixedpoint', 'kinsol']) @@ -235,6 +269,9 @@ def main(argv=None) -> int: if not args.model.is_file(): parser.error(f"File not found: {args.model}") + if (args.calibrate is not None or args.calibrationOutput or args.uqFile) and handler is not _runSSP: + parser.error("--calibrate, --calibrationOutput and --uqFile are only supported for .ssp files") + # A zero step size makes the solver loop forever or crash, so catch it here with a # usage error instead of after the model was already instantiated. for name in ('stepSize', 'initialStepSize', 'tolerance'): diff --git a/src/pip/pyproject.toml b/src/pip/pyproject.toml index 952876aca..1ab130cf3 100644 --- a/src/pip/pyproject.toml +++ b/src/pip/pyproject.toml @@ -21,6 +21,11 @@ dependencies = [ ] [project.optional-dependencies] +# pip install OMSimulator[calibration] +calibration = [ + "numpy", + "scipy" +] # pip install OMSimulator[gui] gui = [ "PySide6", diff --git a/testsuite/tests/CMakeLists.txt b/testsuite/tests/CMakeLists.txt index b2149716c..c11c029b3 100644 --- a/testsuite/tests/CMakeLists.txt +++ b/testsuite/tests/CMakeLists.txt @@ -1,5 +1,6 @@ add_subdirectory(api) add_subdirectory(AircraftVehicleDemonstrator) +add_subdirectory(calibration) add_subdirectory(CompositeModels) add_subdirectory(dcmotor) add_subdirectory(parker) diff --git a/testsuite/tests/calibration/CMakeLists.txt b/testsuite/tests/calibration/CMakeLists.txt new file mode 100644 index 000000000..77446e4c8 --- /dev/null +++ b/testsuite/tests/calibration/CMakeLists.txt @@ -0,0 +1,3 @@ +oms_add_test( + HelloWorldCalibration.py +) diff --git a/testsuite/tests/calibration/HelloWorldCalibration.py b/testsuite/tests/calibration/HelloWorldCalibration.py new file mode 100644 index 000000000..f371d642c --- /dev/null +++ b/testsuite/tests/calibration/HelloWorldCalibration.py @@ -0,0 +1,110 @@ +## status: correct +## linux: yes +## ucrt64: no +## win: no +## mac: no + +## Calibrates a and x_start of HelloWorld (der(x) = a*x, x(0) = x_start) +## against reference data x(t) = exp(-t), i.e., a = -1 and x_start = 1, +## described by Calibration activities of the UQ layered standard: +## cal_01 compares all samples, cal_02 only those from t = 0.5 on, while the +## model still starts at t = 0. + +import math +import zipfile + +from OMSimulator import SSP, CRef, Capi +from OMSimulator.calibration import calibrateSSP, findUQFiles + +Capi.setLogFile("HelloWorldCalibration.log") + +with open("reference.csv", "w") as f: + f.write("time,x|1\n") + for i in range(11): + f.write(f"{0.1*i},{math.exp(-0.1*i)}\n") + +with open("setup.xml", "w") as f: + f.write(''' + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +''') + +model = SSP() +model.addResource("../resources/HelloWorld.fmu", new_name="resources/HelloWorld.fmu") +model.addComponent(CRef("default", "HelloWorld"), "resources/HelloWorld.fmu") +model.getElement(CRef("default", "HelloWorld")).implementation = "me" +model.addResource("setup.xml", new_name="extra/uq/setup.xml") +model.addResource("reference.csv", new_name="extra/uq/reference.csv") +model.addMetaDataReference(None, "extra/uq/setup.xml", "quality", "application/x-uncertainty-quantification") + +for activity in ("cal_01", "cal_02"): + result = calibrateSSP(model, activity, output="HelloWorldCalibration.ssp") + print(f"{activity}: success: {result.success}, samples: {len(result.time)}, from t = {result.time[0]:.1f}") + for name, value in result.parameters.items(): + print(f"{name}: {value:.4f}") + for name, value in result.rmse.items(): + print(f"RMSE {name} < 1e-5: {value < 1e-5}") + +with zipfile.ZipFile("HelloWorldCalibration.ssp") as ssp: + print(sorted(n for n in ssp.namelist() if n.startswith("extra/"))) + with ssp.open("extra/uq/cal_02_calibrated.ssv") as ssv: + print(b'name="HelloWorld.a"' in ssv.read()) + with ssp.open("extra/uq/setup.xml") as uq: + print(b' Date: Tue, 20 Oct 2026 14:40:13 +0200 Subject: [PATCH 2/3] report invalid calibration problems as errors (#1674) A wrong activity id, a missing UQ file or an unknown observed variable made `OMSimulator --calibrate` print a Python traceback. They are user errors, so print them as a single "error:" line on stderr and exit with status 1. OMEdit shows that line in its calibration dialog. Co-Authored-By: Claude Opus 5.5 --- src/OMSimulatorPython/cli.py | 7 ++++++- 1 file changed, 6 insertions(+), 1 deletion(-) diff --git a/src/OMSimulatorPython/cli.py b/src/OMSimulatorPython/cli.py index e4f9de74d..95c6a0dd4 100644 --- a/src/OMSimulatorPython/cli.py +++ b/src/OMSimulatorPython/cli.py @@ -111,7 +111,12 @@ def progress(evaluation: int, parameters: dict, cost: float) -> None: _logInfo(f"evaluation {evaluation}: cost {cost:.6g}") output = args.calibrationOutput or path.with_name(f"{path.stem}.calibrated.ssp") - result = calibrateSSP(ssp, args.calibrate or None, args.uqFile, output, progress) + try: + result = calibrateSSP(ssp, args.calibrate or None, args.uqFile, output, progress) + except (ValueError, KeyError) as e: + # an invalid calibration problem is a user error, not a crash + print(f"error: {e.args[0] if e.args else e}", file=sys.stderr, flush=True) + sys.exit(1) _logInfo( '*** Calibration Info ***\n' f'- status : {"converged" if result.success else "not converged"} ({result.message})\n' From cba3cee94e9833edf205f6bc8fd57990ac98500d Mon Sep 17 00:00:00 2001 From: Adrian Pop Date: Tue, 20 Oct 2026 15:04:36 +0200 Subject: [PATCH 3/3] coverage: run in the build-deps OMSimulator image (#1674) calibration.HelloWorldCalibration failed in the coverage workflow with "No module named 'numpy'": the bare GitHub runner only had what the workflow installed with apt, and the calibration needs the optional OMSimulator[calibration] dependencies. Run the job in docker.openmodelica.org/build-deps:ubuntu-24.04-omsimulator instead, the image the linux64-noble Jenkins stage already uses. It has the compilers, flex, ccache, gcovr and the testsuite's Python packages (numpy and scipy included), so the apt step goes away and the coverage runs see the same environment as Jenkins. Co-Authored-By: Claude Opus 5.5 --- .github/workflows/coverage.yml | 12 +++++------- 1 file changed, 5 insertions(+), 7 deletions(-) diff --git a/.github/workflows/coverage.yml b/.github/workflows/coverage.yml index 75b3c70f2..c4cce430e 100644 --- a/.github/workflows/coverage.yml +++ b/.github/workflows/coverage.yml @@ -18,19 +18,17 @@ jobs: coverage: runs-on: ubuntu-24.04 timeout-minutes: 60 + # The image the linux64-noble Jenkins stage uses: compilers, flex for + # omc-diff, ccache, gcovr and the Python packages of the testsuite. + # Everything else OMSimulator needs is vendored in 3rdParty/. + container: + image: docker.openmodelica.org/build-deps:ubuntu-24.04-omsimulator steps: - uses: actions/checkout@3d3c42e5aac5ba805825da76410c181273ba90b1 # v7.0.1 with: submodules: recursive - # Everything else OMSimulator needs is vendored in 3rdParty/. flex is for - # omc-diff, which the testsuite compares result files with. - - name: Install dependencies - run: | - sudo apt-get update - sudo apt-get install -y cmake g++ flex python3 gcovr ccache - # Only 3rdParty benefits much: the instrumented sources are rebuilt at -O0 # anyway. It still takes the bulk of a cold build off the critical path. - name: Restore ccache