diff --git a/.github/workflows/docs.yml b/.github/workflows/docs.yml index 7ff5123..97db902 100644 --- a/.github/workflows/docs.yml +++ b/.github/workflows/docs.yml @@ -61,7 +61,7 @@ jobs: deploy: runs-on: ubuntu-latest needs: build - if: (github.event_name == 'push' && startsWith(github.ref, 'refs/tags/')) || github.event_name == 'workflow_dispatch' + if: (github.event_name == 'push' && startsWith(github.ref, 'refs/tags/')) permissions: pages: write id-token: write diff --git a/Cargo.toml b/Cargo.toml index e6ef2da..d469a25 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -10,4 +10,40 @@ license = "BSD-3-Clause" [workspace.dependencies] nalgebra = "0.34.1" -pyo3 = { version = "0.27.1", default-features = false, features = ["macros"] } \ No newline at end of file +pyo3 = { version = "0.27.1", default-features = false, features = ["macros"] } + +[workspace.lints.clippy] +# Enable pedantic lints for higher code quality standards +pedantic = { level = "deny", priority = -1 } + +# Specific allows to reduce noise from pedantic +module_name_repetitions = "allow" # Common in large codebases +similar_names = "allow" # Scientific code often has x, y, z, etc. +must_use_candidate = "allow" # Too noisy for mathematical functions +cast_precision_loss = "allow" # Acceptable in numerical code (usize to f64) +cast_possible_truncation = "allow" # Acceptable when converting calculated values +cast_sign_loss = "allow" # Acceptable in constrained numerical contexts +cast_possible_wrap = "allow" # Acceptable in controlled numerical contexts + +# Additional strict lints +missing_errors_doc = "warn" +missing_panics_doc = "warn" +missing_safety_doc = "warn" +undocumented_unsafe_blocks = "warn" + +# Performance and correctness +inefficient_to_string = "warn" +manual_ok_or = "warn" +redundant_closure_for_method_calls = "warn" + +# Complexity and style +cognitive_complexity = "warn" +too_many_lines = "warn" + +# Selected nursery lints +option_if_let_else = "warn" +suboptimal_flops = "allow" # Suggested optimisations often reduce readability in tests + +# Test-specific allows +float_cmp = "allow" # Exact float equality is standard in tests +cast_lossless = "allow" # i32 to f64 conversions in tests/benchmarks \ No newline at end of file diff --git a/docs/tutorials/notebooks/ode_fitting_diffsol.ipynb b/docs/tutorials/notebooks/ode_fitting_diffsol.ipynb index 47d9ae7..04149fe 100644 --- a/docs/tutorials/notebooks/ode_fitting_diffsol.ipynb +++ b/docs/tutorials/notebooks/ode_fitting_diffsol.ipynb @@ -35,7 +35,13 @@ } }, "outputs": [], - "source": "# Import plotting utilities\nimport diffid\nimport matplotlib.pyplot as plt\nimport numpy as np\nfrom diffid.plotting import ode_fit" + "source": [ + "# Import plotting utilities\n", + "import diffid\n", + "import matplotlib.pyplot as plt\n", + "import numpy as np\n", + "from diffid.plotting import ode_fit" + ] }, { "cell_type": "markdown", @@ -114,7 +120,34 @@ } }, "outputs": [], - "source": "# Time points\nt_span = np.linspace(0, 4, 100)\n\n# True logistic growth solution with known parameters\n# y(t) = y0 * exp(r*t) / (1 + y0 * (exp(r*t) - 1) / k)\ny0 = 0.1\nr_true = 1.0\nk_true = 1.0\n\n# Generate clean data\nexp_rt = np.exp(r_true * t_span)\ny_data = y0 * exp_rt / (1 + y0 * (exp_rt - 1) / k_true)\n\n# Add some noise\nnp.random.seed(42)\nnoise_level = 0.01\ny_noisy = y_data + np.random.normal(0, noise_level, size=y_data.shape)\n\n# Diffid expects data as [time, observation] columns\ndata = np.column_stack((t_span, y_noisy))\n\nprint(f\"Generated {len(t_span)} data points\")\nprint(f\"Time range: [{t_span[0]:.2f}, {t_span[-1]:.2f}]\")\nprint(f\"Value range: [{y_noisy.min():.3f}, {y_noisy.max():.3f}]\")\nprint(f\"Noise level: {noise_level}\")\nprint(f\"\\nTrue parameters: r = {r_true}, k = {k_true}\")" + "source": [ + "# Time points\n", + "t_span = np.linspace(0, 4, 100)\n", + "\n", + "# True logistic growth solution with known parameters\n", + "# y(t) = y0 * exp(r*t) / (1 + y0 * (exp(r*t) - 1) / k)\n", + "y0 = 0.1\n", + "r_true = 1.0\n", + "k_true = 1.0\n", + "\n", + "# Generate clean data\n", + "exp_rt = np.exp(r_true * t_span)\n", + "y_data = y0 * exp_rt / (1 + y0 * (exp_rt - 1) / k_true)\n", + "\n", + "# Add some noise\n", + "np.random.seed(42)\n", + "noise_level = 0.01\n", + "y_noisy = y_data + np.random.normal(0, noise_level, size=y_data.shape)\n", + "\n", + "# Diffid expects data as [time, observation] columns\n", + "data = np.column_stack((t_span, y_noisy))\n", + "\n", + "print(f\"Generated {len(t_span)} data points\")\n", + "print(f\"Time range: [{t_span[0]:.2f}, {t_span[-1]:.2f}]\")\n", + "print(f\"Value range: [{y_noisy.min():.3f}, {y_noisy.max():.3f}]\")\n", + "print(f\"Noise level: {noise_level}\")\n", + "print(f\"\\nTrue parameters: r = {r_true}, k = {k_true}\")" + ] }, { "cell_type": "markdown", @@ -184,7 +217,24 @@ } }, "outputs": [], - "source": "# Create builder\nbuilder = (\n diffid.DiffsolBuilder()\n .with_diffsl(dsl_model)\n .with_data(data)\n .with_tolerances(1e-6, 1e-8) # Relative and absolute tolerances\n .with_parameter(\"r\", 10.0) # Initial guess (deliberately wrong)\n .with_parameter(\"k\", 10.0) # Initial guess (deliberately wrong)\n .with_parallel(True) # Enable parallel evaluation\n)\n\nproblem = builder.build()\n\nprint(\"Problem built successfully!\")\nprint(\"\\nInitial parameter guesses: r = 50.0, k = 50.0\")\nprint(f\"(Far from true values: r = {r_true}, k = {k_true})\")" + "source": [ + "# Create builder\n", + "builder = (\n", + " diffid.DiffsolBuilder()\n", + " .with_diffsl(dsl_model)\n", + " .with_data(data)\n", + " .with_tolerances(1e-6, 1e-8) # Relative and absolute tolerances\n", + " .with_parameter(\"r\", 10.0) # Initial guess\n", + " .with_parameter(\"k\", 10.0) # Initial guess\n", + " .with_parallel(True) # Enable parallel evaluation\n", + ")\n", + "\n", + "problem = builder.build()\n", + "\n", + "print(\"Problem built successfully!\")\n", + "print(\"\\nInitial parameter guesses: r = 50.0, k = 50.0\")\n", + "print(f\"(Far from true values: r = {r_true}, k = {k_true})\")" + ] }, { "cell_type": "markdown", @@ -513,4 +563,4 @@ }, "nbformat": 4, "nbformat_minor": 4 -} \ No newline at end of file +} diff --git a/docs/tutorials/notebooks/optimisation_basics.ipynb b/docs/tutorials/notebooks/optimisation_basics.ipynb index b52fc68..1dfccff 100644 --- a/docs/tutorials/notebooks/optimisation_basics.ipynb +++ b/docs/tutorials/notebooks/optimisation_basics.ipynb @@ -249,7 +249,7 @@ "\n", "for name, result in results.items():\n", " print(\n", - " f\"{name:<15} {str(result.success):<10} {result.value:<15.3e} \"\n", + " f\"{name:<15} {result.success!s:<10} {result.value:<15.3e} \"\n", " f\"{result.iterations:<12} {result.evaluations}\"\n", " )\n", "\n", diff --git a/docs/tutorials/notebooks/parallel_optimisation.ipynb b/docs/tutorials/notebooks/parallel_optimisation.ipynb index b56d358..92a80d9 100644 --- a/docs/tutorials/notebooks/parallel_optimisation.ipynb +++ b/docs/tutorials/notebooks/parallel_optimisation.ipynb @@ -26,7 +26,19 @@ } }, "outputs": [], - "source": "import multiprocessing\nimport time\n\nimport diffid\nimport matplotlib.pyplot as plt\nimport numpy as np\nfrom scipy.integrate import solve_ivp\n\n# Detect available cores\nn_cores = multiprocessing.cpu_count()\nprint(f\"Available CPU cores: {n_cores}\")" + "source": [ + "import multiprocessing\n", + "import time\n", + "\n", + "import diffid\n", + "import matplotlib.pyplot as plt\n", + "import numpy as np\n", + "from scipy.integrate import solve_ivp\n", + "\n", + "# Detect available cores\n", + "n_cores = multiprocessing.cpu_count()\n", + "print(f\"Available CPU cores: {n_cores}\")" + ] }, { "cell_type": "markdown", @@ -136,7 +148,39 @@ } }, "outputs": [], - "source": "# Sequential execution (parallel=False)\nprint(\"Running CMA-ES (sequential)...\")\n\nstart = time.time()\n\nresult_seq = (\n diffid.DiffsolBuilder()\n .with_diffsl(model_str)\n .with_data(data)\n .with_parameter(\"alpha\", 0.8)\n .with_parameter(\"beta\", 0.3)\n .with_parameter(\"delta\", 0.05)\n .with_parameter(\"gamma\", 0.3)\n .with_cost(diffid.SSE())\n .with_parallel(False) # Sequential evaluation\n .with_optimiser(diffid.CMAES().with_max_iter(100).with_step_size(0.3))\n .build()\n .optimise()\n)\n\ntime_seq = time.time() - start\n\nprint(\"\\n\" + \"=\" * 60)\nprint(\"CMA-ES (Sequential)\")\nprint(\"=\" * 60)\nprint(f\"Solution: {result_seq.x}\")\nprint(f\"True params: {list(true_params.values())}\")\nprint(f\"Final SSE: {result_seq.value:.3f}\")\nprint(f\"Evaluations: {result_seq.evaluations}\")\nprint(f\"Time: {time_seq:.2f}s\")\nprint(f\"Time per eval: {time_seq / result_seq.evaluations * 1000:.2f}ms\")" + "source": [ + "# Sequential execution (parallel=False)\n", + "print(\"Running CMA-ES (sequential)...\")\n", + "\n", + "start = time.time()\n", + "\n", + "result_seq = (\n", + " diffid.DiffsolBuilder()\n", + " .with_diffsl(model_str)\n", + " .with_data(data)\n", + " .with_parameter(\"alpha\", 0.8)\n", + " .with_parameter(\"beta\", 0.3)\n", + " .with_parameter(\"delta\", 0.05)\n", + " .with_parameter(\"gamma\", 0.3)\n", + " .with_cost(diffid.SSE())\n", + " .with_parallel(False) # Sequential evaluation\n", + " .with_optimiser(diffid.CMAES().with_max_iter(100).with_step_size(0.3))\n", + " .build()\n", + " .optimise()\n", + ")\n", + "\n", + "time_seq = time.time() - start\n", + "\n", + "print(\"\\n\" + \"=\" * 60)\n", + "print(\"CMA-ES (Sequential)\")\n", + "print(\"=\" * 60)\n", + "print(f\"Solution: {result_seq.x}\")\n", + "print(f\"True params: {list(true_params.values())}\")\n", + "print(f\"Final SSE: {result_seq.value:.3f}\")\n", + "print(f\"Evaluations: {result_seq.evaluations}\")\n", + "print(f\"Time: {time_seq:.2f}s\")\n", + "print(f\"Time per eval: {time_seq / result_seq.evaluations * 1000:.2f}ms\")" + ] }, { "cell_type": "markdown", @@ -159,7 +203,42 @@ } }, "outputs": [], - "source": "# Parallel execution (parallel=True)\nprint(\"Running CMA-ES (parallel)...\")\n\nstart = time.time()\n\nresult_par = (\n diffid.DiffsolBuilder()\n .with_diffsl(model_str)\n .with_data(data)\n .with_parameter(\"alpha\", 0.8)\n .with_parameter(\"beta\", 0.3)\n .with_parameter(\"delta\", 0.05)\n .with_parameter(\"gamma\", 0.3)\n .with_cost(diffid.SSE())\n .with_parallel(True) # Parallel evaluation!\n .with_optimiser(diffid.CMAES().with_max_iter(100).with_step_size(0.3))\n .build()\n .optimise()\n)\n\ntime_par = time.time() - start\n\nprint(\"\\n\" + \"=\" * 60)\nprint(\"CMA-ES (Parallel)\")\nprint(\"=\" * 60)\nprint(f\"Solution: {result_par.x}\")\nprint(f\"True params: {list(true_params.values())}\")\nprint(f\"Final SSE: {result_par.value:.3f}\")\nprint(f\"Evaluations: {result_par.evaluations}\")\nprint(f\"Time: {time_par:.2f}s\")\nprint(f\"Time per eval: {time_par / result_par.evaluations * 1000:.2f}ms\")\n\nspeedup = time_seq / time_par\nprint(f\"\\nšŸš€ Speedup: {speedup:.2f}x (with {n_cores} cores)\")" + "source": [ + "# Parallel execution (parallel=True)\n", + "print(\"Running CMA-ES (parallel)...\")\n", + "\n", + "start = time.time()\n", + "\n", + "result_par = (\n", + " diffid.DiffsolBuilder()\n", + " .with_diffsl(model_str)\n", + " .with_data(data)\n", + " .with_parameter(\"alpha\", 0.8)\n", + " .with_parameter(\"beta\", 0.3)\n", + " .with_parameter(\"delta\", 0.05)\n", + " .with_parameter(\"gamma\", 0.3)\n", + " .with_cost(diffid.SSE())\n", + " .with_parallel(True) # Parallel evaluation!\n", + " .with_optimiser(diffid.CMAES().with_max_iter(100).with_step_size(0.3))\n", + " .build()\n", + " .optimise()\n", + ")\n", + "\n", + "time_par = time.time() - start\n", + "\n", + "print(\"\\n\" + \"=\" * 60)\n", + "print(\"CMA-ES (Parallel)\")\n", + "print(\"=\" * 60)\n", + "print(f\"Solution: {result_par.x}\")\n", + "print(f\"True params: {list(true_params.values())}\")\n", + "print(f\"Final SSE: {result_par.value:.3f}\")\n", + "print(f\"Evaluations: {result_par.evaluations}\")\n", + "print(f\"Time: {time_par:.2f}s\")\n", + "print(f\"Time per eval: {time_par / result_par.evaluations * 1000:.2f}ms\")\n", + "\n", + "speedup = time_seq / time_par\n", + "print(f\"\\nšŸš€ Speedup: {speedup:.2f}x (with {n_cores} cores)\")" + ] }, { "cell_type": "code", @@ -177,7 +256,47 @@ } }, "outputs": [], - "source": "# Parallel with larger population (more parallel work per generation)\nprint(\"Running CMA-ES (parallel, large population)...\")\n\nstart = time.time()\n\nresult_par_large = (\n diffid.DiffsolBuilder()\n .with_diffsl(model_str)\n .with_data(data)\n .with_parameter(\"alpha\", 0.8)\n .with_parameter(\"beta\", 0.3)\n .with_parameter(\"delta\", 0.05)\n .with_parameter(\"gamma\", 0.3)\n .with_cost(diffid.SSE())\n .with_parallel(True)\n .with_optimiser(\n diffid.CMAES()\n .with_max_iter(50)\n .with_step_size(0.3)\n .with_population_size(2 * n_cores) # Match population to available cores\n )\n .build()\n .optimise()\n)\n\ntime_par_large = time.time() - start\n\nprint(\"\\n\" + \"=\" * 60)\nprint(f\"CMA-ES (Parallel, Population={2 * n_cores})\")\nprint(\"=\" * 60)\nprint(f\"Solution: {result_par_large.x}\")\nprint(f\"True params: {list(true_params.values())}\")\nprint(f\"Final SSE: {result_par_large.value:.3f}\")\nprint(f\"Evaluations: {result_par_large.evaluations}\")\nprint(f\"Time: {time_par_large:.2f}s\")\nprint(f\"Time per eval: {time_par_large / result_par_large.evaluations * 1000:.2f}ms\")\n\nspeedup_large = time_seq / time_par_large\nprint(f\"\\nšŸš€ Speedup: {speedup_large:.2f}x\")" + "source": [ + "# Parallel with larger population (more parallel work per generation)\n", + "print(\"Running CMA-ES (parallel, large population)...\")\n", + "\n", + "start = time.time()\n", + "\n", + "result_par_large = (\n", + " diffid.DiffsolBuilder()\n", + " .with_diffsl(model_str)\n", + " .with_data(data)\n", + " .with_parameter(\"alpha\", 0.8)\n", + " .with_parameter(\"beta\", 0.3)\n", + " .with_parameter(\"delta\", 0.05)\n", + " .with_parameter(\"gamma\", 0.3)\n", + " .with_cost(diffid.SSE())\n", + " .with_parallel(True)\n", + " .with_optimiser(\n", + " diffid.CMAES()\n", + " .with_max_iter(50)\n", + " .with_step_size(0.3)\n", + " .with_population_size(2 * n_cores) # Match population to available cores\n", + " )\n", + " .build()\n", + " .optimise()\n", + ")\n", + "\n", + "time_par_large = time.time() - start\n", + "\n", + "print(\"\\n\" + \"=\" * 60)\n", + "print(f\"CMA-ES (Parallel, Population={2 * n_cores})\")\n", + "print(\"=\" * 60)\n", + "print(f\"Solution: {result_par_large.x}\")\n", + "print(f\"True params: {list(true_params.values())}\")\n", + "print(f\"Final SSE: {result_par_large.value:.3f}\")\n", + "print(f\"Evaluations: {result_par_large.evaluations}\")\n", + "print(f\"Time: {time_par_large:.2f}s\")\n", + "print(f\"Time per eval: {time_par_large / result_par_large.evaluations * 1000:.2f}ms\")\n", + "\n", + "speedup_large = time_seq / time_par_large\n", + "print(f\"\\nšŸš€ Speedup: {speedup_large:.2f}x\")" + ] }, { "cell_type": "markdown", @@ -200,7 +319,60 @@ } }, "outputs": [], - "source": "# Define an expensive objective function\ndef expensive_objective(params):\n \"\"\"\n An expensive objective function that simulates computation.\n In practice, this could be a complex simulation, ML model, etc.\n \"\"\"\n alpha, beta, delta, gamma = params\n\n # Simulate expensive computation (ODE integration with scipy)\n sol = solve_ivp(\n lotka_volterra,\n [0, 100],\n [10.0, 5.0],\n args=(alpha, beta, delta, gamma),\n t_eval=t_data,\n method=\"RK45\",\n )\n\n if not sol.success:\n return 1e10 # Return large value for failed integrations\n\n # Compute SSE\n y_pred = sol.y.T\n sse = np.sum((y_pred - y_observed) ** 2)\n return sse\n\n\n# Test single evaluation time\nn_evals = 20\ntest_params = [[0.8 + i * 0.02, 0.3, 0.05, 0.3] for i in range(n_evals)]\n\nprint(f\"Testing {n_evals} sequential evaluations...\")\n\nstart = time.time()\nseq_results = [expensive_objective(p) for p in test_params]\ntime_seq_python = time.time() - start\nprint(\n f\"Sequential: {time_seq_python:.2f}s ({time_seq_python / n_evals * 1000:.1f}ms/eval)\"\n)\n\n# Note: Multiprocessing in notebooks has pickling limitations\n# In a regular Python script, you would use:\nprint(\"\"\"\nšŸ’” To parallelise in a script, use multiprocessing:\n\nfrom concurrent.futures import ProcessPoolExecutor\n\ndef evaluate_batch_parallel(params_list, n_workers=8):\n with ProcessPoolExecutor(max_workers=n_workers) as executor:\n return list(executor.map(expensive_objective, params_list))\n\n# This provides near-linear speedup for expensive functions\n\"\"\")" + "source": [ + "# Define an expensive objective function\n", + "def expensive_objective(params):\n", + " \"\"\"\n", + " An expensive objective function that simulates computation.\n", + " In practice, this could be a complex simulation, ML model, etc.\n", + " \"\"\"\n", + " alpha, beta, delta, gamma = params\n", + "\n", + " # Simulate expensive computation (ODE integration with scipy)\n", + " sol = solve_ivp(\n", + " lotka_volterra,\n", + " [0, 100],\n", + " [10.0, 5.0],\n", + " args=(alpha, beta, delta, gamma),\n", + " t_eval=t_data,\n", + " method=\"RK45\",\n", + " )\n", + "\n", + " if not sol.success:\n", + " return 1e10 # Return large value for failed integrations\n", + "\n", + " # Compute SSE\n", + " y_pred = sol.y.T\n", + " return np.sum((y_pred - y_observed) ** 2)\n", + "\n", + "\n", + "# Test single evaluation time\n", + "n_evals = 20\n", + "test_params = [[0.8 + i * 0.02, 0.3, 0.05, 0.3] for i in range(n_evals)]\n", + "\n", + "print(f\"Testing {n_evals} sequential evaluations...\")\n", + "\n", + "start = time.time()\n", + "seq_results = [expensive_objective(p) for p in test_params]\n", + "time_seq_python = time.time() - start\n", + "print(\n", + " f\"Sequential: {time_seq_python:.2f}s ({time_seq_python / n_evals * 1000:.1f}ms/eval)\"\n", + ")\n", + "\n", + "# Note: Multiprocessing in notebooks has pickling limitations\n", + "# In a regular Python script, you would use:\n", + "print(\"\"\"\n", + "šŸ’” To parallelise in a script, use multiprocessing:\n", + "\n", + "from concurrent.futures import ProcessPoolExecutor\n", + "\n", + "def evaluate_batch_parallel(params_list, n_workers=8):\n", + " with ProcessPoolExecutor(max_workers=n_workers) as executor:\n", + " return list(executor.map(expensive_objective, params_list))\n", + "\n", + "# This provides near-linear speedup for expensive functions\n", + "\"\"\")" + ] }, { "cell_type": "code", @@ -315,7 +487,81 @@ } }, "outputs": [], - "source": "# Test DiffsolBuilder scaling with data size\ndata_sizes = [100, 200, 500]\nscaling_results = []\n\nprint(\"Testing DiffsolBuilder scaling with data size...\\n\")\n\nfor n_points in data_sizes:\n # Generate data with different sizes\n t_test = np.linspace(0, 100, n_points)\n sol_test = solve_ivp(\n lotka_volterra,\n [t_test[0], t_test[-1]],\n [10.0, 5.0],\n args=(\n true_params[\"alpha\"],\n true_params[\"beta\"],\n true_params[\"delta\"],\n true_params[\"gamma\"],\n ),\n t_eval=t_test,\n method=\"RK45\",\n )\n y_test = sol_test.y.T + np.random.normal(0, 0.3, (n_points, 2))\n data_test = np.column_stack((t_test, y_test))\n\n # Sequential\n start = time.time()\n _ = (\n diffid.DiffsolBuilder()\n .with_diffsl(model_str)\n .with_data(data_test)\n .with_parameter(\"alpha\", 0.8)\n .with_parameter(\"beta\", 0.3)\n .with_parameter(\"delta\", 0.05)\n .with_parameter(\"gamma\", 0.3)\n .with_cost(diffid.SSE())\n .with_parallel(False)\n .with_optimiser(diffid.CMAES().with_max_iter(30).with_step_size(0.3))\n .build()\n .optimise()\n )\n t_seq_scale = time.time() - start\n\n # Parallel\n start = time.time()\n _ = (\n diffid.DiffsolBuilder()\n .with_diffsl(model_str)\n .with_data(data_test)\n .with_parameter(\"alpha\", 0.8)\n .with_parameter(\"beta\", 0.3)\n .with_parameter(\"delta\", 0.05)\n .with_parameter(\"gamma\", 0.3)\n .with_cost(diffid.SSE())\n .with_parallel(True)\n .with_optimiser(diffid.CMAES().with_max_iter(30).with_step_size(0.3))\n .build()\n .optimise()\n )\n t_par_scale = time.time() - start\n\n sp = t_seq_scale / t_par_scale\n scaling_results.append(\n {\n \"n_points\": n_points,\n \"time_seq\": t_seq_scale,\n \"time_par\": t_par_scale,\n \"speedup\": sp,\n }\n )\n print(\n f\"N={n_points:3d}: Sequential={t_seq_scale:.2f}s, Parallel={t_par_scale:.2f}s, Speedup={sp:.2f}x\"\n )" + "source": [ + "# Test DiffsolBuilder scaling with data size\n", + "data_sizes = [100, 200, 500]\n", + "scaling_results = []\n", + "\n", + "print(\"Testing DiffsolBuilder scaling with data size...\\n\")\n", + "\n", + "for n_points in data_sizes:\n", + " # Generate data with different sizes\n", + " t_test = np.linspace(0, 100, n_points)\n", + " sol_test = solve_ivp(\n", + " lotka_volterra,\n", + " [t_test[0], t_test[-1]],\n", + " [10.0, 5.0],\n", + " args=(\n", + " true_params[\"alpha\"],\n", + " true_params[\"beta\"],\n", + " true_params[\"delta\"],\n", + " true_params[\"gamma\"],\n", + " ),\n", + " t_eval=t_test,\n", + " method=\"RK45\",\n", + " )\n", + " y_test = sol_test.y.T + np.random.normal(0, 0.3, (n_points, 2))\n", + " data_test = np.column_stack((t_test, y_test))\n", + "\n", + " # Sequential\n", + " start = time.time()\n", + " _ = (\n", + " diffid.DiffsolBuilder()\n", + " .with_diffsl(model_str)\n", + " .with_data(data_test)\n", + " .with_parameter(\"alpha\", 0.8)\n", + " .with_parameter(\"beta\", 0.3)\n", + " .with_parameter(\"delta\", 0.05)\n", + " .with_parameter(\"gamma\", 0.3)\n", + " .with_cost(diffid.SSE())\n", + " .with_parallel(False)\n", + " .with_optimiser(diffid.CMAES().with_max_iter(30).with_step_size(0.3))\n", + " .build()\n", + " .optimise()\n", + " )\n", + " t_seq_scale = time.time() - start\n", + "\n", + " # Parallel\n", + " start = time.time()\n", + " _ = (\n", + " diffid.DiffsolBuilder()\n", + " .with_diffsl(model_str)\n", + " .with_data(data_test)\n", + " .with_parameter(\"alpha\", 0.8)\n", + " .with_parameter(\"beta\", 0.3)\n", + " .with_parameter(\"delta\", 0.05)\n", + " .with_parameter(\"gamma\", 0.3)\n", + " .with_cost(diffid.SSE())\n", + " .with_parallel(True)\n", + " .with_optimiser(diffid.CMAES().with_max_iter(30).with_step_size(0.3))\n", + " .build()\n", + " .optimise()\n", + " )\n", + " t_par_scale = time.time() - start\n", + "\n", + " sp = t_seq_scale / t_par_scale\n", + " scaling_results.append(\n", + " {\n", + " \"n_points\": n_points,\n", + " \"time_seq\": t_seq_scale,\n", + " \"time_par\": t_par_scale,\n", + " \"speedup\": sp,\n", + " }\n", + " )\n", + " print(\n", + " f\"N={n_points:3d}: Sequential={t_seq_scale:.2f}s, Parallel={t_par_scale:.2f}s, Speedup={sp:.2f}x\"\n", + " )" + ] }, { "cell_type": "code", diff --git a/examples/predator_prey/generate_data_diffrax.py b/examples/predator_prey/generate_data_diffrax.py index 0334f98..4a0c0f6 100644 --- a/examples/predator_prey/generate_data_diffrax.py +++ b/examples/predator_prey/generate_data_diffrax.py @@ -17,7 +17,7 @@ SEED = 8 -def lotka_volterra(t, state, params): +def lotka_volterra(_t, state, params): x, y = state alpha, beta, delta, gamma = params return jnp.array([alpha * x - beta * x * y, delta * x * y - gamma * y]) diff --git a/examples/predator_prey/predator_prey_diffeqpy.py b/examples/predator_prey/predator_prey_diffeqpy.py index 22b28f7..8d921a1 100644 --- a/examples/predator_prey/predator_prey_diffeqpy.py +++ b/examples/predator_prey/predator_prey_diffeqpy.py @@ -58,7 +58,8 @@ def simulate(params): gen_path = pathlib.Path(__file__).with_name("generate_data_diffrax.py") spec = importlib.util.spec_from_file_location("pp_gen", gen_path) module = importlib.util.module_from_spec(spec) - assert spec is not None and spec.loader is not None + assert spec is not None + assert spec.loader is not None spec.loader.exec_module(module) module.main(data_path) data = np.load(str(data_path)) diff --git a/examples/predator_prey/predator_prey_diffrax.py b/examples/predator_prey/predator_prey_diffrax.py index 53e0d99..2bb1b4b 100644 --- a/examples/predator_prey/predator_prey_diffrax.py +++ b/examples/predator_prey/predator_prey_diffrax.py @@ -29,7 +29,7 @@ TRUE_PARAMS = np.array([1.1, 0.4, 0.1, 0.4]) # [alpha, beta, delta, gamma] -def lotka_volterra(t, state, params): +def lotka_volterra(_t, state, params): """Lotka-Volterra predator-prey dynamics. dx/dt = alpha*x - beta*x*y (prey growth and predation) @@ -81,7 +81,8 @@ def simulate(params): gen_path = pathlib.Path(__file__).with_name("generate_data_diffrax.py") spec = importlib.util.spec_from_file_location("pp_gen", gen_path) module = importlib.util.module_from_spec(spec) - assert spec is not None and spec.loader is not None + assert spec is not None + assert spec.loader is not None spec.loader.exec_module(module) module.main(data_path) data = np.load(str(data_path)) diff --git a/examples/predator_prey/predator_prey_diffsol.py b/examples/predator_prey/predator_prey_diffsol.py index 9c57caf..d671fbd 100644 --- a/examples/predator_prey/predator_prey_diffsol.py +++ b/examples/predator_prey/predator_prey_diffsol.py @@ -24,7 +24,8 @@ gen_path = pathlib.Path(__file__).with_name("generate_data_diffrax.py") spec = importlib.util.spec_from_file_location("pp_gen", gen_path) module = importlib.util.module_from_spec(spec) - assert spec is not None and spec.loader is not None + assert spec is not None + assert spec.loader is not None spec.loader.exec_module(module) module.main(data_path) data = np.load(str(data_path)) diff --git a/pyproject.toml b/pyproject.toml index ecde2dd..7b37a57 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -63,6 +63,7 @@ examples = [ # Non-JAX dependencies "matplotlib>=3.8", + "scipy>=1.16.2", "diffeqpy", ] @@ -74,22 +75,39 @@ fix = true [tool.ruff.lint] select = [ "A", # flake8-builtins: Check for Python builtins being used as variables or parameters + "ARG", # flake8-unused-arguments: Check for unused function arguments "B", # flake8-bugbear: Find likely bugs and design problems + "C90", # mccabe: Check for code complexity "E", # pycodestyle errors - "W", # pycodestyle warnings "F", # pyflakes: Detect various errors by parsing the source file + "FURB", # refurb: Modernization suggestions "I", # isort: Check and enforce import ordering "ISC", # flake8-implicit-str-concat: Check for implicit string concatenation + "PERF", # Perflint: Performance anti-patterns + "PLE", # Pylint errors + "PLR", # Pylint refactor: Check for code smells and complexity + "PLW", # Pylint warnings + "PT", # flake8-pytest-style: Check pytest best practices + "PTH", # flake8-use-pathlib: Prefer pathlib over os.path + "RET", # flake8-return: Check return statement consistency + "RUF", # Ruff-specific rules + "SIM", # flake8-simplify: Suggestions for simplifying code "TID", # flake8-tidy-imports: Validate import hygiene "UP", # pyupgrade: Automatically upgrade syntax for newer versions of Python + "W", # pycodestyle warnings "SLF001", # flake8-string-format: Check for private object name access ] -ignore = ["E501","E741"] +ignore = [ + "E501", # Line too long (handled by formatter) + "E741", # Ambiguous variable name + "PLR0913", # Too many arguments to function call + "PLR2004", # Magic value used in comparison +] [tool.ruff.lint.per-file-ignores] -"tests/*" = ["SLF001"] -"**.ipynb" = ["E402", "E703"] +"tests/*" = ["SLF001", "ARG001", "PLR2004", "RUF001", "RUF002"] # Allow private access, unused fixtures, magic values, and Greek letters in tests +"**.ipynb" = ["E402", "E703", "PLR2004", "ARG001", "RUF001"] # Relax some rules for notebooks [tool.ruff.lint.flake8-tidy-imports] ban-relative-imports = "all" diff --git a/python/Cargo.toml b/python/Cargo.toml index f7f26fb..78e5f85 100644 --- a/python/Cargo.toml +++ b/python/Cargo.toml @@ -5,6 +5,26 @@ edition.workspace = true license.workspace = true readme = "../README.md" +[lints.clippy] +# Inherit workspace pedantic settings +pedantic = { level = "deny", priority = -1 } + +# Python FFI-specific allows (override workspace for FFI needs) +needless_pass_by_value = "allow" # PyO3 functions require owned values +unused_self = "allow" # PyO3 methods need &self even if unused +match_wildcard_for_single_variants = "allow" # FFI enums less likely to change +unnecessary_wraps = "allow" # PyO3 functions often need Result wrapping +missing_panics_doc = "allow" # PyO3 functions can panic on Python errors +missing_errors_doc = "allow" # FFI error handling is implicit +module_name_repetitions = "allow" +similar_names = "allow" +must_use_candidate = "allow" +cast_precision_loss = "allow" +cast_possible_truncation = "allow" +cast_sign_loss = "allow" +cast_possible_wrap = "allow" +float_cmp = "allow" + [package.metadata.maturin] name = "diffid" python-source = "src" diff --git a/python/src/builders.rs b/python/src/builders.rs index b95c9ec..7b2a9a1 100644 --- a/python/src/builders.rs +++ b/python/src/builders.rs @@ -18,6 +18,13 @@ use pyo3_stub_gen::derive::{gen_stub_pyclass, gen_stub_pymethods}; use crate::optimisers::Optimiser; use crate::{DynProblem, PyCostMetric, PyProblem}; +// Type aliases +type BoxedScalarFn = Box f64 + Send + Sync>; +type BoxedGradientFn = Box Vec + Send + Sync>; +type ParameterSpec = (String, f64, Option<(f64, f64)>); +type ScalarBuilderWithFn = ScalarProblemBuilder; +type ScalarBuilderWithGrad = ScalarProblemBuilder; + // Python Objective Function Wrapper pub(crate) struct PyObjectiveFn { callable: Py, @@ -39,8 +46,7 @@ impl PyObjectiveFn { return match array.len() { 1 => Ok(array[0]), n => Err(PyValueError::new_err(format!( - "Objective array must contain exactly one element, got {}", - n + "Objective array must contain exactly one element, got {n}" ))), }; } @@ -49,8 +55,7 @@ impl PyObjectiveFn { return match values.len() { 1 => Ok(values[0]), n => Err(PyValueError::new_err(format!( - "Objective sequence must contain exactly one element, got {}", - n + "Objective sequence must contain exactly one element, got {n}" ))), }; } @@ -62,12 +67,10 @@ impl PyObjectiveFn { let ty_name = result .get_type() .name() - .map(|n| n.to_string()) - .unwrap_or_else(|_| "unknown".to_string()); + .map_or_else(|_| "unknown".to_string(), |n| n.to_string()); Err(PyTypeError::new_err(format!( - "Objective callable must return a float, numpy array, or single-element sequence; got {}", - ty_name + "Objective callable must return a float, numpy array, or single-element sequence; got {ty_name}" ))) }) } @@ -102,14 +105,11 @@ impl PyGradientFn { enum ScalarBuilderState { Empty(ScalarProblemBuilder), WithFunction { - builder: ScalarProblemBuilder f64 + Send + Sync>, NoGradient>, + builder: ScalarBuilderWithFn, py_callable: Arc, }, WithGradient { - builder: ScalarProblemBuilder< - Box f64 + Send + Sync>, - Box Vec + Send + Sync>, - >, + builder: ScalarBuilderWithGrad, py_callable: Arc, py_gradient: Arc, }, @@ -121,10 +121,16 @@ enum ScalarBuilderState { pub struct PyScalarBuilder { state: ScalarBuilderState, pub(crate) default_optimiser: Option, - parameter_specs: Vec<(String, f64, Option<(f64, f64)>)>, + parameter_specs: Vec, config: HashMap, } +impl Default for PyScalarBuilder { + fn default() -> Self { + Self::new() + } +} + #[cfg_attr(feature = "stubgen", gen_stub_pymethods)] #[pymethods] impl PyScalarBuilder { @@ -148,7 +154,7 @@ impl PyScalarBuilder { ScalarBuilderState::WithFunction { py_callable, .. } => { // Recreate the builder with a new closure from the Arc let objective = Arc::clone(py_callable); - let boxed_fn: Box f64 + Send + Sync> = + let boxed_fn: BoxedScalarFn = Box::new(move |x: &[f64]| objective.call(x).unwrap_or(f64::INFINITY)); let mut builder = ScalarProblemBuilder::new().with_function(boxed_fn); @@ -174,14 +180,13 @@ impl PyScalarBuilder { } => { // Recreate both closures let objective = Arc::clone(py_callable); - let boxed_fn: Box f64 + Send + Sync> = + let boxed_fn: BoxedScalarFn = Box::new(move |x: &[f64]| objective.call(x).unwrap_or(f64::INFINITY)); let grad = Arc::clone(py_gradient); - let boxed_grad: Box Vec + Send + Sync> = - Box::new(move |x: &[f64]| { - grad.call(x).unwrap_or_else(|_| vec![f64::NAN; x.len()]) - }); + let boxed_grad: BoxedGradientFn = Box::new(move |x: &[f64]| { + grad.call(x).unwrap_or_else(|_| vec![f64::NAN; x.len()]) + }); let mut builder = ScalarProblemBuilder::new() .with_function(boxed_fn) @@ -255,7 +260,7 @@ impl PyScalarBuilder { let py_fn = Arc::new(PyObjectiveFn::new(obj)); let objective = Arc::clone(&py_fn); - let boxed_fn: Box f64 + Send + Sync> = + let boxed_fn: BoxedScalarFn = Box::new(move |x: &[f64]| objective.call(x).unwrap_or(f64::INFINITY)); slf.state = match std::mem::replace( @@ -283,7 +288,7 @@ impl PyScalarBuilder { let py_grad = Arc::new(PyGradientFn::new(obj)); let grad = Arc::clone(&py_grad); - let boxed_grad: Box Vec + Send + Sync> = + let boxed_grad: BoxedGradientFn = Box::new(move |x: &[f64]| grad.call(x).unwrap_or_else(|_| vec![f64::NAN; x.len()])); slf.state = match std::mem::replace( @@ -316,14 +321,12 @@ impl PyScalarBuilder { if let Some((lower, upper)) = bounds { if lower >= upper { return Err(PyValueError::new_err(format!( - "Invalid bounds for parameter '{}': lower bound ({}) must be less than upper bound ({})", - name, lower, upper + "Invalid bounds for parameter '{name}': lower bound ({lower}) must be less than upper bound ({upper})" ))); } if !initial_value.is_finite() { return Err(PyValueError::new_err(format!( - "Invalid initial value for parameter '{}': must be finite, got {}", - name, initial_value + "Invalid initial value for parameter '{name}': must be finite, got {initial_value}" ))); } } @@ -334,7 +337,7 @@ impl PyScalarBuilder { // Convert Option<(f64, f64)> to ParameterRange let range: diffid_core::problem::ParameterRange = - bounds.map(|b| b.into()).unwrap_or_else(|| Unbounded.into()); + bounds.map_or_else(|| Unbounded.into(), std::convert::Into::into); slf.state = match std::mem::replace( &mut slf.state, @@ -373,7 +376,7 @@ impl PyScalarBuilder { ScalarBuilderState::WithFunction { py_callable, .. } => { // Recreate a fresh builder from the Arc let objective = Arc::clone(py_callable); - let boxed_fn: Box f64 + Send + Sync> = + let boxed_fn: BoxedScalarFn = Box::new(move |x: &[f64]| objective.call(x).unwrap_or(f64::INFINITY)); let mut builder = ScalarProblemBuilder::new().with_function(boxed_fn); @@ -397,14 +400,13 @@ impl PyScalarBuilder { } => { // Recreate fresh builder with both function and gradient let objective = Arc::clone(py_callable); - let boxed_fn: Box f64 + Send + Sync> = + let boxed_fn: BoxedScalarFn = Box::new(move |x: &[f64]| objective.call(x).unwrap_or(f64::INFINITY)); let grad = Arc::clone(py_gradient); - let boxed_grad: Box Vec + Send + Sync> = - Box::new(move |x: &[f64]| { - grad.call(x).unwrap_or_else(|_| vec![f64::NAN; x.len()]) - }); + let boxed_grad: BoxedGradientFn = Box::new(move |x: &[f64]| { + grad.call(x).unwrap_or_else(|_| vec![f64::NAN; x.len()]) + }); let mut builder = ScalarProblemBuilder::new() .with_function(boxed_fn) @@ -463,7 +465,7 @@ fn convert_array_to_dmatrix(data: &PyReadonlyArrayDyn<'_, f64>) -> PyResult, - parameter_specs: Vec<(String, f64, Option<(f64, f64)>)>, + parameter_specs: Vec, config: HashMap, } @@ -494,7 +496,7 @@ impl PyDiffsolBuilder { self.__copy__() } - /// Register the DiffSL program describing the system dynamics. + /// Register the `DiffSL` program describing the system dynamics. fn with_diffsl(mut slf: PyRefMut<'_, Self>, dsl: String) -> PyRefMut<'_, Self> { slf.inner = std::mem::take(&mut slf.inner).with_diffsl(dsl); slf @@ -502,7 +504,7 @@ impl PyDiffsolBuilder { /// Attach observed data used to fit the differential equation. /// - /// The first column must contain the time samples (t_span) and the remaining + /// The first column must contain the time samples (`t_span`) and the remaining /// columns the observed trajectories. fn with_data<'py>( mut slf: PyRefMut<'py, Self>, @@ -531,8 +533,7 @@ impl PyDiffsolBuilder { "sparse" => DiffsolBackend::Sparse, other => { return Err(PyValueError::new_err(format!( - "Unknown backend '{}'. Expected 'dense' or 'sparse'", - other + "Unknown backend '{other}'. Expected 'dense' or 'sparse'" ))) } }; @@ -647,7 +648,7 @@ impl PyDiffsolBuilder { pub struct PyVectorBuilder { inner: VectorProblemBuilder, pub(crate) default_optimiser: Option, - parameter_specs: Vec<(String, f64, Option<(f64, f64)>)>, + parameter_specs: Vec, config: HashMap, } @@ -692,19 +693,15 @@ impl PyVectorBuilder { let params_array = PyArray1::from_slice(py, params); let result = objective.call1(py, (params_array,)).map_err( |e| -> Box { - Box::new(std::io::Error::new( - std::io::ErrorKind::Other, - format!("Objective call failed: {}", e), - )) + Box::new(std::io::Error::other(format!("Objective call failed: {e}"))) }, )?; let array: PyReadonlyArray1 = result.extract(py).map_err( |e| -> Box { - Box::new(std::io::Error::new( - std::io::ErrorKind::Other, - format!("Failed to extract array: {}", e), - )) + Box::new(std::io::Error::other(format!( + "Failed to extract array: {e}" + ))) }, )?; @@ -753,14 +750,12 @@ impl PyVectorBuilder { if let Some((lower, upper)) = bounds { if lower >= upper { return Err(PyValueError::new_err(format!( - "Invalid bounds for parameter '{}': lower bound ({}) must be less than upper bound ({})", - name, lower, upper + "Invalid bounds for parameter '{name}': lower bound ({lower}) must be less than upper bound ({upper})" ))); } if !initial_value.is_finite() { return Err(PyValueError::new_err(format!( - "Invalid initial value for parameter '{}': must be finite, got {}", - name, initial_value + "Invalid initial value for parameter '{name}': must be finite, got {initial_value}" ))); } } diff --git a/python/src/diffid/_diffid.pyi b/python/src/diffid/_diffid.pyi index 99d2b73..5d42420 100644 --- a/python/src/diffid/_diffid.pyi +++ b/python/src/diffid/_diffid.pyi @@ -3,11 +3,9 @@ import builtins import datetime -import typing - import numpy import numpy.typing - +import typing from diffid.sampler import DynamicNestedSampler, MetropolisHastings @typing.final @@ -57,12 +55,14 @@ class Adam: def init( self, initial: typing.Sequence[builtins.float], - bounds: typing.Sequence[tuple[builtins.float, builtins.float]] | None = None, + bounds: typing.Optional[ + typing.Sequence[tuple[builtins.float, builtins.float]] + ] = None, ) -> AdamState: r""" - Initialize ask-tell optimization state. + Initialize ask-tell optimisation state. - Returns an AdamState object that can be used for incremental optimization + Returns an `AdamState` object that can be used for incremental optimisation via the ask-tell interface. Parameters @@ -75,7 +75,7 @@ class Adam: Returns ------- AdamState - State object for ask-tell optimization + State object for ask-tell optimisation Examples -------- @@ -92,9 +92,9 @@ class Adam: @typing.final class AdamState: r""" - Ask-tell state for incremental Adam optimization. + Ask-tell state for incremental Adam optimisation. - This state object allows step-by-step control over the optimization process. + This state object allows step-by-step control over the optimisation process. Use `ask()` to get points to evaluate, and `tell()` to provide results. Examples @@ -112,7 +112,7 @@ class AdamState: """ def ask(self) -> typing.Any: r""" - Get the next action: evaluate points or optimization complete. + Get the next action: evaluate points or optimisation complete. Returns ------- @@ -126,7 +126,7 @@ class AdamState: >>> if isinstance(result, diffid.Evaluate): ... print(f"Need to evaluate {len(result.points)} points") >>> elif isinstance(result, diffid.Done): - ... print(f"Optimisation complete: {result.result}") + ... print(f"optimisation complete: {result.result}") """ def tell( self, result: tuple[builtins.float, typing.Sequence[builtins.float]] @@ -142,9 +142,9 @@ class AdamState: Raises ------ - TellError - If called after optimization has terminated or if result format is invalid - EvaluationError + `TellError` + If called after optimisation has terminated or if result format is invalid + `EvaluationError` If the evaluation failed or contained invalid values Examples @@ -176,14 +176,14 @@ class AdamState: """ def best( self, - ) -> tuple[builtins.list[builtins.float], builtins.float] | None: + ) -> typing.Optional[tuple[builtins.list[builtins.float], builtins.float]]: r""" Get the current best point and value found so far. Returns ------- tuple[list[float], float] | None - (best_point, best_value) or None if no valid evaluations yet + (`best_point`, `best_value`) or None if no valid evaluations yet """ def current_position(self) -> builtins.list[builtins.float]: r""" @@ -244,12 +244,14 @@ class CMAES: def init( self, initial: typing.Sequence[builtins.float], - bounds: typing.Sequence[tuple[builtins.float, builtins.float]] | None = None, + bounds: typing.Optional[ + typing.Sequence[tuple[builtins.float, builtins.float]] + ] = None, ) -> CMAESState: r""" - Initialize ask-tell optimization state. + Initialize ask-tell optimisation state. - Returns a CMAESState object that can be used for incremental optimization + Returns a `CMAESState` object that can be used for incremental optimisation via the ask-tell interface. Parameters @@ -262,7 +264,7 @@ class CMAES: Returns ------- CMAESState - State object for ask-tell optimization + State object for ask-tell optimisation Examples -------- @@ -279,9 +281,9 @@ class CMAES: @typing.final class CMAESState: r""" - Ask-tell state for incremental CMA-ES optimization. + Ask-tell state for incremental CMA-ES optimisation. - This state object allows step-by-step control over the optimization process. + This state object allows step-by-step control over the optimisation process. Use `ask()` to get a population of points to evaluate, and `tell()` to provide results. Examples @@ -298,7 +300,7 @@ class CMAESState: """ def ask(self) -> typing.Any: r""" - Get the next action: evaluate points or optimization complete. + Get the next action: evaluate points or optimisation complete. Returns ------- @@ -309,7 +311,7 @@ class CMAESState: Notes ----- CMA-ES evaluates a population of points each iteration. The number - of points returned depends on the population_size setting. + of points returned depends on the `population_size` setting. """ def tell(self, results: typing.Sequence[builtins.float]) -> None: r""" @@ -319,14 +321,14 @@ class CMAESState: ---------- results : list[float] List of objective function values corresponding to the points - from the last ask() call. Must match the number of points. + from the last `ask()` call. Must match the number of points. Raises ------ - TellError - If called after optimization has terminated or if wrong number + `TellError` + If called after optimisation has terminated or if wrong number of results provided - EvaluationError + `EvaluationError` If evaluations failed or contained invalid values """ def iterations(self) -> builtins.int: @@ -349,14 +351,14 @@ class CMAESState: """ def best( self, - ) -> tuple[builtins.list[builtins.float], builtins.float] | None: + ) -> typing.Optional[tuple[builtins.list[builtins.float], builtins.float]]: r""" Get the current best point and value found so far. Returns ------- tuple[list[float], float] | None - (best_point, best_value) or None if no valid evaluations yet + (`best_point`, `best_value`) or None if no valid evaluations yet """ def mean(self) -> builtins.list[builtins.float]: r""" @@ -401,13 +403,13 @@ class DiffsolBuilder: def __deepcopy__(self, _memo: dict) -> DiffsolBuilder: ... def with_diffsl(self, dsl: builtins.str) -> DiffsolBuilder: r""" - Register the DiffSL program describing the system dynamics. + Register the `DiffSL` program describing the system dynamics. """ def with_data(self, data: numpy.typing.NDArray[numpy.float64]) -> DiffsolBuilder: r""" Attach observed data used to fit the differential equation. - The first column must contain the time samples (t_span) and the remaining + The first column must contain the time samples (`t_span`) and the remaining columns the observed trajectories. """ def remove_data(self) -> DiffsolBuilder: @@ -418,7 +420,9 @@ class DiffsolBuilder: r""" Choose whether to use dense or sparse diffusion solvers. """ - def with_parallel(self, parallel: builtins.bool | None = None) -> DiffsolBuilder: + def with_parallel( + self, parallel: typing.Optional[builtins.bool] = None + ) -> DiffsolBuilder: r""" Opt into parallel proposal generation when supported by the backend. """ @@ -435,7 +439,7 @@ class DiffsolBuilder: self, name: builtins.str, initial_value: builtins.float, - bounds: tuple[builtins.float, builtins.float] | None = None, + bounds: typing.Optional[tuple[builtins.float, builtins.float]] = None, ) -> DiffsolBuilder: r""" Register a named optimisation variable in the order it appears in vectors. @@ -563,12 +567,14 @@ class NelderMead: def init( self, initial: typing.Sequence[builtins.float], - bounds: typing.Sequence[tuple[builtins.float, builtins.float]] | None = None, + bounds: typing.Optional[ + typing.Sequence[tuple[builtins.float, builtins.float]] + ] = None, ) -> NelderMeadState: r""" - Initialize ask-tell optimization state. + Initialize ask-tell optimisation state. - Returns a NelderMeadState object that can be used for incremental optimization + Returns a `NelderMeadState` object that can be used for incremental optimisation via the ask-tell interface. Parameters @@ -581,7 +587,7 @@ class NelderMead: Returns ------- NelderMeadState - State object for ask-tell optimization + State object for ask-tell optimisation Examples -------- @@ -598,9 +604,9 @@ class NelderMead: @typing.final class NelderMeadState: r""" - Ask-tell state for incremental Nelder-Mead optimization. + Ask-tell state for incremental Nelder-Mead optimisation. - This state object allows step-by-step control over the optimization process. + This state object allows step-by-step control over the optimisation process. Use `ask()` to get points to evaluate, and `tell()` to provide results. Examples @@ -617,7 +623,7 @@ class NelderMeadState: """ def ask(self) -> typing.Any: r""" - Get the next action: evaluate points or optimization complete. + Get the next action: evaluate points or optimisation complete. Returns ------- @@ -636,9 +642,9 @@ class NelderMeadState: Raises ------ - TellError - If called after optimization has terminated - EvaluationError + `TellError` + If called after optimisation has terminated + `EvaluationError` If the evaluation failed or contained invalid values """ def iterations(self) -> builtins.int: @@ -661,14 +667,14 @@ class NelderMeadState: """ def best( self, - ) -> tuple[builtins.list[builtins.float], builtins.float] | None: + ) -> typing.Optional[tuple[builtins.list[builtins.float], builtins.float]]: r""" Get the current best point and value from the simplex. Returns ------- tuple[list[float], float] | None - (best_point, best_value) or None if no valid evaluations yet + (`best_point`, `best_value`) or None if no valid evaluations yet """ def __repr__(self) -> builtins.str: ... def __str__(self) -> builtins.str: ... @@ -676,14 +682,14 @@ class NelderMeadState: @typing.final class NestedSamplesIterator: r""" - Iterator for NestedSamples posterior + Iterator for `NestedSamples` posterior """ def __iter__(self) -> NestedSamplesIterator: ... def __next__( self, - ) -> ( - tuple[builtins.list[builtins.float], builtins.float, builtins.float] | None - ): ... + ) -> typing.Optional[ + tuple[builtins.list[builtins.float], builtins.float, builtins.float] + ]: ... @typing.final class OptimisationResults: @@ -698,7 +704,7 @@ class OptimisationResults: Returns ------- numpy.ndarray - Best parameter vector as a NumPy array + Best parameter vector as a `NumPy` array """ @property def value(self) -> builtins.float: @@ -748,7 +754,7 @@ class OptimisationResults: @property def covariance( self, - ) -> builtins.list[builtins.list[builtins.float]] | None: + ) -> typing.Optional[builtins.list[builtins.list[builtins.float]]]: r""" Estimated covariance of the search distribution, if available. """ @@ -762,7 +768,7 @@ class OptimisationResults: """ def __bool__(self) -> builtins.bool: r""" - Return truthiness based on optimization success. + Return truthiness based on optimisation success. Allows using `if result:` instead of `if result.success:`. """ @@ -778,27 +784,27 @@ class Problem: """ def evaluate_gradient( self, x: typing.Sequence[builtins.float] - ) -> builtins.list[builtins.float] | None: + ) -> typing.Optional[builtins.list[builtins.float]]: r""" Evaluate the gradient of the objective function at `x` if available. """ def optimise( self, - initial: typing.Sequence[builtins.float] | None = None, - optimiser: NelderMead | CMAES | Adam | None = None, + initial: typing.Optional[typing.Sequence[builtins.float]] = None, + optimiser: typing.Optional[NelderMead | CMAES | Adam] = None, ) -> OptimisationResults: r""" Solve the problem starting from `initial` using the supplied optimiser. """ def sample( self, - initial: typing.Sequence[builtins.float] | None = None, - sampler: MetropolisHastings | DynamicNestedSampler | None = None, + initial: typing.Optional[typing.Sequence[builtins.float]] = None, + sampler: typing.Optional[MetropolisHastings | DynamicNestedSampler] = None, ) -> typing.Any: r""" Sample from the problem starting from `initial` using the supplied sampler. """ - def get_config(self, _key: builtins.str) -> builtins.float | None: + def get_config(self, _key: builtins.str) -> typing.Optional[builtins.float]: r""" Return the numeric configuration value stored under `key` if present. """ @@ -816,7 +822,7 @@ class Problem: tuple[ builtins.str, builtins.float, - tuple[builtins.float, builtins.float] | None, + typing.Optional[tuple[builtins.float, builtins.float]], ] ]: ... def initial_values(self) -> builtins.list[builtins.float]: ... @@ -851,7 +857,7 @@ class SamplesIterator: def __iter__(self) -> SamplesIterator: ... def __next__( self, - ) -> builtins.list[builtins.list[builtins.float]] | None: ... + ) -> typing.Optional[builtins.list[builtins.list[builtins.float]]]: ... @typing.final class ScalarBuilder: @@ -880,7 +886,7 @@ class ScalarBuilder: self, name: builtins.str, initial_value: builtins.float, - bounds: tuple[builtins.float, builtins.float] | None = None, + bounds: typing.Optional[tuple[builtins.float, builtins.float]] = None, ) -> ScalarBuilder: r""" Register a named optimisation variable in the order it appears in vectors. @@ -919,7 +925,7 @@ class VectorBuilder: self, name: builtins.str, initial_value: builtins.float, - bounds: tuple[builtins.float, builtins.float] | None = None, + bounds: typing.Optional[tuple[builtins.float, builtins.float]] = None, ) -> VectorBuilder: r""" Register a named optimisation variable in the order it appears in vectors. diff --git a/python/src/diffid/errors.py b/python/src/diffid/errors.py index e22c963..187ea05 100644 --- a/python/src/diffid/errors.py +++ b/python/src/diffid/errors.py @@ -161,10 +161,10 @@ def __init__(self): __all__ = [ + "AlreadyTerminated", + "BuildError", "DiffidError", "EvaluationError", - "BuildError", - "TellError", "ResultCountMismatch", - "AlreadyTerminated", + "TellError", ] diff --git a/python/src/diffid/sampler.pyi b/python/src/diffid/sampler.pyi index 6054e40..5b7baa7 100644 --- a/python/src/diffid/sampler.pyi +++ b/python/src/diffid/sampler.pyi @@ -3,11 +3,9 @@ import builtins import datetime -import typing - import numpy import numpy.typing - +import typing from diffid._diffid import NestedSamplesIterator, Problem, SamplesIterator @typing.final @@ -27,17 +25,19 @@ class DynamicNestedSampler: def run( self, problem: Problem, - initial: typing.Sequence[builtins.float] | None = None, + initial: typing.Optional[typing.Sequence[builtins.float]] = None, ) -> NestedSamples: ... def init( self, initial: typing.Sequence[builtins.float], - bounds: typing.Sequence[tuple[builtins.float, builtins.float]] | None = None, + bounds: typing.Optional[ + typing.Sequence[tuple[builtins.float, builtins.float]] + ] = None, ) -> DynamicNestedSamplerState: r""" Initialize ask-tell sampling state. - Returns a DynamicNestedSamplerState object that can be used for incremental + Returns a `DynamicNestedSamplerState` object that can be used for incremental sampling via the ask-tell interface. Parameters @@ -92,7 +92,7 @@ class DynamicNestedSamplerState: ------- Evaluate | Done Either Evaluate(points) requiring function evaluations, - or Done(result) indicating completion with NestedSamples. + or Done(result) indicating completion with `NestedSamples`. """ def tell(self, results: typing.Sequence[builtins.float]) -> None: r""" @@ -105,7 +105,7 @@ class DynamicNestedSamplerState: Raises ------ - TellError + `TellError` If called after sampling has terminated or if wrong number of results provided """ @@ -146,12 +146,14 @@ class MetropolisHastings: def init( self, initial: typing.Sequence[builtins.float], - bounds: typing.Sequence[tuple[builtins.float, builtins.float]] | None = None, + bounds: typing.Optional[ + typing.Sequence[tuple[builtins.float, builtins.float]] + ] = None, ) -> MetropolisHastingsState: r""" Initialize ask-tell sampling state. - Returns a MetropolisHastingsState object that can be used for incremental + Returns a `MetropolisHastingsState` object that can be used for incremental sampling via the ask-tell interface. Parameters @@ -224,7 +226,7 @@ class MetropolisHastingsState: Raises ------ - TellError + `TellError` If called after sampling has terminated or if wrong number of results provided """ @@ -284,7 +286,7 @@ class NestedSamples: r""" Iterate over posterior samples. - Yields tuples of (position, log_likelihood, log_weight). + Yields tuples of (position, `log_likelihood`, `log_weight`). """ def __getitem__( self, idx: builtins.int @@ -295,12 +297,12 @@ class NestedSamples: Parameters ---------- idx : int - Sample index (0 to num_samples - 1) + Sample index (0 to `num_samples` - 1) Returns ------- tuple[list[float], float, float] - Tuple of (position, log_likelihood, log_weight) + Tuple of (position, `log_likelihood`, `log_weight`) """ @typing.final @@ -344,7 +346,7 @@ class Samples: Parameters ---------- idx : int - Chain index (0 to num_chains - 1) + Chain index (0 to `num_chains` - 1) Returns ------- diff --git a/python/src/errors.rs b/python/src/errors.rs index a56ad97..19722ef 100644 --- a/python/src/errors.rs +++ b/python/src/errors.rs @@ -10,41 +10,41 @@ fn get_exception_class<'py>(py: Python<'py>, name: &str) -> PyResult PyErr { Python::attach(|py| { match get_exception_class(py, "EvaluationError") { Ok(exc_class) => { - let message = format!("Evaluation failed: {}", err); + let message = format!("Evaluation failed: {err}"); // Create exception with message match exc_class.call1((message,)) { - Ok(exc_instance) => PyErr::from_value(exc_instance.into()), - Err(_) => PyValueError::new_err(format!("Evaluation failed: {}", err)), + Ok(exc_instance) => PyErr::from_value(exc_instance), + Err(_) => PyValueError::new_err(format!("Evaluation failed: {err}")), } } Err(_) => { // Fallback to ValueError if custom exception not available - PyValueError::new_err(format!("Evaluation failed: {}", err)) + PyValueError::new_err(format!("Evaluation failed: {err}")) } } }) } -/// Convert build errors to Python BuildError +/// Convert build errors to Python `BuildError` pub fn build_error_to_py(err: impl std::fmt::Display) -> PyErr { Python::attach(|py| match get_exception_class(py, "BuildError") { Ok(exc_class) => { - let message = format!("{}", err); + let message = format!("{err}"); match exc_class.call1((message,)) { - Ok(exc_instance) => PyErr::from_value(exc_instance.into()), - Err(_) => PyValueError::new_err(format!("{}", err)), + Ok(exc_instance) => PyErr::from_value(exc_instance), + Err(_) => PyValueError::new_err(format!("{err}")), } } - Err(_) => PyValueError::new_err(format!("{}", err)), + Err(_) => PyValueError::new_err(format!("{err}")), }) } -/// Convert Rust TellError to Python TellError hierarchy +/// Convert Rust `TellError` to Python `TellError` hierarchy pub fn tell_error_to_py(err: CoreTellError) -> PyErr { Python::attach(|py| { match err { @@ -53,23 +53,21 @@ pub fn tell_error_to_py(err: CoreTellError) -> PyErr { Ok(exc_class) => { // Call with expected and received arguments match exc_class.call1((expected, got)) { - Ok(exc_instance) => PyErr::from_value(exc_instance.into()), + Ok(exc_instance) => PyErr::from_value(exc_instance), Err(_) => PyValueError::new_err(format!( - "Expected {} evaluation results, but received {}", - expected, got + "Expected {expected} evaluation results, but received {got}" )), } } Err(_) => PyValueError::new_err(format!( - "Expected {} evaluation results, but received {}", - expected, got + "Expected {expected} evaluation results, but received {got}" )), } } CoreTellError::AlreadyTerminated => { match get_exception_class(py, "AlreadyTerminated") { Ok(exc_class) => match exc_class.call0() { - Ok(exc_instance) => PyErr::from_value(exc_instance.into()), + Ok(exc_instance) => PyErr::from_value(exc_instance), Err(_) => PyValueError::new_err( "Cannot provide results to an already terminated optimisation", ), @@ -83,13 +81,13 @@ pub fn tell_error_to_py(err: CoreTellError) -> PyErr { // For other TellError variants, use base TellError class match get_exception_class(py, "TellError") { Ok(exc_class) => { - let message = format!("Tell error: {}", err); + let message = format!("Tell error: {err}"); match exc_class.call1((message,)) { - Ok(exc_instance) => PyErr::from_value(exc_instance.into()), - Err(_) => PyValueError::new_err(format!("Tell error: {}", err)), + Ok(exc_instance) => PyErr::from_value(exc_instance), + Err(_) => PyValueError::new_err(format!("Tell error: {err}")), } } - Err(_) => PyValueError::new_err(format!("Tell error: {}", err)), + Err(_) => PyValueError::new_err(format!("Tell error: {err}")), } } } diff --git a/python/src/lib.rs b/python/src/lib.rs index d2452f5..f397075 100644 --- a/python/src/lib.rs +++ b/python/src/lib.rs @@ -41,17 +41,16 @@ type ParameterSpecEntry = (String, f64, Option<(f64, f64)>); // Import objective types for the problem enum use diffid_core::problem::{DiffsolObjective, ScalarObjective, VectorObjective}; +// Type aliases to reduce complexity warnings +type BoxedScalarFn = Box f64 + Send + Sync>; +type BoxedGradientFn = Box Vec + Send + Sync>; +type ScalarProblemType = Problem>; +type ScalarGradientProblemType = Problem>; + // Enum to hold different Problem types internally pub(crate) enum DynProblem { - Scalar(Problem f64 + Send + Sync>>>), - ScalarWithGradient( - Problem< - ScalarObjective< - Box f64 + Send + Sync>, - Box Vec + Send + Sync>, - >, - >, - ), + Scalar(ScalarProblemType), + ScalarWithGradient(ScalarGradientProblemType), Vector(Problem), Diffsol(Problem), } @@ -228,7 +227,7 @@ impl PyProblem { fn evaluate(&self, x: Vec) -> PyResult { self.inner.evaluate(&x).map_err(|e| { crate::errors::evaluation_error_to_py(diffid_core::errors::EvaluationError::message( - format!("{}", e), + format!("{e}"), )) }) } @@ -239,7 +238,7 @@ impl PyProblem { DynProblem::ScalarWithGradient(p) => match p.evaluate_with_gradient(&x) { Ok((_val, grad_opt)) => Ok(grad_opt), Err(e) => Err(crate::errors::evaluation_error_to_py( - diffid_core::errors::EvaluationError::message(format!("{}", e)), + diffid_core::errors::EvaluationError::message(format!("{e}")), )), }, _ => Ok(None), @@ -346,13 +345,13 @@ impl PyProblem { /// Return a detailed string representation of the problem. fn __repr__(&self) -> String { let dim = self.inner.dimension(); - format!("Problem(dimension={})", dim) + format!("Problem(dimension={dim})") } /// Return a concise string representation of the problem. fn __str__(&self) -> String { let dim = self.inner.dimension(); - format!("{}-dimensional problem", dim) + format!("{dim}-dimensional problem") } } diff --git a/python/src/optimisers.rs b/python/src/optimisers.rs index fdffb45..dd2feab 100644 --- a/python/src/optimisers.rs +++ b/python/src/optimisers.rs @@ -161,7 +161,7 @@ impl PyNelderMead { /// Initialize ask-tell optimisation state. /// - /// Returns a NelderMeadState object that can be used for incremental optimisation + /// Returns a `NelderMeadState` object that can be used for incremental optimisation /// via the ask-tell interface. /// /// Parameters @@ -192,9 +192,7 @@ impl PyNelderMead { initial: Vec, bounds: Option>, ) -> PyResult { - let bounds = bounds - .map(Bounds::new) - .unwrap_or_else(|| Bounds::unbounded_like(&initial)); + let bounds = bounds.map_or_else(|| Bounds::unbounded_like(&initial), Bounds::new); let (state, _initial_point) = self.inner.init(initial, bounds); Ok(PyNelderMeadState { inner: state }) @@ -288,7 +286,7 @@ impl PyCMAES { /// Initialize ask-tell optimisation state. /// - /// Returns a CMAESState object that can be used for incremental optimisation + /// Returns a `CMAESState` object that can be used for incremental optimisation /// via the ask-tell interface. /// /// Parameters @@ -315,9 +313,7 @@ impl PyCMAES { /// ... state.tell(values) #[pyo3(signature = (initial, bounds=None))] fn init(&self, initial: Vec, bounds: Option>) -> PyResult { - let bounds = bounds - .map(Bounds::new) - .unwrap_or_else(|| Bounds::unbounded_like(&initial)); + let bounds = bounds.map_or_else(|| Bounds::unbounded_like(&initial), Bounds::new); let (state, _initial_point) = self.inner.init(initial, bounds); Ok(PyCMAESState { inner: state }) @@ -406,7 +402,7 @@ impl PyAdam { /// Initialize ask-tell optimisation state. /// - /// Returns an AdamState object that can be used for incremental optimisation + /// Returns an `AdamState` object that can be used for incremental optimisation /// via the ask-tell interface. /// /// Parameters @@ -433,9 +429,7 @@ impl PyAdam { /// ... state.tell(values) #[pyo3(signature = (initial, bounds=None))] fn init(&self, initial: Vec, bounds: Option>) -> PyResult { - let bounds = bounds - .map(Bounds::new) - .unwrap_or_else(|| Bounds::unbounded_like(&initial)); + let bounds = bounds.map_or_else(|| Bounds::unbounded_like(&initial), Bounds::new); let (state, _initial_point) = self.inner.init(initial, bounds); Ok(PyAdamState { inner: state }) @@ -503,9 +497,9 @@ impl PyAdamState { /// /// Raises /// ------ - /// TellError + /// `TellError` /// If called after optimisation has terminated or if result format is invalid - /// EvaluationError + /// `EvaluationError` /// If the evaluation failed or contained invalid values /// /// Examples @@ -546,7 +540,7 @@ impl PyAdamState { /// Returns /// ------- /// tuple[list[float], float] | None - /// (best_point, best_value) or None if no valid evaluations yet + /// (`best_point`, `best_value`) or None if no valid evaluations yet fn best(&self) -> Option<(Vec, f64)> { self.inner .best() @@ -569,7 +563,7 @@ impl PyAdamState { self.inner.iterations(), self.inner.evaluations(), match self.inner.best() { - Some((_, value)) => format!("{:.6}", value), + Some((_, value)) => format!("{value:.6}"), None => "None".to_string(), } ) @@ -631,9 +625,9 @@ impl PyNelderMeadState { /// /// Raises /// ------ - /// TellError + /// `TellError` /// If called after optimisation has terminated - /// EvaluationError + /// `EvaluationError` /// If the evaluation failed or contained invalid values fn tell(&mut self, result: f64) -> PyResult<()> { self.inner.tell(result).map_err(tell_error_to_py) @@ -664,7 +658,7 @@ impl PyNelderMeadState { /// Returns /// ------- /// tuple[list[float], float] | None - /// (best_point, best_value) or None if no valid evaluations yet + /// (`best_point`, `best_value`) or None if no valid evaluations yet fn best(&self) -> Option<(Vec, f64)> { self.inner .best() @@ -677,7 +671,7 @@ impl PyNelderMeadState { self.inner.iterations(), self.inner.evaluations(), match self.inner.best() { - Some((_, value)) => format!("{:.6}", value), + Some((_, value)) => format!("{value:.6}"), None => "None".to_string(), } ) @@ -728,7 +722,7 @@ impl PyCMAESState { /// Notes /// ----- /// CMA-ES evaluates a population of points each iteration. The number - /// of points returned depends on the population_size setting. + /// of points returned depends on the `population_size` setting. fn ask(&self, py: Python<'_>) -> Py { match self.inner.ask() { AskResult::Evaluate(points) => Py::new(py, PyEvaluate { points }).unwrap().into_any(), @@ -744,14 +738,14 @@ impl PyCMAESState { /// ---------- /// results : list[float] /// List of objective function values corresponding to the points - /// from the last ask() call. Must match the number of points. + /// from the last `ask()` call. Must match the number of points. /// /// Raises /// ------ - /// TellError + /// `TellError` /// If called after optimisation has terminated or if wrong number /// of results provided - /// EvaluationError + /// `EvaluationError` /// If evaluations failed or contained invalid values fn tell(&mut self, results: Vec) -> PyResult<()> { self.inner.tell(results).map_err(tell_error_to_py) @@ -782,7 +776,7 @@ impl PyCMAESState { /// Returns /// ------- /// tuple[list[float], float] | None - /// (best_point, best_value) or None if no valid evaluations yet + /// (`best_point`, `best_value`) or None if no valid evaluations yet fn best(&self) -> Option<(Vec, f64)> { self.inner .best() @@ -815,7 +809,7 @@ impl PyCMAESState { self.inner.iterations(), self.inner.evaluations(), match self.inner.best() { - Some((_, value)) => format!("{:.6}", value), + Some((_, value)) => format!("{value:.6}"), None => "None".to_string(), }, self.inner.sigma() diff --git a/python/src/results.rs b/python/src/results.rs index 38deada..6b005a8 100644 --- a/python/src/results.rs +++ b/python/src/results.rs @@ -70,7 +70,7 @@ pub struct PyDone { impl PyDone { fn __repr__(&self, py: Python<'_>) -> PyResult { let result_repr = self.result.bind(py).repr()?.to_string(); - Ok(format!("Done(result={})", result_repr)) + Ok(format!("Done(result={result_repr})")) } fn __str__(&self) -> String { @@ -79,7 +79,7 @@ impl PyDone { } impl PyDone { - /// Create a Done variant with OptimisationResults + /// Create a Done variant with `OptimisationResults` pub fn with_optimisation_results(py: Python<'_>, results: OptimisationResults) -> Self { let py_results = PyOptimisationResults { inner: results }; Self { @@ -94,7 +94,7 @@ impl PyDone { } } - /// Create a Done variant with NestedSamples + /// Create a Done variant with `NestedSamples` pub fn with_nested_samples(py: Python<'_>, samples: PyNestedSamples) -> Self { Self { result: Py::new(py, samples).unwrap().into_any(), @@ -117,7 +117,7 @@ impl PyOptimisationResults { /// Returns /// ------- /// numpy.ndarray - /// Best parameter vector as a NumPy array + /// Best parameter vector as a `NumPy` array #[getter] fn x<'py>(&self, py: Python<'py>) -> Bound<'py, PyArray1> { self.inner.x.to_pyarray(py) diff --git a/python/src/samplers.rs b/python/src/samplers.rs index 1bfc918..91f3273 100644 --- a/python/src/samplers.rs +++ b/python/src/samplers.rs @@ -171,7 +171,7 @@ impl PySamples { /// Parameters /// ---------- /// idx : int - /// Chain index (0 to num_chains - 1) + /// Chain index (0 to `num_chains` - 1) /// /// Returns /// ------- @@ -305,7 +305,7 @@ impl PyNestedSamples { /// Iterate over posterior samples. /// - /// Yields tuples of (position, log_likelihood, log_weight). + /// Yields tuples of (position, `log_likelihood`, `log_weight`). fn __iter__(slf: PyRef<'_, Self>) -> PyResult { let posterior: Vec<_> = slf .inner @@ -331,12 +331,12 @@ impl PyNestedSamples { /// Parameters /// ---------- /// idx : int - /// Sample index (0 to num_samples - 1) + /// Sample index (0 to `num_samples` - 1) /// /// Returns /// ------- /// tuple[list[float], float, float] - /// Tuple of (position, log_likelihood, log_weight) + /// Tuple of (position, `log_likelihood`, `log_weight`) fn __getitem__(&self, idx: isize) -> PyResult<(Vec, f64, f64)> { let posterior = self.inner.posterior(); let len = posterior.len() as isize; @@ -366,7 +366,7 @@ impl PyNestedSamples { } } -/// Iterator for NestedSamples posterior +/// Iterator for `NestedSamples` posterior #[cfg_attr(feature = "stubgen", gen_stub_pyclass)] #[pyclass] struct NestedSamplesIterator { @@ -440,7 +440,7 @@ impl PyMetropolisHastings { /// Initialize ask-tell sampling state. /// - /// Returns a MetropolisHastingsState object that can be used for incremental + /// Returns a `MetropolisHastingsState` object that can be used for incremental /// sampling via the ask-tell interface. /// /// Parameters @@ -471,9 +471,7 @@ impl PyMetropolisHastings { initial: Vec, bounds: Option>, ) -> PyResult { - let bounds = bounds - .map(Bounds::new) - .unwrap_or_else(|| Bounds::unbounded_like(&initial)); + let bounds = bounds.map_or_else(|| Bounds::unbounded_like(&initial), Bounds::new); let state = self.inner.init(initial, bounds); Ok(PyMetropolisHastingsState { inner: state }) @@ -536,7 +534,7 @@ impl PyDynamicNestedSampler { /// Initialize ask-tell sampling state. /// - /// Returns a DynamicNestedSamplerState object that can be used for incremental + /// Returns a `DynamicNestedSamplerState` object that can be used for incremental /// sampling via the ask-tell interface. /// /// Parameters @@ -567,9 +565,7 @@ impl PyDynamicNestedSampler { initial: Vec, bounds: Option>, ) -> PyResult { - let bounds = bounds - .map(Bounds::new) - .unwrap_or_else(|| Bounds::unbounded_like(&initial)); + let bounds = bounds.map_or_else(|| Bounds::unbounded_like(&initial), Bounds::new); let (state, _initial_points) = self.inner.init(initial, bounds); Ok(PyDynamicNestedSamplerState { inner: state }) @@ -634,7 +630,7 @@ impl PyMetropolisHastingsState { /// /// Raises /// ------ - /// TellError + /// `TellError` /// If called after sampling has terminated or if wrong number /// of results provided fn tell(&mut self, results: Vec) -> PyResult<()> { @@ -709,7 +705,7 @@ impl PyDynamicNestedSamplerState { /// ------- /// Evaluate | Done /// Either Evaluate(points) requiring function evaluations, - /// or Done(result) indicating completion with NestedSamples. + /// or Done(result) indicating completion with `NestedSamples`. fn ask(&self, py: Python<'_>) -> Py { match self.inner.ask() { AskResult::Evaluate(points) => Py::new(py, PyEvaluate { points }).unwrap().into_any(), @@ -732,7 +728,7 @@ impl PyDynamicNestedSamplerState { /// /// Raises /// ------ - /// TellError + /// `TellError` /// If called after sampling has terminated or if wrong number /// of results provided fn tell(&mut self, results: Vec) -> PyResult<()> { diff --git a/rust/Cargo.toml b/rust/Cargo.toml index c27fb3a..1fce866 100644 --- a/rust/Cargo.toml +++ b/rust/Cargo.toml @@ -12,6 +12,9 @@ keywords = ["optimisation", "time-series", "differential-equations", "sampler", categories = ["algorithms", "science", "simulation"] include = ["src/**", "benches/**", "Cargo.toml", "../README.md", "../LICENSE"] +[lints] +workspace = true + [dependencies] diffsol = { version = "0.10.1" } nalgebra.workspace = true diff --git a/rust/benches/diffsol_benches.rs b/rust/benches/diffsol_benches.rs index 7a8b94f..046716b 100644 --- a/rust/benches/diffsol_benches.rs +++ b/rust/benches/diffsol_benches.rs @@ -6,13 +6,13 @@ use std::time::Duration; macro_rules! build_logistic_problem { ($backend:expr, $parallel:expr) => {{ - let dsl = r#" + let dsl = r" in_i { r = 1, k = 1 } r { 1 } k { 1 } u_i { y = 0.1 } F_i { (r * y) * (1 - (y / k)) } -"#; +"; let t_span: Vec = (0..20).map(|i| i as f64 * 0.1).collect(); let data_values: Vec = t_span.iter().map(|t| 0.1 * (*t).exp()).collect(); @@ -41,7 +41,7 @@ fn bench_diffsol_single_eval(c: &mut Criterion) { for backend in [DiffsolBackend::Dense, DiffsolBackend::Sparse] { let problem = build_logistic_problem!(backend, false); group.bench_with_input( - BenchmarkId::new("evaluate", format!("{:?}", backend)), + BenchmarkId::new("evaluate", format!("{backend:?}")), &problem, |b, problem| { b.iter(|| { @@ -62,14 +62,14 @@ fn bench_diffsol_population_eval(c: &mut Criterion) { let population: Vec> = (0..64) .map(|i| { - let scale = 0.8 + (i as f64) * 0.01; + let scale = f64::from(i).mul_add(0.01, 0.8); vec![1.0 * scale, 1.0 / scale] }) .collect(); for backend in [DiffsolBackend::Dense, DiffsolBackend::Sparse] { for ¶llel in &[false, true] { - let label = format!("{:?}_parallel={}", backend, parallel); + let label = format!("{backend:?}_parallel={parallel}"); let problem = build_logistic_problem!(backend, parallel); group.bench_with_input( diff --git a/rust/src/builders/diffsol.rs b/rust/src/builders/diffsol.rs index 5fcdff8..10abd21 100644 --- a/rust/src/builders/diffsol.rs +++ b/rust/src/builders/diffsol.rs @@ -16,7 +16,8 @@ pub enum DiffsolBackend { Sparse, } -#[derive(Debug, Clone)] +#[must_use] +#[derive(Debug, Clone, Copy)] pub struct DiffsolConfig { pub rtol: f64, pub atol: f64, @@ -76,6 +77,7 @@ impl DiffsolConfig { } } +#[must_use] #[derive(Clone)] pub struct DiffsolProblemBuilder { equations: Option, @@ -109,7 +111,7 @@ impl DiffsolProblemBuilder { self } - /// Registers the DiffSL differential equation system. + /// Registers the `DiffSL` differential equation system. pub fn with_diffsl(mut self, equations: String) -> Self { self.equations = Some(equations); self @@ -193,6 +195,13 @@ impl DiffsolProblemBuilder { } /// Build the problem + /// + /// # Errors + /// + /// Returns an error if: + /// - No data has been provided via `with_data` + /// - Data has fewer than 2 columns (needs time + at least one observation) + /// - No differential equations have been provided via `with_diffsl` pub fn build(self) -> Result, ProblemBuilderError> { // Unpack data and verify let data_with_t = self.data.as_ref().ok_or(ProblemBuilderError::MissingData)?; @@ -202,7 +211,7 @@ impl DiffsolProblemBuilder { got: data_with_t.ncols(), }); } - let t_span: Vec = data_with_t.column(0).iter().cloned().collect(); + let t_span: Vec = data_with_t.column(0).iter().copied().collect(); let data = data_with_t.columns(1, data_with_t.ncols() - 1).into_owned(); // Check costs and provide default if empty diff --git a/rust/src/builders/mod.rs b/rust/src/builders/mod.rs index b80faab..93fa505 100644 --- a/rust/src/builders/mod.rs +++ b/rust/src/builders/mod.rs @@ -20,11 +20,11 @@ impl std::fmt::Display for ProblemBuilderError { fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result { match self { Self::MissingData => write!(f, "Missing data"), - Self::BuildFailed(msg) => write!(f, "build failed: {}", msg), + Self::BuildFailed(msg) => write!(f, "build failed: {msg}"), Self::MissingSystem => write!(f, "Missing system"), Self::MissingVectorFn => write!(f, "Missing vector function"), Self::DimensionMismatch { expected, got } => { - write!(f, "expected {} elements, got {}", expected, got) + write!(f, "expected {expected} elements, got {got}") } } } @@ -36,11 +36,11 @@ impl From for ProblemBuilderError { ProblemError::DimensionMismatch { expected, got } => { ProblemBuilderError::DimensionMismatch { expected, got } } - ProblemError::BuildFailed(msg) => ProblemBuilderError::BuildFailed(msg), + ProblemError::BuildFailed(msg) + | ProblemError::EvaluationFailed(msg) + | ProblemError::SolverError(msg) => ProblemBuilderError::BuildFailed(msg), // For other variants, wrap them in BuildFailed ProblemError::External(e) => ProblemBuilderError::BuildFailed(e.to_string()), - ProblemError::EvaluationFailed(msg) => ProblemBuilderError::BuildFailed(msg), - ProblemError::SolverError(msg) => ProblemBuilderError::BuildFailed(msg), } } } @@ -56,11 +56,11 @@ mod tests { #[test] fn test_diffsol_builder() { - let dsl = r#" + let dsl = r" in_i { r = 1, k = 1} u_i { y = 0.1 } F_i { (r * y) * (1 - (y / k)) } -"#; +"; let t_span: Vec = (0..5).map(|i| i as f64 * 0.1).collect(); let data_values = vec![0.1, 0.2, 0.3, 0.4, 0.5]; @@ -130,10 +130,10 @@ F_i { (r * y) * (1 - (y / k)) } // Perfect fit should have zero cost (a=1, b=1 gives [1,2,3,4,5]) let cost = problem.evaluate(&[1.0, 1.0]).expect("evaluation failed"); - assert!(cost.abs() < 1e-10, "expected near-zero cost, got {}", cost); + assert!(cost.abs() < 1e-10, "expected near-zero cost, got {cost}"); // Non-perfect fit should have positive cost let cost = problem.evaluate(&[0.5, 0.5]).expect("evaluation failed"); - assert!(cost > 0.0, "expected positive cost, got {}", cost); + assert!(cost > 0.0, "expected positive cost, got {cost}"); } } diff --git a/rust/src/builders/scalar.rs b/rust/src/builders/scalar.rs index 702ef98..99f228f 100644 --- a/rust/src/builders/scalar.rs +++ b/rust/src/builders/scalar.rs @@ -3,6 +3,7 @@ use crate::optimisers::Optimiser; use crate::prelude::{ParameterSpec, Problem}; use crate::problem::{NoFunction, NoGradient, ParameterRange, ScalarObjective}; +#[must_use] #[derive(Clone)] pub struct ScalarProblemBuilder { f: F, @@ -87,6 +88,11 @@ where F: Fn(&[f64]) -> f64 + Send + Sync + 'static, { /// Build the problem + /// + /// # Errors + /// + /// Currently this method always succeeds, but returns `Result` for consistency + /// with other builder `build` methods that may fail. pub fn build(self) -> Result>, ProblemBuilderError> { // Build objective let objective = ScalarObjective::new(self.f); @@ -102,6 +108,12 @@ where F: Fn(&[f64]) -> f64 + Send + Sync + 'static, G: Fn(&[f64]) -> Vec + Send + Sync + 'static, { + /// Build the problem + /// + /// # Errors + /// + /// Currently this method always succeeds, but returns `Result` for consistency + /// with other builder `build` methods that may fail. pub fn build(self) -> Result>, ProblemBuilderError> { // Build objective let objective = ScalarObjective::with_gradient(self.f, self.gradient); diff --git a/rust/src/builders/vector.rs b/rust/src/builders/vector.rs index 3b475da..731b3a4 100644 --- a/rust/src/builders/vector.rs +++ b/rust/src/builders/vector.rs @@ -5,6 +5,7 @@ use crate::prelude::{ParameterSpec, Problem}; use crate::problem::{ParameterRange, VectorFn, VectorObjective}; use std::sync::Arc; +#[must_use] #[derive(Clone)] pub struct VectorProblemBuilder { function: Option, @@ -96,6 +97,13 @@ impl VectorProblemBuilder { } /// Build the problem + /// + /// # Errors + /// + /// Returns an error if: + /// - No function has been provided via `with_function` + /// - No data has been provided via `with_data` + /// - The objective fails to build (e.g., dimension mismatches) pub fn build(self) -> Result, ProblemBuilderError> { // Default costs let mut costs = self.costs; diff --git a/rust/src/common.rs b/rust/src/common.rs index 874fc89..f9864a8 100644 --- a/rust/src/common.rs +++ b/rust/src/common.rs @@ -180,7 +180,7 @@ impl Bounds { let width = hi - lo; let sigma = width * scale; let draw = rng.sample::(StandardNormal); - (base + draw * sigma).clamp(lo, hi) + draw.mul_add(sigma, base).clamp(lo, hi) } else { // At least one bound is non-finite let base = if lo.is_finite() { @@ -192,7 +192,7 @@ impl Bounds { }; let sigma = scale; let offset = rng.sample::(StandardNormal); - base + offset * sigma + offset.mul_add(sigma, base) } }) .collect() diff --git a/rust/src/cost/mod.rs b/rust/src/cost/mod.rs index e38c7d1..6266d98 100644 --- a/rust/src/cost/mod.rs +++ b/rust/src/cost/mod.rs @@ -197,8 +197,7 @@ impl GaussianNll { pub fn new(weight: Option, variance: f64) -> Self { assert!( variance > 0.0 && variance.is_finite(), - "Variance must be positive and finite, got {}", - variance + "Variance must be positive and finite, got {variance}" ); let log_term = (2.0 * PI * variance).ln(); @@ -227,7 +226,7 @@ impl CostMetric for GaussianNll { let sse: f64 = residuals.iter().map(|&r| r * r).sum(); // NLL = (n/2) * ln(2Ļ€ĻƒĀ²) + (1/2ĻƒĀ²) * Ī£rĀ² - (0.5 * n * self.log_term + 0.5 * sse / self.variance) * self.weight + (0.5 * n).mul_add(self.log_term, 0.5 * sse / self.variance) * self.weight } fn name(&self) -> &'static str { @@ -246,7 +245,7 @@ impl CostMetric for GaussianNll { let n = residuals.len() as f64; let sse: f64 = residuals.iter().map(|&r| r * r).sum(); - let cost = (0.5 * n * self.log_term + 0.5 * sse / self.variance) * self.weight; + let cost = (0.5 * n).mul_add(self.log_term, 0.5 * sse / self.variance) * self.weight; if sensitivities.is_empty() { return Some((cost, Vec::new())); @@ -313,7 +312,7 @@ mod tests { // Create a 2x1 sensitivity matrix with values [0.5, 0.5] let triplets = vec![(0, 0, 0.5), (1, 0, 0.5)]; let sens_matrix: NalgebraMat = - Matrix::try_from_triplets(2, 1, triplets, Default::default()).unwrap(); + Matrix::try_from_triplets(2, 1, triplets, diffsol::NalgebraContext).unwrap(); let (cost, grad) = metric .evaluate_with_sensitivities(&residuals, &[sens_matrix]) @@ -347,7 +346,7 @@ mod tests { // Create a 2x1 sensitivity matrix with values [0.5, 0.5] let triplets = vec![(0, 0, 0.5), (1, 0, 0.5)]; let sens_matrix: NalgebraMat = - Matrix::try_from_triplets(2, 1, triplets, Default::default()).unwrap(); + Matrix::try_from_triplets(2, 1, triplets, diffsol::NalgebraContext).unwrap(); let (cost, grad) = metric .evaluate_with_sensitivities(&residuals, &[sens_matrix]) @@ -395,7 +394,7 @@ mod tests { // Create a 2x1 sensitivity matrix with values [0.5, 0.5] let triplets = vec![(0, 0, 0.5), (1, 0, 0.5)]; let sens_matrix: NalgebraMat = - Matrix::try_from_triplets(2, 1, triplets, Default::default()).unwrap(); + Matrix::try_from_triplets(2, 1, triplets, diffsol::NalgebraContext).unwrap(); let (cost, grad) = metric .evaluate_with_sensitivities(&residuals, &[sens_matrix]) diff --git a/rust/src/errors.rs b/rust/src/errors.rs index be24fb1..1ae16ed 100644 --- a/rust/src/errors.rs +++ b/rust/src/errors.rs @@ -34,7 +34,7 @@ impl EvaluationError { impl fmt::Display for EvaluationError { fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result { match self { - Self::User(e) => write!(f, "Evaluation failed:: {}", e), + Self::User(e) => write!(f, "Evaluation failed:: {e}"), Self::NonFiniteValue => write!(f, "Evaluation failed::NonFiniteValue"), Self::NonFiniteGradient => write!(f, "Evaluation failed::NonFiniteGradient"), } @@ -45,8 +45,7 @@ impl std::error::Error for EvaluationError { fn source(&self) -> Option<&(dyn std::error::Error + 'static)> { match self { Self::User(e) => Some(e.as_ref()), - Self::NonFiniteValue => None, - Self::NonFiniteGradient => None, + Self::NonFiniteValue | Self::NonFiniteGradient => None, } } } @@ -78,7 +77,7 @@ pub enum TellError { ResultCountMismatch { expected: usize, got: usize }, /// Gradient dimension doesn't match point dimension GradientDimensionMismatch { expected: usize, got: usize }, - /// Wrapper for EvaluationError + /// Wrapper for `EvaluationError` EvaluationFailed(EvaluationError), } @@ -89,12 +88,12 @@ impl fmt::Display for TellError { match self { TellError::AlreadyTerminated => write!(f, "Algorithm already terminated"), TellError::ResultCountMismatch { expected, got } => { - write!(f, "Expected {} results, got {}", expected, got) + write!(f, "Expected {expected} results, got {got}") } TellError::GradientDimensionMismatch { expected, got } => { - write!(f, "Expected gradient dimension {}, got {}", expected, got) + write!(f, "Expected gradient dimension {expected}, got {got}") } - TellError::EvaluationFailed(e) => write!(f, "Evaluation failed {:?}", e), + TellError::EvaluationFailed(e) => write!(f, "Evaluation failed {e:?}"), } } } diff --git a/rust/src/optimisers/adam.rs b/rust/src/optimisers/adam.rs index 27812a4..666d654 100644 --- a/rust/src/optimisers/adam.rs +++ b/rust/src/optimisers/adam.rs @@ -7,6 +7,7 @@ use std::error::Error as StdError; use std::time::{Duration, Instant}; /// Configuration for the Adam optimiser +#[must_use] #[derive(Clone, Debug)] pub struct Adam { max_iter: usize, @@ -152,9 +153,9 @@ impl MomentumState { for (i, g) in gradient.iter().enumerate() { // Update biased first moment estimate - self.m[i] = beta1 * self.m[i] + (1.0 - beta1) * g; + self.m[i] = beta1.mul_add(self.m[i], (1.0 - beta1) * g); // Update biased second moment estimate - self.v[i] = beta2 * self.v[i] + (1.0 - beta2) * g * g; + self.v[i] = beta2.mul_add(self.v[i], (1.0 - beta2) * g * g); // Compute bias-corrected estimates let m_hat = self.m[i] / bias_correction1; @@ -224,6 +225,12 @@ impl AdamState { /// Report the evaluation result (value and gradient) for the last point from `ask()` /// /// Pass `Err` if the objective function failed to evaluate + /// + /// # Errors + /// + /// Returns an error if: + /// - The optimizer has already terminated + /// - The gradient dimension doesn't match the problem dimension pub fn tell(&mut self, result: T) -> Result<(), TellError> where T: TryInto, @@ -241,7 +248,7 @@ impl AdamState { self.history .push(EvaluatedPoint::new(self.x.clone(), f64::NAN)); self.phase = AdamPhase::Terminated(TerminationReason::FunctionEvaluationFailed( - format!("{}", err), + format!("{err}"), )); return Ok(()); } @@ -284,6 +291,12 @@ impl AdamState { } /// Get the current best point and value + /// + /// # Panics + /// + /// This method will not panic in practice, as it only compares finite values. + /// The `unwrap()` is safe because `partial_cmp` only returns `None` for NaN comparisons, + /// which are filtered out by the `is_finite()` check. pub fn best(&self) -> Option<(&[f64], f64)> { self.history .iter() @@ -403,6 +416,11 @@ impl Adam { /// Run optimisation using a closure for evaluation /// /// The closure should return `(value, gradient)` for a given point + /// + /// # Panics + /// + /// Panics if the optimisation state machine enters an unexpected state after + /// a tell error. This should not occur under normal operation. pub fn run( &self, mut objective: F, @@ -434,7 +452,7 @@ impl Adam { match state.ask() { AskResult::Done(results) => results, - _ => panic!("Unexpected state after tell error"), + AskResult::Evaluate(_) => panic!("Unexpected state after tell error"), } } @@ -509,6 +527,7 @@ mod tests { } /// fallible sphere function + #[allow(clippy::unnecessary_wraps)] fn sphere_fallible(x: &[f64]) -> Result<(f64, Vec), std::io::Error> { let value: f64 = x.iter().map(|xi| xi * xi).sum(); let grad: Vec = x.iter().map(|xi| 2.0 * xi).collect(); @@ -856,7 +875,7 @@ mod tests { "Tell after termination should return AlreadyTerminated error" ); } - _ => panic!("Should have terminated after max_iter=1"), + AskResult::Evaluate(_) => panic!("Should have terminated after max_iter=1"), } } @@ -886,7 +905,9 @@ mod tests { TerminationReason::FunctionEvaluationFailed(_) )); } - _ => panic!("Should have terminated due to non-finite gradient"), + AskResult::Evaluate(_) => { + panic!("Should have terminated due to non-finite gradient") + } } } AskResult::Done(_) => panic!("Should not terminate after first evaluation"), @@ -1039,8 +1060,7 @@ mod tests { for val in &points[0] { assert!( *val >= -2.0 && *val <= 2.0, - "Point {} violates bounds [-2.0, 2.0]", - val + "Point {val} violates bounds [-2.0, 2.0]" ); } result = sphere_infallible(&points[0]); @@ -1109,7 +1129,7 @@ mod tests { TerminationReason::FunctionEvaluationFailed(_) )); } - _ => panic!("Should terminate after evaluation error"), + AskResult::Evaluate(_) => panic!("Should terminate after evaluation error"), } } AskResult::Done(_) => panic!("Should not terminate after first evaluation"), diff --git a/rust/src/optimisers/cmaes.rs b/rust/src/optimisers/cmaes.rs index c3ae1c0..99559bb 100644 --- a/rust/src/optimisers/cmaes.rs +++ b/rust/src/optimisers/cmaes.rs @@ -1,5 +1,3 @@ -#![allow(unexpected_cfgs)] - use crate::common::{AskResult, Bounds, Point}; use crate::errors::{EvaluationError, TellError}; use crate::optimisers::{ @@ -12,6 +10,7 @@ use std::cmp::Ordering; use std::time::{Duration, Instant}; /// Configuration for the CMA-ES optimiser +#[must_use] #[derive(Clone, Debug)] pub struct CMAES { max_iter: usize, @@ -69,7 +68,7 @@ impl CMAES { fn compute_population_size(&self, dim: usize) -> usize { self.population_size.unwrap_or_else(|| { if dim > 0 { - let suggested = (4.0 + (3.0 * (dim as f64).ln())).floor() as usize; + let suggested = 3.0f64.mul_add((dim as f64).ln(), 4.0).floor() as usize; suggested.max(4).max(2 * dim) } else { 4 @@ -146,9 +145,9 @@ impl StrategyParameters { let c_sigma = (mu_eff + 2.0) / (dim_f + mu_eff + 5.0); let d_sigma = Self::compute_d_sigma(mu_eff, dim_f, c_sigma); let c_c = (4.0 + mu_eff / dim_f) / (dim_f + 4.0 + 2.0 * mu_eff / dim_f); - let c1 = 2.0 / ((dim_f + 1.3).powi(2) + mu_eff); + let c1 = 2.0 / (dim_f + 1.3).mul_add(dim_f + 1.3, mu_eff); let c_mu = ((1.0 - c1) - .min(2.0 * (mu_eff - 2.0 + 1.0 / mu_eff) / ((dim_f + 2.0).powi(2) + mu_eff))) + .min(2.0 * (mu_eff - 2.0 + 1.0 / mu_eff) / (dim_f + 2.0).mul_add(dim_f + 2.0, mu_eff))) .max(0.0); let chi_n = dim_f.sqrt() * (1.0 - 1.0 / (4.0 * dim_f) + 1.0 / (21.0 * dim_f.powi(2))); @@ -168,7 +167,7 @@ impl StrategyParameters { pub fn compute_d_sigma(mu_eff: f64, dim_f: f64, c_sigma: f64) -> f64 { let sqrt_term = ((mu_eff - 1.0) / (dim_f + 1.0)).max(0.0).sqrt(); - 1.0 + c_sigma + 2.0 * (sqrt_term - 1.0).max(0.0) + 2.0f64.mul_add((sqrt_term - 1.0).max(0.0), 1.0 + c_sigma) } } @@ -253,10 +252,9 @@ impl CMAESState { let mut distribution = DistributionState::new(&initial_point); distribution.sigma = config.step_size.max(1e-12); - let rng = match config.seed { - Some(seed) => StdRng::seed_from_u64(seed), - None => StdRng::from_os_rng(), - }; + let rng = config + .seed + .map_or_else(StdRng::from_os_rng, StdRng::seed_from_u64); Self { config, @@ -288,6 +286,12 @@ impl CMAESState { } /// Report the evaluation results for the points from `ask()` + /// + /// # Errors + /// + /// Returns an error if: + /// - The optimiser has already terminated + /// - The number of results doesn't match the expected population size pub fn tell(&mut self, results: I) -> Result<(), TellError> where I: IntoIterator, @@ -301,13 +305,8 @@ impl CMAESState { // We collect into a Result> first to handle errors early let values: Vec = results .into_iter() - .map(|r| r.try_into()) - .map(|res| { - match res { - Ok(eval) => eval.value(), - Err(_) => f64::INFINITY, // Error as infinite cost - } - }) + .map(std::convert::TryInto::try_into) + .map(|res| res.map_or(f64::INFINITY, |eval| eval.value())) .collect(); // Take ownership of current phase @@ -320,14 +319,14 @@ impl CMAESState { match phase { CMAESPhase::EvaluatingInitial { initial_point } => { - self.handle_initial_evaluated(initial_point, values)?; + self.handle_initial_evaluated(&initial_point, &values)?; } CMAESPhase::AwaitingPopulation { candidates, z_vectors, old_mean, } => { - self.handle_population_evaluated(candidates, z_vectors, old_mean, values)?; + self.handle_population_evaluated(candidates, z_vectors, &old_mean, values)?; } CMAESPhase::Terminated(_) => unreachable!(), } @@ -370,8 +369,8 @@ impl CMAESState { // Phase Handlers fn handle_initial_evaluated( &mut self, - initial_point: Point, - results: Vec, + initial_point: &Point, + results: &[f64], ) -> Result<(), TellError> { if results.len() != 1 { return Err(TellError::ResultCountMismatch { @@ -400,7 +399,7 @@ impl CMAESState { &mut self, candidates: Vec, z_vectors: Vec>, - old_mean: DVector, + old_mean: &DVector, results: Vec, ) -> Result<(), TellError> { let expected = candidates.len(); @@ -441,7 +440,7 @@ impl CMAESState { } // Update CMA-ES state - let termination_reason = self.update_distribution(&population, &old_mean); + let termination_reason = self.update_distribution(&population, old_mean); // Update final population self.final_population = population.iter().map(|(pt, _)| pt.clone()).collect(); @@ -505,7 +504,7 @@ impl CMAESState { let step = &self.distribution.eigenvectors * (&step_matrix * &z); let candidate_vec = &self.distribution.mean + step * self.distribution.sigma; - let mut candidate: Point = candidate_vec.iter().cloned().collect(); + let mut candidate: Point = candidate_vec.iter().copied().collect(); // Apply bounds self.bounds.clamp(&mut candidate); @@ -666,6 +665,11 @@ impl CMAES { /// Run optimisation using a closure for evaluation /// /// This is a convenience wrapper around the ask/tell interface + /// + /// # Panics + /// + /// Panics if the optimisation state machine enters an unexpected state after + /// a tell error. This should not occur under normal operation. pub fn run( &self, mut objective: F, @@ -698,10 +702,18 @@ impl CMAES { match state.ask() { AskResult::Done(opt_results) => opt_results, - _ => panic!("Unexpected state after tell error"), + AskResult::Evaluate(_) => panic!("Unexpected state after tell error"), } } + /// Run optimisation using batch evaluation + /// + /// Similar to `run`, but evaluates all population points in a single call + /// + /// # Panics + /// + /// Panics if the optimisation state machine enters an unexpected state after + /// a tell error. This should not occur under normal operation. pub fn run_batch( &self, objective: F, @@ -715,7 +727,7 @@ impl CMAES { { let (mut state, first_point) = self.init(initial, bounds); - let mut results = objective(&vec![first_point]); + let mut results = objective(&[first_point]); loop { if state.tell(results).is_err() { @@ -734,7 +746,7 @@ impl CMAES { match state.ask() { AskResult::Done(opt_results) => opt_results, - _ => panic!("Unexpected state after tell error"), + AskResult::Evaluate(_) => panic!("Unexpected state after tell error"), } } } @@ -745,6 +757,7 @@ mod tests { use crate::builders::ScalarProblemBuilder; use std::convert::Infallible; + #[allow(clippy::unnecessary_wraps)] fn sphere(x: &[f64]) -> Result { Ok(x.iter().map(|xi| xi * xi).sum()) } @@ -762,7 +775,7 @@ mod tests { match state.tell(current_results) { Ok(()) => {} Err(TellError::AlreadyTerminated) => break, - Err(e) => panic!("Unexpected error: {:?}", e), + Err(e) => panic!("Unexpected error: {e:?}"), } match state.ask() { @@ -869,7 +882,7 @@ mod tests { .zip(&covariance) .for_each(|(row_i, row_j)| { row_i.iter().zip(row_j).for_each(|(a, b)| { - assert!((a - b).abs() < 1e-12, "covariance matrix must be symmetric") + assert!((a - b).abs() < 1e-12, "covariance matrix must be symmetric"); }); }); @@ -882,8 +895,7 @@ mod tests { assert!( eigenvalues.iter().all(|&eig| eig >= -1e-10), - "covariance must be positive semi-definite: {:?}", - eigenvalues + "covariance must be positive semi-definite: {eigenvalues:?}" ); } @@ -904,9 +916,7 @@ mod tests { assert!( (computed - expected).abs() < 1e-12, - "d_sigma mismatch: expected {}, got {}", - expected, - computed + "d_sigma mismatch: expected {expected}, got {computed}" ); // For this case, the sqrt term is less than 1, so it should clamp to 0 @@ -946,7 +956,7 @@ mod tests { let updated = CMAESState::update_covariance(&cov, c1, c_mu, &p_c, h_sigma, c_c, &rank_mu); for (exp, got) in expected.iter().zip(updated.iter()) { - assert!((exp - got).abs() < 1e-12, "expected {} got {}", exp, got); + assert!((exp - got).abs() < 1e-12, "expected {exp} got {got}"); } } @@ -968,10 +978,11 @@ mod tests { let updated = CMAESState::update_covariance(&cov, c1, c_mu, &p_c, h_sigma, c_c, &rank_mu); for (exp, got) in expected.iter().zip(updated.iter()) { - assert!((exp - got).abs() < 1e-12, "expected {} got {}", exp, got); + assert!((exp - got).abs() < 1e-12, "expected {exp} got {got}"); } } + #[allow(clippy::unnecessary_wraps)] fn sphere_cmaes(x: &[f64]) -> Result { Ok(x.iter().map(|xi| xi * xi).sum()) } @@ -1183,8 +1194,7 @@ mod tests { for &val in &first_point[..] { assert!( (-2.0..=2.0).contains(&val), - "Initial point {:?} violates bounds", - first_point + "Initial point {first_point:?} violates bounds" ); } @@ -1200,8 +1210,7 @@ mod tests { for &val in point { assert!( (-2.0..=2.0).contains(&val), - "Point {:?} violates bounds", - point + "Point {point:?} violates bounds" ); } } @@ -1375,6 +1384,7 @@ mod tests { } } + #[allow(clippy::unnecessary_wraps)] fn rosenbrock_cmaes(x: &[f64]) -> Result { let a = 1.0; let b = 100.0; diff --git a/rust/src/optimisers/mod.rs b/rust/src/optimisers/mod.rs index a071ae5..54badde 100644 --- a/rust/src/optimisers/mod.rs +++ b/rust/src/optimisers/mod.rs @@ -81,6 +81,11 @@ impl ScalarOptimiser { /// * `initial` - Initial point (will be auto-expanded if needed) /// * `bounds` - Optional parameter bounds /// + /// # Panics + /// + /// May panic when using Nelder-Mead if the objective function returns an empty vector. + /// This should not occur under normal operation. + /// /// # Returns /// Optimisation results including best point, value, and diagnostics /// @@ -111,7 +116,7 @@ impl ScalarOptimiser { ScalarOptimiser::CMAES(cm) => cm.run_batch(objective, initial, bounds), ScalarOptimiser::NelderMead(nm) => nm.run( |x| { - let result = objective(&vec![x.to_vec()]); + let result = objective(&[x.to_vec()]); result.into_iter().next().unwrap() // ToDO: This needs proper error integration }, initial, @@ -203,7 +208,7 @@ impl Optimiser { pub fn as_scalar(&self) -> Option<&ScalarOptimiser> { match self { Optimiser::Scalar(opt) => Some(opt), - _ => None, + Optimiser::Gradient(_) => None, } } @@ -211,7 +216,7 @@ impl Optimiser { pub fn as_scalar_mut(&mut self) -> Option<&mut ScalarOptimiser> { match self { Optimiser::Scalar(opt) => Some(opt), - _ => None, + Optimiser::Gradient(_) => None, } } @@ -221,7 +226,7 @@ impl Optimiser { pub fn as_gradient(&self) -> Option<&GradientOptimiser> { match self { Optimiser::Gradient(opt) => Some(opt), - _ => None, + Optimiser::Scalar(_) => None, } } @@ -229,7 +234,7 @@ impl Optimiser { pub fn as_gradient_mut(&mut self) -> Option<&mut GradientOptimiser> { match self { Optimiser::Gradient(opt) => Some(opt), - _ => None, + Optimiser::Scalar(_) => None, } } @@ -237,10 +242,14 @@ impl Optimiser { /// /// Returns `Ok(ScalarOptimiser)` on success, or `Err(self)` if this is not /// a scalar optimiser. + /// + /// # Errors + /// + /// Returns `Err(self)` if this optimiser is a gradient optimiser, not a scalar optimiser. pub fn into_scalar(self) -> Result { match self { Optimiser::Scalar(opt) => Ok(opt), - other => Err(other), + other @ Optimiser::Gradient(_) => Err(other), } } @@ -248,10 +257,14 @@ impl Optimiser { /// /// Returns `Ok(GradientOptimiser)` on success, or `Err(self)` if this is not /// a gradient optimiser. + /// + /// # Errors + /// + /// Returns `Err(self)` if this optimiser is a scalar optimiser, not a gradient optimiser. pub fn into_gradient(self) -> Result { match self { Optimiser::Gradient(opt) => Ok(opt), - other => Err(other), + other @ Optimiser::Scalar(_) => Err(other), } } @@ -450,7 +463,7 @@ impl fmt::Display for TerminationReason { write!(f, "Patience elapsed") } TerminationReason::FunctionEvaluationFailed(msg) => { - write!(f, "Function evaluation failed: {}", msg) + write!(f, "Function evaluation failed: {msg}") } } } @@ -565,12 +578,12 @@ mod tests { match nm { ScalarOptimiser::NelderMead(_) => (), - _ => panic!("Expected NelderMead variant"), + ScalarOptimiser::CMAES(_) => panic!("Expected NelderMead variant"), } match cm { ScalarOptimiser::CMAES(_) => (), - _ => panic!("Expected CMAES variant"), + ScalarOptimiser::NelderMead(_) => panic!("Expected CMAES variant"), } } @@ -890,12 +903,14 @@ mod tests { let result = optimiser.run( |x| -> Result<(f64, Vec), ProblemError> { let (value, gradient_opt) = problem.evaluate_with_gradient(x)?; - match gradient_opt { - Some(g) => Ok((value, g)), - None => Err(ProblemError::EvaluationFailed( - "Adam optimiser requires an available gradient".to_string(), - )), - } + gradient_opt.map_or_else( + || { + Err(ProblemError::EvaluationFailed( + "Adam optimiser requires an available gradient".to_string(), + )) + }, + |g| Ok((value, g)), + ) }, vec![1.0, 2.0], Bounds::unbounded(2), @@ -906,7 +921,7 @@ mod tests { TerminationReason::FunctionEvaluationFailed(ref msg) => { assert!(msg.contains("requires an available gradient")); } - other => panic!("expected FunctionEvaluationFailed, got {:?}", other), + other => panic!("expected FunctionEvaluationFailed, got {other:?}"), } } @@ -1023,9 +1038,7 @@ mod tests { for (x1, x2) in result_one.x.iter().zip(result_two.x.iter()) { assert!( (x1 - x2).abs() < 1e-10, - "expected identical optima: {} vs {}", - x1, - x2 + "expected identical optima: {x1} vs {x2}" ); } assert_eq!(result_one.covariance, result_two.covariance); @@ -1121,7 +1134,7 @@ mod tests { .zip(&covariance) .for_each(|(row_i, row_j)| { row_i.iter().zip(row_j).for_each(|(a, b)| { - assert!((a - b).abs() < 1e-12, "covariance matrix must be symmetric") + assert!((a - b).abs() < 1e-12, "covariance matrix must be symmetric"); }); }); @@ -1134,8 +1147,7 @@ mod tests { assert!( eigenvalues.iter().all(|&eig| eig >= -1e-10), - "covariance must be positive semi-definite: {:?}", - eigenvalues + "covariance must be positive semi-definite: {eigenvalues:?}" ); } @@ -1156,9 +1168,7 @@ mod tests { assert!( (computed - expected).abs() < 1e-12, - "d_sigma mismatch: expected {}, got {}", - expected, - computed + "d_sigma mismatch: expected {expected}, got {computed}" ); // For this case, the sqrt term is less than 1, so it should clamp to 0 @@ -1198,7 +1208,7 @@ mod tests { let updated = CMAESState::update_covariance(&cov, c1, c_mu, &p_c, h_sigma, c_c, &rank_mu); for (exp, got) in expected.iter().zip(updated.iter()) { - assert!((exp - got).abs() < 1e-12, "expected {} got {}", exp, got); + assert!((exp - got).abs() < 1e-12, "expected {exp} got {got}"); } } @@ -1220,7 +1230,7 @@ mod tests { let updated = CMAESState::update_covariance(&cov, c1, c_mu, &p_c, h_sigma, c_c, &rank_mu); for (exp, got) in expected.iter().zip(updated.iter()) { - assert!((exp - got).abs() < 1e-12, "expected {} got {}", exp, got); + assert!((exp - got).abs() < 1e-12, "expected {exp} got {got}"); } } } diff --git a/rust/src/optimisers/nelder_mead.rs b/rust/src/optimisers/nelder_mead.rs index ce52862..a846220 100644 --- a/rust/src/optimisers/nelder_mead.rs +++ b/rust/src/optimisers/nelder_mead.rs @@ -7,6 +7,7 @@ use std::cmp::Ordering; use std::time::{Duration, Instant}; /// Configuration for the Nelder-Mead optimiser +#[must_use] #[derive(Clone, Debug)] pub struct NelderMead { max_iter: usize, @@ -101,6 +102,11 @@ impl NelderMead { /// Run optimisation using a closure for evaluation /// /// This is a convenience wrapper around the ask/tell interface + /// + /// # Panics + /// + /// Panics if the optimisation state machine enters an unexpected state after + /// a tell error. This should not occur under normal operation. pub fn run( &self, mut objective: F, @@ -133,7 +139,7 @@ impl NelderMead { // Should not reach here normally match state.ask() { AskResult::Done(results) => results, - _ => panic!("Unexpected state after tell error"), + AskResult::Evaluate(_) => panic!("Unexpected state after tell error"), } } } @@ -211,7 +217,8 @@ impl NelderMeadState { NelderMeadPhase::EvaluatingInitial => { AskResult::Evaluate(vec![self.initial_point.clone()]) } - NelderMeadPhase::BuildingSimplex { pending_point, .. } => { + NelderMeadPhase::BuildingSimplex { pending_point, .. } + | NelderMeadPhase::Shrinking { pending_point, .. } => { AskResult::Evaluate(vec![pending_point.clone()]) } NelderMeadPhase::AwaitingReflection { @@ -223,15 +230,16 @@ impl NelderMeadState { NelderMeadPhase::AwaitingContraction { contract_point, .. } => { AskResult::Evaluate(vec![contract_point.clone()]) } - NelderMeadPhase::Shrinking { pending_point, .. } => { - AskResult::Evaluate(vec![pending_point.clone()]) - } } } /// Report the evaluation result for the last point from `ask()` /// /// Pass `Err` if the objective function failed to evaluate + /// + /// # Errors + /// + /// Returns an error if the optimiser has already terminated. pub fn tell(&mut self, result: T) -> Result<(), TellError> where T: TryInto, @@ -246,7 +254,7 @@ impl NelderMeadState { Err(e) => { let err: EvaluationError = e.into(); self.phase = NelderMeadPhase::Terminated( - TerminationReason::FunctionEvaluationFailed(format!("{}", err)), + TerminationReason::FunctionEvaluationFailed(format!("{err}")), ); return Ok(()); } @@ -543,10 +551,10 @@ impl NelderMeadState { fn compute_simplex_vertex(&self, dim: usize) -> Point { let mut point = self.initial_point.clone(); - if point[dim] != 0.0 { - point[dim] *= 1.0 + self.config.step_size; - } else { + if point[dim] == 0.0 { point[dim] = self.config.step_size; + } else { + point[dim] *= 1.0 + self.config.step_size; } // Ensure the point differs from the initial point @@ -691,6 +699,7 @@ mod tests { use super::*; use std::convert::Infallible; + #[allow(clippy::unnecessary_wraps)] fn rosenbrock(x: &[f64]) -> Result { let a = 1.0; let b = 100.0; @@ -712,8 +721,7 @@ mod tests { loop { match state.tell(current_value) { Ok(()) => {} - Err(TellError::AlreadyTerminated) => break, - _ => break, + Err(_) => break, } match state.ask() { @@ -740,6 +748,7 @@ mod tests { assert!(results.value < 1e-6); } + #[allow(clippy::unnecessary_wraps)] fn sphere(x: &[f64]) -> Result { Ok(x.iter().map(|xi| xi * xi).sum()) } @@ -928,8 +937,7 @@ mod tests { for &val in point { assert!( (-1.0..=1.0).contains(&val), - "Point {:?} violates bounds", - point + "Point {point:?} violates bounds" ); } } diff --git a/rust/src/problem/diffsol_problem.rs b/rust/src/problem/diffsol_problem.rs index 8339c43..130d14d 100644 --- a/rust/src/problem/diffsol_problem.rs +++ b/rust/src/problem/diffsol_problem.rs @@ -89,7 +89,7 @@ impl DiffsolObjective { .rtol(self.config.rtol) .build_from_diffsl(&self.dsl) .map_err(|e| { - ProblemError::BuildFailed(format!("Failed to build ODE model: {}", e)) + ProblemError::BuildFailed(format!("Failed to build ODE model: {e}")) })?; Ok(DiffsolSimulator::Dense(Box::new(diff_system))) } @@ -99,7 +99,7 @@ impl DiffsolObjective { .rtol(self.config.rtol) .build_from_diffsl(&self.dsl) .map_err(|e| { - ProblemError::BuildFailed(format!("Failed to build ODE model: {}", e)) + ProblemError::BuildFailed(format!("Failed to build ODE model: {e}")) })?; Ok(DiffsolSimulator::Sparse(Box::new(diff_system))) } @@ -228,27 +228,25 @@ impl Objective for DiffsolObjective { self.with_simulator_cached(x, |problem| match problem { DiffsolSimulator::Dense(p) => { let solver_result = p.bdf::(); - let mut solver = match solver_result { - Ok(s) => s, - Err(_) => return Ok(FAILED_SOLVE_PENALTY), + let Ok(mut solver) = solver_result else { + return Ok(FAILED_SOLVE_PENALTY); }; - match Self::solve_safely(|| solver.solve_dense(&self.t_span)) { - Ok(solution) => self.calculate_cost(&solution), - Err(_) => Ok(FAILED_SOLVE_PENALTY), - } + Self::solve_safely(|| solver.solve_dense(&self.t_span)).map_or_else( + |_| Ok(FAILED_SOLVE_PENALTY), + |solution| self.calculate_cost(&solution), + ) } DiffsolSimulator::Sparse(p) => { let solver_result = p.bdf::(); - let mut solver = match solver_result { - Ok(s) => s, - Err(_) => return Ok(FAILED_SOLVE_PENALTY), + let Ok(mut solver) = solver_result else { + return Ok(FAILED_SOLVE_PENALTY); }; - match Self::solve_safely(|| solver.solve_dense(&self.t_span)) { - Ok(solution) => self.calculate_cost(&solution), - Err(_) => Ok(FAILED_SOLVE_PENALTY), - } + Self::solve_safely(|| solver.solve_dense(&self.t_span)).map_or_else( + |_| Ok(FAILED_SOLVE_PENALTY), + |solution| self.calculate_cost(&solution), + ) } }) } @@ -305,13 +303,12 @@ mod tests { use super::*; use crate::cost::{GaussianNll, SumSquaredError}; - #[allow(dead_code)] fn build_logistic_problem(backend: DiffsolBackend) -> DiffsolObjective { - let dsl = r#" + let dsl = r" in_i {r = 1, k = 1 } u_i { y = 0.1 } F_i { (r * y) * (1 - (y / k)) } -"#; +"; let t_span: Vec = (0..6).map(|i| i as f64 * 0.2).collect(); let data_values: Vec = t_span.iter().map(|t| 0.1 * (*t).exp()).collect(); @@ -364,7 +361,7 @@ F_i { (r * y) * (1 - (y / k)) } // Build a 2x1 sensitivity matrix (2 elements) which mismatches the 3 residuals let triplets = vec![(0, 0, 0.5), (1, 0, 0.5)]; let wrong_size_sens: NalgebraMat = - Matrix::try_from_triplets(2, 1, triplets, Default::default()).unwrap(); + Matrix::try_from_triplets(2, 1, triplets, diffsol::NalgebraContext).unwrap(); let result = std::panic::catch_unwind(|| { metric @@ -384,7 +381,7 @@ F_i { (r * y) * (1 - (y / k)) } let sensitivities: Vec> = (0..3) .map(|param_idx| { let triplets = vec![(param_idx, 0, 1.0)]; - Matrix::try_from_triplets(3, 1, triplets, Default::default()).unwrap() + Matrix::try_from_triplets(3, 1, triplets, diffsol::NalgebraContext).unwrap() }) .collect(); @@ -415,6 +412,7 @@ F_i { (r * y) * (1 - (y / k)) } let eps = 1e-5_f64; // Compare against finite-difference approximation of problem.evaluate + #[allow(clippy::needless_range_loop)] for i in 0..params.len() { let mut params_fd = params; let fd = finite_difference(&mut params_fd, i, eps, |p| { @@ -427,11 +425,7 @@ F_i { (r * y) * (1 - (y / k)) } let diff = (fd - g).abs(); assert!( diff < 1e-6, - "gradient mismatch for param {}: fd={} grad={} diff={}", - i, - fd, - g, - diff + "gradient mismatch for param {i}: fd={fd} grad={g} diff={diff}" ); } } diff --git a/rust/src/problem/mod.rs b/rust/src/problem/mod.rs index 2e74c8b..607ae83 100644 --- a/rust/src/problem/mod.rs +++ b/rust/src/problem/mod.rs @@ -27,7 +27,7 @@ pub enum ProblemError { impl From for ProblemError { fn from(e: DiffsolError) -> Self { - ProblemError::BuildFailed(format!("{}", e)) + ProblemError::BuildFailed(format!("{e}")) } } @@ -40,13 +40,13 @@ impl From> for ProblemError { impl std::fmt::Display for ProblemError { fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result { match self { - Self::EvaluationFailed(msg) => write!(f, "evaluation failed: {}", msg), - Self::SolverError(msg) => write!(f, "solver failed: {}", msg), + Self::EvaluationFailed(msg) => write!(f, "evaluation failed: {msg}"), + Self::SolverError(msg) => write!(f, "solver failed: {msg}"), Self::DimensionMismatch { expected, got } => { - write!(f, "expected {} elements, got {}", expected, got) + write!(f, "expected {expected} elements, got {got}") } - Self::BuildFailed(msg) => write!(f, "build failed: {}", msg), - Self::External(err) => write!(f, "external error: {}", err), + Self::BuildFailed(msg) => write!(f, "build failed: {msg}"), + Self::External(err) => write!(f, "external error: {err}"), } } } @@ -111,11 +111,11 @@ impl ParameterSet { } pub fn push(&mut self, spec: ParameterSpec) { - self.0.push(spec) + self.0.push(spec); } pub fn clear(&mut self) { - self.0.clear() + self.0.clear(); } pub fn take(&mut self) -> Vec { @@ -148,6 +148,15 @@ impl ParameterSet { } } +impl<'a> IntoIterator for &'a ParameterSet { + type Item = &'a ParameterSpec; + type IntoIter = std::slice::Iter<'a, ParameterSpec>; + + fn into_iter(self) -> Self::IntoIter { + self.iter() + } +} + #[derive(Debug, Clone, PartialEq)] pub struct ParameterRange(RangeInclusive); @@ -191,7 +200,7 @@ impl From for ParameterRange { } } -/// Implement for RangeInclusive +/// Implement for `RangeInclusive` impl From> for ParameterRange { fn from(range: RangeInclusive) -> Self { Self(range) @@ -201,6 +210,15 @@ impl From> for ParameterRange { /// An Objective trait, used to define the core /// evaluation of a `problem`. pub trait Objective: Send + Sync { + /// Evaluate the objective function at point `x` + /// + /// # Errors + /// + /// Returns [`ProblemError`] if evaluation fails: + /// - [`ProblemError::EvaluationFailed`]: Objective function computation fails + /// - [`ProblemError::SolverError`]: ODE/DAE solver fails to converge + /// - [`ProblemError::DimensionMismatch`]: Input or output dimensions are incorrect + /// - [`ProblemError::External`]: External function propagates an error fn evaluate(&self, x: &[f64]) -> Result; fn gradient(&self, _x: &[f64]) -> Option> { @@ -212,6 +230,15 @@ pub trait Objective: Send + Sync { false } + /// Evaluate the objective function and its gradient at point `x` + /// + /// # Errors + /// + /// Returns [`ProblemError`] if evaluation fails: + /// - [`ProblemError::EvaluationFailed`]: Objective or gradient computation fails + /// - [`ProblemError::SolverError`]: ODE/DAE solver fails to converge + /// - [`ProblemError::DimensionMismatch`]: Input or output dimensions are incorrect + /// - [`ProblemError::External`]: External function propagates an error fn evaluate_with_gradient(&self, x: &[f64]) -> Result<(f64, Option>), ProblemError> { Ok((self.evaluate(x)?, self.gradient(x))) } @@ -296,6 +323,11 @@ pub struct VectorObjective { } impl VectorObjective { + /// Create a new vector objective + /// + /// # Errors + /// + /// Returns [`ProblemError::BuildFailed`] if the data vector is empty. pub fn new( f: VectorFn, data: Vec, @@ -352,10 +384,20 @@ impl Problem { self.objective.has_gradient() } + /// Evaluate the problem at point `x` + /// + /// # Errors + /// + /// Returns [`ProblemError`] if evaluation fails. See [`Objective::evaluate`] for details. pub fn evaluate(&self, x: &[f64]) -> Result { self.objective.evaluate(x) } + /// Evaluate the problem and its gradient at point `x` + /// + /// # Errors + /// + /// Returns [`ProblemError`] if evaluation fails. See [`Objective::evaluate_with_gradient`] for details. pub fn evaluate_with_gradient( &self, x: &[f64], @@ -418,12 +460,14 @@ impl Problem { grad_opt.run( |x| { let (value, grad) = self.evaluate_with_gradient(x)?; - match grad { - Some(grad) => Ok((value, grad)), - None => Err(ProblemError::EvaluationFailed( - "Gradient optimiser requires gradient".to_string(), - )), - } + grad.map_or_else( + || { + Err(ProblemError::EvaluationFailed( + "Gradient optimiser requires gradient".to_string(), + )) + }, + |grad| Ok((value, grad)), + ) }, x0, self.parameters.bounds(), @@ -510,8 +554,7 @@ mod tests { .expect("evaluation failed"); assert!( cost.abs() < 1e-10, - "expected near-zero cost with true params, got {}", - cost + "expected near-zero cost with true params, got {cost}" ); // Test with wrong parameters @@ -613,8 +656,7 @@ mod tests { let cost = problem.evaluate(&[1.0]).expect("evaluation failed"); assert!( (cost - 1.0).abs() < 1e-10, - "expected RMSE of 1.0, got {}", - cost + "expected RMSE of 1.0, got {cost}" ); } diff --git a/rust/src/sampler/dynamic_nested/mcmc_proposal.rs b/rust/src/sampler/dynamic_nested/mcmc_proposal.rs index e872afb..353af92 100644 --- a/rust/src/sampler/dynamic_nested/mcmc_proposal.rs +++ b/rust/src/sampler/dynamic_nested/mcmc_proposal.rs @@ -86,7 +86,7 @@ impl MCMCProposalGenerator { /// Compute per-dimension step sizes based on bounds width /// /// For each dimension i, calculates: - /// scale[i] = step_size_factor Ɨ (upper[i] - lower[i]) + /// scale[i] = `step_size_factor` Ɨ (upper[i] - lower[i]) /// /// For unbounded dimensions, uses a default scale of 1.0. /// This makes step sizes adaptive to the natural scale of each parameter. diff --git a/rust/src/sampler/dynamic_nested/mod.rs b/rust/src/sampler/dynamic_nested/mod.rs index 845a5ac..b595d3e 100644 --- a/rust/src/sampler/dynamic_nested/mod.rs +++ b/rust/src/sampler/dynamic_nested/mod.rs @@ -65,6 +65,7 @@ pub struct DynamicNestedSamplerState { } /// Configurable Dynamic Nested Sampling engine +#[must_use] #[derive(Clone, Debug)] pub struct DynamicNestedSampler { live_points: usize, @@ -151,12 +152,11 @@ impl DynamicNestedSampler { /// * `bounds` - Parameter bounds for sampling /// /// # Returns - /// Tuple of (state, initial_candidates) where candidates should be evaluated + /// Tuple of (state, `initial_candidates`) where candidates should be evaluated pub fn init(&self, _initial: Point, bounds: Bounds) -> (DynamicNestedSamplerState, Vec) { - let mut rng = match self.seed { - Some(seed) => StdRng::seed_from_u64(seed), - None => StdRng::from_rng(&mut rand::rng()), - }; + let mut rng = self + .seed + .map_or_else(|| StdRng::from_rng(&mut rand::rng()), StdRng::seed_from_u64); let dimension = bounds.dimension(); @@ -230,8 +230,8 @@ impl DynamicNestedSamplerState { } DNSPhase::AwaitingReplacementBatch { pending_positions, .. - } => AskResult::Evaluate(pending_positions.clone()), - DNSPhase::AwaitingExpansion { + } + | DNSPhase::AwaitingExpansion { pending_positions, .. } => AskResult::Evaluate(pending_positions.clone()), } @@ -243,7 +243,13 @@ impl DynamicNestedSamplerState { /// Errors are treated as infinite values (rejected). /// /// # Arguments - /// * `results` - Evaluation results matching the last ask() request + /// * `results` - Evaluation results matching the last `ask()` request + /// + /// # Errors + /// + /// Returns an error if: + /// - The sampler has already terminated + /// - The number of results doesn't match the expected count /// /// # Returns /// `Ok(())` on success, or `TellError` if already terminated or result count mismatch @@ -267,16 +273,13 @@ impl DynamicNestedSamplerState { // Convert results to Vec, treating errors as INFINITY let values: Vec = results .into_iter() - .map(|r| match r.try_into() { - Ok(eval) => eval.value(), - Err(_) => f64::INFINITY, - }) + .map(|r| r.try_into().map_or(f64::INFINITY, |eval| eval.value())) .collect(); // Dispatch to phase handler match &self.phase { DNSPhase::InitialisingLivePoints { .. } => self.handle_initialisation(values), - DNSPhase::AwaitingReplacementBatch { .. } => self.handle_replacement_batch(values), + DNSPhase::AwaitingReplacementBatch { .. } => self.handle_replacement_batch(&values), DNSPhase::AwaitingExpansion { .. } => self.handle_expansion(values), DNSPhase::Terminated(_) => unreachable!("Already checked above"), } @@ -358,7 +361,7 @@ impl DynamicNestedSamplerState { /// /// Evaluates a batch of MCMC proposals and selects the best valid one (likelihood > threshold). /// If all proposals are rejected, the removed point is restored. - fn handle_replacement_batch(&mut self, values: Vec) -> Result<(), TellError> { + fn handle_replacement_batch(&mut self, values: &[f64]) -> Result<(), TellError> { // Extract phase data let (positions, removed, threshold) = match &self.phase { DNSPhase::AwaitingReplacementBatch { @@ -439,6 +442,11 @@ impl DynamicNestedSamplerState { } /// Start next iteration of the main sampling loop + /// + /// Note: Returns Result for API consistency across sampler methods, + /// though this function currently never errors. Maintains uniform + /// interface for potential future error cases. + #[allow(clippy::unnecessary_wraps)] fn start_next_iteration(&mut self) -> Result<(), TellError> { // Check termination conditions if self.sampler_state.live_points().is_empty() { @@ -568,6 +576,10 @@ impl DynamicNestedSampler { /// nested sampling calculations. /// * `initial` - Initial point (currently unused, reserved for future) /// * `bounds` - Parameter bounds + /// + /// # Panics + /// + /// Panics if `tell()` fails during the sampling loop, which indicates a bug in the sampler. pub fn run(&self, mut objective: F, initial: Point, bounds: Bounds) -> NestedSamples where F: FnMut(&[f64]) -> R, @@ -586,7 +598,9 @@ impl DynamicNestedSampler { results = points.iter().map(|p| objective(p)).collect(); } AskResult::Done(SamplingResults::Nested(samples)) => return samples, - _ => unreachable!("DynamicNestedSampler always returns Nested results"), + AskResult::Done(_) => { + unreachable!("DynamicNestedSampler always returns Nested results") + } } } } @@ -602,6 +616,11 @@ impl DynamicNestedSampler { /// nested sampling calculations. /// * `initial` - Initial point (currently unused, reserved for future) /// * `bounds` - Parameter bounds + /// + /// # Panics + /// + /// Panics if the sampling state machine enters an unexpected state after a tell error. + /// This should not occur under normal operation. pub fn run_batch(&self, objective: F, initial: Point, bounds: Bounds) -> NestedSamples where F: Fn(&[Vec]) -> Vec, @@ -622,7 +641,9 @@ impl DynamicNestedSampler { results = objective(&points); } AskResult::Done(SamplingResults::Nested(samples)) => return samples, - _ => unreachable!("DynamicNestedSampler always returns Nested results"), + AskResult::Done(_) => { + unreachable!("DynamicNestedSampler always returns Nested results") + } } } @@ -702,8 +723,7 @@ mod tests { let mean = nested.mean()[0]; assert!( mean.is_finite(), - "posterior mean must be finite, got {:.4}", - mean + "posterior mean must be finite, got {mean:.4}" ); // The posterior should be concentrated within the prior bounds supplied in the builder. @@ -718,7 +738,7 @@ mod tests { let evidence_sum: f64 = nested .posterior() .iter() - .map(|sample| sample.evidence_weight()) + .map(super::results::NestedSample::evidence_weight) .sum(); assert!(evidence_sum.is_finite() && evidence_sum > 0.0); } @@ -922,7 +942,7 @@ mod tests { state.tell(results).unwrap(); let mut iterations = 0; - let max_allowed = 100000; // Safety limit to prevent infinite loop + let max_allowed = 100_000; // Safety limit to prevent infinite loop loop { match state.ask() { @@ -934,9 +954,10 @@ mod tests { state.tell(results).unwrap(); iterations += 1; - if iterations > max_allowed { - panic!("Exceeded safety limit - sampler not terminating"); - } + assert!( + iterations <= max_allowed, + "Exceeded safety limit - sampler not terminating" + ); } AskResult::Done(SamplingResults::Nested(samples)) => { // Should terminate eventually @@ -953,15 +974,15 @@ mod tests { | SamplerTermination::InsufficientLivePoints => { // Expected termination reasons } - other => { + other @ SamplerTermination::EvaluationFailed(_) => { // Other reasons are also ok, just document what we see - eprintln!("Terminated with reason: {:?}", other); + eprintln!("Terminated with reason: {other:?}"); } } } break; } - _ => panic!("Unexpected result"), + AskResult::Done(_) => panic!("Unexpected result"), } } } diff --git a/rust/src/sampler/dynamic_nested/results.rs b/rust/src/sampler/dynamic_nested/results.rs index 87f6620..f56f567 100644 --- a/rust/src/sampler/dynamic_nested/results.rs +++ b/rust/src/sampler/dynamic_nested/results.rs @@ -239,7 +239,7 @@ fn information_from_samples(log_z: f64, samples: &[NestedSample]) -> f64 { } // Shift by the maximum log-weight to maintain numerical stability. - let max_delta = deltas.iter().cloned().fold(f64::NEG_INFINITY, f64::max); + let max_delta = deltas.iter().copied().fold(f64::NEG_INFINITY, f64::max); if !max_delta.is_finite() { return 0.0; diff --git a/rust/src/sampler/dynamic_nested/scheduler.rs b/rust/src/sampler/dynamic_nested/scheduler.rs index 18ffe3c..a38d160 100644 --- a/rust/src/sampler/dynamic_nested/scheduler.rs +++ b/rust/src/sampler/dynamic_nested/scheduler.rs @@ -20,7 +20,7 @@ impl Scheduler { let baseline = baseline_live.max(MIN_LIVE_POINTS); let expansion = expansion_factor.max(0.0); let tol = termination_tol.abs().max(1e-10); - let max_live = ((baseline as f64) * (1.0 + 4.0 * expansion)).ceil() as usize; + let max_live = ((baseline as f64) * 4.0f64.mul_add(expansion, 1.0)).ceil() as usize; Self { baseline_live: baseline, @@ -35,7 +35,7 @@ impl Scheduler { /// Compute the desired live-set size given the estimated information gain. pub fn target(&mut self, information: f64, current_live: usize) -> usize { let info = information.max(0.0); - let scale = 1.0 + self.expansion_factor * info.sqrt(); + let scale = self.expansion_factor.mul_add(info.sqrt(), 1.0); let mut desired = (self.baseline_live as f64 * scale).round() as usize; desired = desired.clamp(MIN_LIVE_POINTS, self.max_live); diff --git a/rust/src/sampler/dynamic_nested/state.rs b/rust/src/sampler/dynamic_nested/state.rs index c129b9d..14af038 100644 --- a/rust/src/sampler/dynamic_nested/state.rs +++ b/rust/src/sampler/dynamic_nested/state.rs @@ -69,10 +69,7 @@ impl RemovedPoint { impl SamplerState { /// Initialise state from an initial live set, inferring problem dimension. pub fn new(live_points: Vec) -> Self { - let dimension = live_points - .first() - .map(|point| point.position.len()) - .unwrap_or(0); + let dimension = live_points.first().map_or(0, |point| point.position.len()); Self { live_points, posterior: Vec::new(), @@ -114,15 +111,6 @@ impl SamplerState { .fold(f64::NEG_INFINITY, f64::max) } - /// Worst log-likelihood among the current live points. - #[allow(dead_code)] - pub fn min_log_likelihood(&self) -> f64 { - self.live_points - .iter() - .map(|p| p.log_likelihood) - .fold(f64::INFINITY, f64::min) - } - /// Index of the live point with the lowest likelihood. pub fn worst_index(&self) -> Option { self.live_points diff --git a/rust/src/sampler/metropolis_hastings.rs b/rust/src/sampler/metropolis_hastings.rs index eda9acb..36a4850 100644 --- a/rust/src/sampler/metropolis_hastings.rs +++ b/rust/src/sampler/metropolis_hastings.rs @@ -7,6 +7,7 @@ use rand::{Rng, SeedableRng}; use rand_distr::StandardNormal; use std::time::Instant; +#[must_use] #[derive(Clone, Debug)] pub struct MetropolisHastings { num_chains: usize, @@ -60,6 +61,10 @@ impl MetropolisHastings { /// /// Internally uses the ask/tell interface for consistency. For external /// control of the evaluation loop, use `init()`, `ask()`, and `tell()` directly. + /// + /// # Panics + /// + /// Panics if `tell()` fails during the sampling loop, which indicates a bug in the sampler. pub fn run(&self, mut objective: F, initial: Point, bounds: Bounds) -> Samples where F: FnMut(&[f64]) -> R, @@ -80,7 +85,9 @@ impl MetropolisHastings { AskResult::Done(SamplingResults::MCMC(samples)) => { return samples; } - _ => unreachable!("MetropolisHastings always returns MCMC results"), + AskResult::Done(_) => { + unreachable!("MetropolisHastings always returns MCMC results") + } } } } @@ -92,6 +99,11 @@ impl MetropolisHastings { /// /// Internally uses the ask/tell interface for consistency. For external /// control of the evaluation loop, use `init()`, `ask()`, and `tell()` directly. + /// + /// # Panics + /// + /// Panics if the sampling state machine enters an unexpected state after a tell error. + /// This should not occur under normal operation. pub fn run_batch(&self, objective: F, initial: Point, bounds: Bounds) -> Samples where F: Fn(&[Vec]) -> Vec, @@ -101,7 +113,7 @@ impl MetropolisHastings { let initial_point = initial; let mut state = self.init(initial_point.clone(), bounds); // ToDo: performance improvement, remove clone - let mut results = objective(&vec![initial_point]); + let mut results = objective(&[initial_point]); loop { // Call ask and break if an error is encountered @@ -116,7 +128,9 @@ impl MetropolisHastings { AskResult::Done(SamplingResults::MCMC(samples)) => { return samples; } - _ => unreachable!("MetropolisHastings always returns MCMC results"), + AskResult::Done(_) => { + unreachable!("MetropolisHastings always returns MCMC results") + } } } @@ -185,10 +199,9 @@ impl MetropolisHastings { /// } /// ``` pub fn init(&self, initial: Vec, bounds: Bounds) -> MetropolisHastingsState { - let mut seed_rng = match self.seed { - Some(s) => StdRng::seed_from_u64(s), - None => StdRng::from_os_rng(), - }; + let mut seed_rng = self + .seed + .map_or_else(StdRng::from_os_rng, StdRng::seed_from_u64); // Move initial and clamp it let mut initial_point = initial; @@ -277,6 +290,12 @@ impl MetropolisHastingsState { /// AskResult::Done(results) => { /* done */ } /// } /// ``` + /// + /// # Errors + /// + /// Returns an error if: + /// - The sampler has already terminated + /// - The number of results doesn't match the expected count pub fn tell(&mut self, results: I) -> Result<(), TellError> where I: IntoIterator, @@ -289,10 +308,7 @@ impl MetropolisHastingsState { let values: Vec = results .into_iter() - .map(|r| match r.try_into() { - Ok(eval) => eval.value(), - Err(_) => f64::INFINITY, - }) + .map(|r| r.try_into().map_or(f64::INFINITY, |eval| eval.value())) .collect(); if values.len() != self.chains.len() { diff --git a/rust/src/sampler/mod.rs b/rust/src/sampler/mod.rs index b9e40b2..871be4b 100644 --- a/rust/src/sampler/mod.rs +++ b/rust/src/sampler/mod.rs @@ -136,7 +136,7 @@ impl Sampler { pub fn as_scalar(&self) -> Option<&ScalarSampler> { match self { Sampler::Scalar(s) => Some(s), - _ => None, + Sampler::Gradient(_) => None, } } @@ -144,7 +144,7 @@ impl Sampler { pub fn as_scalar_mut(&mut self) -> Option<&mut ScalarSampler> { match self { Sampler::Scalar(s) => Some(s), - _ => None, + Sampler::Gradient(_) => None, } } @@ -152,7 +152,7 @@ impl Sampler { pub fn as_gradient(&self) -> Option<&GradientSampler> { match self { Sampler::Gradient(g) => Some(g), - _ => None, + Sampler::Scalar(_) => None, } } @@ -160,25 +160,33 @@ impl Sampler { pub fn as_gradient_mut(&mut self) -> Option<&mut GradientSampler> { match self { Sampler::Gradient(g) => Some(g), - _ => None, + Sampler::Scalar(_) => None, } } // ā”€ā”€ā”€ Type extraction - consuming ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ā”€ /// Try to convert into a scalar sampler + /// + /// # Errors + /// + /// Returns `Err(self)` if this sampler is a gradient sampler, not a scalar sampler. pub fn into_scalar(self) -> Result { match self { Sampler::Scalar(s) => Ok(s), - other => Err(other), + other @ Sampler::Gradient(_) => Err(other), } } /// Try to convert into a gradient sampler + /// + /// # Errors + /// + /// Returns `Err(self)` if this sampler is a scalar sampler, not a gradient sampler. pub fn into_gradient(self) -> Result { match self { Sampler::Gradient(g) => Ok(g), - other => Err(other), + other @ Sampler::Scalar(_) => Err(other), } } @@ -367,7 +375,7 @@ impl SamplingResults { pub fn as_mcmc(&self) -> Option<&Samples> { match self { SamplingResults::MCMC(s) => Some(s), - _ => None, + SamplingResults::Nested(_) => None, } } @@ -375,7 +383,7 @@ impl SamplingResults { pub fn as_mcmc_mut(&mut self) -> Option<&mut Samples> { match self { SamplingResults::MCMC(s) => Some(s), - _ => None, + SamplingResults::Nested(_) => None, } } @@ -383,7 +391,7 @@ impl SamplingResults { pub fn as_nested(&self) -> Option<&NestedSamples> { match self { SamplingResults::Nested(s) => Some(s), - _ => None, + SamplingResults::MCMC(_) => None, } } @@ -391,7 +399,7 @@ impl SamplingResults { pub fn as_nested_mut(&mut self) -> Option<&mut NestedSamples> { match self { SamplingResults::Nested(s) => Some(s), - _ => None, + SamplingResults::MCMC(_) => None, } } @@ -569,7 +577,7 @@ mod tests { } break; } - _ => panic!("Expected MCMC results"), + AskResult::Done(_) => panic!("Expected MCMC results"), } } } @@ -596,7 +604,7 @@ mod tests { state.tell(results).unwrap(); } AskResult::Done(SamplingResults::MCMC(samples)) => break samples, - _ => panic!("Expected MCMC results"), + AskResult::Done(_) => panic!("Expected MCMC results"), } }; @@ -624,7 +632,7 @@ mod tests { let err = state.tell(results).unwrap_err(); assert!(matches!(err, TellError::ResultCountMismatch { .. })); } - _ => panic!("Expected Evaluate"), + AskResult::Done(_) => panic!("Expected Evaluate"), } } @@ -644,13 +652,13 @@ mod tests { .tell(points.into_iter().map(|_| Ok::(1.0))) .expect("tell should succeed"); } - _ => panic!("Expected Evaluate on first ask"), + AskResult::Done(_) => panic!("Expected Evaluate on first ask"), } // Now should be Done match state.ask() { AskResult::Done(_) => (), - _ => panic!("Expected Done after iteration 0"), + AskResult::Evaluate(_) => panic!("Expected Done after iteration 0"), } // Trying to tell again should fail @@ -693,7 +701,7 @@ mod tests { ); break; } - _ => panic!("Expected MCMC results"), + AskResult::Done(_) => panic!("Expected MCMC results"), } } } @@ -732,7 +740,7 @@ mod tests { assert_eq!(iteration, 5); break; } - _ => panic!("Expected MCMC results"), + AskResult::Done(_) => panic!("Expected MCMC results"), } } } @@ -774,7 +782,7 @@ mod tests { assert!(eval_count > 0, "Should have evaluated some points"); break; } - _ => panic!("Expected Nested results"), + AskResult::Done(_) => panic!("Expected Nested results"), } } } @@ -804,7 +812,7 @@ mod tests { state.tell(results).unwrap(); } AskResult::Done(SamplingResults::Nested(samples)) => break samples, - _ => panic!("Expected Nested results"), + AskResult::Done(_) => panic!("Expected Nested results"), } }; @@ -866,14 +874,13 @@ mod tests { .collect(); state.tell(results).unwrap(); return; // Test complete - } else { - // Still in initialization or expansion phase, provide correct results - let results: Vec<_> = points - .iter() - .map(|x| Ok::(0.5 * x[0].powi(2))) - .collect(); - state.tell(results).unwrap(); } + // Still in initialisation or expansion phase, provide correct results + let results: Vec<_> = points + .iter() + .map(|x| Ok::(0.5 * x[0].powi(2))) + .collect(); + state.tell(results).unwrap(); } AskResult::Done(_) => { panic!("Should not complete before testing error handling"); @@ -923,7 +930,7 @@ mod tests { assert!(samples.draws() > 0); break; } - _ => panic!("Expected Nested results"), + AskResult::Done(_) => panic!("Expected Nested results"), } } } @@ -958,7 +965,7 @@ mod tests { DNSPhase::AwaitingReplacementBatch { .. } => seen_single_replacement = true, DNSPhase::AwaitingExpansion { .. } => seen_expansion = true, DNSPhase::Terminated(_) => break, - _ => {} + DNSPhase::InitialisingLivePoints { .. } => {} } match state.ask() { diff --git a/rust/src/types.rs b/rust/src/types.rs index f2a0884..7a474c9 100644 --- a/rust/src/types.rs +++ b/rust/src/types.rs @@ -19,7 +19,7 @@ impl ScalarEvaluation { } } -/// TryFrom f64 with validation +/// `TryFrom` f64 with validation impl TryFrom for ScalarEvaluation { type Error = EvaluationError; @@ -60,7 +60,7 @@ impl GradientEvaluation { } } -/// TryFrom (f64, Vec) with validation +/// `TryFrom` (f64, Vec) with validation impl TryFrom<(f64, Vec)> for GradientEvaluation { type Error = EvaluationError; diff --git a/rust/tests/diffsol_optimisation.rs b/rust/tests/diffsol_optimisation.rs index 6003559..d418743 100644 --- a/rust/tests/diffsol_optimisation.rs +++ b/rust/tests/diffsol_optimisation.rs @@ -3,16 +3,16 @@ use nalgebra::DMatrix; #[test] fn diffsol_builder_supports_end_to_end_optimisation() { - let dsl = r#" + let dsl = r" in_i { a = 1 } u_i { y = 0.1 } F_i { a * y } -"#; +"; let true_param = 1.2_f64; let y0 = 0.1_f64; - let t_span: Vec = (0..20).map(|i| i as f64 * 0.05).collect(); + let t_span: Vec = (0..20).map(|i| f64::from(i) * 0.05).collect(); let data_values: Vec = t_span.iter().map(|t| y0 * (true_param * t).exp()).collect(); let mut data = Vec::with_capacity(t_span.len() * 2); diff --git a/rust/tests/dynamic_nested.rs b/rust/tests/dynamic_nested.rs index e584042..2f352f0 100644 --- a/rust/tests/dynamic_nested.rs +++ b/rust/tests/dynamic_nested.rs @@ -28,8 +28,7 @@ fn dynamic_nested_sampler_integration() { // Relax tolerance - nested sampling with limited live points has variability assert!( (posterior_mean - 0.5).abs() < 0.2, - "Posterior mean {} should be close to 0.5", - posterior_mean + "Posterior mean {posterior_mean} should be close to 0.5" ); // Verify posterior samples are valid @@ -47,13 +46,13 @@ fn build_logistic_objective( backend: DiffsolBackend, parallel: bool, ) -> (impl Fn(&[f64]) -> f64, Bounds) { - let dsl = r#" + let dsl = r" in_i { r = 1, k = 1 } u_i { y = 0.1 } F_i { (r * y) * (1 - (y / k)) } -"#; +"; - let t_span: Vec = (0..20).map(|i| i as f64 * 0.1).collect(); + let t_span: Vec = (0..20).map(|i| f64::from(i) * 0.1).collect(); let data_values: Vec = t_span.iter().map(|t| 0.1 * (*t).exp()).collect(); // Data matrix: column-major format with t_span first, then data_values @@ -142,11 +141,7 @@ fn dynamic_nested_sampler_parallel_vs_sequential_consistency() { let mean_diff = (p - s).abs(); assert!( mean_diff < 0.5, - "Parameter {} means should be similar: parallel={}, sequential={}, diff={}", - i, - p, - s, - mean_diff + "Parameter {i} means should be similar: parallel={p}, sequential={s}, diff={mean_diff}" ); } } diff --git a/rust/tests/dynamic_nested_advanced.rs b/rust/tests/dynamic_nested_advanced.rs index d130439..a4d0a4a 100644 --- a/rust/tests/dynamic_nested_advanced.rs +++ b/rust/tests/dynamic_nested_advanced.rs @@ -16,8 +16,8 @@ fn gaussian_evidence_accuracy() { let problem = ScalarProblemBuilder::new() .with_function(move |x: &[f64]| { - let log_norm = sigma.ln() + 0.5 * (2.0 * std::f64::consts::PI).ln(); - 0.5 * (x[0] / sigma).powi(2) + log_norm + let log_norm = 0.5f64.mul_add((2.0 * std::f64::consts::PI).ln(), sigma.ln()); + 0.5f64.mul_add((x[0] / sigma).powi(2), log_norm) }) .with_parameter("x", 0.0, (prior_lower, prior_upper)) .build() @@ -67,7 +67,10 @@ fn bimodal_distribution_samples_both_modes() { let max_log = log_p1.max(log_p2); let log_sum = max_log + ((log_p1 - max_log).exp() + (log_p2 - max_log).exp()).ln(); // Return negative log likelihood - -(log_sum - 2.0_f64.ln() - sigma.ln() - 0.5 * (2.0 * std::f64::consts::PI).ln()) + -0.5f64.mul_add( + -(2.0 * std::f64::consts::PI).ln(), + log_sum - 2.0_f64.ln() - sigma.ln(), + ) }) .with_parameter("x", 0.0, (-10.0, 10.0)) .build() @@ -97,12 +100,8 @@ fn bimodal_distribution_samples_both_modes() { } } - assert!(near_mode1 > 0, "should sample from first mode at x={}", mu1); - assert!( - near_mode2 > 0, - "should sample from second mode at x={}", - mu2 - ); + assert!(near_mode1 > 0, "should sample from first mode at x={mu1}"); + assert!(near_mode2 > 0, "should sample from second mode at x={mu2}"); } /// Test that infinite likelihoods are handled gracefully. @@ -185,8 +184,8 @@ fn high_dimensional_problem() { ); for (i, &m) in nested.mean().iter().enumerate() { - assert!(m.is_finite(), "mean[{}] should be finite", i); - assert!(m.abs() < 0.5, "mean[{}] = {} should be near origin", i, m); + assert!(m.is_finite(), "mean[{i}] should be finite"); + assert!(m.abs() < 0.5, "mean[{i}] = {m} should be near origin"); } } @@ -231,7 +230,7 @@ fn very_small_evidence() { let offset: f64 = 100.0; let problem = ScalarProblemBuilder::new() - .with_function(move |x: &[f64]| offset + 0.5 * x[0].powi(2)) + .with_function(move |x: &[f64]| 0.5f64.mul_add(x[0].powi(2), offset)) .with_parameter("x", 0.0, (-5.0, 5.0)) .build() .expect("failed to build problem"); @@ -289,8 +288,7 @@ fn posterior_weights_normalize() { assert!( (weight_sum - 1.0).abs() < 0.05, - "weights should sum to ~1.0, got {}", - weight_sum + "weights should sum to ~1.0, got {weight_sum}" ); } @@ -330,9 +328,7 @@ fn mean_matches_weighted_average() { assert!( (manual_mean - reported_mean).abs() < 1e-6, - "manual mean {} should match reported mean {}", - manual_mean, - reported_mean + "manual mean {manual_mean} should match reported mean {reported_mean}" ); } @@ -400,10 +396,7 @@ fn respects_parameter_bounds() { let x = sample.position[0]; assert!( x >= lower && x <= upper, - "sample {} should be within bounds [{}, {}]", - x, - lower, - upper + "sample {x} should be within bounds [{lower}, {upper}]" ); } } @@ -445,8 +438,8 @@ fn information_increases_with_constraint() { let make_problem = |sigma: f64| { ScalarProblemBuilder::new() .with_function(move |x: &[f64]| { - let log_norm = sigma.ln() + 0.5 * (2.0 * std::f64::consts::PI).ln(); - 0.5 * (x[0] / sigma).powi(2) + log_norm + let log_norm = 0.5f64.mul_add((2.0 * std::f64::consts::PI).ln(), sigma.ln()); + 0.5f64.mul_add((x[0] / sigma).powi(2), log_norm) }) .with_parameter("x", 0.0, (-10.0, 10.0)) .build() @@ -490,8 +483,8 @@ fn terminates_with_high_information() { let problem = ScalarProblemBuilder::new() .with_function(move |x: &[f64]| { - let log_norm = sigma.ln() + 0.5 * (2.0 * std::f64::consts::PI).ln(); - 0.5 * (x[0] / sigma).powi(2) + log_norm + let log_norm = 0.5f64.mul_add((2.0 * std::f64::consts::PI).ln(), sigma.ln()); + 0.5f64.mul_add((x[0] / sigma).powi(2), log_norm) }) .with_parameter("x", 0.0, ParameterRange::from((-5.0, 5.0))) .build() diff --git a/tests/integration/test_mathematical_suite.py b/tests/integration/test_mathematical_suite.py index 515d924..54c0d43 100644 --- a/tests/integration/test_mathematical_suite.py +++ b/tests/integration/test_mathematical_suite.py @@ -88,7 +88,7 @@ def make_cmaes() -> diffid.CMAES: ], ) @pytest.mark.parametrize( - "objective, dimension, initial, expected, position_tol, fun_tol", + ("objective", "dimension", "initial", "expected", "position_tol", "fun_tol"), [ pytest.param( sphere, diff --git a/tests/test_docs.py b/tests/test_docs.py index 3ac02b0..e981144 100644 --- a/tests/test_docs.py +++ b/tests/test_docs.py @@ -30,6 +30,7 @@ def test_mkdocs_build_succeeds(self): """ result = subprocess.run( ["mkdocs", "build", "--strict"], + check=False, cwd=ROOT, capture_output=True, text=True, @@ -50,7 +51,7 @@ def test_all_notebooks_are_valid_json(self): assert len(notebooks) > 0, "No notebooks found" for notebook_path in notebooks: - with open(notebook_path, encoding="utf-8") as f: + with notebook_path.open(encoding="utf-8") as f: nb_data = json.load(f) # Validate basic notebook structure diff --git a/tests/unit/test_diffsol.py b/tests/unit/test_diffsol.py index 404162d..095ec32 100644 --- a/tests/unit/test_diffsol.py +++ b/tests/unit/test_diffsol.py @@ -188,7 +188,7 @@ def build_problem(cost_metric=None): ) assert pytest.approx(expected_gaussian, rel=1e-6, abs=1e-9) == gaussian_cost - with pytest.raises(ValueError): + with pytest.raises(ValueError, match="variance must be positive"): diffid.GaussianNLL(0.0) diff --git a/tests/unit/test_python_autodiff.py b/tests/unit/test_python_autodiff.py index fb05e15..ff832f2 100644 --- a/tests/unit/test_python_autodiff.py +++ b/tests/unit/test_python_autodiff.py @@ -49,7 +49,7 @@ def test_python_builder_gradient_with_jax(): @pytest.mark.parametrize( - "test_point,description", + ("test_point", "description"), [ (np.array([0.0, 0.0, 0.0], dtype=np.float32), "origin"), (np.array([1.0, 1.0, 1.0], dtype=np.float32), "all ones"), diff --git a/uv.lock b/uv.lock index de4eaaf..5ea6896 100644 --- a/uv.lock +++ b/uv.lock @@ -484,6 +484,7 @@ examples = [ { name = "jax", marker = "python_full_version < '3.14' and platform_machine == 'x86_64' and sys_platform == 'darwin'" }, { name = "jaxlib", marker = "python_full_version < '3.14' and platform_machine == 'x86_64' and sys_platform == 'darwin'" }, { name = "matplotlib" }, + { name = "scipy" }, ] [package.metadata] @@ -520,6 +521,7 @@ examples = [ { name = "jax", marker = "python_full_version < '3.14' and platform_machine == 'x86_64' and sys_platform == 'darwin'", specifier = "==0.4.38" }, { name = "jaxlib", marker = "python_full_version < '3.14' and platform_machine == 'x86_64' and sys_platform == 'darwin'", specifier = "==0.4.38" }, { name = "matplotlib", specifier = ">=3.8" }, + { name = "scipy", specifier = ">=1.16.2" }, ] [[package]]