"Fossies" - the Fresh Open Source Software Archive 
Member "mapm_4.9.5a/mapmasn0.c" (21 Feb 2010, 4798 Bytes) of package /linux/misc/old/mapm-4.9.5a.tar.gz:
As a special service "Fossies" has tried to format the requested source page into HTML format using (guessed) C and C++ source code syntax highlighting (style:
standard) with prefixed line numbers and
code folding option.
Alternatively you can here
view or
download the uninterpreted source code file.
For more information about "mapmasn0.c" see the
Fossies "Dox" file reference documentation.
1
2 /*
3 * M_APM - mapmasn0.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: mapmasn0.c,v 1.8 2007/12/03 01:49:49 mike Exp $
23 *
24 * This file contains the 'ARC' family of functions; ARC-SIN,
25 * ARC-COS, ARC-TAN when the input arg is very close to 0 (zero).
26 *
27 * $Log: mapmasn0.c,v $
28 * Revision 1.8 2007/12/03 01:49:49 mike
29 * Update license
30 *
31 * Revision 1.7 2003/06/02 16:51:13 mike
32 * *** empty log message ***
33 *
34 * Revision 1.6 2003/06/02 16:49:48 mike
35 * tweak the decimal places
36 *
37 * Revision 1.5 2003/06/02 16:47:39 mike
38 * tweak arctan algorithm some more
39 *
40 * Revision 1.4 2003/05/31 22:38:07 mike
41 * optimize arctan by using fewer digits as subsequent
42 * terms get smaller
43 *
44 * Revision 1.3 2002/11/03 21:36:43 mike
45 * Updated function parameters to use the modern style
46 *
47 * Revision 1.2 2000/12/02 20:11:37 mike
48 * add comments
49 *
50 * Revision 1.1 2000/12/02 20:08:27 mike
51 * Initial revision
52 */
53
54 #include "m_apm_lc.h"
55
56 /****************************************************************************/
57 /*
58 Calculate arcsin using the identity :
59
60 x
61 arcsin (x) == arctan [ --------------- ]
62 sqrt(1 - x^2)
63
64 */
65 void M_arcsin_near_0(M_APM rr, int places, M_APM aa)
66 {
67 M_APM tmp5, tmp6;
68
69 tmp5 = M_get_stack_var();
70 tmp6 = M_get_stack_var();
71
72 M_cos_to_sin(tmp5, (places + 8), aa);
73 m_apm_divide(tmp6, (places + 8), aa, tmp5);
74 M_arctan_near_0(rr, places, tmp6);
75
76 M_restore_stack(2);
77 }
78 /****************************************************************************/
79 /*
80 Calculate arccos using the identity :
81
82 arccos (x) == PI / 2 - arcsin (x)
83
84 */
85 void M_arccos_near_0(M_APM rr, int places, M_APM aa)
86 {
87 M_APM tmp1, tmp2;
88
89 tmp1 = M_get_stack_var();
90 tmp2 = M_get_stack_var();
91
92 M_check_PI_places(places);
93 M_arcsin_near_0(tmp1, (places + 4), aa);
94 m_apm_subtract(tmp2, MM_lc_HALF_PI, tmp1);
95 m_apm_round(rr, places, tmp2);
96
97 M_restore_stack(2);
98 }
99 /****************************************************************************/
100 /*
101 calculate arctan (x) with the following series:
102
103 x^3 x^5 x^7 x^9
104 arctan (x) = x - --- + --- - --- + --- ...
105 3 5 7 9
106
107 */
108 void M_arctan_near_0(M_APM rr, int places, M_APM aa)
109 {
110 M_APM tmp0, tmp2, tmpR, tmpS, digit, term;
111 int tolerance, dplaces, local_precision;
112 long m1;
113
114 tmp0 = M_get_stack_var();
115 tmp2 = M_get_stack_var();
116 tmpR = M_get_stack_var();
117 tmpS = M_get_stack_var();
118 term = M_get_stack_var();
119 digit = M_get_stack_var();
120
121 tolerance = aa->m_apm_exponent - (places + 4);
122 dplaces = (places + 8) - aa->m_apm_exponent;
123
124 m_apm_copy(term, aa);
125 m_apm_copy(tmpS, aa);
126 m_apm_multiply(tmp0, aa, aa);
127 m_apm_round(tmp2, (dplaces + 8), tmp0);
128
129 m1 = 1L;
130
131 while (TRUE)
132 {
133 /*
134 * do the subtraction term
135 */
136
137 m_apm_multiply(tmp0, term, tmp2);
138
139 if ((tmp0->m_apm_exponent < tolerance) || (tmp0->m_apm_sign == 0))
140 {
141 m_apm_round(rr, places, tmpS);
142 break;
143 }
144
145 local_precision = dplaces + tmp0->m_apm_exponent;
146
147 if (local_precision < 20)
148 local_precision = 20;
149
150 m1 += 2;
151 m_apm_set_long(digit, m1);
152 m_apm_round(term, local_precision, tmp0);
153 m_apm_divide(tmp0, local_precision, term, digit);
154 m_apm_subtract(tmpR, tmpS, tmp0);
155
156 /*
157 * do the addition term
158 */
159
160 m_apm_multiply(tmp0, term, tmp2);
161
162 if ((tmp0->m_apm_exponent < tolerance) || (tmp0->m_apm_sign == 0))
163 {
164 m_apm_round(rr, places, tmpR);
165 break;
166 }
167
168 local_precision = dplaces + tmp0->m_apm_exponent;
169
170 if (local_precision < 20)
171 local_precision = 20;
172
173 m1 += 2;
174 m_apm_set_long(digit, m1);
175 m_apm_round(term, local_precision, tmp0);
176 m_apm_divide(tmp0, local_precision, term, digit);
177 m_apm_add(tmpS, tmpR, tmp0);
178 }
179
180 M_restore_stack(6); /* restore the 6 locals we used here */
181 }
182 /****************************************************************************/