From e15f0fd082d61aeb90427fe7c297f9008235abbd Mon Sep 17 00:00:00 2001 From: lachlangrose Date: Thu, 29 Oct 2026 14:50:02 +1030 Subject: [PATCH 1/3] test: pin the recipe and inequality-constraint problems (strict xfail) Recipe: to_recipe_json / from_recipe_json do not store a feature's build parameters (nelements, interpolatortype, regularisation), so a model built with anything but the defaults rebuilds with the defaults and gives a different field (differences of 12 to 60 on a field of range about 560). The default-parameter case passes; the three others are strict xfails. Inequality: ADMM enforces the bounds with the arguments the existing builder test uses, but with its default arguments returns essentially the unconstrained solution (4 of 24 bounds violated, worst 0.056), and with solver cg the inequality data are dropped without a warning. A control (the unconstrained solution does violate bounds) and the working ADMM configuration pass; the default-ADMM case and the cg-warning expectation are strict xfails. Co-Authored-By: Claude Sonnet 5.5 --- .../tests/test_inequality_enforcement.py | 122 ++++++++++++++++++ .../modelling/test_recipe_reproduces_field.py | 64 +++++++++ 2 files changed, 186 insertions(+) create mode 100644 packages/loop_interpolation/tests/test_inequality_enforcement.py create mode 100644 tests/unit/modelling/test_recipe_reproduces_field.py diff --git a/packages/loop_interpolation/tests/test_inequality_enforcement.py b/packages/loop_interpolation/tests/test_inequality_enforcement.py new file mode 100644 index 00000000..b0c4a018 --- /dev/null +++ b/packages/loop_interpolation/tests/test_inequality_enforcement.py @@ -0,0 +1,122 @@ +"""Inequality constraints must be enforced, or at least not silently ignored. + +``InterpolatorBuilder.add_inequality_constraints`` takes rows ``x, y, z, lower, upper`` and the +documented meaning is ``lower <= f(x, y, z) <= upper``. Inequality constraints are kept apart +from the equality system and are only used by the ADMM solver. + +The scenario is a sparse model on a gaussian dome whose true field is known: one borehole (three +contact points on the 0.3, 0.5 and 0.7 horizons) and two orientations. Thirty random points whose +true value is clearly above or below the 0.5 horizon become one-sided bounds on it: +``0.5 <= f`` above, ``f <= 0.5`` below. Without the bounds this data leaves 4 of the 24 bounded +points on the wrong side of 0.5. + +What was observed (not diagnosed): + +* ADMM with ``nmajor=10, admm_weight=0.01, maxiter=50`` (the arguments the existing builder test + uses) enforces every bound, and the error against the truth falls. +* ADMM with its default arguments, or with only some of those arguments, returns essentially the + unconstrained solution: 4 of 24 bounds still violated, worst by 0.056. +* With the ``cg`` solver the inequality data have no effect, and nothing is logged. + +The last two are strict xfails. When they start to pass they XPASS, which fails the run: remove +the markers then. +""" +import logging + +import numpy as np +import pytest +from loop_common.geometry import BoundingBox +from loop_interpolation import InterpolatorBuilder +from scipy.optimize import brentq + +TOLERANCE = 1e-3 +LEVEL = 0.5 + + +def truth(xyz): + r2 = (xyz[:, 0] - 0.5) ** 2 + (xyz[:, 1] - 0.5) ** 2 + return xyz[:, 2] - 0.2 * np.exp(-r2 / (2 * 0.2**2)) + + +def truth_gradient(xyz): + r2 = (xyz[:, 0] - 0.5) ** 2 + (xyz[:, 1] - 0.5) ** 2 + e = 0.2 * np.exp(-r2 / (2 * 0.2**2)) / 0.2**2 + return np.c_[e * (xyz[:, 0] - 0.5), e * (xyz[:, 1] - 0.5), np.ones(len(xyz))] + + +def sparse_data(): + contacts = [] + for level in (0.3, 0.5, 0.7): + z = brentq(lambda z, level=level: truth(np.array([[0.5, 0.5, z]]))[0] - level, 0.0, 1.0) + contacts.append([0.5, 0.5, z, level]) + orientation_xyz = np.array([[0.4, 0.4, 0.5], [0.6, 0.6, 0.5]]) + normals = np.c_[orientation_xyz, truth_gradient(orientation_xyz)] + + rng = np.random.default_rng(9) + candidates = rng.uniform(0.1, 0.9, (30, 3)) + rows = [] + for point, value in zip(candidates, truth(candidates)): + if value > LEVEL + 0.05: + rows.append([*point, LEVEL, 10.0]) # f >= 0.5 + elif value < LEVEL - 0.05: + rows.append([*point, -10.0, LEVEL]) # f <= 0.5 + return np.array(contacts), normals, np.array(rows) + + +def builder(with_inequality=True): + contacts, normals, inequality = sparse_data() + bounding_box = BoundingBox(np.array([0.0, 0.0, 0.0]), np.array([1.0, 1.0, 1.0])) + built = InterpolatorBuilder( + "FDI", bounding_box=bounding_box, nelements=4000, buffer=0.1 + ) + built.add_value_constraints(contacts).add_normal_constraints(normals) + if with_inequality: + built.add_inequality_constraints(inequality) + return built + + +def violations(interpolator): + """How far each bounded point is outside its bounds (0 where satisfied).""" + _, _, inequality = sparse_data() + value = interpolator.evaluate_value(inequality[:, :3]) + below = np.maximum(inequality[:, 3] - value, 0.0) + above = np.maximum(value - inequality[:, 4], 0.0) + return below + above + + +def test_the_unconstrained_solution_violates_some_bounds(): + """The scenario is only meaningful if the data alone get some bounded points wrong.""" + interpolator = builder(with_inequality=False).setup_interpolator().solve().build() + assert np.count_nonzero(violations(interpolator) > TOLERANCE) >= 3 + + +def test_admm_with_the_arguments_of_the_existing_builder_test_enforces_the_bounds(): + interpolator = ( + builder() + .setup_interpolator() + .use_solver("admm", nmajor=10, admm_weight=0.01, maxiter=50) + .solve() + .build() + ) + assert violations(interpolator).max() <= TOLERANCE + + +@pytest.mark.xfail( + strict=True, + reason="ADMM with its default arguments returns essentially the unconstrained solution " + "here (4 of 24 bounds violated, worst 0.056); cause not diagnosed", +) +def test_admm_with_default_arguments_enforces_the_bounds(): + interpolator = builder().setup_interpolator().use_solver("admm").solve().build() + assert violations(interpolator).max() <= TOLERANCE + + +@pytest.mark.xfail( + strict=True, + reason="expectation, not a documented contract: with solver 'cg' the inequality data are " + "dropped without a warning", +) +def test_cg_with_inequality_constraints_warns_that_they_are_not_used(caplog): + with caplog.at_level(logging.WARNING): + builder().setup_interpolator().use_solver("cg").solve().build() + assert any("inequality" in record.getMessage().lower() for record in caplog.records) diff --git a/tests/unit/modelling/test_recipe_reproduces_field.py b/tests/unit/modelling/test_recipe_reproduces_field.py new file mode 100644 index 00000000..27ecac30 --- /dev/null +++ b/tests/unit/modelling/test_recipe_reproduces_field.py @@ -0,0 +1,64 @@ +"""A model recipe should rebuild a model that gives the same scalar field. + +``GeologicalModel.to_recipe_json`` / ``from_recipe_json`` are documented as a way to describe and +rebuild a model. The existing recipe tests check the names of features, the data and the +stratigraphic column, but not the field, and the field depends on how each feature was built: +the number of elements, the interpolator type and the regularisation weight. + +Those build parameters are not stored: a feature entry in the recipe has only ``name``, ``kind``, +``faults`` and ``regions``. A model built with the defaults rebuilds exactly, because the rebuild +uses the defaults. A model built with anything else rebuilds with the defaults and gives a +different field (on the Claudius data, differences of 12 to 62 on a field of range about 560). + +The non-default cases are strict xfails. When the recipe stores the build parameters they will +XPASS, which fails the run: remove the markers then. +""" +import numpy as np +import pytest + +from LoopStructural import GeologicalModel +from LoopStructural.datasets import load_claudius + +RECIPE_DROPS_BUILD_PARAMETERS = pytest.mark.xfail( + strict=True, + reason="to_recipe_json stores only name, kind, faults and regions for a feature, not the " + "build parameters (nelements, interpolatortype, regularisation, ...), so the rebuilt " + "model uses the defaults", +) + +CASES = [ + pytest.param({}, id="defaults"), + pytest.param( + {"nelements": 2000}, id="nelements", marks=RECIPE_DROPS_BUILD_PARAMETERS + ), + pytest.param( + {"interpolatortype": "PLI"}, id="interpolatortype", marks=RECIPE_DROPS_BUILD_PARAMETERS + ), + pytest.param( + {"regularisation": 10.0}, id="regularisation", marks=RECIPE_DROPS_BUILD_PARAMETERS + ), +] + + +def build_model(**build_parameters): + data, bb = load_claudius() + model = GeologicalModel(bb[0, :], bb[1, :]) + model.data = data + model.stratigraphic_column.clear(basement=False) + model.create_and_add_foliation("strati", **build_parameters) + model.update() + return model, bb + + +@pytest.mark.parametrize("build_parameters", CASES) +def test_rebuilt_model_gives_the_same_field(build_parameters): + model, bb = build_model(**build_parameters) + + rebuilt = GeologicalModel.from_recipe_json(model.to_recipe_json()) + rebuilt.update() + + points = np.random.default_rng(0).uniform(bb[0], bb[1], (300, 3)) + original = model.evaluate_feature_value("strati", points) + again = rebuilt.evaluate_feature_value("strati", points) + np.testing.assert_allclose(again, original, atol=1e-6) + From 150c145b283cd3409759e0977794b9c89a52752f Mon Sep 17 00:00:00 2001 From: lachlangrose Date: Thu, 29 Oct 2026 15:23:38 +1030 Subject: [PATCH 2/3] fix: store build parameters in model recipes, and refuse what they cannot rebuild A recipe rebuilt every feature with the default build parameters and as a foliation, so a model built with non-default nelements, interpolatortype or regularisation gave a different field (12 to 60 on a field of range ~560), a fault came back as a foliation, and a model with an unconformity could not be rebuilt at all. - create_and_add_feature records the parameters it was *requested* with on the feature (the element count realised by the interpolator would not reproduce the grid) - the recipe stores them per foliation (build.feature_type, build.parameters) and the reader replays them through create_and_add_feature - the stratigraphic column's group_mapping is stored - feature_recipe_kind tests subclasses first, so faults and unconformities are no longer labelled structural_frame / foliation - to_recipe_* raise LoopValueError for a feature it cannot rebuild yet (non-foliation, custom data=, a parameter JSON cannot hold, no creation record) instead of writing a wrong recipe - recipes written before this (no build entry) still load, with a warning; an older non-foliation entry is refused The three strict xfails from the recipe test now pass and their markers are removed. Co-Authored-By: Claude Sonnet 5.5 --- COMPAT.md | 1 + .../modelling/core/_model_serializer.py | 112 ++++++++++++- .../modelling/core/geological_model.py | 5 +- .../modelling/core/stratigraphic_column.py | 3 + tests/unit/modelling/test_recipe_limits.py | 148 ++++++++++++++++++ .../modelling/test_recipe_reproduces_field.py | 34 +--- 6 files changed, 268 insertions(+), 35 deletions(-) create mode 100644 tests/unit/modelling/test_recipe_limits.py diff --git a/COMPAT.md b/COMPAT.md index 11ee6262..5075c172 100644 --- a/COMPAT.md +++ b/COMPAT.md @@ -37,6 +37,7 @@ unchanged, but the new path is preferred going forward. |---|---|---| | `LoopStructural.modelling.features._feature_converters.add_fold_to_feature` (private module, imported directly by the QGIS plugin) | `GeologicalModel.add_fold_to_feature` | 2026-07-24 | | `LoopStructural.modelling.features._feature_converters.convert_feature_to_structural_frame` | `GeologicalModel.convert_feature_to_structural_frame` | 2026-07-24 | +| `GeologicalModel.to_recipe_json()` / `to_recipe_dict()` / `save_recipe()` silently wrote a lossy recipe for any model: only `name`, `kind`, `faults` for each feature (build parameters dropped, kinds mislabelled), and a model with a fault, unconformity or custom `data=` rebuilt wrongly or not at all | Each foliation entry now stores the parameters it was requested with (`build`), `group_mapping` is stored, and `from_recipe_*` replays them, so a model built with non-default `nelements`, `interpolatortype` or `regularisation` rebuilds to the same field. **A model containing a fault, unconformity or any other non-foliation feature, a feature built with its own `data=`, or a parameter JSON cannot hold, now raises `LoopValueError` from `to_recipe_*` instead of writing a recipe that rebuilds wrongly.** Older recipes (no `build`) still load, with a warning, but an older entry whose kind is not `foliation` is refused. Faults and unconformities are not yet supported | 2026-10-29 | | `GeologicalModel.create_and_add_fault(..., faults=[])` | `GeologicalModel.create_and_add_fault(..., faults=None)` (list built internally, same effective default) | 2026-07-27 | | `GeologicalModel.create_and_add_intrusion(..., intrusion_frame_parameters={}, geometric_scaling_parameters={})` | `GeologicalModel.create_and_add_intrusion(..., intrusion_frame_parameters=None, geometric_scaling_parameters=None)` (dicts built internally, same effective default) | 2026-07-27 | | `GeologicalModel.get_fault_surfaces(faults=[])` | `GeologicalModel.get_fault_surfaces(faults=None)` (list built internally, same effective default) | 2026-07-27 | diff --git a/LoopStructural/modelling/core/_model_serializer.py b/LoopStructural/modelling/core/_model_serializer.py index cbf4d112..70dbe545 100644 --- a/LoopStructural/modelling/core/_model_serializer.py +++ b/LoopStructural/modelling/core/_model_serializer.py @@ -13,6 +13,7 @@ import json import pathlib +import numpy as np import pandas as pd from ...geometry import BoundingBox @@ -24,19 +25,96 @@ logger = getLogger(__name__) +# Parameters of create_and_add_feature that are not stored with a feature's creation record: +# ``index`` is the position (the recipe keeps the order), ``faults`` are other features (the recipe +# keeps their names), and ``data`` is a data frame (see ``custom_data``). +_NOT_RECORDED = ("index", "faults", "data") + +# Feature types a recipe can rebuild from the parameters stored with it. +_SUPPORTED_FEATURE_TYPES = ("foliation",) + + +def _json_safe(value: object, where: str) -> object: + """Return ``value`` as something ``json`` can write, or raise ``LoopValueError``.""" + if value is None or isinstance(value, (bool, int, float, str)): + return value + if isinstance(value, np.generic): + return value.item() + if isinstance(value, np.ndarray): + return value.tolist() + if isinstance(value, (list, tuple)): + return [_json_safe(v, where) for v in value] + if isinstance(value, dict): + return {str(k): _json_safe(v, f"{where}[{k!r}]") for k, v in value.items()} + raise LoopValueError( + f"{where} is a {type(value).__name__}, which cannot be written to a recipe. " + "Only numbers, strings, booleans, None, lists and dictionaries of these can be stored." + ) + + class ModelSerializer: @staticmethod def feature_recipe_kind(feature): + # subclasses first: FaultSegment is a StructuralFrame and UnconformityFeature a + # GeologicalFeature + if isinstance(feature, FaultSegment): + return "fault" + if isinstance(feature, UnconformityFeature): + return "unconformity" if isinstance(feature, GeologicalFeature): return "foliation" if isinstance(feature, StructuralFrame): return "structural_frame" - if isinstance(feature, UnconformityFeature): - return "unconformity" - if isinstance(feature, FaultSegment): - return "fault" return feature.__class__.__name__.lower() + @staticmethod + def record_creation(feature: object, feature_type: str, params: dict) -> None: + """Remember how ``feature`` was created, for the recipe. + + Called by ``GeologicalModel.create_and_add_feature`` with the parameters that were + *requested* (for example ``nelements``, not the element count that was realised). + """ + feature._recipe_creation = { + "feature_type": feature_type, + "parameters": {k: v for k, v in params.items() if k not in _NOT_RECORDED}, + "custom_data": params.get("data") is not None, + } + + @staticmethod + def feature_build_entry(feature: object) -> dict: + """Return the ``build`` part of a feature's recipe entry. + + Raises + ------ + LoopValueError + if the feature cannot be rebuilt from a recipe yet + """ + kind = ModelSerializer.feature_recipe_kind(feature) + creation = getattr(feature, "_recipe_creation", None) + if creation is None: + raise LoopValueError( + f"Feature '{feature.name}' ({kind}) cannot be stored in a recipe: there is no " + "record of how it was created. Recipes store features created with " + "GeologicalModel.create_and_add_feature (or a create_and_add_* method)." + ) + feature_type = creation["feature_type"] + if kind != "foliation" or feature_type not in _SUPPORTED_FEATURE_TYPES: + raise LoopValueError( + f"Feature '{feature.name}' ({kind}, created as '{feature_type}') cannot be " + "stored in a recipe yet: recipes currently support foliations only." + ) + if creation["custom_data"]: + raise LoopValueError( + f"Feature '{feature.name}' was created with its own data=, which a recipe cannot " + "store: a recipe rebuilds a feature from the rows of the model data whose " + "feature_name matches. Put the data in model.data instead." + ) + parameters = { + key: _json_safe(value, f"parameter '{key}' of feature '{feature.name}'") + for key, value in creation["parameters"].items() + } + return {"feature_type": feature_type, "parameters": parameters} + @staticmethod def to_dict(model): result = {} @@ -68,6 +146,7 @@ def to_recipe_dict(model, data_reference=None): if getattr(fault, "name", None) ], "regions": [], + "build": ModelSerializer.feature_build_entry(feature), } recipe["model"]["features"].append(feature_entry) if data_reference is not None: @@ -132,10 +211,27 @@ def from_recipe_dict(cls, recipe): feature_name = feature_entry.get("name") if not isinstance(feature_name, str): raise TypeError("each feature entry must include a string name") - feature_data = model.data.loc[model.data["feature_name"] == feature_name].copy() - if feature_data.empty: - feature_data = None - feature = model.create_and_add_foliation(feature_name, data=feature_data) + build = feature_entry.get("build") + if build is None: + # a recipe written before build parameters were stored + kind = feature_entry.get("kind", "foliation") + if kind != "foliation": + raise LoopValueError( + f"Cannot rebuild feature '{feature_name}' of kind '{kind}' from this " + "recipe: recipes currently support foliations only." + ) + logger.warning( + f"The recipe has no build parameters for feature '{feature_name}' (it was " + "written by an older version); building it with the default parameters" + ) + build = {"feature_type": "foliation", "parameters": {}} + if not isinstance(build, dict) or not isinstance(build.get("parameters", {}), dict): + raise TypeError("a feature's build entry must be a dictionary with parameters") + feature = model.create_and_add_feature( + build.get("feature_type", "foliation"), + feature_name, + **build.get("parameters", {}), + ) if feature is None: raise ValueError(f"Could not recreate feature '{feature_name}' from recipe") feature_map[feature_name] = feature diff --git a/LoopStructural/modelling/core/geological_model.py b/LoopStructural/modelling/core/geological_model.py index cabaff0c..bf844e12 100644 --- a/LoopStructural/modelling/core/geological_model.py +++ b/LoopStructural/modelling/core/geological_model.py @@ -708,7 +708,10 @@ def create_and_add_feature(self, feature_type: str, name: str, **params): feature : BaseFeature the created geological feature, or None if it could not be built """ - return FeatureBuilderRegistry.create(feature_type, self, name, **params) + feature = FeatureBuilderRegistry.create(feature_type, self, name, **params) + if feature is not None: + ModelSerializer.record_creation(feature, feature_type, params) + return feature @public_api(tier="stable") def create_and_add_foliation( diff --git a/LoopStructural/modelling/core/stratigraphic_column.py b/LoopStructural/modelling/core/stratigraphic_column.py index a2d1971f..b1f83c04 100644 --- a/LoopStructural/modelling/core/stratigraphic_column.py +++ b/LoopStructural/modelling/core/stratigraphic_column.py @@ -641,6 +641,7 @@ def to_dict(self): """ return { "elements": [element.to_dict() for element in self.order], + "group_mapping": dict(self.group_mapping), } def update_from_dict(self, data): @@ -658,6 +659,7 @@ def update_from_dict(self, data): else: element = StratigraphicUnit.from_dict(element_data) self.add_element(element) + self.group_mapping = dict(data.get("group_mapping", {})) @classmethod def from_dict(cls, data): @@ -675,6 +677,7 @@ def from_dict(cls, data): else: element = StratigraphicUnit.from_dict(element_data) column.add_element(element) + column.group_mapping = dict(data.get("group_mapping", {})) return column def get_isovalues(self, where: str = 'bottom') -> dict[str, float]: diff --git a/tests/unit/modelling/test_recipe_limits.py b/tests/unit/modelling/test_recipe_limits.py new file mode 100644 index 00000000..741693cf --- /dev/null +++ b/tests/unit/modelling/test_recipe_limits.py @@ -0,0 +1,148 @@ +"""What a recipe can and cannot store, and that it says so instead of rebuilding something wrong.""" +import json + +import numpy as np +import pandas as pd +import pytest + +from LoopStructural import GeologicalModel +from LoopStructural.datasets import load_claudius +from LoopStructural.modelling.core._model_serializer import ModelSerializer +from LoopStructural.utils import LoopValueError + + +def claudius_model(**build_parameters): + data, bb = load_claudius() + model = GeologicalModel(bb[0, :], bb[1, :]) + model.data = data + model.stratigraphic_column.clear(basement=False) + model.create_and_add_foliation("strati", **{"nelements": 1000, **build_parameters}) + return model, data, bb + + +def fault_model(): + data = pd.DataFrame( + [ + [5, 5, 5, 0, 0.70710678, 0.0, 0.70710678, 0, "fault"], + [5, 5, 5, 0, -0.70710678, 0.0, 0.70710678, 1, "fault"], + [8, 5, 5, 0, 0, 0, 1, np.nan, "strati"], + ], + columns=["X", "Y", "Z", "val", "nx", "ny", "nz", "coord", "feature_name"], + ) + model = GeologicalModel([0, 0, 0], [10, 10, 10]) + model.data = data + model.create_and_add_fault( + "fault", 1, nelements=1000, buffer=0.5, major_axis=10, minor_axis=3, intermediate_axis=10 + ) + model.create_and_add_foliation("strati", nelements=1000, faults=[model["fault"]]) + return model + + +# ------------------------------------------------------------------ what is stored +def test_the_requested_parameters_are_stored_not_the_realised_ones(): + model, _, _ = claudius_model(nelements=2000, regularisation=10.0) + entry = json.loads(model.to_recipe_json())["model"]["features"][0] + parameters = entry["build"]["parameters"] + assert entry["build"]["feature_type"] == "foliation" + assert parameters["nelements"] == 2000 + assert parameters["regularisation"] == 10.0 + assert parameters["interpolatortype"] == "FDI" + # the realised element count is different, and storing it would not reproduce the grid + model.update() + assert model["strati"].builder.build_arguments["nelements"] != 2000 + + +def test_data_and_other_features_are_not_stored_as_parameters(): + model, _, _ = claudius_model() + parameters = json.loads(model.to_recipe_json())["model"]["features"][0]["build"]["parameters"] + assert not {"data", "faults", "index"} & set(parameters) + + +def test_group_mapping_survives_a_round_trip(): + model, _, _ = claudius_model() + model.stratigraphic_column.add_unit("unit_0", thickness=60, id=0) + model.stratigraphic_column.group_mapping["Group_0"] = "strati" + rebuilt = GeologicalModel.from_recipe_json(model.to_recipe_json()) + assert rebuilt.stratigraphic_column.group_mapping == {"Group_0": "strati"} + + +def test_a_column_without_a_group_mapping_still_loads(): + """Recipes written before group_mapping was stored have no such key.""" + model, _, _ = claudius_model() + model.stratigraphic_column.add_unit("unit_0", thickness=60, id=0) + recipe = model.to_recipe_dict() + del recipe["model"]["stratigraphic_column"]["group_mapping"] + rebuilt = GeologicalModel.from_recipe_dict(recipe) + assert rebuilt.stratigraphic_column.group_mapping == {} + + +# ------------------------------------------------------------------ kind detection +def test_a_fault_is_recognised_as_a_fault(): + model = fault_model() + assert ModelSerializer.feature_recipe_kind(model["fault"]) == "fault" + assert ModelSerializer.feature_recipe_kind(model["strati"]) == "foliation" + + +def test_an_unconformity_is_recognised_as_an_unconformity(): + model, _, _ = claudius_model() + model.add_unconformity(model["strati"], 0) + kinds = {f.name: ModelSerializer.feature_recipe_kind(f) for f in model.features} + assert kinds["__strati_unconformity"] == "unconformity" + + +# ------------------------------------------------------------------ what is refused +def test_a_model_with_a_fault_is_refused_not_rebuilt_without_the_fault(): + with pytest.raises(LoopValueError, match=r"fault.*foliations only"): + fault_model().to_recipe_json() + + +def test_a_model_with_an_unconformity_is_refused(): + model, _, _ = claudius_model() + model.add_unconformity(model["strati"], 0) + with pytest.raises(LoopValueError, match=r"__strati_unconformity.*unconformity"): + model.to_recipe_json() + + +def test_a_feature_built_with_its_own_data_is_refused(): + model, data, _ = claudius_model() + model.create_and_add_foliation( + "strati_2", data=data[data["feature_name"] == "strati"].iloc[::2], nelements=1000 + ) + with pytest.raises(LoopValueError, match=r"strati_2.*its own data"): + model.to_recipe_json() + + +def test_a_parameter_json_cannot_hold_is_refused(): + model, _, _ = claudius_model(solver=lambda a, b: None) + with pytest.raises(LoopValueError, match=r"parameter 'solver' of feature 'strati'"): + model.to_recipe_json() + + +def test_a_feature_with_no_creation_record_is_refused(): + model, _, _ = claudius_model() + del model["strati"]._recipe_creation + with pytest.raises(LoopValueError, match="no record of how it was created"): + model.to_recipe_json() + + +# ------------------------------------------------------------------ older recipes +def test_an_older_recipe_without_build_parameters_loads_with_a_warning(caplog): + import logging + + model, _, _ = claudius_model() + recipe = model.to_recipe_dict() + del recipe["model"]["features"][0]["build"] + with caplog.at_level(logging.WARNING): + rebuilt = GeologicalModel.from_recipe_dict(recipe) + assert [f.name for f in rebuilt.features] == ["strati"] + assert any("no build parameters" in record.getMessage() for record in caplog.records) + + +def test_an_older_recipe_entry_that_is_not_a_foliation_is_refused(): + model, _, _ = claudius_model() + recipe = model.to_recipe_dict() + entry = recipe["model"]["features"][0] + del entry["build"] + entry["kind"] = "structural_frame" + with pytest.raises(LoopValueError, match="foliations only"): + GeologicalModel.from_recipe_dict(recipe) diff --git a/tests/unit/modelling/test_recipe_reproduces_field.py b/tests/unit/modelling/test_recipe_reproduces_field.py index 27ecac30..a6ac2bb7 100644 --- a/tests/unit/modelling/test_recipe_reproduces_field.py +++ b/tests/unit/modelling/test_recipe_reproduces_field.py @@ -1,17 +1,9 @@ """A model recipe should rebuild a model that gives the same scalar field. -``GeologicalModel.to_recipe_json`` / ``from_recipe_json`` are documented as a way to describe and -rebuild a model. The existing recipe tests check the names of features, the data and the -stratigraphic column, but not the field, and the field depends on how each feature was built: -the number of elements, the interpolator type and the regularisation weight. - -Those build parameters are not stored: a feature entry in the recipe has only ``name``, ``kind``, -``faults`` and ``regions``. A model built with the defaults rebuilds exactly, because the rebuild -uses the defaults. A model built with anything else rebuilds with the defaults and gives a -different field (on the Claudius data, differences of 12 to 62 on a field of range about 560). - -The non-default cases are strict xfails. When the recipe stores the build parameters they will -XPASS, which fails the run: remove the markers then. +The field depends on how each feature was built: the number of elements, the interpolator type +and the regularisation weight. The recipe stores the parameters each feature was *requested* +with (not the element count that was realised, which would not reproduce the grid), and the +rebuilt model repeats the request. """ import numpy as np import pytest @@ -19,23 +11,13 @@ from LoopStructural import GeologicalModel from LoopStructural.datasets import load_claudius -RECIPE_DROPS_BUILD_PARAMETERS = pytest.mark.xfail( - strict=True, - reason="to_recipe_json stores only name, kind, faults and regions for a feature, not the " - "build parameters (nelements, interpolatortype, regularisation, ...), so the rebuilt " - "model uses the defaults", -) - CASES = [ pytest.param({}, id="defaults"), + pytest.param({"nelements": 2000}, id="nelements"), + pytest.param({"interpolatortype": "PLI"}, id="interpolatortype"), + pytest.param({"regularisation": 10.0}, id="regularisation"), pytest.param( - {"nelements": 2000}, id="nelements", marks=RECIPE_DROPS_BUILD_PARAMETERS - ), - pytest.param( - {"interpolatortype": "PLI"}, id="interpolatortype", marks=RECIPE_DROPS_BUILD_PARAMETERS - ), - pytest.param( - {"regularisation": 10.0}, id="regularisation", marks=RECIPE_DROPS_BUILD_PARAMETERS + {"nelements": 2000, "regularisation": 10.0, "damp": False}, id="several-parameters" ), ] From ee1e5c34a1a032ebe4b08c487f8bb09d956df089 Mon Sep 17 00:00:00 2001 From: lachlangrose Date: Thu, 29 Oct 2026 18:02:55 +1030 Subject: [PATCH 3/3] fix: deep-copy recipe parameters at creation; test mutation and array parameters Co-Authored-By: Claude Sonnet 5.5 --- .../modelling/core/_model_serializer.py | 6 ++++- tests/unit/modelling/test_recipe_limits.py | 27 +++++++++++++++++++ 2 files changed, 32 insertions(+), 1 deletion(-) diff --git a/LoopStructural/modelling/core/_model_serializer.py b/LoopStructural/modelling/core/_model_serializer.py index 70dbe545..5696c5ab 100644 --- a/LoopStructural/modelling/core/_model_serializer.py +++ b/LoopStructural/modelling/core/_model_serializer.py @@ -10,6 +10,7 @@ delegate to the staticmethods here. """ +import copy import json import pathlib @@ -76,7 +77,10 @@ def record_creation(feature: object, feature_type: str, params: dict) -> None: """ feature._recipe_creation = { "feature_type": feature_type, - "parameters": {k: v for k, v in params.items() if k not in _NOT_RECORDED}, + # deep copy: a list/dict/array the caller changes later must not change the recipe + "parameters": { + k: copy.deepcopy(v) for k, v in params.items() if k not in _NOT_RECORDED + }, "custom_data": params.get("data") is not None, } diff --git a/tests/unit/modelling/test_recipe_limits.py b/tests/unit/modelling/test_recipe_limits.py index 741693cf..fe9c673b 100644 --- a/tests/unit/modelling/test_recipe_limits.py +++ b/tests/unit/modelling/test_recipe_limits.py @@ -146,3 +146,30 @@ def test_an_older_recipe_entry_that_is_not_a_foliation_is_refused(): entry["kind"] = "structural_frame" with pytest.raises(LoopValueError, match="foliations only"): GeologicalModel.from_recipe_dict(recipe) + + +# ------------------------------------------------------- parameters are copied, not aliased +def test_changing_a_parameter_after_creation_does_not_change_the_recipe(): + model, _, _ = claudius_model() + values = [1, 2, 3] + weights = {"a": [1.0, 2.0]} + ModelSerializer.record_creation( + model["strati"], "foliation", {"nelements": 1000, "v": values, "w": weights} + ) + before = model.to_recipe_json() + values.append(4) + weights["a"].append(3.0) + weights["b"] = 1 + assert model.to_recipe_json() == before + parameters = json.loads(before)["model"]["features"][0]["build"]["parameters"] + assert parameters["v"] == [1, 2, 3] + assert parameters["w"] == {"a": [1.0, 2.0]} + + +def test_a_numpy_array_parameter_is_stored_as_a_list_and_copied(): + model, _, _ = claudius_model() + vector = np.array([0.0, 1.0, 0.5]) + ModelSerializer.record_creation(model["strati"], "foliation", {"nelements": 1000, "vec": vector}) + vector[0] = 99.0 + parameters = json.loads(model.to_recipe_json())["model"]["features"][0]["build"]["parameters"] + assert parameters["vec"] == [0.0, 1.0, 0.5]