Skip to content

Commit 52ee9d6

Browse files
committed
Use Eigen::Ref for input
1 parent 632a734 commit 52ee9d6

File tree

2 files changed

+22
-22
lines changed

2 files changed

+22
-22
lines changed

src/finitediff.cpp

Lines changed: 11 additions & 11 deletions
Original file line numberDiff line numberDiff line change
@@ -65,7 +65,7 @@ double get_denominator(const AccuracyOrder accuracy)
6565

6666
// Compute the gradient of a function at a point using finite differences.
6767
void finite_gradient(
68-
const Eigen::VectorXd& x,
68+
const Eigen::Ref<const Eigen::VectorXd>& x,
6969
const std::function<double(const Eigen::VectorXd&)>& f,
7070
Eigen::VectorXd& grad,
7171
const AccuracyOrder accuracy,
@@ -93,7 +93,7 @@ void finite_gradient(
9393
}
9494

9595
void finite_jacobian(
96-
const Eigen::VectorXd& x,
96+
const Eigen::Ref<const Eigen::VectorXd>& x,
9797
const std::function<Eigen::VectorXd(const Eigen::VectorXd&)>& f,
9898
Eigen::MatrixXd& jac,
9999
const AccuracyOrder accuracy,
@@ -121,7 +121,7 @@ void finite_jacobian(
121121
}
122122

123123
void finite_hessian(
124-
const Eigen::VectorXd& x,
124+
const Eigen::Ref<const Eigen::VectorXd>& x,
125125
const std::function<double(const Eigen::VectorXd&)>& f,
126126
Eigen::MatrixXd& hess,
127127
const AccuracyOrder accuracy,
@@ -159,8 +159,8 @@ void finite_hessian(
159159

160160
// Compare if two gradients are close enough.
161161
bool compare_gradient(
162-
const Eigen::VectorXd& x,
163-
const Eigen::VectorXd& y,
162+
const Eigen::Ref<const Eigen::VectorXd>& x,
163+
const Eigen::Ref<const Eigen::VectorXd>& y,
164164
const double test_eps,
165165
const std::string& msg)
166166
{
@@ -185,8 +185,8 @@ bool compare_gradient(
185185

186186
// Compare if two jacobians are close enough.
187187
bool compare_jacobian(
188-
const Eigen::MatrixXd& x,
189-
const Eigen::MatrixXd& y,
188+
const Eigen::Ref<const Eigen::MatrixXd>& x,
189+
const Eigen::Ref<const Eigen::MatrixXd>& y,
190190
const double test_eps,
191191
const std::string& msg)
192192
{
@@ -216,16 +216,16 @@ bool compare_jacobian(
216216

217217
// Compare if two hessians are close enough.
218218
bool compare_hessian(
219-
const Eigen::MatrixXd& x,
220-
const Eigen::MatrixXd& y,
219+
const Eigen::Ref<const Eigen::MatrixXd>& x,
220+
const Eigen::Ref<const Eigen::MatrixXd>& y,
221221
const double test_eps,
222222
const std::string& msg)
223223
{
224224
return compare_jacobian(x, y, test_eps, msg);
225225
}
226226

227227
// Flatten the matrix rowwise
228-
Eigen::VectorXd flatten(const Eigen::MatrixXd& X)
228+
Eigen::VectorXd flatten(const Eigen::Ref<const Eigen::MatrixXd>& X)
229229
{
230230
Eigen::VectorXd x(X.size());
231231
for (int i = 0; i < X.rows(); i++) {
@@ -237,7 +237,7 @@ Eigen::VectorXd flatten(const Eigen::MatrixXd& X)
237237
}
238238

239239
// Unflatten rowwise
240-
Eigen::MatrixXd unflatten(const Eigen::VectorXd& x, int dim)
240+
Eigen::MatrixXd unflatten(const Eigen::Ref<const Eigen::VectorXd>& x, int dim)
241241
{
242242
assert(x.size() % dim == 0);
243243
Eigen::MatrixXd X(x.size() / dim, dim);

src/finitediff.hpp

Lines changed: 11 additions & 11 deletions
Original file line numberDiff line numberDiff line change
@@ -32,7 +32,7 @@ enum AccuracyOrder {
3232
* @param[in] eps Value of the finite difference step.
3333
*/
3434
void finite_gradient(
35-
const Eigen::VectorXd& x,
35+
const Eigen::Ref<const Eigen::VectorXd>& x,
3636
const std::function<double(const Eigen::VectorXd&)>& f,
3737
Eigen::VectorXd& grad,
3838
const AccuracyOrder accuracy = SECOND,
@@ -48,7 +48,7 @@ void finite_gradient(
4848
* @param[in] eps Value of the finite difference step.
4949
*/
5050
void finite_jacobian(
51-
const Eigen::VectorXd& x,
51+
const Eigen::Ref<const Eigen::VectorXd>& x,
5252
const std::function<Eigen::VectorXd(const Eigen::VectorXd&)>& f,
5353
Eigen::MatrixXd& jac,
5454
const AccuracyOrder accuracy = SECOND,
@@ -64,7 +64,7 @@ void finite_jacobian(
6464
* @param[in] eps Value of the finite difference step.
6565
*/
6666
void finite_hessian(
67-
const Eigen::VectorXd& x,
67+
const Eigen::Ref<const Eigen::VectorXd>& x,
6868
const std::function<double(const Eigen::VectorXd&)>& f,
6969
Eigen::MatrixXd& hess,
7070
const AccuracyOrder accuracy = SECOND,
@@ -81,8 +81,8 @@ void finite_hessian(
8181
* @return A boolean for if x and y are close to the same value.
8282
*/
8383
bool compare_gradient(
84-
const Eigen::VectorXd& x,
85-
const Eigen::VectorXd& y,
84+
const Eigen::Ref<const Eigen::VectorXd>& x,
85+
const Eigen::Ref<const Eigen::VectorXd>& y,
8686
const double test_eps = 1e-4,
8787
const std::string& msg = "compare_gradient ");
8888

@@ -97,8 +97,8 @@ bool compare_gradient(
9797
* @return A boolean for if x and y are close to the same value.
9898
*/
9999
bool compare_jacobian(
100-
const Eigen::MatrixXd& x,
101-
const Eigen::MatrixXd& y,
100+
const Eigen::Ref<const Eigen::MatrixXd>& x,
101+
const Eigen::Ref<const Eigen::MatrixXd>& y,
102102
const double test_eps = 1e-4,
103103
const std::string& msg = "compare_jacobian ");
104104

@@ -113,15 +113,15 @@ bool compare_jacobian(
113113
* @return A boolean for if x and y are close to the same value.
114114
*/
115115
bool compare_hessian(
116-
const Eigen::MatrixXd& x,
117-
const Eigen::MatrixXd& y,
116+
const Eigen::Ref<const Eigen::MatrixXd>& x,
117+
const Eigen::Ref<const Eigen::MatrixXd>& y,
118118
const double test_eps = 1e-4,
119119
const std::string& msg = "compare_hessian ");
120120

121121
/// @brief Flatten the matrix rowwise
122-
Eigen::VectorXd flatten(const Eigen::MatrixXd& X);
122+
Eigen::VectorXd flatten(const Eigen::Ref<const Eigen::MatrixXd>& X);
123123

124124
/// @brief Unflatten rowwise
125-
Eigen::MatrixXd unflatten(const Eigen::VectorXd& x, int dim);
125+
Eigen::MatrixXd unflatten(const Eigen::Ref<const Eigen::VectorXd>& x, int dim);
126126

127127
} // namespace fd

0 commit comments

Comments
 (0)