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 */
m_apm_arcsinh(M_APM rr,int places,M_APM aa)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 */
m_apm_arccosh(M_APM rr,int places,M_APM aa)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 */
m_apm_arctanh(M_APM rr,int places,M_APM aa)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 /****************************************************************************/
156