Repository: fitbenchmarking/fitbenchmarking

Issue:
Add Jacobian check
We should check Jacobians before fitting (outside of timed loop) and use scipy if not any good


Discussion:
github-actions[bot]: ## Unit & system test results
  3 files   -   1    3 suites   - 1   12m 58s ⏱️ -53s
601 tests  -   8  588 ✅  -  10  10 💤 ± 0  3 ❌ +2 
606 runs   - 540  588 ✅  - 483  15 💤  - 59  3 ❌ +2 

For more details on these failures, see [this check](https://github.com/fitbenchmarking/fitbenchmarking/runs/35443041508).

Results for commit f283d0fb. ± Comparison against base commit d926ea32.

<details>
  <summary>This pull request <b>removes</b> 8 tests.</summary>

```
fitbenchmarking.parsing.tests.test_parsers.TestParsers ‑ test_factory[nist-/home/runner/work/fitbenchmarking/fitbenchmarking/fitbenchmarking/parsing/tests/nist/basic.dat]
fitbenchmarking.parsing.tests.test_parsers.TestParsers ‑ test_factory[nist-/home/runner/work/fitbenchmarking/fitbenchmarking/fitbenchmarking/parsing/tests/nist/start_end_x.dat]
fitbenchmarking.parsing.tests.test_parsers.TestParsers ‑ test_function_evaluation[/home/runner/work/fitbenchmarking/fitbenchmarking/fitbenchmarking/parsing/tests/nist/function_evaluations.json-nist]
fitbenchmarking.parsing.tests.test_parsers.TestParsers ‑ test_hessian_evaluation[/home/runner/work/fitbenchmarking/fitbenchmarking/fitbenchmarking/parsing/tests/nist/hessian_evaluations.json-nist]
fitbenchmarking.parsing.tests.test_parsers.TestParsers ‑ test_jacobian_evaluation[/home/runner/work/fitbenchmarking/fitbenchmarking/fitbenchmarking/parsing/tests/nist/jacobian_evaluations.json-nist]
fitbenchmarking.parsing.tests.test_parsers.TestParsers ‑ test_parsers[expected0-nist-/home/runner/work/fitbenchmarking/fitbenchmarking/fitbenchmarking/parsing/tests/nist/start_end_x.dat]
fitbenchmarking.parsing.tests.test_parsers.TestParsers ‑ test_parsers[expected1-nist-/home/runner/work/fitbenchmarking/fitbenchmarking/fitbenchmarking/parsing/tests/nist/basic.dat]
fitbenchmarking.parsing.tests.test_parsers.TestParsers ‑ test_sparsej_evaluation[/home/runner/work/fitbenchmarking/fitbenchmarking/fitbenchmarking/parsing/tests/nist/sparse_jacobian_evaluations.json-nist]
```
</details>

[test-results]:data:application/gzip;base64,H4sIAPlTgWcC/1WMyw6DIBBFf8Ww7mJEQejPNAhMQuqjQViZ/nvHitYuz7k3Z2UYBr+we9XcKrbkkE5wOZoU5omw6xQJmtI2SqgPeizZWlJCqZ96hhepGk6DJgwlugsf4xzJbJeYpz0qC/w3d1OS4hCX4pevQTuPY0gEDLlqHGDPQfG+h9ahBqt5bbhswSgtOHrvG2TvD11LWU8HAQAA

coveralls: 
[![Coverage Status](https://coveralls.io/builds/71679415/badge)](https://coveralls.io/builds/71679415)

coverage: 91.424% (-0.04%) from 91.462%
when pulling **f283d0fb2082bb04df90c921a2640a8952feee3f on 1360-add-jacobian-check**
into **d926ea32e1e444689ea75160d9be2ec3b4f580d1 on master**.


Patch:
diff --git a/docs/source/users/running.rst b/docs/source/users/running.rst
index bdb319d9..ed6f6530 100644
--- a/docs/source/users/running.rst
+++ b/docs/source/users/running.rst
@@ -192,6 +192,16 @@ They can be changed using the arguments in the table below.
      - ``--dont_run_dash``
      - | Use this option if you have decided
        | not to run the Dash app.
+   * - *Check jacobian*
+     -
+     - ``--check_jacobian``
+     - | Use this option to check the jacobian 
+       | against a finite difference approximation.
+   * - *Don't check jacobian*
+     -
+     - ``--dont_check_jacobian``
+     - | Use this option if you have decided
+       | not to check the jacobian.
 
 **For example, to change the results directory:**
 
diff --git a/fitbenchmarking/cli/main.py b/fitbenchmarking/cli/main.py
index 81a08b33..905b00ba 100755
--- a/fitbenchmarking/cli/main.py
+++ b/fitbenchmarking/cli/main.py
@@ -178,8 +178,8 @@ of the Fitbenchmarking docs. """
         help="Set the port for Dash.",
     )
 
-    group1 = parser.add_mutually_exclusive_group()
-    group1.add_argument(
+    make_plots_group = parser.add_mutually_exclusive_group()
+    make_plots_group.add_argument(
         "--make_plots",
         action="store_true",
         help=(
@@ -187,7 +187,7 @@ of the Fitbenchmarking docs. """
             "during runtime."
         ),
     )
-    group1.add_argument(
+    make_plots_group.add_argument(
         "--dont_make_plots",
         action="store_true",
         help=(
@@ -196,8 +196,8 @@ of the Fitbenchmarking docs. """
         ),
     )
 
-    group2 = parser.add_mutually_exclusive_group()
-    group2.add_argument(
+    results_browser_group = parser.add_mutually_exclusive_group()
+    results_browser_group.add_argument(
         "--results_browser",
         action="store_true",
         help=(
@@ -206,7 +206,7 @@ of the Fitbenchmarking docs. """
             "of a fit benchmark."
         ),
     )
-    group2.add_argument(
+    results_browser_group.add_argument(
         "--no_results_browser",
         action="store_true",
         help=(
@@ -216,8 +216,8 @@ of the Fitbenchmarking docs. """
         ),
     )
 
-    group3 = parser.add_mutually_exclusive_group()
-    group3.add_argument(
+    pbar_group = parser.add_mutually_exclusive_group()
+    pbar_group.add_argument(
         "--pbar",
         action="store_true",
         help=(
@@ -225,7 +225,7 @@ of the Fitbenchmarking docs. """
             "see the progress bar during runtime."
         ),
     )
-    group3.add_argument(
+    pbar_group.add_argument(
         "--no_pbar",
         action="store_true",
         help=(
@@ -268,8 +268,8 @@ of the Fitbenchmarking docs. """
         help="Specify the file path to write the logs to.",
     )
 
-    group2 = parser.add_mutually_exclusive_group()
-    group2.add_argument(
+    append_log_group = parser.add_mutually_exclusive_group()
+    append_log_group.add_argument(
         "--append_log",
         action="store_true",
         help=(
@@ -279,7 +279,7 @@ of the Fitbenchmarking docs. """
             "subsequent run."
         ),
     )
-    group2.add_argument(
+    append_log_group.add_argument(
         "--overwrite_log",
         action="store_true",
         help=(
@@ -317,8 +317,8 @@ of the Fitbenchmarking docs. """
         ),
     )
 
-    group4 = parser.add_mutually_exclusive_group()
-    group4.add_argument(
+    run_dash_group = parser.add_mutually_exclusive_group()
+    run_dash_group.add_argument(
         "--run_dash",
         action="store_true",
         help=(
@@ -326,7 +326,7 @@ of the Fitbenchmarking docs. """
             "interactive plots."
         ),
     )
-    group4.add_argument(
+    run_dash_group.add_argument(
         "--dont_run_dash",
         action="store_true",
         help=(
@@ -334,6 +334,23 @@ of the Fitbenchmarking docs. """
             "run dash for interactive plots."
         ),
     )
+
+    jacobian_check_group = parser.add_mutually_exclusive_group()
+    jacobian_check_group.add_argument(
+        "--check_jacobian",
+        action="store_true",
+        help=(
+            "Use this option to check the jacobian against"
+            "a finite difference approximation."
+        ),
+    )
+    jacobian_check_group.add_argument(
+        "--dont_check_jacobian",
+        action="store_true",
+        help=(
+            "Use this option if you have decided not to check the jacobian."
+        ),
+    )
     return parser
 
 
@@ -547,6 +564,13 @@ def main():
     elif args.dont_run_dash:
         options_dictionary["run_dash"] = False
 
+    # Check if check_jacobian in options.py should be overridden, and if so,
+    # add to options_dictionary
+    if args.check_jacobian:
+        options_dictionary["check_jacobian"] = True
+    elif args.dont_check_jacobian:
+        options_dictionary["check_jacobian"] = False
+
     # Check if benchmark in options.py should be overridden, and if so,
     # add to options_dictionary
     if args.pbar:
diff --git a/fitbenchmarking/core/fitting_benchmarking.py b/fitbenchmarking/core/fitting_benchmarking.py
index 3fd867b3..87a5b962 100644
--- a/fitbenchmarking/core/fitting_benchmarking.py
+++ b/fitbenchmarking/core/fitting_benchmarking.py
@@ -9,6 +9,7 @@ import timeit
 
 import numpy as np
 from codecarbon import EmissionsTracker
+from scipy.optimize import approx_fprime, check_grad
 from tqdm import tqdm, trange
 from tqdm.contrib.logging import logging_redirect_tqdm
 
@@ -409,6 +410,21 @@ class Fit:
                             jacobian.name() if jacobian.name() else "default",
                         )
 
+                    if (
+                        self._options.check_jacobian
+                        and not controller.problem.multifit
+                    ):
+                        params = list(
+                            controller.starting_values[
+                                controller.parameter_set
+                            ].values()
+                        )
+                        self._check_jacobian(
+                            func=cost_func.eval_cost,
+                            jac=cost_func.jac_cost,
+                            params=params,
+                        )
+
                     #######################
                     # Loops over Hessians #
                     #######################
@@ -629,3 +645,37 @@ class Fit:
             controller.check_bounds_respected()
 
         return accuracy, runtimes, energy
+
+    def _check_jacobian(self, func, jac, params):
+        """
+        Check how similar the jacobian is to a finite difference
+        approximation of the function
+
+        :param func: The function for which to compute the jacobian
+        :type func: function
+
+        :param jac: The jacobian of the function being used in Fitbenchmarking
+        :type jac: function
+
+        :param params: The values at which to evaluate the jacobian matrix
+        :type params: list
+
+        """
+        analytical_grad = jac(params)
+        finite_diff_grad = approx_fprime(params, func)
+
+        max_range = (np.ptp(analytical_grad) + np.ptp(finite_diff_grad)) / 2
+
+        error_from_check_grad = check_grad(func, jac, params)
+        normalized_error = error_from_check_grad / max_range
+
+        if normalized_error > 10**-3:
+            LOGGER.warning(
+                "A relative error of %.6f was detected between "
+                "the jacobian computed by Fitbenchmarking and the one "
+                "obtained through a finite difference approximation. "
+                "This might depend on either the initial parameters "
+                "provided or the jacobian function, if this has also "
+                "been provided by the user.",
+                normalized_error,
+            )
diff --git a/fitbenchmarking/utils/options.py b/fitbenchmarking/utils/options.py
index 8d16f386..e6bd8e90 100644
--- a/fitbenchmarking/utils/options.py
+++ b/fitbenchmarking/utils/options.py
@@ -292,6 +292,7 @@ class Options:
         ],
         "results_browser": [True, False],
         "run_dash": [True, False],
+        "check_jacobian": [True, False],
         "colour_map": plt.colormaps(),
     }
     VALID_LOGGING = {
@@ -469,6 +470,7 @@ class Options:
         "comparison_mode": "both",
         "results_browser": True,
         "run_dash": True,
+        "check_jacobian": True,
         "table_type": [
             "acc",
             "runtime",
@@ -628,6 +630,13 @@ class Options:
                 output.getboolean, "run_dash", additional_options
             )
 
+        if "check_jacobian" in additional_options:
+            self.check_jacobian = additional_options["check_jacobian"]
+        else:
+            self.check_jacobian = self.read_value(
+                output.getboolean, "check_jacobian", additional_options
+            )
+
         if "pbar" in additional_options:
             self.pbar = additional_options["pbar"]
         else:
@@ -805,6 +814,7 @@ class Options:
             "make_plots": self.make_plots,
             "results_browser": self.results_browser,
             "run_dash": self.run_dash,
+            "check_jacobian": self.check_jacobian,
             "pbar": self.pbar,
             "table_type": list_to_string(self.table_type),
             "run_name": self.run_name,


Test patch:
diff --git a/fitbenchmarking/core/tests/test_fitting_benchmarking.py b/fitbenchmarking/core/tests/test_fitting_benchmarking.py
index fd038090..5f7327ac 100644
--- a/fitbenchmarking/core/tests/test_fitting_benchmarking.py
+++ b/fitbenchmarking/core/tests/test_fitting_benchmarking.py
@@ -519,7 +519,10 @@ class JacobianTests(unittest.TestCase):
         assert len(results) == 2
 
     @patch(f"{FITTING_DIR}.Fit._loop_over_hessians")
-    def test_loop_over_jacobians_sparsity_check_true(self, loop_over_hessians):
+    @patch(f"{FITTING_DIR}.Fit._check_jacobian")
+    def test_loop_over_jacobians_sparsity_check_true(
+        self, check_jacobian, loop_over_hessians
+    ):
         """
         The test checks _loop_over_jacobians method
         handles the check for sparsity correctly
@@ -530,7 +533,37 @@ class JacobianTests(unittest.TestCase):
         self.controller.jacobian_enabled_solvers = ["Newton-CG"]
         self.controller.sparsity_enabled_solvers = ["Newton-CG"]
         results = self.fit._loop_over_jacobians(self.controller)
-        assert len(results) == 3
+        assert len(results) == 3 == check_jacobian.call_count
+
+    @patch(f"{FITTING_DIR}.check_grad")
+    @patch(f"{FITTING_DIR}.np.ptp")
+    def test_check_jacobian_raises_warning(self, mock_ptp, mock_check_grad):
+        """
+        The test checks _check_jacobian method raises
+        a warning when the jacobian check fails
+        """
+        with self.assertLogs(LOGGER, level="WARNING") as log:
+            mock_check_grad.return_value = 0.01
+            mock_ptp.return_value = 1
+            self.controller.cost_func.jacobian = Analytic(
+                self.controller.cost_func.problem
+            )
+            self.fit._check_jacobian(
+                func=self.controller.cost_func.eval_cost,
+                jac=self.controller.cost_func.jac_cost,
+                params=[94.9, 0.009, 90.1, 113.0, 20.0, 73.8, 140.0, 20.0],
+            )
+            self.assertTrue(
+                (
+                    "A relative error of 0.010000 was detected between "
+                    "the jacobian computed by Fitbenchmarking and the one "
+                    "obtained through a finite difference approximation. "
+                    "This might depend on either the initial parameters "
+                    "provided or the jacobian function, if this has also "
+                    "been provided by the user."
+                )
+                in log.output[0]
+            )
 
 
 class MinimizersTests(unittest.TestCase):