Stan  1.0
probability, sampling & optimization
matrix_error_handling.hpp
Go to the documentation of this file.
1 #ifndef __STAN__MATH__MATRIX_ERROR_HANDLING_HPP__
2 #define __STAN__MATH__MATRIX_ERROR_HANDLING_HPP__
3 
4 #include <stan/math/matrix.hpp>
5 #include <stan/meta/traits.hpp>
7 #include <boost/type_traits/common_type.hpp>
8 
9 namespace stan {
10 
11  namespace math {
12 
13  // FIXME: update warnings
14  template <typename T_size1, typename T_size2, typename T_result,
15  class Policy>
16  inline bool check_size_match(const char* function,
17  T_size1 i,
18  T_size2 j,
19  T_result* result,
20  const Policy&) {
22  typedef typename boost::common_type<T_size1,T_size2>::type common_type;
23  if (static_cast<common_type>(i) != static_cast<common_type>(j)) {
24  std::ostringstream msg;
25  msg << "i and j must be same. Found i=%1%, j=" << j;
26  T_result tmp = raise_domain_error<T_result,T_size1>(function,
27  msg.str().c_str(),
28  i,
29  Policy());
30  if (result != 0)
31  *result = tmp;
32  return false;
33  }
34  return true;
35  }
36 
37  template <typename T_size1, typename T_size2, typename T_result>
38  inline
39  bool check_size_match(const char* function,
40  T_size1 i,
41  T_size2 j,
42  T_result* result) {
43  return check_size_match(function,i,j,result,default_policy());
44  }
45 
46  template <typename T_size1, typename T_size2>
47  inline
48  bool check_size_match(const char* function,
49  T_size1 i,
50  T_size2 j,
51  T_size1* result = 0) {
52  return check_size_match(function,i,j,result,default_policy());
53  }
54 
55 
56 
69  template <typename T_y, typename T_result, class Policy>
70  inline bool check_symmetric(const char* function,
71  const Eigen::Matrix<T_y,Eigen::Dynamic,Eigen::Dynamic>& y,
72  const char* name,
73  T_result* result,
74  const Policy&) {
75  typedef
77  size_type;
78  size_type k = y.rows();
79  if (k == 1)
80  return true;
81  for (size_type m = 0; m < k; ++m) {
82  for (size_type n = m + 1; n < k; ++n) {
83  if (fabs(y(m,n) - y(n,m)) > CONSTRAINT_TOLERANCE) {
84  std::ostringstream message;
85  message << name << " is not symmetric. "
86  << name << "[" << m << "," << n << "] is %1%, but "
87  << name << "[" << n << "," << m
88  << "] element is " << y(n,m);
89  T_result tmp
90  = policies::raise_domain_error<T_y>(function,
91  message.str().c_str(),
92  y(m,n), Policy());
93  if (result != 0)
94  *result = tmp;
95  return false;
96  }
97  }
98  }
99  return true;
100  }
101 
102 
103  template <typename T_y, typename T_result>
104  inline bool check_symmetric(const char* function,
105  const Eigen::Matrix<T_y,Eigen::Dynamic,Eigen::Dynamic>& y,
106  const char* name,
107  T_result* result) {
108  return check_symmetric(function,y,name,result,default_policy());
109  }
110 
111  template <typename T>
112  inline bool check_symmetric(const char* function,
113  const Eigen::Matrix<T,Eigen::Dynamic,Eigen::Dynamic>& y,
114  const char* name,
115  T* result = 0) {
116  return check_symmetric(function,y,name,result,default_policy());
117  }
118 
119 
120 
121 
134  // FIXME: update warnings (message has (0,0) item)
135  template <typename T_y, typename T_result, class Policy>
136  inline bool check_pos_definite(const char* function,
137  const Eigen::Matrix<T_y,Eigen::Dynamic,Eigen::Dynamic>& y,
138  const char* name,
139  T_result* result,
140  const Policy&) {
141  typedef
143  size_type;
144  if (y.rows() == 1 && y(0,0) <= CONSTRAINT_TOLERANCE) {
145  std::ostringstream message;
146  message << name << " is not positive definite. "
147  << name << "(0,0) is %1%.";
148  T_result tmp = policies::raise_domain_error<T_y>(function,
149  message.str().c_str(),
150  y(0,0), Policy());
151  if (result != 0)
152  *result = tmp;
153  return false;
154  }
155  Eigen::LDLT< Eigen::Matrix<T_y,Eigen::Dynamic,Eigen::Dynamic> > cholesky
156  = y.ldlt();
157  if((cholesky.vectorD().array() <= CONSTRAINT_TOLERANCE).any()) {
158  std::ostringstream message;
159  message << name << " is not positive definite. "
160  << name << "(0,0) is %1%.";
161  T_result tmp = policies::raise_domain_error<T_y>(function,
162  message.str().c_str(),
163  y(0,0), Policy());
164  if (result != 0)
165  *result = tmp;
166  return false;
167  }
168  return true;
169  }
170 
171 
172  template <typename T_y, typename T_result>
173  inline bool check_pos_definite(const char* function,
174  const Eigen::Matrix<T_y,Eigen::Dynamic,Eigen::Dynamic>& y,
175  const char* name,
176  T_result* result) {
177  return check_pos_definite(function,y,name,result,default_policy());
178  }
179 
180 
181  template <typename T>
182  inline bool check_pos_definite(const char* function,
183  const Eigen::Matrix<T,Eigen::Dynamic,Eigen::Dynamic>& y,
184  const char* name,
185  T* result = 0) {
186  return check_pos_definite(function,y,name,result,default_policy());
187  }
188 
189 
190 
191 
204  // FIXME: update warnings
205  template <typename T_y, typename T_result, class Policy>
206  inline bool check_cov_matrix(const char* function,
207  const Eigen::Matrix<T_y,Eigen::Dynamic,Eigen::Dynamic>& y,
208  const char* name,
209  T_result* result,
210  const Policy&) {
211  if (!check_size_match(function, y.rows(), y.cols(), result, Policy()))
212  return false;
213  if (!check_positive(function, y.rows(), "rows", result, Policy()))
214  return false;
215  if (!check_symmetric(function, y, name, result, Policy()))
216  return false;
217  if (!check_pos_definite(function, y, name, result, Policy()))
218  return false;
219  return true;
220  }
221 
222  template <typename T_y, typename T_result>
223  inline bool check_cov_matrix(const char* function,
224  const Eigen::Matrix<T_y,Eigen::Dynamic,Eigen::Dynamic>& y,
225  const char* name,
226  T_result* result) {
227  return check_cov_matrix(function,y,name,result,default_policy());
228  }
229 
230 
231  template <typename T>
232  inline bool check_cov_matrix(const char* function,
233  const Eigen::Matrix<T,Eigen::Dynamic,Eigen::Dynamic>& y,
234  const char* name,
235  T* result = 0) {
236  return check_cov_matrix(function,y,name,result,default_policy());
237  }
238 
239 
240 
256  // FIXME: update warnings
257  template <typename T_y, typename T_result, class Policy>
258  inline bool check_corr_matrix(const char* function,
259  const Eigen::Matrix<T_y,Eigen::Dynamic,Eigen::Dynamic>& y,
260  const char* name,
261  T_result* result,
262  const Policy&) {
263  if (!check_size_match(function, y.rows(), y.cols(), result, Policy()))
264  return false;
265  if (!check_positive(function, y.rows(), "rows", result, Policy()))
266  return false;
267  if (!check_symmetric(function, y, "y", result, Policy()))
268  return false;
270  k = 0; k < y.rows(); ++k) {
271  if (fabs(y(k,k) - 1.0) > CONSTRAINT_TOLERANCE) {
272  std::ostringstream message;
273  message << name << " is not a valid correlation matrix. "
274  << name << "(" << k << "," << k
275  << ") is %1%, but should be near 1.0";
276  T_result tmp
277  = policies::raise_domain_error<T_y>(function,
278  message.str().c_str(),
279  y(k,k), Policy());
280  if (result != 0)
281  *result = tmp;
282  return false;
283  }
284  }
285  if (!check_pos_definite(function, y, "y", result, Policy()))
286  return false;
287  return true;
288  }
289 
290 
291  template <typename T_y, typename T_result>
292  inline bool check_corr_matrix(const char* function,
293  const Eigen::Matrix<T_y,Eigen::Dynamic,Eigen::Dynamic>& y,
294  const char* name,
295  T_result* result) {
296  return check_corr_matrix(function,y,name,result,default_policy());
297  }
298 
299  template <typename T>
300  inline bool check_corr_matrix(const char* function,
301  const Eigen::Matrix<T,Eigen::Dynamic,Eigen::Dynamic>& y,
302  const char* name,
303  T* result = 0) {
304  return check_corr_matrix(function,y,name,result,default_policy());
305  }
306 
318  // FIXME: update warnings
319  template <typename T_covar, typename T_result, class Policy>
320  inline bool check_cov_matrix(const char* function,
321  const Eigen::Matrix<T_covar,Eigen::Dynamic,Eigen::Dynamic>& Sigma,
322  T_result* result,
323  const Policy&) {
324  if (!check_size_match(function, Sigma.rows(), Sigma.cols(),
325  result, Policy()))
326  return false;
327  if (!check_positive(function, Sigma.rows(), "rows", result, Policy()))
328  return false;
329  if (!check_symmetric(function, Sigma, "Sigma", result, Policy()))
330  return false;
331  return true;
332  }
333  template <typename T_covar, typename T_result>
334  inline bool check_cov_matrix(const char* function,
335  const Eigen::Matrix<T_covar,Eigen::Dynamic,Eigen::Dynamic>& Sigma,
336  T_result* result) {
337  return check_cov_matrix(function,Sigma,result,default_policy());
338 
339  }
340  template <typename T>
341  inline bool check_cov_matrix(const char* function,
342  const Eigen::Matrix<T,Eigen::Dynamic,Eigen::Dynamic>& Sigma,
343  T* result = 0) {
344  return check_cov_matrix(function,Sigma,result,default_policy());
345 
346  }
347 
348 
349 
350 
365  template <typename T_prob,
366  typename T_result,
367  class Policy>
368  bool check_simplex(const char* function,
369  const Eigen::Matrix<T_prob,Eigen::Dynamic,1>& theta,
370  const char* name,
371  T_result* result,
372  const Policy&) {
375  if (theta.size() == 0) {
376  std::string message(name);
377  message += " is not a valid simplex. %1% elements in the vector.";
378  T_result tmp = raise_domain_error<size_t, size_t>(function,
379  message.c_str(),
380  0,
381  Policy());
382  if (result != 0)
383  *result = tmp;
384  return false;
385  }
386  if (fabs(1.0 - theta.sum()) > CONSTRAINT_TOLERANCE) {
387  std::stringstream msg;
388  T_prob sum = theta.sum();
389  msg << "in function check_simplex(%1%), ";
390  msg << name << " is not a valid simplex.";
391  msg << " The sum of the elements should be 1, but is " << sum;
392  T_result tmp = raise_domain_error<T_result,T_prob>(function,
393  msg.str().c_str(),
394  sum,
395  Policy());
396  if (result != 0)
397  *result = tmp;
398  return false;
399  }
400  for (size_t n = 0; n < theta.size(); n++) {
401  if (!(theta[n] >= 0)) {
402  std::ostringstream stream;
403  stream << name << " is not a valid simplex."
404  << " The element at " << n
405  << " is %1%, but should be greater than or equal to 0";
406  T_result tmp
407  = raise_domain_error<T_result,T_prob>(function,
408  stream.str().c_str(),
409  theta[n],
410  Policy());
411  if (result != 0)
412  *result = tmp;
413  return false;
414  }
415  }
416  return true;
417  }
418  template <typename T_y,
419  typename T_result> // = typename T_prob_vector::value_type,
420  inline bool check_simplex(const char* function,
421  const Eigen::Matrix<T_y,Eigen::Dynamic,1>& theta,
422  const char* name,
423  T_result* result) {
424  return check_simplex(function,theta,name,result,default_policy());
425  }
426  template <typename T>
427  inline bool check_simplex(const char* function,
428  const Eigen::Matrix<T,Eigen::Dynamic,1>& theta,
429  const char* name,
430  T* result = 0) {
431  return check_simplex(function,theta,name,result,default_policy());
432  }
433 
434 
435 
436 
451  template <typename T_y, typename T_result, class Policy>
452  bool check_ordered(const char* function,
453  const Eigen::Matrix<T_y,Eigen::Dynamic,1>& y,
454  const char* name,
455  T_result* result,
456  const Policy&) {
458  typedef typename Eigen::Matrix<T_y,Eigen::Dynamic,1>::size_type size_t;
459  if (y.size() == 0) {
460  return true;
461  }
462  for (size_t n = 1; n < y.size(); n++) {
463  if (!(y[n] > y[n-1])) {
464  std::ostringstream stream;
465  stream << name << " is not a valid ordered vector."
466  << " The element at " << n
467  << " is %1%, but should be greater than the previous element, "
468  << y[n-1];
469  T_result tmp = raise_domain_error<T_result,T_y>(function,
470  stream.str().c_str(),
471  y[n],
472  Policy());
473  if (result != 0)
474  *result = tmp;
475  return false;
476  }
477  }
478  return true;
479  }
480  template <typename T_y, typename T_result>
481  bool check_ordered(const char* function,
482  const Eigen::Matrix<T_y,Eigen::Dynamic,1>& y,
483  const char* name,
484  T_result* result) {
485  return check_ordered(function,y,name,result,default_policy());
486  }
487  template <typename T>
488  bool check_ordered(const char* function,
489  const Eigen::Matrix<T,Eigen::Dynamic,1>& y,
490  const char* name,
491  T* result = 0) {
492  return check_ordered(function,y,name,result,default_policy());
493  }
494 
509  template <typename T_y, typename T_result, class Policy>
510  bool check_positive_ordered(const char* function,
511  const Eigen::Matrix<T_y,Eigen::Dynamic,1>& y,
512  const char* name,
513  T_result* result,
514  const Policy&) {
516  typedef typename Eigen::Matrix<T_y,Eigen::Dynamic,1>::size_type size_t;
517  if (y.size() == 0) {
518  return true;
519  }
520  if (y[0] < 0) {
521  std::ostringstream stream;
522  stream << name << " is not a valid positive_ordered vector."
523  << " The element at 0 is %1%, but should be postive.";
524  T_result tmp = raise_domain_error<T_result,T_y>(function,
525  stream.str().c_str(),
526  y[0],
527  Policy());
528  if (result != 0)
529  *result = tmp;
530  return false;
531  }
532  for (size_t n = 1; n < y.size(); n++) {
533  if (!(y[n] > y[n-1])) {
534  std::ostringstream stream;
535  stream << name << " is not a valid ordered vector."
536  << " The element at " << n
537  << " is %1%, but should be greater than the previous element, "
538  << y[n-1];
539  T_result tmp = raise_domain_error<T_result,T_y>(function,
540  stream.str().c_str(),
541  y[n],
542  Policy());
543  if (result != 0)
544  *result = tmp;
545  return false;
546  }
547  }
548  return true;
549  }
550  template <typename T_y, typename T_result>
551  bool check_positive_ordered(const char* function,
552  const Eigen::Matrix<T_y,Eigen::Dynamic,1>& y,
553  const char* name,
554  T_result* result) {
555  return check_positive_ordered(function,y,name,result,default_policy());
556  }
557  template <typename T>
558  bool check_positive_ordered(const char* function,
559  const Eigen::Matrix<T,Eigen::Dynamic,1>& y,
560  const char* name,
561  T* result = 0) {
562  return check_positive_ordered(function,y,name,result,default_policy());
563  }
564 
565 
566 
567  // error_handling functions for Eigen matrix
568 
569  template <typename T_y, typename T_result, class Policy>
570  inline bool check_not_nan(const char* function,
571  const Eigen::Matrix<T_y,Eigen::Dynamic,Eigen::Dynamic>& y,
572  const char* name,
573  T_result* result,
574  const Policy&) {
575  for (int i = 0; i < y.rows(); i++) {
576  for (int j = 0; j < y.cols(); j++) {
577  if (boost::math::isnan(y(i,j))) {
578  std::ostringstream message;
579  message << name << "[" << i << "," << j
580  << "] is %1%, but must not be nan!";
581  T_result tmp
582  = policies::raise_domain_error<T_y>(function,
583  message.str().c_str(),
584  y(i,j), Policy());
585  if (result != 0)
586  *result = tmp;
587  return false;
588  }
589  }
590  }
591  return true;
592  }
593  template <typename T_y, typename T_result>
594  inline bool check_not_nan(const char* function,
595  const Eigen::Matrix<T_y,Eigen::Dynamic,Eigen::Dynamic>& y,
596  const char* name,
597  T_result* result) {
598  return check_not_nan(function,y,name,result,default_policy());
599  }
600  template <typename T>
601  inline bool check_not_nan(const char* function,
602  const Eigen::Matrix<T,Eigen::Dynamic,Eigen::Dynamic>& y,
603  const char* name,
604  T* result = 0) {
605  return check_not_nan(function,y,name,result,default_policy());
606  }
607 
608  }
609 }
610 
611 #endif
internal::traits< Derived >::Index size_type
bool isnan(const stan::agrad::var v)
Checks if the given number is NaN.
var fabs(const var &a)
Return the absolute value of the variable (cmath).
Definition: agrad.hpp:2023
T_result raise_domain_error(const char *function, const char *message, const T_val &val, const Policy &)
bool check_corr_matrix(const char *function, const Eigen::Matrix< T_y, Eigen::Dynamic, Eigen::Dynamic > &y, const char *name, T_result *result, const Policy &)
Return true if the specified matrix is a valid correlation matrix.
bool check_simplex(const char *function, const Eigen::Matrix< T_prob, Eigen::Dynamic, 1 > &theta, const char *name, T_result *result, const Policy &)
Return true if the specified vector is simplex.
bool check_ordered(const char *function, const Eigen::Matrix< T_y, Eigen::Dynamic, 1 > &y, const char *name, T_result *result, const Policy &)
Return true if the specified vector is sorted into increasing order.
bool check_symmetric(const char *function, const Eigen::Matrix< T_y, Eigen::Dynamic, Eigen::Dynamic > &y, const char *name, T_result *result, const Policy &)
Return true if the specified matrix is symmetric.
T sum(const std::vector< T > &xs)
Return the sum of the values in the specified standard vector.
Definition: matrix.hpp:1131
bool check_pos_definite(const char *function, const Eigen::Matrix< T_y, Eigen::Dynamic, Eigen::Dynamic > &y, const char *name, T_result *result, const Policy &)
Return true if the specified matrix is positive definite.
bool check_positive_ordered(const char *function, const Eigen::Matrix< T_y, Eigen::Dynamic, 1 > &y, const char *name, T_result *result, const Policy &)
Return true if the specified vector contains only non-negative values and is sorted into increasing o...
bool check_not_nan(const char *function, const T_y &y, const char *name, T_result *result, const Policy &)
Checks if the variable y is nan.
bool check_cov_matrix(const char *function, const Eigen::Matrix< T_y, Eigen::Dynamic, Eigen::Dynamic > &y, const char *name, T_result *result, const Policy &)
Return true if the specified matrix is a valid covariance matrix.
const double CONSTRAINT_TOLERANCE
The tolerance for checking arithmetic bounds In rank and in simplexes.
bool check_positive(const char *function, const T_y &y, const char *name, T_result *result, const Policy &)
boost::math::policies::policy default_policy
Default error-handling policy from Boost.
bool check_size_match(const char *function, T_size1 i, T_size2 j, T_result *result, const Policy &)
Probability, optimization and sampling library.
Definition: agrad.cpp:6

     [ Stan Home Page ] © 2011–2012, Stan Development Team.