Feature/bump version (#36)

* remove branch from badges

* fix badges with main branch

* bump version

* increase levenberg and nalgebra dependencies to current ones

* update changelog

* fix clippy lints

* update readme

* fix docs some more

* appease clippy, improve docs

* more clippy

* bump crate version

* support link

CHANGES.md

@@ -3,6 +3,11 @@
 This is the changelog for the `varpro` library.
 See also here for a [version history](https://crates.io/crates/varpro/versions).
 
+## 0.10.0
+
+- Updated dependencies to current versions.
+- Updates to documentation
+
 ## 0.9.0
 
 - We can now calculate confidence bands using via `FitStatistics`. Only available

Cargo.toml

@@ -1,6 +1,6 @@
 [package]
 name = "varpro"
-version = "0.9.0"
+version = "0.10.0"
 authors = ["geo-ant"]
 edition = "2021"
 license = "MIT"
@@ -16,8 +16,8 @@ members = ["shared_test_code"]
 
 [dependencies]
 thiserror = "1"
-levenberg-marquardt = "0.13"
-nalgebra = { version = "0.32" } #, features = ["rayon"]}
+levenberg-marquardt = "0.14"
+nalgebra = { version = "0.33" } #, features = ["rayon"]}
 num-traits = "0.2"
 distrs = "0.2"
 # rayon = "1.6"

README.md

@@ -11,6 +11,11 @@ Nonlinear function fitting made simple. This library provides robust and fast
 least-squares fitting of a wide class of model functions to data.
 It uses the VarPro algorithm to achieve this, hence the name.
 
+## Support
+
+If you like this crate, recommend it to others. You can also [buy me a coffee](https://ko-fi.com/geoant)
+if you feel like it.
+
 ## Introduction
 
 This crate implements a powerful algorithm to fit model functions to data, 
@@ -28,8 +33,8 @@ which is a notoriously ill-conditioned problem.
 
 ### What is VarPro?
 
-Variable Projection (VarPro) is an algorithm that exploits that its fitting
-problem can be separated into linear and nonlinear parameters.
+Variable Projection (VarPro) is an algorithm that tatkes advantage of
+the fact that its fitting problem can be separated into linear and nonlinear parameters.
 First, the linear parameters are eliminated using some clever linear algebra. Then,
 the fitting problem is rewritten so that it depends only on the nonlinear parameters.
 Finally, this reduced problem is solved by using a general purpose nonlinear minimization algorithm,
@@ -145,7 +150,7 @@ than just using a nonlinear solver without the magic of varpro.
 But this crate offers an additional way to eek out the last bits of performance.
 
 The `SeparableNonlinearModel` trait can be manually implemented to describe a
-model function. This often allows us to shave of the last hundreds of microseconds
+model function. This often allows us to shave off the last hundreds of microseconds
 from the computation, e.g. by caching intermediate calculations. The crate documentation
 contains detailed examples.
 

benches/double_exponential_without_noise.rs

@@ -1,23 +1,17 @@
 use criterion::{criterion_group, criterion_main, Criterion};
 use levenberg_marquardt::LeastSquaresProblem;
 use levenberg_marquardt::LevenbergMarquardt;
-use nalgebra::ComplexField;
-
 use nalgebra::DVector;
 use nalgebra::DefaultAllocator;
-
 use nalgebra::Dyn;
 use nalgebra::OVector;
 use nalgebra::RawStorageMut;
-
 use nalgebra::Storage;
 use nalgebra::U1;
-
 use pprof::criterion::{Output, PProfProfiler};
 use shared_test_code::models::DoubleExpModelWithConstantOffsetSepModel;
 use shared_test_code::models::DoubleExponentialDecayFittingWithOffsetLevmar;
 use shared_test_code::*;
-
 use varpro::prelude::SeparableNonlinearModel;
 use varpro::solvers::levmar::LevMarProblem;
 use varpro::solvers::levmar::LevMarProblemBuilder;
@@ -39,12 +33,10 @@ fn build_problem<Model>(
 ) -> LevMarProblem<Model, false>
 where
     Model: SeparableNonlinearModel<ScalarType = f64>,
-    DefaultAllocator: nalgebra::allocator::Allocator<f64, Dyn>,
-    DefaultAllocator: nalgebra::allocator::Allocator<f64, Dyn, Dyn>,
-    <DefaultAllocator as nalgebra::allocator::Allocator<f64, Dyn>>::Buffer: Storage<f64, Dyn>,
-    <DefaultAllocator as nalgebra::allocator::Allocator<f64, Dyn>>::Buffer: RawStorageMut<f64, Dyn>,
-    DefaultAllocator:
-        nalgebra::allocator::Allocator<(<f64 as ComplexField>::RealField, usize), Dyn>,
+    DefaultAllocator: nalgebra::allocator::Allocator<Dyn>,
+    DefaultAllocator: nalgebra::allocator::Allocator<Dyn, Dyn>,
+    <DefaultAllocator as nalgebra::allocator::Allocator<Dyn>>::Buffer<f64>: Storage<f64, Dyn>,
+    <DefaultAllocator as nalgebra::allocator::Allocator<Dyn>>::Buffer<f64>: RawStorageMut<f64, Dyn>,
 {
     let DoubleExponentialParameters {
         tau1,

shared_test_code/Cargo.toml

@@ -7,8 +7,8 @@ edition = "2021"
 
 [dependencies]
 varpro = {path = "../"}
-nalgebra = "0.32"
+nalgebra = "0.33"
 num-traits = "0.2"
-levenberg-marquardt = "0.13"
+levenberg-marquardt = "0.14"
 approx = "0.5"
 rand = "0.8"

shared_test_code/src/lib.rs

@@ -54,8 +54,8 @@ pub fn evaluate_complete_model_at_params<Model>(
 where
     Model::ScalarType: Scalar + ComplexField,
     Model: SeparableNonlinearModel,
-    DefaultAllocator: nalgebra::allocator::Allocator<Model::ScalarType, Dyn, Dyn>,
-    DefaultAllocator: nalgebra::allocator::Allocator<Model::ScalarType, Dyn>,
+    DefaultAllocator: nalgebra::allocator::Allocator<Dyn, Dyn>,
+    DefaultAllocator: nalgebra::allocator::Allocator<Dyn>,
 {
     let original_params = model.params();
     model
@@ -80,8 +80,8 @@ pub fn evaluate_complete_model_at_params_mrhs<Model>(
 where
     Model::ScalarType: Scalar + ComplexField,
     Model: SeparableNonlinearModel,
-    DefaultAllocator: nalgebra::allocator::Allocator<Model::ScalarType, Dyn, Dyn>,
-    DefaultAllocator: nalgebra::allocator::Allocator<Model::ScalarType, Dyn>,
+    DefaultAllocator: nalgebra::allocator::Allocator<Dyn, Dyn>,
+    DefaultAllocator: nalgebra::allocator::Allocator<Dyn>,
 {
     let original_params = model.params();
     model

shared_test_code/src/models.rs

@@ -2,11 +2,6 @@ use nalgebra::{DVector, Dyn, OMatrix, OVector, U1};
 use varpro::model::SeparableModel;
 use varpro::{model::errors::ModelError, prelude::*};
 
-/// contains a double exponential model with constant offset that does not
-/// perform caching, as preliminary results indicate that using caching can
-/// lead to performance penalties for large models
-pub mod uncached;
-
 /// A separable model for double exponential decay
 /// with a constant offset
 /// f_j = c1*exp(-x_j/tau1) + c2*exp(-x_j/tau2) + c3

src/basis_function/mod.rs

@@ -9,6 +9,7 @@ use nalgebra::{DVector, Scalar};
 ///
 /// * `$\vec{x}$` is an independent variable, like time, spatial coordinates and so on,
 /// * `$\alpha_j$` are scalar parameter of this function
+///
 /// The functions must have at least one parameter argument additional ot the location parameter, i.e. a function
 /// `$f(\vec{x})$` does not satisfy the trait, but e.g. `$f(\vec{x},\alpha_1)$` does. The function
 /// must be vector valued and should produce a vector of the same size as `\vec{x}`. (If you want to
@@ -18,15 +19,20 @@ use nalgebra::{DVector, Scalar};
 /// If it is violated, then calculations using the basis function will return errors in the [SeparableModel](crate::model::SeparableModel).
 ///
 /// # Variadic Functions
+///
 /// Since Rust does not have variadic functions or generics, this trait is implemented for all functions up to
 /// a maximum number of arguments. Right now this is implemented for up to 10 arguments.
 ///
 /// # Generic Parameters
+///
 /// ## Scalar Type
+///
 /// The functions must have an interface `Fn(&DVector<ScalarType>,ScalarType)-> DVector<ScalarType>`,
 /// `Fn(&DVector<ScalarType>,ScalarType,ScalarType)-> DVector<ScalarType>` and so on.
 /// All numeric types, including the return value must be of the same scalar type.
+///
 /// ## ArgList : The argument list
+///
 /// This type is of no consequence for the user because it will be correctly inferred when
 /// passing a function. It is a [nifty trick](https://geo-ant.github.io/blog/2021/rust-traits-and-variadic-functions/)
 /// that allows us to implement this trait for functions taking different arguments. Just FYI: The
@@ -42,9 +48,10 @@ where
     /// # Arguments
     /// * `x`: the vector `$\vec{x}$`
     /// * `params`: The parameters `$(\alpha_1,...,\alpha_N)$` as a slice. The slice must have
-    /// the correct number of arguments for calling the function (no more, no less).
+    ///    the correct number of arguments for calling the function (no more, no less).
     ///
     /// # Effect and Panics
+    ///
     /// If the slice has fewer elements than the parameter argument list. If the slice has more element,
     /// only the first `$N$` elements are used and dispatched to the parameters in order, i.e.
     /// `$\alpha_1$`=`param[0]`, ..., `$\alpha_N$`=`param[N-1]`. Calling eval will result in evaluating

src/lib.rs

@@ -47,10 +47,11 @@
 //! Lets look at the components of this equation in more detail. The vector valued function
 //! `$\vec{f}(\vec{\alpha},\vec{c}) \in \mathbb{C}^{N_{data}}$` is the model we want to fit. It depends on
 //! two vector valued parameters:
+//!
 //! * `$\vec{\alpha}=(\alpha_1,\dots,\alpha_{N_{params}})^T$` is the vector of nonlinear model parameters.
-//! We will get back to these later.
+//!    We will get back to these later.
 //! * `$\vec{c}=(c_1,\dots,c_{N_{basis}})^T$` is the vector of coefficients for the basis functions.
-//! Those are the linear model parameters.
+//!    Those are the linear model parameters.
 //!
 //! Note that we call `$\vec{\alpha}$` the _nonlinear parameters_ and `$\vec{c}$` the _coefficients_ of the model
 //! just to make the distinction between the two types of parameters clear. The coefficients are
@@ -200,7 +201,7 @@
 //! ```
 //!
 //! We'll see in the example how the [function](crate::model::builder::SeparableModelBuilder::function) method
-//! and the [partial_deriv](crate::model::builder::SeparableModelBuilderProxyWithDerivatives::partial_deriv)
+//! and the [partial_deriv](crate::model::builder::SeparableModelBuilder::partial_deriv)
 //! methods let us add the function and the derivative as base functions.
 //!
 //! There is a second type of basis function, which corresponds to coefficient `$c_3$`. This is a constant
@@ -349,11 +350,11 @@
 //! statement is the following:
 //!
 //!  * We have not only one observation but a set `$\{\vec{y}_s\}$`, `$s=1,...,S$` of
-//! observations.
+//!    observations.
 //!  * We want to fit the separable nonlinear function `$\vec{f}(\vec{\alpha},\vec{c})$`
-//! to all vectors of observations, but in such a way that the linear parameters
-//! are allowed to vary with `$s$`, but the nonlinear parameters
-//! are the same for the whole dataset.
+//!    to all vectors of observations, but in such a way that the linear parameters
+//!    are allowed to vary with `$s$`, but the nonlinear parameters
+//!    are the same for the whole dataset.
 //!
 //! This is called global fitting with multiple right hand sides,
 //! because the nonlinear parameters are not

src/model/builder/mod.rs

@@ -54,9 +54,9 @@ pub mod error;
 /// **Function Arguments and Output**
 ///
 /// * The first argument of the function must be a reference to a `&DVector` type
-/// that accepts the independent variable (the `$\vec{x}$` values) and the other
-/// parameters must be scalars that are the nonlinear parameters that the basis
-/// function depends on.
+///   that accepts the independent variable (the `$\vec{x}$` values) and the other
+///   parameters must be scalars that are the nonlinear parameters that the basis
+///   function depends on.
 ///
 /// So if we want to model a basis function `$\vec{f_1}(\vec{x},\vec{\alpha})$`
 /// where `$\vec{\alpha}=(\alpha_1,\alpha_2)$` we would write the function in Rust as
@@ -129,14 +129,14 @@ pub mod error;
 /// ** Rules You Must Abide By **
 ///
 /// * Basis functions must be **nonlinear** in the parameters they take. If they aren't, you can always
-/// rewrite the problem so that the linear parameters go in the coefficient vector `$\vec{c}$`. This
-/// means that each partial derivative also depend on all the parameters that the basis function depends
-/// on.
+///   rewrite the problem so that the linear parameters go in the coefficient vector `$\vec{c}$`. This
+///   means that each partial derivative also depend on all the parameters that the basis function depends
+///   on.
 ///
 /// * Derivatives must take the same parameter arguments *and in the same order* as the original
-/// basis function. This means if basis function `$\vec{f}_j$` is given as `$\vec{f}_j(\vec{x},a,b)$`,
-/// then the derivatives must also be given with the parameters `$a,b$` in the same order, i.e.
-/// `$\partial/\partial a \vec{f}_j(\vec{x},a,b)$`, `$\partial/\partial b \vec{f}_j(\vec{x},a,b)$`.
+///   basis function. This means if basis function `$\vec{f}_j$` is given as `$\vec{f}_j(\vec{x},a,b)$`,
+///   then the derivatives must also be given with the parameters `$a,b$` in the same order, i.e.
+///   `$\partial/\partial a \vec{f}_j(\vec{x},a,b)$`, `$\partial/\partial b \vec{f}_j(\vec{x},a,b)$`.
 ///
 /// **Rules Enforced at Compile Time**
 ///
@@ -333,8 +333,8 @@ where
     /// * The list of parameters must only contain unique names
     /// * The list of parameter names must not be empty
     /// * Parameter names must not contain a comma. This is a precaution because
-    /// `&["alpha,beta"]` most likely indicates a typo for `&["alpha","beta"]`. Any other form
-    /// of punctuation is allowed.
+    ///   `&["alpha,beta"]` most likely indicates a typo for `&["alpha","beta"]`. Any other form
+    ///   of punctuation is allowed.
     pub fn new<StrCollection>(parameter_names: StrCollection) -> Self
     where
         StrCollection: IntoIterator,

src/model/builder/modelfunction_builder/mod.rs

@@ -16,7 +16,7 @@ use crate::model::model_basis_function::ModelBasisFunction;
 /// * the model parameters are unique and non-empty
 /// * the function parameters are unique, non-empty and a subset of the model params
 /// * a derivative is provided for each parameter that the function depends on
-/// (all other derivatives are zero, because the function does not depends on other params)
+///   (all other derivatives are zero, because the function does not depend on other params)
 ///
 #[doc(hidden)]
 pub struct ModelBasisFunctionBuilder<ScalarType>
@@ -38,13 +38,17 @@ where
     /// begin constructing a modelfunction for a specific model. The modelfunction must take
     /// a subset of the model parameters. This is the first step in creating a function, because
     /// the modelfunction must have partial derivatives specified for each parameter it takes.
+    ///
     /// # Arguments
+    ///
     /// * `model_parameters`: the model parameters of the model to which this function belongs. This is important
-    /// so the builder understands how the parameters of the function relate to the parameters of the model.
+    ///    so the builder understands how the parameters of the function relate to the parameters of the model.
     /// * `function_parameters`: the parameters that the function takes. Those must be in the order
-    /// of the parameter vector. The paramters must not be empty, nor may they contain duplicates
+    ///    of the parameter vector. The paramters must not be empty, nor may they contain duplicates
     /// * `function`: the actual function.
+    ///
     /// # Result
+    ///
     /// A model builder that can be used to add derivatives.
     pub fn new<FuncType, StrCollection, ArgList>(
         model_parameters: Vec<String>,
@@ -88,9 +92,9 @@ where
     /// Add a derivative for the function with respect to the given parameter.
     /// # Arguments
     /// * `parameter`: the parameter with respect to which the derivative is taken.
-    /// The parameter must be inside the set of model parameters. Furthermore the
+    ///    The parameter must be inside the set of model parameters. Furthermore the
     /// * `derivative`: the partial derivative of the function with which the
-    /// builder was created.
+    ///    builder was created.
     pub fn partial_deriv<FuncType, ArgList>(mut self, parameter: &str, derivative: FuncType) -> Self
     where
         FuncType: BasisFunction<ScalarType, ArgList> + 'static,

src/model/detail.rs

@@ -9,7 +9,9 @@ use crate::model::builder::error::ModelBuildError;
 /// * the set of parameters is not empty
 /// * the set of parameters contains only unique elements
 /// * any of the parameter names contains a comma. This indicates most likely a typo when giving the parameter list
+///
 /// # Returns
+///
 /// Ok if the conditions hold, otherwise an error variant.
 pub fn check_parameter_names<StrType>(param_names: &[StrType]) -> Result<(), ModelBuildError>
 where
@@ -80,11 +82,13 @@ where
 /// of model parameters and dispatch them to our model function
 /// # Arguments
 /// * `model_parameters`: the parameters that the complete model that this function belongs to
-/// depends on
+///    depends on
 /// * `function_parameters`: the parameters (in right to left order) that the basisfunction depends
-/// on. Must be a subset of the model parameters.
+///    on. Must be a subset of the model parameters.
 /// * `function` a basis function that depends on a number of parameters
+///
 /// # Result
+///
 /// Say our model depends on parameters `$\vec{p}=(\alpha,\beta,\gamma)$` and the `function` argument
 /// is a basis function `$f(\vec{x},\gamma,\beta)$`. Then calling `create_wrapped_basis_function(&["alpha","beta","gamma"],&["gamma","alpha"],f)`
 /// creates a wrapped function `$\tilde{f}(\vec{x},\vec{p})$` which can be called with `$\tilde{f}(\vec{x},\vec{p})=f(\vec{x},\gamma,\beta$`.

src/model/mod.rs

@@ -237,8 +237,8 @@ pub mod test;
 #[allow(clippy::type_complexity)]
 pub trait SeparableNonlinearModel
 where
-    DefaultAllocator: nalgebra::allocator::Allocator<Self::ScalarType, Dyn>,
-    DefaultAllocator: nalgebra::allocator::Allocator<Self::ScalarType, Dyn, Dyn>,
+    DefaultAllocator: nalgebra::allocator::Allocator<Dyn>,
+    DefaultAllocator: nalgebra::allocator::Allocator<Dyn, Dyn>,
 {
     /// the scalar number type for this model, which should be
     /// a real or complex number type, commonly either `f64` or `f32`.
@@ -312,11 +312,11 @@ where
     /// # Arguments
     ///
     /// * `derivative_index`: The index of the nonlinear parameter with respect to which
-    /// partial derivative should be evaluated. We use _zero based indexing_
-    /// here! Put in more simple terms, say your model has three nonlinear parameters
-    /// `a,b,c`, so your vector of nonlinear parameters is `$\vec{\alpha} = (a,b,c)$`.
-    /// Then index 0 requests `$\partial/\partial_a$`, index 1 requests `$\partial/\partial_b$`
-    /// and index 2 requests `$\partial/\partial_c$`.
+    ///    partial derivative should be evaluated. We use _zero based indexing_
+    ///    here! Put in more simple terms, say your model has three nonlinear parameters
+    ///    `a,b,c`, so your vector of nonlinear parameters is `$\vec{\alpha} = (a,b,c)$`.
+    ///    Then index 0 requests `$\partial/\partial_a$`, index 1 requests `$\partial/\partial_b$`
+    ///    and index 2 requests `$\partial/\partial_c$`.
     ///
     /// # Result
     ///