spdcalc 2.0.1

SPDCalc, the Spontaneous Parametric Downconversion Calculator
Documentation 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132

extern crate ndarray;
// extern crate optimize;
// use nelder_mead::{*, params::*, bounds::*};
// use argmin::prelude::*;
// use argmin::solver::neldermead::NelderMead;
// use serde::{Deserialize, Serialize};

use argmin::core::{CostFunction, Executor};
// use ndarray::{Array, ArrayView1};
use argmin::solver::neldermead::NelderMead;
// use optimize::*;

// #[derive(Clone, Default, Serialize, Deserialize)]
// struct NelderMead1d<F: FnMut(f64) -> f64> {
//   func: F,
// }
//
// impl<F : FnMut(f64) -> f64> ArgminOp for NelderMead1d<F>
// where
//   F: serde::Serialize + std::clone::Clone + std::marker::Sync +
// std::marker::Send {
//   /// Type of the parameter vector
//   type Param = Vec<f64>;
//   /// Type of the return value computed by the cost function
//   type Output = f64;
//   /// Type of the Hessian. Can be `()` if not needed.
//   type Hessian = ();
//   /// Type of the Jacobian. Can be `()` if not needed.
//   type Jacobian = ();
//
//   /// Apply the cost function to a parameter `p`
//   fn apply(&self, p: &Self::Param) -> Result<Self::Output, Error> {
//     Ok((self.func)(p[0]))
//   }
// }

struct Cost1d<F>
where
  F: Fn(f64) -> f64,
{
  func: F,
  min: f64,
  max: f64,
}
impl<F> CostFunction for Cost1d<F>
where
  F: Fn(f64) -> f64,
{
  type Param = f64;
  type Output = f64;
  fn cost(&self, x: &Self::Param) -> Result<Self::Output, argmin::core::Error> {
    if x > &self.max || x < &self.min {
      Ok(f64::INFINITY)
    } else {
      let f = &self.func;
      Ok(f(*x))
    }
  }
}

/// nelder mead optimization. Returns x
pub fn nelder_mead_1d<F>(
  func: F,
  guess: (f64, f64),
  max_iter: u64,
  min: f64,
  max: f64,
  tolerance: f64,
) -> f64
where
  F: Fn(f64) -> f64,
{
  let cost = Cost1d { func, min, max };
  let solver = NelderMead::new(vec![guess.0, guess.1])
    .with_sd_tolerance(tolerance)
    .unwrap();

  let res = Executor::new(cost, solver)
    .configure(|state| state.max_iters(max_iter))
    .run()
    .unwrap();

  use argmin::core::State;
  *res.state().get_best_param().unwrap()

  // let minimizer = NelderMeadBuilder::default()
  //   .xtol(tolerance)
  //   .ftol(tolerance)
  //   .maxiter(max_iter as usize)
  //   .build()
  //   .unwrap();

  // let cost = |args : ArrayView1<f64>| {
  //   let x = args[0];
  //   if x > max || x < min {
  //     std::f64::INFINITY
  //   }
  //   // high cost if x outside bounds
  //   else {
  //     func(x)
  //   }
  // };

  // let ans = minimizer.minimize(&cost, Array::from_vec(vec![guess]).view());

  // ans[0]

  // let (x, _fx) = minimize(
  //   |args| func(args[0]),
  //   vec![guess],
  //   (max-min)/100.,
  //   Params::default(),
  //   Bounds { min: vec![min], max: vec![max] },
  //   max_iter
  // );
  //
  // x[0]

  // let cost = NelderMead1d { func };
  //
  // // Set up solver -- note that the proper choice of the vertices is very
  // important! let solver = NelderMead::new()
  //   .with_initial_params(vec![guess])
  //   .sd_tolerance(0.0001);
  //
  //   // Run solver
  // let res = Executor::new(cost, solver, vec![])
  //   .max_iters(100)
  //   .run()?;
  //
  // res
}