From f0500e2deb91c83d9f8e052d77c8a93272ca4cdd Mon Sep 17 00:00:00 2001 From: Timothy Nunn Date: Fri, 3 Oct 2025 13:02:00 +0100 Subject: [PATCH 1/4] Scale linear system of equations in plane_stress and skip test if incompatible --- process/tf_coil.py | 17 ++++++++++++++--- tests/conftest.py | 35 +++++++++++++++++++---------------- tests/unit/test_tfcoil.py | 2 +- 3 files changed, 34 insertions(+), 20 deletions(-) diff --git a/process/tf_coil.py b/process/tf_coil.py index 9927b2e541..153c8bd102 100644 --- a/process/tf_coil.py +++ b/process/tf_coil.py @@ -5066,12 +5066,23 @@ def plane_stress(nu, rad, ey, j, nlayers, n_radial_array): # out of numba, running the code, and then copying the result # back in. This means that the linear algebra solve is not compiled and runs # as if it were written in normal Python. - # This is necessary because numba compiles against the SciPy algebra library, - # not the Numpy one. There are differences in the Scipy solvers depending on - # operating system/architecture/Scipy version that cause tests to fail. + # This is necessary because it uses mpmath solvers that cannot be numba compiled. + # We do this because scipy and numpy linear algebra solvers produce + # slightly different results depending on the system (e.g. Mac and Linux) which mean + # PROCESS can produce different results depending on where it is run! # https://github.com/ukaea/PROCESS/issues/3027 # https://github.com/scipy/scipy/issues/23639 with numba.objmode(cc="float64[:]"): + # These matrices can often end up being very ill-conditioned which can lead to numerical + # instability when solving below. + # Scaling the matrix can help reduce numerical instability by reducing the condition of matrix. + # Here, we scale aa such that the largest element on a given row is 1.0. This does not + # change the solution provided each element of a given row is scaled by the same scalar + # and the corresponding entry in bb is also scaled the same amount. + row_scale = np.max(np.abs(aa), axis=1) + aa /= np.repeat(row_scale[:, np.newaxis], aa.shape[0], axis=1) + bb /= row_scale + cc = np.linalg.solve(aa, bb) # Multiply c by (-1) (John Last, internal CCFE memorandum, 21/05/2013) diff --git a/tests/conftest.py b/tests/conftest.py index 2664e16d38..c711855233 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -92,29 +92,32 @@ def opt_params_only(request: SubRequest) -> bool: return request.config.getoption("--opt-params-only") -@pytest.fixture(scope="session", autouse=True) -def precondition(request): - """Check a user for an outdated system - Warn a user if their system is outdated and could have - regression test floating point issues. +@pytest.fixture +def skip_if_incompatible_system(): + """Skip the test using this fixture if it is detcted that their system is incompatible + and may raise errors because of floating-point rounding error. + """ + if not system_compatible(): + pytest.skip( + "This test could fail on your system due to differences caused by " + "floating-point rounding error" + ) - Exits the test suite if trying to overwrite tests - to stop inaccurate test assets being written. - e.g. "pytest --overwrite" returns True here. - :param request: request fixture to access CLI args - :type request: SubRequest +@pytest.fixture(scope="session", autouse=True) +def running_on_compatible_system_warning(): + """Check for an outdated system + + Warn the user if their system is outdated and could have test floating point issues. """ - compatible = system_compatible() - if compatible: + if system_compatible(): return warnings.warn( """ - \u001b[33m\033[1mYou are running the PROCESS test suite on an outdated system.\033[0m - This can cause floating point rounding errors in regression tests. + \u001b[33m\033[1mYou are running the PROCESS test suite on an incompatible system.\033[0m + This can cause floating point rounding errors in tests. - Please see documentation for information on running PROCESS (and tests) - using a Docker/Singularity container. + Some unit tests may be skipped! """, UserWarning, stacklevel=2, diff --git a/tests/unit/test_tfcoil.py b/tests/unit/test_tfcoil.py index 41e79f2482..007b6060df 100644 --- a/tests/unit/test_tfcoil.py +++ b/tests/unit/test_tfcoil.py @@ -5128,7 +5128,7 @@ class PlaneStressParam(NamedTuple): ), ), ) -def test_plane_stress(planestressparam, monkeypatch): +def test_plane_stress(planestressparam, skip_if_incompatible_system): """ Automatically generated Regression Unit Test for plane_stress. From 58a7e1d635a5e130e8df28df845f897f46a83555 Mon Sep 17 00:00:00 2001 From: Timothy Nunn Date: Fri, 3 Oct 2025 16:02:19 +0100 Subject: [PATCH 2/4] Add warning to docs about 'incompatible' systems --- documentation/installation/installation.md | 5 +++++ process/tf_coil.py | 7 +++---- 2 files changed, 8 insertions(+), 4 deletions(-) diff --git a/documentation/installation/installation.md b/documentation/installation/installation.md index e795bf7c98..ba2410f396 100644 --- a/documentation/installation/installation.md +++ b/documentation/installation/installation.md @@ -14,6 +14,11 @@ PROCESS is a command-line Python program. As such, it should be supported by a w wsl --install ``` +!!! Warning "Users of 'incompatible' systems" + PROCESS is a complex numerical code that performs thousands of calculations and, on certain systems (e.g. Mac and some non-Ubuntu Linux systems), we have observed that floating-point rounding error can slightly alter the results of individual models--they are not wrong, just different. However, in a code like PROCESS, this error _could_ accumulate over many iterations of the optimiser to cause non-convergence/convergence to a very different design point. To be clear, this does not happen at the moment in any of our regression tests, but is theoretically possible. + + In our test suite, we skip a couple of tests if the system is not Linux or has an `ldd` version less than `2.31` because of this floating-point rounding difference: you can look for tests using the `skip_if_incompatible_system` fixture. While we encourage users of all systems to use PROCESS, we highly recommend you check your results on our officially support systems (listed above) and Python versions. + Start by downloading the PROCESS source code using `git`: diff --git a/process/tf_coil.py b/process/tf_coil.py index 153c8bd102..861346ab51 100644 --- a/process/tf_coil.py +++ b/process/tf_coil.py @@ -5066,10 +5066,9 @@ def plane_stress(nu, rad, ey, j, nlayers, n_radial_array): # out of numba, running the code, and then copying the result # back in. This means that the linear algebra solve is not compiled and runs # as if it were written in normal Python. - # This is necessary because it uses mpmath solvers that cannot be numba compiled. - # We do this because scipy and numpy linear algebra solvers produce - # slightly different results depending on the system (e.g. Mac and Linux) which mean - # PROCESS can produce different results depending on where it is run! + # This is necessary because numba compiles against the SciPy algebra library, + # not the Numpy one. There are differences in the Scipy solvers depending on + # operating system/architecture/Scipy version that cause tests to fail. # https://github.com/ukaea/PROCESS/issues/3027 # https://github.com/scipy/scipy/issues/23639 with numba.objmode(cc="float64[:]"): From 1394fab22e0d0f1ff4108827e5df1ffd8e6b44e0 Mon Sep 17 00:00:00 2001 From: Timothy Nunn Date: Tue, 18 Nov 2025 09:11:47 +0000 Subject: [PATCH 3/4] Skip test_solv PF coil integration test on incompatible systems --- tests/integration/test_pfcoil_int.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/integration/test_pfcoil_int.py b/tests/integration/test_pfcoil_int.py index 7180864642..c8b35f64b4 100644 --- a/tests/integration/test_pfcoil_int.py +++ b/tests/integration/test_pfcoil_int.py @@ -1662,7 +1662,7 @@ def test_mtrx(pfcoil: PFCoil): assert_array_almost_equal(bvec, bvec_exp) -def test_solv(pfcoil: PFCoil): +def test_solv(pfcoil: PFCoil, skip_if_incompatible_system): """Test solv() with simple arguments. Running baseline_2019 results in 2D array args with 740 elements: unfeasible From 6ef31be17e004d959cab045395163ed61a9a5210 Mon Sep 17 00:00:00 2001 From: Timothy Nunn Date: Thu, 27 Nov 2025 10:18:17 +0000 Subject: [PATCH 4/4] Update comment and docs --- documentation/installation/installation.md | 6 +++--- process/tf_coil.py | 11 +++++++---- 2 files changed, 10 insertions(+), 7 deletions(-) diff --git a/documentation/installation/installation.md b/documentation/installation/installation.md index ba2410f396..2b377fc245 100644 --- a/documentation/installation/installation.md +++ b/documentation/installation/installation.md @@ -14,10 +14,10 @@ PROCESS is a command-line Python program. As such, it should be supported by a w wsl --install ``` -!!! Warning "Users of 'incompatible' systems" - PROCESS is a complex numerical code that performs thousands of calculations and, on certain systems (e.g. Mac and some non-Ubuntu Linux systems), we have observed that floating-point rounding error can slightly alter the results of individual models--they are not wrong, just different. However, in a code like PROCESS, this error _could_ accumulate over many iterations of the optimiser to cause non-convergence/convergence to a very different design point. To be clear, this does not happen at the moment in any of our regression tests, but is theoretically possible. +!!! Warning "Differences between systems" + PROCESS is a complex numerical code that performs thousands of calculations and, on certain systems (e.g. Mac and some non-Ubuntu Linux systems), we have observed that floating-point rounding error can slightly alter the results of individual models--they are not wrong, just different. However, in a code like PROCESS, this error _could_ accumulate over many iterations of the optimiser to cause non-convergence/convergence to a very different design point. [Issue #3972](https://github.com/ukaea/PROCESS/issues/3972) shows that these differences and shows how these small errors can accumulate to slightly change a solution. - In our test suite, we skip a couple of tests if the system is not Linux or has an `ldd` version less than `2.31` because of this floating-point rounding difference: you can look for tests using the `skip_if_incompatible_system` fixture. While we encourage users of all systems to use PROCESS, we highly recommend you check your results on our officially support systems (listed above) and Python versions. + In some sections of the code, these differences between systems are more pronounced (e.g. where operations occur on ill-conditioned matrices). Therefore, in our test suite, we skip a couple of tests if the system is not Linux or has an `ldd` version less than `2.31` because of this floating-point rounding difference: you can look for tests using the `skip_if_incompatible_system` fixture. While we encourage users of all systems to use PROCESS, we highly recommend you check your results on our officially support systems (listed above) and Python versions. Start by downloading the PROCESS source code using `git`: diff --git a/process/tf_coil.py b/process/tf_coil.py index 861346ab51..4de2820044 100644 --- a/process/tf_coil.py +++ b/process/tf_coil.py @@ -5066,9 +5066,9 @@ def plane_stress(nu, rad, ey, j, nlayers, n_radial_array): # out of numba, running the code, and then copying the result # back in. This means that the linear algebra solve is not compiled and runs # as if it were written in normal Python. - # This is necessary because numba compiles against the SciPy algebra library, - # not the Numpy one. There are differences in the Scipy solvers depending on - # operating system/architecture/Scipy version that cause tests to fail. + # This is done because numba compiles against the SciPy algebra library, + # not the Numpy one. We have observed some odd behaviour when using the SciPy library + # for this specific problem and so opt for using the numpy library instead. # https://github.com/ukaea/PROCESS/issues/3027 # https://github.com/scipy/scipy/issues/23639 with numba.objmode(cc="float64[:]"): @@ -5078,8 +5078,11 @@ def plane_stress(nu, rad, ey, j, nlayers, n_radial_array): # Here, we scale aa such that the largest element on a given row is 1.0. This does not # change the solution provided each element of a given row is scaled by the same scalar # and the corresponding entry in bb is also scaled the same amount. + # NOTE: this does not entirely solve the numerical instability and you can get above-floating point + # differences in the result of this function depending on system. row_scale = np.max(np.abs(aa), axis=1) - aa /= np.repeat(row_scale[:, np.newaxis], aa.shape[0], axis=1) + # The transpose below ensures the scale is repeated along the row, not the column + aa /= np.broadcast_to(row_scale, aa.shape).T bb /= row_scale cc = np.linalg.solve(aa, bb)