Skip to content

Commit

Permalink
update subtree
Browse files Browse the repository at this point in the history
  • Loading branch information
Yi Zhang committed Nov 14, 2024
1 parent 9425768 commit 6fb6f1a
Show file tree
Hide file tree
Showing 19 changed files with 84 additions and 4,166 deletions.
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
#ifndef STAN_MATH_TORSTEN_LINEAR_INTERPOLATION_HPP
#define STAN_MATH_TORSTEN_LINEAR_INTERPOLATION_HPP
#ifndef STAN_MATH_TORSTEN_PMX_LN_INTERPOLATE_HPP
#define STAN_MATH_TORSTEN_PMX_LN_INTERPOLATE_HPP

#include <stan/math/rev.hpp>
#include <stan/math/prim.hpp>
Expand All @@ -12,7 +12,7 @@
namespace torsten {

/**
* Return the values of a piecewise linear function at specifed values of the
* Return the values of a piecewise linear function at specifed values of the
* function argument. The function is specified in terms of a set of x y pairs.
*
* @tparam T0 Type of value to be interpolated
Expand All @@ -26,16 +26,16 @@ namespace torsten {

template <typename T0, typename T1, typename T2>
typename stan::return_type_t<T0, T1, T2>
inline linear_interpolation(const T0& xout,
inline pmx_ln_interpolate(const T0& xout,
const std::vector<T1>& x,
const std::vector<T2>& y) {
stan::math::check_finite("linear_interpolation", "xout", xout);
stan::math::check_finite("linear_interpolation", "x", x);
stan::math::check_finite("linear_interpolation", "y", y);
stan::math::check_nonzero_size("linear_interpolation", "x", x);
stan::math::check_nonzero_size("linear_interpolation", "y", y);
stan::math::check_ordered("linear_interpolation", "x", x);
stan::math::check_matching_sizes("linear_interpolation", "x", x, "y", y);
stan::math::check_finite("pmx_ln_interpolate", "xout", xout);
stan::math::check_finite("pmx_ln_interpolate", "x", x);
stan::math::check_finite("pmx_ln_interpolate", "y", y);
stan::math::check_nonzero_size("pmx_ln_interpolate", "x", x);
stan::math::check_nonzero_size("pmx_ln_interpolate", "y", y);
stan::math::check_ordered("pmx_ln_interpolate", "x", x);
stan::math::check_matching_sizes("pmx_ln_interpolate", "x", x, "y", y);

typename stan::return_type_t<T0, T1, T2> yout;
if (xout < x.front()) {
Expand All @@ -51,15 +51,15 @@ namespace torsten {

template <typename T0, typename T1, typename T2>
std::vector<typename stan::return_type_t<T0, T1, T2>>
inline linear_interpolation(const std::vector<T0>& xout,
inline pmx_ln_interpolate(const std::vector<T0>& xout,
const std::vector<T1>& x,
const std::vector<T2>& y) {
stan::math::check_nonzero_size("linear_interpolation", "xout", xout);
stan::math::check_finite("linear_interpolation", "xout", xout);
stan::math::check_nonzero_size("pmx_ln_interpolate", "xout", xout);
stan::math::check_finite("pmx_ln_interpolate", "xout", xout);

std::vector<typename stan::return_type_t<T0, T1, T2>> yout(xout.size());
std::transform(xout.begin(), xout.end(), yout.begin(),
[&x, &y](const T0& xi) { return linear_interpolation(xi, x, y); });
[&x, &y](const T0& xi) { return pmx_ln_interpolate(xi, x, y); });
return yout;
}
}
Expand Down
Original file line number Diff line number Diff line change
@@ -1,13 +1,13 @@
#include <stan/math/torsten/linear_interpolation.hpp>
#include <stan/math/torsten/pmx_ln_interpolate.hpp>
#include <stan/math/prim.hpp>
#include <gtest/gtest.h>
#include <iostream>
#include <sstream>
#include <vector>
#include <test/unit/util.hpp>

TEST(linear_interpolation, xdbl){
using torsten::linear_interpolation;
TEST(pmx_ln_interpolate, xdbl){
using torsten::pmx_ln_interpolate;
using stan::math::var;
int nx = 5;
std::vector<double> x(nx), y(nx);
Expand All @@ -19,12 +19,12 @@ TEST(linear_interpolation, xdbl){
}

double xout = 3.4;
double yout = linear_interpolation(xout, x, y);
double yout = pmx_ln_interpolate(xout, x, y);
EXPECT_FLOAT_EQ(yout, k * xout);
}

TEST(linear_interpolation, in_range_gradient_x_y){
using torsten::linear_interpolation;
TEST(pmx_ln_interpolate, in_range_gradient_x_y){
using torsten::pmx_ln_interpolate;
using stan::math::var;
using stan::math::value_of;
int nx = 5;
Expand All @@ -38,7 +38,7 @@ TEST(linear_interpolation, in_range_gradient_x_y){

{ // xout is double
double xout = 0.4;
var yout = linear_interpolation(xout, x, y);
var yout = pmx_ln_interpolate(xout, x, y);
int i = 0;
var y0 = y.at(i) + (y.at(i+1) - y.at(i)) / (x.at(i+1) - x.at(i)) * (xout - x.at(i));
EXPECT_FLOAT_EQ(value_of(yout), k * xout);
Expand All @@ -63,7 +63,7 @@ TEST(linear_interpolation, in_range_gradient_x_y){

{ // xout is double
double xout = 2.4;
var yout = linear_interpolation(xout, x, y);
var yout = pmx_ln_interpolate(xout, x, y);
int i = 2;
var y0 = y.at(i) + (y.at(i+1) - y.at(i)) / (x.at(i+1) - x.at(i)) * (xout - x.at(i));
EXPECT_FLOAT_EQ(value_of(yout), k * xout);
Expand All @@ -88,7 +88,7 @@ TEST(linear_interpolation, in_range_gradient_x_y){

{ // xout is var
var xout = 1.4;
var yout = linear_interpolation(xout, x, y);
var yout = pmx_ln_interpolate(xout, x, y);
int i = 1;
var y0 = y.at(i) + (y.at(i+1) - y.at(i)) / (x.at(i+1) - x.at(i)) * (xout - x.at(i));
EXPECT_FLOAT_EQ(value_of(yout), value_of(k * xout));
Expand Down Expand Up @@ -121,8 +121,8 @@ TEST(linear_interpolation, in_range_gradient_x_y){
}
}

TEST(linear_interpolation, out_range_gradient_x_y){
using torsten::linear_interpolation;
TEST(pmx_ln_interpolate, out_range_gradient_x_y){
using torsten::pmx_ln_interpolate;
using stan::math::var;
using stan::math::value_of;
int nx = 5;
Expand All @@ -136,7 +136,7 @@ TEST(linear_interpolation, out_range_gradient_x_y){

{ // xout is double
double xout = -1.2;
var yout = linear_interpolation(xout, x, y);
var yout = pmx_ln_interpolate(xout, x, y);
int i = 0;
var y0 = y.front();
EXPECT_FLOAT_EQ(value_of(yout), value_of(y.front()));
Expand All @@ -161,7 +161,7 @@ TEST(linear_interpolation, out_range_gradient_x_y){

{ // xout is double
double xout = 7.5;
var yout = linear_interpolation(xout, x, y);
var yout = pmx_ln_interpolate(xout, x, y);
int i = 2;
var y0 = y.back();
EXPECT_FLOAT_EQ(value_of(yout), value_of(y.back()));
Expand All @@ -185,7 +185,7 @@ TEST(linear_interpolation, out_range_gradient_x_y){
}
}

TEST(linear_interpolation, linear_example) {
TEST(pmx_ln_interpolate, linear_example) {
int nx = 5, nout = 3;
std::vector<double> x(nx), y(nx), xout(nout), yout;
double youtTrue;
Expand All @@ -199,7 +199,7 @@ TEST(linear_interpolation, linear_example) {
xout[1] = 2.5;
xout[2] = 4.5;

yout = torsten::linear_interpolation(xout, x, y);
yout = torsten::pmx_ln_interpolate(xout, x, y);
for(int i = 0; i < nout; i++){
if(xout[i] <= x[0]){
youtTrue = x[0];
Expand All @@ -212,11 +212,11 @@ TEST(linear_interpolation, linear_example) {
}
}

TEST(linear_interpolation, xgradient){
TEST(pmx_ln_interpolate, xgradient){
using stan::math::var;
int nx = 5, nout = 3;
std::vector<double> xdbl(nx), ydbl(nx), xoutdbl(nout), thisGrad(nx);
Eigen::MatrixXd trueJac = Eigen::MatrixXd::Zero(3, 5);
Eigen::MatrixXd trueJac = Eigen::MatrixXd::Zero(3, 5);

for(int i = 0; i < nx; i++){
xdbl[i] = i;
Expand All @@ -242,7 +242,7 @@ TEST(linear_interpolation, xgradient){
for(int k = 0; k < nout; k++){
xout[k] = xoutdbl[k];
}
yout = torsten::linear_interpolation(xout, x, y);
yout = torsten::pmx_ln_interpolate(xout, x, y);

yout[i].grad(x, thisGrad);

Expand All @@ -252,11 +252,11 @@ TEST(linear_interpolation, xgradient){
}
}

TEST(linear_interpolation, ygradient){
TEST(pmx_ln_interpolate, ygradient){
using stan::math::var;
int nx = 5, nout = 3;
std::vector<double> xdbl(nx), ydbl(nx), xoutdbl(nout), thisGrad(nx);
Eigen::MatrixXd trueJac = Eigen::MatrixXd::Zero(3, 5);
Eigen::MatrixXd trueJac = Eigen::MatrixXd::Zero(3, 5);

for(int i = 0; i < nx; i++){
xdbl[i] = i;
Expand Down Expand Up @@ -284,21 +284,21 @@ TEST(linear_interpolation, ygradient){
for(int k = 0; k < nout; k++){
xout[k] = xoutdbl[k];
}
yout = torsten::linear_interpolation(xout, x, y);
yout = torsten::pmx_ln_interpolate(xout, x, y);

yout[i].grad(y, thisGrad);

for(int j = 0; j < nx; j++){
EXPECT_EQ(trueJac(i, j), thisGrad[j]);
}
}
}

TEST(linear_interpolation, xoutgradient){
TEST(pmx_ln_interpolate, xoutgradient){
using stan::math::var;
int nx = 5, nout = 3;
std::vector<double> xdbl(nx), ydbl(nx), xoutdbl(nout), thisGrad(nx);
Eigen::MatrixXd trueJac = Eigen::MatrixXd::Zero(3, 3);
Eigen::MatrixXd trueJac = Eigen::MatrixXd::Zero(3, 3);

for(int i = 0; i < nx; i++){
xdbl[i] = i;
Expand All @@ -321,7 +321,7 @@ TEST(linear_interpolation, xoutgradient){
for(int k = 0; k < nout; k++){
xout[k] = xoutdbl[k];
}
yout = torsten::linear_interpolation(xout, x, y);
yout = torsten::pmx_ln_interpolate(xout, x, y);

yout[i].grad(xout, thisGrad);

Expand All @@ -331,7 +331,7 @@ TEST(linear_interpolation, xoutgradient){
}
}

TEST(linear_interpolation, error_conditions) {
TEST(pmx_ln_interpolate, error_conditions) {
using stan::math::var;
int nx = 5, nout = 3;
std::vector<double> x(nx), y(nx), xout(nout), yout;
Expand All @@ -346,34 +346,34 @@ TEST(linear_interpolation, error_conditions) {
xout[2] = 4.5;

std::vector<double> xout_bad;
EXPECT_THROW_MSG(torsten::linear_interpolation(xout_bad, x, y),
EXPECT_THROW_MSG(torsten::pmx_ln_interpolate(xout_bad, x, y),
std::invalid_argument,
"xout has size 0");

std::vector<double> x_bad;
EXPECT_THROW_MSG(torsten::linear_interpolation(xout, x_bad, y),
EXPECT_THROW_MSG(torsten::pmx_ln_interpolate(xout, x_bad, y),
std::invalid_argument,
"x has size 0");

std::vector<double> y_bad;
EXPECT_THROW_MSG(torsten::linear_interpolation(xout, x, y_bad),
EXPECT_THROW_MSG(torsten::pmx_ln_interpolate(xout, x, y_bad),
std::invalid_argument,
"y has size 0");

std::vector<double> x3_bad = x;
x3_bad[2] = 0.0;
EXPECT_THROW_MSG(torsten::linear_interpolation(xout, x3_bad, y),
EXPECT_THROW_MSG(torsten::pmx_ln_interpolate(xout, x3_bad, y),
std::domain_error,
"x is not a valid ordered vector");

std::vector<double> x2_bad(nx - 1);
std::vector<double> x2_bad(nx - 1);
for(int i = 0; i < (nx - 1); i++) x2_bad[i] = x[i];
EXPECT_THROW_MSG(torsten::linear_interpolation(xout, x2_bad, y),
EXPECT_THROW_MSG(torsten::pmx_ln_interpolate(xout, x2_bad, y),
std::invalid_argument,
"size of x (4) and size of y (5) must match in size");
}

TEST(linear_interpolation, error_conditions_inf) {
TEST(pmx_ln_interpolate, error_conditions_inf) {
using stan::math::var;
std::stringstream expected_is_inf;
expected_is_inf << "is " << std::numeric_limits<double>::infinity();
Expand All @@ -394,37 +394,37 @@ TEST(linear_interpolation, error_conditions_inf) {

std::vector<double> xout_bad = xout;
xout_bad[0] = inf;
EXPECT_THROW_MSG(torsten::linear_interpolation(xout_bad, x, y),
EXPECT_THROW_MSG(torsten::pmx_ln_interpolate(xout_bad, x, y),
std::domain_error,
"xout");
EXPECT_THROW_MSG(torsten::linear_interpolation(xout_bad, x, y),
EXPECT_THROW_MSG(torsten::pmx_ln_interpolate(xout_bad, x, y),
std::domain_error,
expected_is_inf.str());

xout_bad = xout;
xout_bad[0] = -inf;
EXPECT_THROW_MSG(torsten::linear_interpolation(xout_bad, x, y),
EXPECT_THROW_MSG(torsten::pmx_ln_interpolate(xout_bad, x, y),
std::domain_error,
"xout");
EXPECT_THROW_MSG(torsten::linear_interpolation(xout_bad, x, y),
EXPECT_THROW_MSG(torsten::pmx_ln_interpolate(xout_bad, x, y),
std::domain_error,
expected_is_neg_inf.str());

std::vector<double> x_bad = x;
x_bad[0] = inf;
EXPECT_THROW_MSG(torsten::linear_interpolation(xout, x_bad, y),
EXPECT_THROW_MSG(torsten::pmx_ln_interpolate(xout, x_bad, y),
std::domain_error,
"x");
EXPECT_THROW_MSG(torsten::linear_interpolation(xout, x_bad, y),
"x[1]");
EXPECT_THROW_MSG(torsten::pmx_ln_interpolate(xout, x_bad, y),
std::domain_error,
expected_is_inf.str());

std::vector<double> y_bad = y;
y_bad[0] = -inf;
EXPECT_THROW_MSG(torsten::linear_interpolation(xout, x, y_bad),
EXPECT_THROW_MSG(torsten::pmx_ln_interpolate(xout, x, y_bad),
std::domain_error,
"y");
EXPECT_THROW_MSG(torsten::linear_interpolation(xout, x, y_bad),
EXPECT_THROW_MSG(torsten::pmx_ln_interpolate(xout, x, y_bad),
std::domain_error,
expected_is_neg_inf.str());
}
2 changes: 1 addition & 1 deletion cmdstan/stan/lib/stan_math/stan/math/torsten/torsten.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -16,7 +16,7 @@
#include <stan/math/torsten/pmx_solve_twocpt_rk45.hpp>
#include <stan/math/torsten/pmx_solve_twocpt_bdf.hpp>

#include <stan/math/torsten/linear_interpolation.hpp>
#include <stan/math/torsten/pmx_ln_interpolate.hpp>
#include <stan/math/torsten/univariate_integral.hpp>

#include <stan/math/torsten/dsolve.hpp>
Expand Down

This file was deleted.

Loading

0 comments on commit 6fb6f1a

Please sign in to comment.