Velocity Reviews - Computer Hardware Reviews

Velocity Reviews > Newsgroups > Programming > C Programming > sqrtl(), etc. for Cygwin

Reply
Thread Tools

sqrtl(), etc. for Cygwin

 
 
Fred K
Guest
Posts: n/a
 
      11-07-2012
On Wednesday, November 7, 2012 2:29:23 PM UTC-8, Bill Cunningham wrote:
> Keith Thompson wrote: [snip] > gcc doesn't have sqrtl() on any system. It's part of the runtime > library, not the compiler. That's not just a quibble; most or all > Linux-based systems use the GNU C library, but Cygwin apparently uses > something else. [...] Just for my own information and for future reference. Are functions like this sqrtl and other GNU C function OT here? Isn't this group restricted to ANSI and ISO C ? Bill


Read Keith Thompson's first reply:
"Both sqrtf and sqrtl were added to the ISO C standard in 1999"
 
Reply With Quote
 
 
 
 
Anand Hariharan
Guest
Posts: n/a
 
      11-07-2012
On Nov 7, 4:29*pm, "Bill Cunningham" <(E-Mail Removed)> wrote:
> Keith Thompson wrote:
>
> [snip]
>
> > gcc doesn't have sqrtl() on any system. *It's part of the runtime
> > library, not the compiler. *That's not just a quibble; most or all
> > Linux-based systems use the GNU C library, but Cygwin apparently uses
> > something else.

>
> [...]
>
> Just for my own information and for future reference. Are functions like
> this sqrtl and other GNU C function OT here? Isn't this group restricted to
> ANSI and ISO C ?
>
> Bill


Couple of lines below the point where you snipped his post, Keith
states "Both sqrtf and sqrtl were added to the ISO C standard in
1999".

- Anand
 
Reply With Quote
 
 
 
 
Bill Cunningham
Guest
Posts: n/a
 
      11-07-2012
Fred K wrote:
> On Wednesday, November 7, 2012 2:29:23 PM UTC-8, Bill Cunningham
> wrote:
>> Keith Thompson wrote: [snip] > gcc doesn't have sqrtl() on any
>> system. It's part of the runtime > library, not the compiler. That's
>> not just a quibble; most or all > Linux-based systems use the GNU C
>> library, but Cygwin apparently uses > something else. [...] Just for
>> my own information and for future reference. Are functions like this
>> sqrtl and other GNU C function OT here? Isn't this group restricted
>> to ANSI and ISO C ? Bill

>
> Read Keith Thompson's first reply:
> "Both sqrtf and sqrtl were added to the ISO C standard in 1999"


Oh I see I'm sorry. Missed that. So it's on topic since it's ISO C. But
GNU extensions are not?

Bill


 
Reply With Quote
 
Keith Thompson
Guest
Posts: n/a
 
      11-07-2012
"Steven G. Kargl" <(E-Mail Removed)> writes:
> On Wed, 07 Nov 2012 12:15:58 -0800, Keith Thompson wrote:

[...]
>> Of course you can always use sqrt() instead of sqrtl(), at the cost of
>> some loss of precision (but if you didn't care about that you probably
>> wouldn't be using long double in the first place).

>
> As long as one is disinteresting in x = inf, NaN, and possibly -0. You
> can do
>
> #include <math.h>
> long double
> foo_sqrtl(long double x)
> {
> long double r;
> r = sqrt((double)x);
> return ((r + (x / r))/2);
> }
>
> The tricky part is making the above compliant with IEEE 754.


Cool.

Experiment shows that it works correctly for 1.0 and 2.0. Therefore, by
mathematical induction, it works perfectly for all argument values.

}

(Incidentally, the cast is unnecessary but harmless.)

--
Keith Thompson (The_Other_Keith) http://www.velocityreviews.com/forums/(E-Mail Removed) <http://www.ghoti.net/~kst>
Will write code for food.
"We must do something. This is something. Therefore, we must do this."
-- Antony Jay and Jonathan Lynn, "Yes Minister"
 
Reply With Quote
 
James Dow Allen
Guest
Posts: n/a
 
      11-08-2012
On Nov 8, 6:40*am, pete <(E-Mail Removed)> wrote:

