use ffi::{self, FFI};
use std::os::raw::{c_int, c_void};
use VectorF64;
ffi_wrapper!(MultiFitFSolverType, *mut sys::gsl_multifit_fsolver_type);
pub struct MultiFitFunction(pub sys::gsl_multifit_function);
ffi_wrapper!(
MultiFitFSolver,
*mut sys::gsl_multifit_fsolver,
gsl_multifit_fsolver_free
);
impl MultiFitFSolver {
#[doc(alias = "gsl_multifit_fsolver_alloc")]
pub fn new(t: &MultiFitFSolverType, n: usize, p: usize) -> Option<MultiFitFSolver> {
let tmp = unsafe { sys::gsl_multifit_fsolver_alloc(t.unwrap_shared(), n, p) };
if tmp.is_null() {
None
} else {
Some(MultiFitFSolver::wrap(tmp))
}
}
#[doc(alias = "gsl_multifit_fsolver_set")]
pub fn set(&mut self, f: &mut MultiFitFunction, x: &mut VectorF64) -> ::Value {
::Value::from(unsafe {
sys::gsl_multifit_fsolver_set(self.unwrap_unique(), &mut f.0, x.unwrap_shared())
})
}
#[doc(alias = "gsl_multifit_fsolver_iterate")]
pub fn iterate(&mut self) -> ::Value {
::Value::from(unsafe { sys::gsl_multifit_fsolver_iterate(self.unwrap_unique()) })
}
#[doc(alias = "gsl_multifit_fsolver_name")]
pub fn name(&self) -> String {
unsafe {
let tmp = sys::gsl_multifit_fsolver_name(self.unwrap_shared());
String::from_utf8_lossy(::std::ffi::CStr::from_ptr(tmp).to_bytes()).to_string()
}
}
#[doc(alias = "gsl_multifit_fsolver_position")]
pub fn position(&self) -> VectorF64 {
unsafe { ffi::FFI::wrap(sys::gsl_multifit_fsolver_position(self.unwrap_shared())) }
}
}
ffi_wrapper!(
MultiFitFdfSolver,
*mut sys::gsl_multifit_fdfsolver,
gsl_multifit_fdfsolver_free
);
impl MultiFitFdfSolver {
#[doc(alias = "gsl_multifit_fdfsolver_alloc")]
pub fn new(_type: &MultiFitFdfSolverType, n: usize, p: usize) -> Option<MultiFitFdfSolver> {
let s = unsafe { sys::gsl_multifit_fdfsolver_alloc(_type.unwrap_shared(), n, p) };
if s.is_null() {
None
} else {
Some(MultiFitFdfSolver::wrap(s))
}
}
#[doc(alias = "gsl_multifit_fdfsolver_set")]
pub fn set(&mut self, f: &mut MultiFitFunctionFdf, x: &::VectorF64) -> ::Value {
::Value::from(unsafe {
sys::gsl_multifit_fdfsolver_set(self.unwrap_unique(), f.to_raw(), x.unwrap_shared())
})
}
pub fn x(&self) -> ::VectorF64 {
unsafe { ffi::FFI::soft_wrap((*self.unwrap_shared()).x) }
}
pub fn f(&self) -> ::VectorF64 {
unsafe { ffi::FFI::soft_wrap((*self.unwrap_shared()).f) }
}
pub fn dx(&self) -> ::VectorF64 {
unsafe { ffi::FFI::soft_wrap((*self.unwrap_shared()).dx) }
}
pub fn g(&self) -> ::VectorF64 {
unsafe { ffi::FFI::soft_wrap((*self.unwrap_shared()).g) }
}
pub fn sqrt_wts(&self) -> ::VectorF64 {
unsafe { ffi::FFI::soft_wrap((*self.unwrap_shared()).sqrt_wts) }
}
#[doc(alias = "gsl_multifit_fdfsolver_name")]
pub fn name(&self) -> String {
unsafe {
let tmp = sys::gsl_multifit_fdfsolver_name(self.unwrap_shared());
String::from_utf8_lossy(::std::ffi::CStr::from_ptr(tmp).to_bytes()).to_string()
}
}
#[doc(alias = "gsl_multifit_fdfsolver_iterate")]
pub fn iterate(&mut self) -> ::Value {
::Value::from(unsafe { sys::gsl_multifit_fdfsolver_iterate(self.unwrap_unique()) })
}
#[doc(alias = "gsl_multifit_fdfsolver_position")]
pub fn position(&self) -> ::VectorF64 {
unsafe { ffi::FFI::wrap(sys::gsl_multifit_fdfsolver_position(self.unwrap_shared())) }
}
#[allow(unused_assignments)]
#[doc(alias = "gsl_multifit_test_delta")]
pub fn driver(&mut self, max_iter: usize, epsabs: f64, epsrel: f64) -> ::Value {
let mut status = ::Value::Failure;
let ptr = self.unwrap_shared();
if !ptr.is_null() {
let mut iter = 0usize;
loop {
status = self.iterate();
if status != ::Value::Success {
break;
}
status = ::Value::from(unsafe {
sys::gsl_multifit_test_delta((*ptr).dx, (*ptr).x, epsabs, epsrel)
});
iter += 1;
if status != ::Value::Continue || iter >= max_iter {
break;
}
}
}
status
}
}
ffi_wrapper!(
MultiFitFdfSolverType,
*const sys::gsl_multifit_fdfsolver_type
);
impl MultiFitFdfSolverType {
#[doc(alias = "gsl_multifit_fdfsolver_lmder")]
pub fn lmder() -> MultiFitFdfSolverType {
ffi_wrap!(gsl_multifit_fdfsolver_lmder)
}
#[doc(alias = "gsl_multifit_fdfsolver_lmsder")]
pub fn lmsder() -> MultiFitFdfSolverType {
ffi_wrap!(gsl_multifit_fdfsolver_lmsder)
}
}
pub struct MultiFitFunctionFdf {
pub f: Option<Box<dyn Fn(::VectorF64, ::VectorF64) -> ::Value>>,
pub df: Option<Box<dyn Fn(::VectorF64, ::MatrixF64) -> ::Value>>,
pub fdf: Option<Box<dyn Fn(::VectorF64, ::VectorF64, ::MatrixF64) -> ::Value>>,
pub n: usize,
pub p: usize,
intern: sys::gsl_multifit_function_fdf,
}
impl MultiFitFunctionFdf {
#[doc(alias = "gsl_multifit_function_fdf")]
pub fn new(n: usize, p: usize, nevalf: usize, nevaldf: usize) -> MultiFitFunctionFdf {
MultiFitFunctionFdf {
f: None,
df: None,
fdf: None,
n,
p,
intern: sys::gsl_multifit_function_fdf {
f: Some(f),
df: Some(df),
fdf: Some(fdf),
n,
p,
params: ::std::ptr::null_mut(),
nevalf,
nevaldf,
},
}
}
#[allow(clippy::wrong_self_convention)]
fn to_raw(&mut self) -> *mut sys::gsl_multifit_function_fdf {
self.intern.n = self.n;
self.intern.p = self.p;
self.intern.params = self as *mut MultiFitFunctionFdf as *mut c_void;
&mut self.intern
}
}
unsafe extern "C" fn f(
x: *const sys::gsl_vector,
params: *mut c_void,
pf: *mut sys::gsl_vector,
) -> c_int {
let t = params as *mut MultiFitFunctionFdf;
if let Some(ref i_f) = (*t).f {
i_f(
ffi::FFI::soft_wrap(x as usize as *mut _),
ffi::FFI::soft_wrap(pf),
)
.into()
} else {
::Value::Success.into()
}
}
unsafe extern "C" fn df(
x: *const sys::gsl_vector,
params: *mut c_void,
pdf: *mut sys::gsl_matrix,
) -> c_int {
let t = params as *mut MultiFitFunctionFdf;
if let Some(ref i_df) = (*t).df {
i_df(
ffi::FFI::soft_wrap(x as usize as *mut _),
ffi::FFI::soft_wrap(pdf),
)
.into()
} else {
::Value::Success.into()
}
}
unsafe extern "C" fn fdf(
x: *const sys::gsl_vector,
params: *mut c_void,
pf: *mut sys::gsl_vector,
pdf: *mut sys::gsl_matrix,
) -> c_int {
let t = params as *mut MultiFitFunctionFdf;
if let Some(ref i_fdf) = (*t).fdf {
i_fdf(
ffi::FFI::soft_wrap(x as usize as *mut _),
ffi::FFI::soft_wrap(pf),
ffi::FFI::soft_wrap(pdf),
)
.into()
} else {
::Value::Success.into()
}
}