"Fossies" - the Fresh Open Source Software Archive 
Member "mapm_4.9.5a/mapmgues.c" (21 Feb 2010, 5906 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 "mapmgues.c" see the
Fossies "Dox" file reference documentation.
1
2 /*
3 * M_APM - mapmgues.c
4 *
5 * Copyright (C) 1999 - 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: mapmgues.c,v 1.20 2007/12/03 01:52:55 mike Exp $
23 *
24 * This file contains the functions that generate the initial
25 * 'guesses' for the sqrt, cbrt, log, arcsin, and arccos functions.
26 *
27 * $Log: mapmgues.c,v $
28 * Revision 1.20 2007/12/03 01:52:55 mike
29 * Update license
30 *
31 * Revision 1.19 2003/07/21 20:03:49 mike
32 * check for invalid inputs to _set_double
33 *
34 * Revision 1.18 2003/05/01 21:58:45 mike
35 * remove math.h
36 *
37 * Revision 1.17 2003/04/16 16:52:47 mike
38 * change cbrt guess to use reciprocal value for new cbrt algorithm
39 *
40 * Revision 1.16 2003/04/11 14:18:13 mike
41 * add comments
42 *
43 * Revision 1.15 2003/04/10 22:28:35 mike
44 * .
45 *
46 * Revision 1.14 2003/04/09 21:33:19 mike
47 * induce same error for asin and acos
48 *
49 * Revision 1.13 2003/04/09 20:11:38 mike
50 * induce error of 10 ^ -5 in log guess for known starting
51 * point in the iterative algorithm
52 *
53 * Revision 1.12 2003/03/27 19:32:59 mike
54 * simplify log guess since caller guarantee's limited input magnitude
55 *
56 * Revision 1.11 2002/11/03 21:45:53 mike
57 * Updated function parameters to use the modern style
58 *
59 * Revision 1.10 2001/03/20 22:08:27 mike
60 * delete unneeded logic in asin guess
61 *
62 * Revision 1.9 2000/09/26 17:05:11 mike
63 * guess 1/sqrt instead of sqrt due to new sqrt algorithm
64 *
65 * Revision 1.8 2000/04/10 21:13:13 mike
66 * minor tweaks to log_guess
67 *
68 * Revision 1.7 2000/04/03 17:25:45 mike
69 * added function to estimate the cube root (cbrt)
70 *
71 * Revision 1.6 1999/07/18 01:57:35 mike
72 * adjust arc-sin guess for small exponents
73 *
74 * Revision 1.5 1999/07/09 22:32:50 mike
75 * optimize some functions
76 *
77 * Revision 1.4 1999/05/12 21:22:27 mike
78 * add more comments
79 *
80 * Revision 1.3 1999/05/12 21:00:51 mike
81 * added new sqrt guess function
82 *
83 * Revision 1.2 1999/05/11 02:10:12 mike
84 * added some comments
85 *
86 * Revision 1.1 1999/05/10 20:56:31 mike
87 * Initial revision
88 */
89
90 #include "m_apm_lc.h"
91
92 /****************************************************************************/
93 void M_get_sqrt_guess(M_APM r, M_APM a)
94 {
95 char buf[48];
96 double dd;
97
98 m_apm_to_string(buf, 15, a);
99 dd = atof(buf); /* sqrt algorithm actually finds 1/sqrt */
100 m_apm_set_double(r, (1.0 / sqrt(dd)));
101 }
102 /****************************************************************************/
103 /*
104 * for cbrt, log, asin, and acos we induce an error of 10 ^ -5.
105 * this enables the iterative routine to be more efficient
106 * by knowing exactly how accurate the initial guess is.
107 *
108 * this also prevents some corner conditions where the iterative
109 * functions may terminate too soon.
110 */
111 /****************************************************************************/
112 void M_get_cbrt_guess(M_APM r, M_APM a)
113 {
114 char buf[48];
115 double dd;
116
117 m_apm_to_string(buf, 15, a);
118 dd = atof(buf);
119 dd = log(dd) / 3.0; /* cbrt algorithm actually finds 1/cbrt */
120 m_apm_set_double(r, (1.00001 / exp(dd)));
121 }
122 /****************************************************************************/
123 void M_get_log_guess(M_APM r, M_APM a)
124 {
125 char buf[48];
126 double dd;
127
128 m_apm_to_string(buf, 15, a);
129 dd = atof(buf);
130 m_apm_set_double(r, (1.00001 * log(dd))); /* induce error of 10 ^ -5 */
131 }
132 /****************************************************************************/
133 /*
134 * the implementation of the asin & acos functions
135 * guarantee that 'a' is always < 0.85, so it is
136 * safe to multiply by a number > 1
137 */
138 void M_get_asin_guess(M_APM r, M_APM a)
139 {
140 char buf[48];
141 double dd;
142
143 m_apm_to_string(buf, 15, a);
144 dd = atof(buf);
145 m_apm_set_double(r, (1.00001 * asin(dd))); /* induce error of 10 ^ -5 */
146 }
147 /****************************************************************************/
148 void M_get_acos_guess(M_APM r, M_APM a)
149 {
150 char buf[48];
151 double dd;
152
153 m_apm_to_string(buf, 15, a);
154 dd = atof(buf);
155 m_apm_set_double(r, (1.00001 * acos(dd))); /* induce error of 10 ^ -5 */
156 }
157 /****************************************************************************/
158 /*
159 convert a C 'double' into an M_APM value.
160 */
161 void m_apm_set_double(M_APM atmp, double dd)
162 {
163 char *cp, *p, *ps, buf[64];
164
165 if (dd == 0.0) /* special case for 0 exactly */
166 M_set_to_zero(atmp);
167 else
168 {
169 sprintf(buf, "%.14E", dd);
170
171 if ((cp = strstr(buf, "E")) == NULL)
172 {
173 M_apm_log_error_msg(M_APM_RETURN,
174 "\'m_apm_set_double\', Invalid double input (likely a NAN or +/- INF)");
175
176 M_set_to_zero(atmp);
177 return;
178 }
179
180 if (atoi(cp + sizeof(char)) == 0)
181 *cp = '\0';
182
183 p = cp;
184
185 while (TRUE)
186 {
187 p--;
188 if (*p == '0' || *p == '.')
189 *p = ' ';
190 else
191 break;
192 }
193
194 ps = buf;
195 p = buf;
196
197 while (TRUE)
198 {
199 if ((*p = *ps) == '\0')
200 break;
201
202 if (*ps++ != ' ')
203 p++;
204 }
205
206 m_apm_set_string(atmp, buf);
207 }
208 }
209 /****************************************************************************/