diff --git a/documentation/installation/installation.md b/documentation/installation/installation.md index e795bf7c98..2b377fc245 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 "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 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 9927b2e541..4de2820044 100644 --- a/process/tf_coil.py +++ b/process/tf_coil.py @@ -5066,12 +5066,25 @@ 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[:]"): + # 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. + # 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) + # 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) # 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/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 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.