real.cpp

Go to the documentation of this file.

/*
 * real.cpp - some real valued function implementations
 *
 * Copyright (C) 2008 Stefan Jahn <stefan@lkcc.org>
 *
 * This is free software; you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation; either version 2, or (at your option)
 * any later version.
 *
 * This software is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with this package; see the file COPYING.  If not, write to
 * the Free Software Foundation, Inc., 51 Franklin Street - Fifth Floor,
 * Boston, MA 02110-1301, USA.
 *
 * $Id: real.cpp 1825 2011-03-11 20:42:14Z ela $
 *
 */

#if HAVE_CONFIG_H
# include <config.h>
#endif

#include <math.h>
#include <assert.h>

#include "consts.h"
#include "real.h"

#ifndef HAVE_ROUND
nr_double_t round (const nr_double_t arg) {
  return (arg > 0) ? floor (arg + 0.5) : ceil (arg - 0.5);
}
#endif /* HAVE_ROUND */

#ifndef HAVE_TRUNC
nr_double_t trunc (const double arg) {
  return arg > 0 ? floor (arg) : floor (arg + 1);
}
#endif /* HAVE_TRUNC */

#ifndef HAVE_ACOSH
nr_double_t acosh (const double arg) {
  return log (arg + sqrt (arg * arg - 1.0));
}
#endif /* HAVE_ACOSH */

#ifndef HAVE_ASINH
nr_double_t asinh (const double arg) {
  return log (arg + sqrt (arg * arg + 1.0));
}
#endif /* HAVE_ASINH */

unsigned int
factorial (unsigned int n) {
  unsigned int result = 1;

  /* 13! > 2^32 */
  assert (n < 13);

  if (n == 0)
    return 1;

  for (; n > 1; n--)
    result = result * n;

  return result;
}

nr_double_t real (const nr_double_t r) {
  return r;
}

nr_double_t imag (const nr_double_t r) {
  return 0.0;
}

nr_double_t norm (const nr_double_t r) {
  return r * r;
}

nr_double_t abs (const nr_double_t r) {      
  return fabs (r);   
}
     
nr_double_t conj (const nr_double_t r) {
  return r;
}

nr_double_t limexp (const nr_double_t r) {
  return r < M_LIMEXP ? exp (r) : exp (M_LIMEXP) * (1.0 + (r - M_LIMEXP));
}

nr_double_t signum (const nr_double_t d) {
  if (d == 0) return 0;
  return d < 0 ? -1 : 1;
}

nr_double_t sign (const nr_double_t d) {
  return d < 0 ? -1 : 1;
}

nr_double_t xhypot (const nr_double_t a, const nr_double_t b) {
  nr_double_t c = fabs (a);
  nr_double_t d = fabs (b);
  if (c > d) {
    nr_double_t e = d / c;
    return c * sqrt (1 + e * e);
  }
  else if (d == 0)
    return 0;
  else {
    nr_double_t e = c / d;
    return d * sqrt (1 + e * e);
  }
}

nr_double_t sinc (const nr_double_t d) {
  if (d == 0) return 1;
  return sin (d) / d;
}

nr_double_t fix (const nr_double_t d) {
  return (d > 0) ? floor (d) : ceil (d);
}

nr_double_t step (const nr_double_t d) {
  nr_double_t x = d;
  if (x < 0.0)
    x = 0.0;
  else if (x > 0.0)
    x = 1.0;
  else
    x = 0.5;
  return x;
}