1 #ifndef GENESIS_UTILS_MATH_DISTANCE_H_
2 #define GENESIS_UTILS_MATH_DISTANCE_H_
69 template <
class ForwardIterator>
70 double p_norm( ForwardIterator first, ForwardIterator last,
double p = 2.0 )
73 if( p < 1.0 || ( ! std::isfinite( p ) && ! std::isinf( p ))) {
74 throw std::runtime_error(
"Cannot calculate p-norm with p < 1.0" );
77 assert( std::isfinite( p ) || std::isinf( p ));
85 if( std::isfinite( *it ) ) {
86 if( std::isfinite( p )) {
87 sum += std::pow( std::abs( *it ), p );
89 sum = std::max(
sum, std::abs( *it ));
103 if( std::isfinite( p )) {
104 return std::pow(
sum, 1.0 / p );
119 inline double p_norm( std::vector<double>
const& vec,
double p = 2.0 )
121 return p_norm( vec.begin(), vec.end(), p );
130 template <
class ForwardIterator>
133 return p_norm( first, last, 1.0 );
144 return p_norm( vec.begin(), vec.end(), 1.0 );
153 template <
class ForwardIterator>
156 return p_norm( first, last, 2.0 );
167 return p_norm( vec.begin(), vec.end(), 2.0 );
177 template <
class ForwardIterator>
180 return p_norm( first, last, std::numeric_limits<double>::infinity() );
192 return p_norm( vec.begin(), vec.end(), std::numeric_limits<double>::infinity() );
215 template <
class ForwardIterator>
223 while( it_out != last ) {
224 if( std::isfinite( *it_out ) ) {
226 if( *it_out <= 0.0 ) {
227 throw std::invalid_argument(
228 "Cannot calculate Aitchison norm of non-positive values."
234 while( it_in != last ) {
235 if( std::isfinite( *it_in ) ) {
236 auto const ln = std::log( *it_out / *it_in );
254 return std::sqrt(
sum / ( 2.0 *
static_cast<double>( cnt )));
286 template <
class ForwardIteratorA,
class ForwardIteratorB>
288 ForwardIteratorA first_a, ForwardIteratorA last_a,
289 ForwardIteratorB first_b, ForwardIteratorB last_b,
293 if( p < 1.0 || ( ! std::isfinite( p ) && ! std::isinf( p ))) {
294 throw std::runtime_error(
"Cannot calculate p-norm distance with p < 1.0" );
297 assert( std::isfinite( p ) || std::isinf( p ));
303 if( std::isfinite( p )) {
306 sum += std::pow( std::abs( val_a - val_b ), p );
310 assert( std::isinf( p ));
312 sum = std::max(
sum, std::abs( val_a - val_b ) );
324 if( std::isfinite( p )) {
325 return std::pow(
sum, 1.0 / p );
340 std::vector<double>
const& vec_a, std::vector<double>
const& vec_b,
double p = 2.0
342 return p_norm_distance( vec_a.begin(), vec_a.end(), vec_b.begin(), vec_b.end(), p );
352 template <
class ForwardIteratorA,
class ForwardIteratorB>
354 ForwardIteratorA first_a, ForwardIteratorA last_a,
355 ForwardIteratorB first_b, ForwardIteratorB last_b
368 std::vector<double>
const& vec_a, std::vector<double>
const& vec_b
370 return p_norm_distance( vec_a.begin(), vec_a.end(), vec_b.begin(), vec_b.end(), 1.0 );
380 template <
class ForwardIteratorA,
class ForwardIteratorB>
382 ForwardIteratorA first_a, ForwardIteratorA last_a,
383 ForwardIteratorB first_b, ForwardIteratorB last_b
396 std::vector<double>
const& vec_a, std::vector<double>
const& vec_b
398 return p_norm_distance( vec_a.begin(), vec_a.end(), vec_b.begin(), vec_b.end(), 2.0 );
409 template <
class ForwardIteratorA,
class ForwardIteratorB>
411 ForwardIteratorA first_a, ForwardIteratorA last_a,
412 ForwardIteratorB first_b, ForwardIteratorB last_b
415 first_a, last_a, first_b, last_b, std::numeric_limits<double>::infinity()
429 std::vector<double>
const& vec_a, std::vector<double>
const& vec_b
432 vec_a.begin(), vec_a.end(), vec_b.begin(), vec_b.end(),
433 std::numeric_limits<double>::infinity()
482 #endif // include guard