1 #ifndef STAN_MATH_PRIM_FUN_LB_CONSTRAIN_HPP
2 #define STAN_MATH_PRIM_FUN_LB_CONSTRAIN_HPP
3
4 #include <stan/math/prim/meta.hpp>
5 #include <stan/math/prim/err.hpp>
6 #include <stan/math/prim/fun/constants.hpp>
7 #include <stan/math/prim/fun/identity_constrain.hpp>
8 #include <stan/math/prim/fun/identity_free.hpp>
9 #include <stan/math/prim/fun/add.hpp>
10 #include <stan/math/prim/fun/exp.hpp>
11 #include <stan/math/prim/fun/eval.hpp>
12 #include <stan/math/prim/fun/sum.hpp>
13 #include <stan/math/prim/fun/value_of.hpp>
14 #include <cmath>
15
16 namespace stan {
17 namespace math {
18
19 /**
20 * Return the lower-bounded value for the specified unconstrained input
21 * and specified lower bound.
22 *
23 * <p>The transform applied is
24 *
25 * <p>\f$f(x) = \exp(x) + L\f$
26 *
27 * <p>where \f$L\f$ is the constant lower bound.
28 *
29 * @tparam T Scalar.
30 * @tparam L Scalar.
31 * @param[in] x Unconstrained input
32 * @param[in] lb Lower bound
33 * @return Constrained matrix
34 */
35 template <typename T, typename L, require_all_stan_scalar_t<T, L>* = nullptr,
36 require_all_not_st_var<T, L>* = nullptr>
lb_constrain(const T & x,const L & lb)37 inline auto lb_constrain(const T& x, const L& lb) {
38 if (unlikely(value_of_rec(lb) == NEGATIVE_INFTY)) {
39 return identity_constrain(x, lb);
40 } else {
41 return add(exp(x), lb);
42 }
43 }
44
45 /**
46 * Return the lower-bounded value for the specified unconstrained
47 * input and specified lower bound, incrementing the specified
48 * reference with the log absolute Jacobian determinant of the
49 * transform.
50 *
51 * @tparam T Scalar.
52 * @tparam L Scalar.
53 * @param[in] x unconstrained input
54 * @param[in] lb lower bound on output
55 * @param[in,out] lp reference to log probability to increment
56 * @return lower-bound constrained value corresponding to inputs
57 */
58 template <typename T, typename L, require_all_stan_scalar_t<T, L>* = nullptr,
59 require_all_not_st_var<T, L>* = nullptr>
lb_constrain(const T & x,const L & lb,return_type_t<T,L> & lp)60 inline auto lb_constrain(const T& x, const L& lb, return_type_t<T, L>& lp) {
61 if (value_of_rec(lb) == NEGATIVE_INFTY) {
62 return identity_constrain(x, lb);
63 } else {
64 lp += x;
65 return add(exp(x), lb);
66 }
67 }
68
69 /**
70 * Specialization of `lb_constrain` to apply a scalar lower bound elementwise
71 * to each input.
72 *
73 * @tparam T A type inheriting from `EigenBase`.
74 * @tparam L Scalar.
75 * @param[in] x unconstrained input
76 * @param[in] lb lower bound on output
77 * @return lower-bound constrained value corresponding to inputs
78 */
79 template <typename T, typename L, require_eigen_t<T>* = nullptr,
80 require_stan_scalar_t<L>* = nullptr,
81 require_all_not_st_var<T, L>* = nullptr>
lb_constrain(T && x,L && lb)82 inline auto lb_constrain(T&& x, L&& lb) {
83 return eval(x.unaryExpr([lb](auto&& x) { return lb_constrain(x, lb); }));
84 }
85
86 /**
87 * Specialization of `lb_constrain` to apply a scalar lower bound elementwise
88 * to each input.
89 *
90 * @tparam T A type inheriting from `EigenBase`.
91 * @tparam L Scalar.
92 * @param[in] x unconstrained input
93 * @param[in] lb lower bound on output
94 * @param[in,out] lp reference to log probability to increment
95 * @return lower-bound constrained value corresponding to inputs
96 */
97 template <typename T, typename L, require_eigen_t<T>* = nullptr,
98 require_stan_scalar_t<L>* = nullptr,
99 require_all_not_st_var<T, L>* = nullptr>
lb_constrain(const T & x,const L & lb,return_type_t<T,L> & lp)100 inline auto lb_constrain(const T& x, const L& lb, return_type_t<T, L>& lp) {
101 return eval(
102 x.unaryExpr([lb, &lp](auto&& xx) { return lb_constrain(xx, lb, lp); }));
103 }
104
105 /**
106 * Specialization of `lb_constrain` to apply a matrix of lower bounds
107 * elementwise to each input element.
108 *
109 * @tparam T A type inheriting from `EigenBase`.
110 * @tparam L A type inheriting from `EigenBase`.
111 * @param[in] x unconstrained input
112 * @param[in] lb lower bound on output
113 * @return lower-bound constrained value corresponding to inputs
114 */
115 template <typename T, typename L, require_all_eigen_t<T, L>* = nullptr,
116 require_all_not_st_var<T, L>* = nullptr>
lb_constrain(T && x,L && lb)117 inline auto lb_constrain(T&& x, L&& lb) {
118 check_matching_dims("lb_constrain", "x", x, "lb", lb);
119 return eval(x.binaryExpr(
120 lb, [](auto&& x, auto&& lb) { return lb_constrain(x, lb); }));
121 }
122
123 /**
124 * Specialization of `lb_constrain` to apply a matrix of lower bounds
125 * elementwise to each input element.
126 *
127 * @tparam T A type inheriting from `EigenBase`.
128 * @tparam L A type inheriting from `EigenBase`.
129 * @param[in] x unconstrained input
130 * @param[in] lb lower bound on output
131 * @param[in,out] lp reference to log probability to increment
132 * @return lower-bound constrained value corresponding to inputs
133 */
134 template <typename T, typename L, require_all_eigen_t<T, L>* = nullptr,
135 require_all_not_st_var<T, L>* = nullptr>
lb_constrain(const T & x,const L & lb,return_type_t<T,L> & lp)136 inline auto lb_constrain(const T& x, const L& lb, return_type_t<T, L>& lp) {
137 check_matching_dims("lb_constrain", "x", x, "lb", lb);
138 return eval(x.binaryExpr(
139 lb, [&lp](auto&& xx, auto&& lbb) { return lb_constrain(xx, lbb, lp); }));
140 }
141
142 /**
143 * Specialization of `lb_constrain` to apply a container of lower bounds
144 * elementwise to each input element.
145 *
146 * @tparam T A Any type with a Scalar `scalar_type`.
147 * @tparam L A type inheriting from `EigenBase` or a scalar.
148 * @param[in] x unconstrained input
149 * @param[in] lb lower bound on output
150 * @return lower-bound constrained value corresponding to inputs
151 */
152 template <typename T, typename L, require_not_std_vector_t<L>* = nullptr>
lb_constrain(const std::vector<T> & x,const L & lb)153 inline auto lb_constrain(const std::vector<T>& x, const L& lb) {
154 std::vector<plain_type_t<decltype(lb_constrain(x[0], lb))>> ret(x.size());
155 for (size_t i = 0; i < x.size(); ++i) {
156 ret[i] = lb_constrain(x[i], lb);
157 }
158 return ret;
159 }
160
161 /**
162 * Specialization of `lb_constrain` to apply a container of lower bounds
163 * elementwise to each input element.
164 *
165 * @tparam T A Any type with a Scalar `scalar_type`.
166 * @tparam L A type inheriting from `EigenBase` or a standard vector.
167 * @param[in] x unconstrained input
168 * @param[in] lb lower bound on output
169 * @param[in,out] lp reference to log probability to increment
170 * @return lower-bound constrained value corresponding to inputs
171 */
172 template <typename T, typename L, require_not_std_vector_t<L>* = nullptr>
lb_constrain(const std::vector<T> & x,const L & lb,return_type_t<T,L> & lp)173 inline auto lb_constrain(const std::vector<T>& x, const L& lb,
174 return_type_t<T, L>& lp) {
175 std::vector<plain_type_t<decltype(lb_constrain(x[0], lb))>> ret(x.size());
176 for (size_t i = 0; i < x.size(); ++i) {
177 ret[i] = lb_constrain(x[i], lb, lp);
178 }
179 return ret;
180 }
181
182 /**
183 * Specialization of `lb_constrain` to apply a container of lower bounds
184 * elementwise to each input element.
185 *
186 * @tparam T A Any type with a Scalar `scalar_type`.
187 * @tparam L A type inheriting from `EigenBase` or a standard vector.
188 * @param[in] x unconstrained input
189 * @param[in] lb lower bound on output
190 * @return lower-bound constrained value corresponding to inputs
191 */
192 template <typename T, typename L>
lb_constrain(const std::vector<T> & x,const std::vector<L> & lb)193 inline auto lb_constrain(const std::vector<T>& x, const std::vector<L>& lb) {
194 check_matching_dims("lb_constrain", "x", x, "lb", lb);
195 std::vector<plain_type_t<decltype(lb_constrain(x[0], lb[0]))>> ret(x.size());
196 for (size_t i = 0; i < x.size(); ++i) {
197 ret[i] = lb_constrain(x[i], lb[i]);
198 }
199 return ret;
200 }
201
202 /**
203 * Specialization of `lb_constrain` to apply a container of lower bounds
204 * elementwise to each input element.
205 *
206 * @tparam T A Any type with a Scalar `scalar_type`.
207 * @tparam L A type inheriting from `EigenBase` or a standard vector.
208 * @param[in] x unconstrained input
209 * @param[in] lb lower bound on output
210 * @param[in,out] lp reference to log probability to increment
211 * @return lower-bound constrained value corresponding to inputs
212 */
213 template <typename T, typename L>
lb_constrain(const std::vector<T> & x,const std::vector<L> & lb,return_type_t<T,L> & lp)214 inline auto lb_constrain(const std::vector<T>& x, const std::vector<L>& lb,
215 return_type_t<T, L>& lp) {
216 check_matching_dims("lb_constrain", "x", x, "lb", lb);
217 std::vector<plain_type_t<decltype(lb_constrain(x[0], lb[0]))>> ret(x.size());
218 for (size_t i = 0; i < x.size(); ++i) {
219 ret[i] = lb_constrain(x[i], lb[i], lp);
220 }
221 return ret;
222 }
223
224 /**
225 * Specialization of `lb_constrain` to apply a container of lower bounds
226 * elementwise to each input element. If the `Jacobian` parameter is `true`, the
227 * log density accumulator is incremented with the log absolute Jacobian
228 * determinant of the transform. All of the transforms are specified with their
229 * Jacobians in the *Stan Reference Manual* chapter Constraint Transforms.
230 *
231 * @tparam Jacobian if `true`, increment log density accumulator with log
232 * absolute Jacobian determinant of constraining transform
233 * @tparam T A type inheriting from `Eigen::EigenBase`, a `var_value` with inner
234 * type inheriting from `Eigen::EigenBase`, a standard vector, or a scalar
235 * @tparam L A type inheriting from `Eigen::EigenBase`, a `var_value` with inner
236 * type inheriting from `Eigen::EigenBase`, a standard vector, or a scalar
237 * @param[in] x unconstrained input
238 * @param[in] lb lower bound on output
239 * @param[in, out] lp log density accumulator
240 * @return lower-bound constrained value corresponding to inputs
241 */
242 template <bool Jacobian, typename T, typename L>
lb_constrain(const T & x,const L & lb,return_type_t<T,L> & lp)243 inline auto lb_constrain(const T& x, const L& lb, return_type_t<T, L>& lp) {
244 if (Jacobian) {
245 return lb_constrain(x, lb, lp);
246 } else {
247 return lb_constrain(x, lb);
248 }
249 }
250
251 } // namespace math
252 } // namespace stan
253
254 #endif
255