src/libs/mapm/mapmhasn.c

   1
   2 /*
   3  *  M_APM  -  mapmhasn.c
   4  *
   5  *  Copyright (C) 2000 - 2007   Michael C. Ring
   6  *
   7  *  Permission to use, copy, and distribute this software and its
   8  *  documentation for any purpose with or without fee is hereby granted,
   9  *  provided that the above copyright notice appear in all copies and
  10  *  that both that copyright notice and this permission notice appear
  11  *  in supporting documentation.
  12  *
  13  *  Permission to modify the software is granted. Permission to distribute
  14  *  the modified code is granted. Modifications are to be distributed by
  15  *  using the file 'license.txt' as a template to modify the file header.
  16  *  'license.txt' is available in the official MAPM distribution.
  17  *
  18  *  This software is provided "as is" without express or implied warranty.
  19  */
  20
  21 /*
  22  *      $Id: mapmhasn.c,v 1.7 2007/12/03 01:53:33 mike Exp $
  23  *
  24  *      This file contains the Inverse Hyperbolic SIN, COS, & TAN functions.
  25  *
  26  *      $Log: mapmhasn.c,v $
  27  *      Revision 1.7  2007/12/03 01:53:33  mike
  28  *      Update license
  29  *
  30  *      Revision 1.6  2003/07/24 16:28:50  mike
  31  *      update arcsinh
  32  *
  33  *      Revision 1.5  2003/07/23 23:08:27  mike
  34  *      fix problem with arcsinh when input is a very large
  35  *      negative number.
  36  *
  37  *      Revision 1.4  2003/07/21 20:36:33  mike
  38  *      Modify error messages to be in a consistent format.
  39  *
  40  *      Revision 1.3  2003/03/31 21:53:21  mike
  41  *      call generic error handling function
  42  *
  43  *      Revision 1.2  2002/11/03 21:25:03  mike
  44  *      Updated function parameters to use the modern style
  45  *
  46  *      Revision 1.1  2000/04/03 18:16:29  mike
  47  *      Initial revision
  48  */
  49
  50 #include "m_apm_lc.h"
  51
  52 /****************************************************************************/
  53 /*
  54  *      arcsinh(x) == log [ x + sqrt(x^2 + 1) ]
  55  *
  56  *      also, use arcsinh(-x) == -arcsinh(x)
  57  */
  58 void    m_apm_arcsinh(M_APM rr, int places, M_APM aa)
  59 {
  60 M_APM   tmp0, tmp1, tmp2;
  61
  62 /* result is 0 if input is 0 */
  63
  64 if (aa->m_apm_sign == 0)
  65   {
  66    M_set_to_zero(rr);
  67    return;
  68   }
  69
  70 tmp0 = M_get_stack_var();
  71 tmp1 = M_get_stack_var();
  72 tmp2 = M_get_stack_var();
  73
  74 m_apm_absolute_value(tmp0, aa);
  75 m_apm_multiply(tmp1, tmp0, tmp0);
  76 m_apm_add(tmp2, tmp1, MM_One);
  77 m_apm_sqrt(tmp1, (places + 6), tmp2);
  78 m_apm_add(tmp2, tmp0, tmp1);
  79 m_apm_log(rr, places, tmp2);
  80
  81 rr->m_apm_sign = aa->m_apm_sign;                          /* fix final sign */
  82
  83 M_restore_stack(3);
  84 }
  85 /****************************************************************************/
  86 /*
  87  *      arccosh(x) == log [ x + sqrt(x^2 - 1) ]
  88  *
  89  *      x >= 1.0
  90  */
  91 void    m_apm_arccosh(M_APM rr, int places, M_APM aa)
  92 {
  93 M_APM   tmp1, tmp2;
  94 int     ii;
  95
  96 ii = m_apm_compare(aa, MM_One);
  97
  98 if (ii == -1)       /* x < 1 */
  99   {
 100    M_apm_log_error_msg(M_APM_RETURN, "\'m_apm_arccosh\', Argument < 1");
 101    M_set_to_zero(rr);
 102    return;
 103   }
 104
 105 tmp1 = M_get_stack_var();
 106 tmp2 = M_get_stack_var();
 107
 108 m_apm_multiply(tmp1, aa, aa);
 109 m_apm_subtract(tmp2, tmp1, MM_One);
 110 m_apm_sqrt(tmp1, (places + 6), tmp2);
 111 m_apm_add(tmp2, aa, tmp1);
 112 m_apm_log(rr, places, tmp2);
 113
 114 M_restore_stack(2);
 115 }
 116 /****************************************************************************/
 117 /*
 118  *      arctanh(x) == 0.5 * log [ (1 + x) / (1 - x) ]
 119  *
 120  *      |x| < 1.0
 121  */
 122 void    m_apm_arctanh(M_APM rr, int places, M_APM aa)
 123 {
 124 M_APM   tmp1, tmp2, tmp3;
 125 int     ii, local_precision;
 126
 127 tmp1 = M_get_stack_var();
 128
 129 m_apm_absolute_value(tmp1, aa);
 130
 131 ii = m_apm_compare(tmp1, MM_One);
 132
 133 if (ii >= 0)       /* |x| >= 1.0 */
 134   {
 135    M_apm_log_error_msg(M_APM_RETURN, "\'m_apm_arctanh\', |Argument| >= 1");
 136    M_set_to_zero(rr);
 137    M_restore_stack(1);
 138    return;
 139   }
 140
 141 tmp2 = M_get_stack_var();
 142 tmp3 = M_get_stack_var();
 143
 144 local_precision = places + 8;
 145
 146 m_apm_add(tmp1, MM_One, aa);
 147 m_apm_subtract(tmp2, MM_One, aa);
 148 m_apm_divide(tmp3, local_precision, tmp1, tmp2);
 149 m_apm_log(tmp2, local_precision, tmp3);
 150 m_apm_multiply(tmp1, tmp2, MM_0_5);
 151 m_apm_round(rr, places, tmp1);
 152
 153 M_restore_stack(3);
 154 }
 155 /****************************************************************************/