> I have some portable freestanding implementations of
> *double fs_log2(double x);
> *[et cetera]
> * *http://www.mindspring.com/~pfilandr/C/fs_math/fs_math.h
> * *http://www.mindspring.com/~pfilandr/C/fs_math/fs_math.c


Thank you! I've downloaded these and will start using them as
needed. I particularly like that they're "freestanding"
and in just two files. (Often I've given up trying to use
on-line code, rather than face the annoyance of tracking down
prerequisites.)

James
 
Reply With Quote
 
glen herrmannsfeldt
Guest
Posts: n/a
 
      11-08-2012
Keith Thompson <(E-Mail Removed)> wrote:

(snip)
>> #include <math.h>
>> long double
>> foo_sqrtl(long double x)
>> {
>> long double r;
>> r = sqrt((double)x);
>> return ((r + (x / r))/2);
>> }


(snip)
> Experiment shows that it works correctly for 1.0 and 2.0. Therefore, by
> mathematical induction, it works perfectly for all argument values.


> }


> (Incidentally, the cast is unnecessary but harmless.)


It reminds me of the narrowing casts required by Java.
(Except I think Java doesn't have long double yet, but when it does...)

-- glen
 
Reply With Quote
 
ImpalerCore
Guest
Posts: n/a
 
      11-08-2012
On Nov 7, 6:40*pm, pete <(E-Mail Removed)> wrote:
> ImpalerCore wrote:
>
> > On Nov 7, 1:40 pm, James Dow Allen <(E-Mail Removed)> wrote:
> > > I've been using Cygwin lately, and am largely satisfied.
> > > It has gcc, but lacks sqrtl(), atan2l(), etc.

>
> > > What's the recommended remedy?
> > > I Googled "source sqrtl" and got to codeforge.com
> > > but only AFTER passing one painful Gotcha test
> > > was I told I'd need another password.

>
> > > I'd bother with making a free(?) account if I were sure
> > > codeforge.com is the best piece of cheese in the maze,
> > > but it was just a random Google hit and
> > > it seemed wiser to ask the experts here.
> > > (I downloaded Cygwin just a few months ago --
> > > why doesn't its gcc have sqrtl?)

>
> > \code
> > #include <stdio.h>
> > #define c_sizeof_array(arr) (sizeof (arr) / sizeof ((arr)[0]))

>
> > const double tolerance = 0.0000001;
> > const int max_iterations = 10;

>
> > double heron_sqrt( double x );

>
> > int main( void )
> > {
> > * size_t idx;
> > * char tmp_str[40];

>
> > * double input[] = {
> > * * 0, 1, 10, 100, 1000, 10000, 0.1, 0.25, -1
> > * };

>
> > * for ( idx = 0; idx < c_sizeof_array (input); ++idx )
> > * {
> > * * snprintf( tmp_str, sizeof (tmp_str), "sqrt( %.2f )", input[idx]);
> > * * printf( "%-20s = %.10f\n", tmp_str, heron_sqrt( input[idx] ) );
> > * }

>
> > * return 0;
> > }

>
> > double heron_sqrt( double x )
> > {
> > * /*
> > * ** The algorithm used to calculate the square root 'x' of number 'S'
> > * ** is based on Heron's method.
> > * **
> > * ** Let x(0) be an arbitrary positive number.
> > * ** Let x(n+1) be the average of x(n) and 'S / x(n)'
> > * ** Repeat until the desired accuracy is achieved.
> > * **
> > * ** The closer x(0) is to the square root value, the fewer number
> > * ** of iterations are needed to converge to the accuracy required.
> > * **
> > * ** The global variable 'max_iterations' allows the person to
> > * ** limit the accuracy of the result.
> > * **/

>
> > * double sqrt_x = 0;
> > * double guess = 0;
> > * double error = 0;
> > * int itr = 0;

>
> > * if ( x == 0.0 ) {
> > * * sqrt_x = 0;
> > * }
> > * else
> > * {
> > * * sqrt_x = x * 0.5;

>
> > * * do
> > * * {
> > * * * guess = sqrt_x;
> > * * * sqrt_x = 0.5 * ( guess + x / guess );
> > * * * error = sqrt_x * sqrt_x - x;
> > * * * error = error < 0 ? -error : error;
> > * * * ++itr;
> > * * } while ( error > tolerance && itr < max_iterations );
> > * }

>
> > * return sqrt_x;
> > }
> > \endcode

>
> > There are some issues, like how you want to deal with negative
> > numbers, but you should be able to patch the algorithm up to create a
> > workable substitute for 'sqrtl' by replacing the type with 'long
> > double'.

>
> I prefer this implementation of Hero's method
> for finding square roots:
>
> double
> square_root(double x)
> {
> * * double root = x;
>
> * * if (root > 0) {
> * * * * double next = root / 2 + 0.5;
>
> * * * * do {
> * * * * * * root = next;
> * * * * * * next = (x / root + root) / 2;
> * * * * } while (root > next);
> * * }
> * * return root;
>
> }
>
> (square_root(x)) returns the value of the largest double
> which is less than or equal to the square root of (x).


Thanks for the alternate implementation. The version I gave was from
an example from when I was mucking around contract programming, but
with the contracts stripped out.

\code
#include <stdio.h>
#include <clover/macros.h>
#include <clover/contract.h>

const double tolerance = 0.0000001;
const int max_iterations = 10;

double heron_sqrt( double x );

int main( void )
{
size_t idx;
char tmp_str[40];

double input[] = {
0, 1, 10, 100, 1000, 10000, 0.1, 0.25, -1
};

for ( idx = 0; idx < c_sizeof_array (input); ++idx )
{
snprintf( tmp_str, sizeof (tmp_str), "sqrt( %.2f )", input[idx] );
printf( "%-20s = %.10f\n", tmp_str, heron_sqrt( input[idx] ) );
}

return 0;
}

double heron_sqrt( double x )
{
/*
* The algorithm used to calculate the square root 'x' of number 'S'
* is based on Heron's method.
*
* Let x(0) be an arbitrary positive number.
* Let x(n+1) be the average of x(n) and 'S / x(n)'
* Repeat until the desired accuracy is achieved.
*
* The closer x(0) is to the square root value, the fewer number
* of iterations are needed to converge to the accuracy required.
*
* The global variable 'max_iterations' allows the person to
* limit the accuracy of the result to cause postcondition
* violations.
*/

double sqrt_x = 0;
double guess = 0;
double error = 0;
int itr = 0;
C_POST( double x_sq; )

C_REQUIRE( x >= 0 );

if ( x == 0.0 ) {
sqrt_x = 0;
}
else
{
sqrt_x = x * 0.5;

do
{
guess = sqrt_x;
sqrt_x = 0.5 * ( guess + x / guess );
error = sqrt_x * sqrt_x - x;
error = error < 0 ? -error : error;
++itr;
} while ( error > tolerance && itr < max_iterations );
}

C_POST( x_sq = sqrt_x * sqrt_x; )
C_ENSURE( x_sq - tolerance <= x && x <= x_sq + tolerance );

return sqrt_x;
}
\endcode

C_REQUIRE and C_ENSURE are basically glorified asserts that separate
pre and post conditions (which are either elided or evaluated
depending on the level of contract evaluation), and C_POST is a macro
to hide or declare the variables and statements when evaluating
postconditions. I used 'max_iterations' as a driver to test
postcondition evaluation.

> I have some portable freestanding implementations of
> *double fs_log2(double x);
> *double fs_exp2(double x);
> *long double fs_powl(long double x, long double y);
> *long double fs_sqrtl(long double x);
> *long double fs_logl(long double x);
> *long double fs_expl(long double x);
> *long double fs_cosl(long double x);
> *long double fs_fmodl(long double x, long double y);
> declared here:
> * *http://www.mindspring.com/~pfilandr/C/fs_math/fs_math.h
> and defined here:
> * *http://www.mindspring.com/~pfilandr/C/fs_math/fs_math.c


I've visited your site several times in the past and appreciate your
generosity of sharing your code. I hope to do the same eventually
once I get my code to a respectable state.

> In terms of accuracy, my square root function is as good as any.
> The pow, log, and exp functions are not quite as good.
>
> I like the performance of my cosine function.
> Given:
> * * #include <math.h>
> * * #define DEGREES_TO_RADIANS(A) * ((A) * atan(1) / 45)
> on the implementation that I am using,
> using my fs_cos,
> the value of (fs_cos(DEGREES_TO_RADIANS(-3540)) - 0.5))
> is 9.992007e-016.
> The value obtatined by using the supplied standard library
> cos function in (cos(DEGREES_TO_RADIANS(-3540)) - 0.5)
> is -1.060133e-015.
> Mathematically, the value should be zero.


One of these days I need to take a deeper look at your math and sort
function implementations.

Best regards,
John D.
 
Reply With Quote
 
ImpalerCore
Guest
Posts: n/a
 
      11-08-2012
On Nov 7, 6:40*pm, Keith Thompson <(E-Mail Removed)> wrote:
> "Steven G. Kargl" <(E-Mail Removed)> writes:
> > On Wed, 07 Nov 2012 12:15:58 -0800, Keith Thompson wrote:

> [...]
> >> Of course you can always use sqrt() instead of sqrtl(), at the cost of
> >> some loss of precision (but if you didn't care about that you probably
> >> wouldn't be using long double in the first place).

>
> > As long as one is disinteresting in x = inf, NaN, and possibly -0. *You
> > can do

>
> > #include <math.h>
> > long double
> > foo_sqrtl(long double x)
> > {
> > * *long double r;
> > * *r = sqrt((double)x);
> > * *return ((r + (x / r))/2);
> > }

>
> > The tricky part is making the above compliant with IEEE 754.

>
> Cool.
>
> Experiment shows that it works correctly for 1.0 and 2.0. *Therefore, by
> mathematical induction, it works perfectly for all argument values.
>
> }
>
> (Incidentally, the cast is unnecessary but harmless.)


I believe that the cast is there to make the loss of precision more
explicit, e.g. for documentation rather than semantic purposes, and to
avoid potential slaps on the wrist from the compiler. I tend to favor
that casting style as well for any type conversion that loses
precision for statements that do not include a declaration of the
variable type.

\code snippet
double x;
....
int y = x; /* ok, type explicit in statement */
\endcode

\code snippet
double a;
int b;

....

b = (int)a;
\endcode

Best regards,
John D.
 
Reply With Quote
 
Noob
Guest
Posts: n/a
 
      11-08-2012
pete wrote:

> Steven G. Kargl wrote:
>
>> As long as one is disinteresting in x = inf, NaN, and possibly -0.

>
> Is that a hyphenated zero?


https://en.wikipedia.org/wiki/Signed_zero

 
Reply With Quote
 
Philip Lantz
Guest
Posts: n/a
 
      11-08-2012
Steven G. Kargl wrote:
> glen herrmannsfeldt wrote:
> > Steven G. Kargl wrote:
> >> #include <math.h>
> >> long double
> >> foo_sqrtl(long double x)
> >> {
> >> long double r;
> >> r = sqrt((double)x);
> >> return ((r + (x / r))/2);
> >> }

>
> > ... the last one should be
> >
> > return(x/r+(r-x/r)/2);

>
> I beg to differ.


Mathematically, ((r + x/r)/2) and (x/r + (r - x/r)/2) are equivalent. Is
there a numerical reason to prefer one over the other? (Obviously, if
there are no such concerns, one might as well use the simpler one.)

 
Reply With Quote
 
 
 
Reply

Thread Tools

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off
Trackbacks are On
Pingbacks are On
Refbacks are Off


Similar Threads
Thread Thread Starter Forum Replies Last Post
diff between pc-msvc pc-cygwin and svg-GDI - which is fastest vc Firefox 10 02-28-2004 01:05 AM
confirm unsubscribe from cygwin@cygwin.com cygwin-help@cygwin.com Python 0 09-05-2003 04:42 PM
WELCOME to cygwin@cygwin.com cygwin-help@cygwin.com Python 1 09-05-2003 07:46 AM
confirm unsubscribe from cygwin-announce@cygwin.com cygwin-announce-help@cygwin.com Python 0 09-05-2003 01:29 AM
confirm unsubscribe from cygwin@cygwin.com cygwin-help@cygwin.com Python 0 09-04-2003 06:34 PM



Advertisments