glibc/sysdeps/ieee754/dbl-64/mpsqrt.c
Ulrich Drepper bb3f4825c4 Update.
2003-11-28  Ulrich Drepper  <drepper@redhat.com>

	* sysdeps/x86_64/fpu/libm-test-ulps: Add some more minor changes
	to compensate other setup.

2003-11-27  Andreas Jaeger  <aj@suse.de>

	* sysdeps/x86_64/fpu/libm-test-ulps: Add ulps for new atan2 test.

	* math/libm-test.inc (atan2_test): Add test that run infinitly.
	Reported by "Willus" <etc231etc231@willus.com>.

2003-11-27  Michael Matz  <matz@suse.de>

	* sysdeps/ieee754/dbl-64/mpsqrt.c (fastiroot): Fix 64-bit problem
	with wrong types.

2003-11-28  Jakub Jelinek  <jakub@redhat.com>

	* posix/regexec.c (acquire_init_state_context): Make inline.
	Add always_inline attribute.
	(check_matching): Add BE macro.  Move if (cur_state->has_backref)
	into if (dfa->nbackref).
	(sift_states_backward): Fix comment.
	(transit_state): Add BE macro.  Move if (next_state->has_backref)
	into if (dfa->nbackref && next_state).  Don't check for next_state
	!= NULL twice.
	* posix/regcomp.c (peek_token): Use opr.ctx_type instead of opr.idx
	for ANCHOR.
	(parse_expression): Only call init_word_char if word context will be
	needed.

	* posix/bug-regex11.c (tests): Add new tests.

	* posix/tst-regex.c: Include getopt.h.
	(timing): New variable.
	(main): Set timing to 1 if --timing argument is present.
	Add 2 new tests.
	(run_test, run_test_backwards): Handle timing.

2003-11-27  Jakub Jelinek  <jakub@redhat.com>

	* posix/regex_internal.h (re_string_t): Remove mbs_case field.
	Add offsets, valid_raw_len, raw_len, raw_stop, mbs_allocated and
	offsets_needed fields.  Change icase, is_utf8 and map_notascii
	type from int bitfield to unsigned char.
	(MBS_ALLOCATED, MBS_CASE_ALLOCATED): Remove.
	(build_wcs_upper_buffer): Change prototype to return int.
	(re_string_peek_byte_case, re_string_fetch_byte_case): Remove
	defines, add prototypes.
	* posix/regex_internal.c (re_string_allocate): Don't initialize
	stop here.  Don't initialize mbs_case.  Set valid_raw_len.
	Use mbs_allocated instead of MBS_* macros.
	(re_string_construct): Don't initialize stop and valid_len here.
	Don't initialize mbs_case.  Use mbs_allocated instead of MBS_*
	macros.  Reallocate buffers if build_wcs_upper_buffer converted
	too few bytes.  Set valid_len to bufs_len only for single byte
	no translation and set in that case valid_raw_len as well.
	(re_string_realloc_buffers): Reallocate offsets if not NULL.
	Use mbs_allocated instead of MBS_ALLOCATED.  Don't reallocate
	mbs_case.
	(re_string_construct_common): Initialize raw_len, mbs_allocated,
	stop and raw_stop.
	(build_wcs_buffer): Apply pstr->trans before mbrtowc instead of
	after it.  Set valid_raw_len.  Don't set mbs_case.
	(build_wcs_upper_buffer): Return REG_NOERROR or REG_ESPACE.
	Only use the fast path if !pstr->offsets_needed.  Apply pstr->trans
	before mbrtowc instead of after it.  If upper case character
	uses different number of bytes than lower case, goto to the
	slow path.  Don't call towupper unnecessarily twice.  Set
	valid_raw_len as well.  Handle in the slow path the case if
	lower and upper case use different number of characters.
	Don't set mbs_case.
	(re_string_skip_chars): Use valid_raw_len instead of valid_len.
	(build_upper_buffer): Don't set mbs_case.  Add BE macro.  Set
	valid_raw_len.
	(re_string_translate_buffer): Set mbs instead of mbs_case.  Set
	valid_raw_len.
	(re_string_reconstruct): Use raw_len/raw_stop to initialize
	len/stop.  Clear valid_raw_len and offsets_needed when clearing
	valid_len.  Use mbs_allocated instead of MBS_* macros.
	Check original offset against valid_raw_len instead of valid_len.
	Remove mbs_case handling.  Adjust valid_raw_len together with
	valid_len.  If is_utf8 and looking for tip context, apply
	pstr->trans first.  If buffers start with partial multi-byte
	character, initialize mbs array as well if mbs_allocated.
	Check return value of build_wcs_upper_buffer.
	(re_string_peek_byte_case): New function.
	(re_string_fetch_byte_case): New function.
	(re_string_destruct): Use mbs_allocated instead of MBS_ALLOCATED.
	Don't free mbs_case.  Free offsets.
	* posix/regcomp.c (init_dfa): Only check if charset name is UTF-8
	if mb_cur_max == 6.
	* posix/regexec.c (re_search_internal): Initialize input.raw_stop
	as well.  Use valid_raw_len instead of valid_len when looking
	through fastmap.  Adjust registers through input.offsets.
	(extend_buffers): Allow build_wcs_upper_buffer to fail.
	* posix/bug-regex18.c (tests): Enable #ifdefed out tests.  Add new
	tests.
2003-11-29 06:13:09 +00:00

104 lines
4.4 KiB
C

/*
* IBM Accurate Mathematical Library
* written by International Business Machines Corp.
* Copyright (C) 2001 Free Software Foundation
*
* This program is free software; you can redistribute it and/or modify
* it under the terms of the GNU Lesser General Public License as published by
* the Free Software Foundation; either version 2.1 of the License, or
* (at your option) any later version.
*
* This program is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU Lesser General Public License for more details.
*
* You should have received a copy of the GNU Lesser General Public License
* along with this program; if not, write to the Free Software
* Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
*/
/****************************************************************************/
/* MODULE_NAME:mpsqrt.c */
/* */
/* FUNCTION:mpsqrt */
/* fastiroot */
/* */
/* FILES NEEDED:endian.h mpa.h mpsqrt.h */
/* mpa.c */
/* Multi-Precision square root function subroutine for precision p >= 4. */
/* The relative error is bounded by 3.501*r**(1-p), where r=2**24. */
/* */
/****************************************************************************/
#include "endian.h"
#include "mpa.h"
/****************************************************************************/
/* Multi-Precision square root function subroutine for precision p >= 4. */
/* The relative error is bounded by 3.501*r**(1-p), where r=2**24. */
/* Routine receives two pointers to Multi Precision numbers: */
/* x (left argument) and y (next argument). Routine also receives precision */
/* p as integer. Routine computes sqrt(*x) and stores result in *y */
/****************************************************************************/
double fastiroot(double);
void __mpsqrt(mp_no *x, mp_no *y, int p) {
#include "mpsqrt.h"
int i,m,ex,ey;
double dx,dy;
mp_no
mphalf = {0,{0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,
0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,
0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0}},
mp3halfs = {0,{0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,
0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,
0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0}};
mp_no mpxn,mpz,mpu,mpt1,mpt2;
/* Prepare multi-precision 1/2 and 3/2 */
mphalf.e =0; mphalf.d[0] =ONE; mphalf.d[1] =HALFRAD;
mp3halfs.e=1; mp3halfs.d[0]=ONE; mp3halfs.d[1]=ONE; mp3halfs.d[2]=HALFRAD;
ex=EX; ey=EX/2; __cpy(x,&mpxn,p); mpxn.e -= (ey+ey);
__mp_dbl(&mpxn,&dx,p); dy=fastiroot(dx); __dbl_mp(dy,&mpu,p);
__mul(&mpxn,&mphalf,&mpz,p);
m=mp[p];
for (i=0; i<m; i++) {
__mul(&mpu,&mpu,&mpt1,p);
__mul(&mpt1,&mpz,&mpt2,p);
__sub(&mp3halfs,&mpt2,&mpt1,p);
__mul(&mpu,&mpt1,&mpt2,p);
__cpy(&mpt2,&mpu,p);
}
__mul(&mpxn,&mpu,y,p); EY += ey;
return;
}
/***********************************************************/
/* Compute a double precision approximation for 1/sqrt(x) */
/* with the relative error bounded by 2**-51. */
/***********************************************************/
double fastiroot(double x) {
union {int i[2]; double d;} p,q;
double y,z, t;
int n;
static const double c0 = 0.99674, c1 = -0.53380, c2 = 0.45472, c3 = -0.21553;
p.d = x;
p.i[HIGH_HALF] = (p.i[HIGH_HALF] & 0x3FFFFFFF ) | 0x3FE00000 ;
q.d = x;
y = p.d;
z = y -1.0;
n = (q.i[HIGH_HALF] - p.i[HIGH_HALF])>>1;
z = ((c3*z + c2)*z + c1)*z + c0; /* 2**-7 */
z = z*(1.5 - 0.5*y*z*z); /* 2**-14 */
p.d = z*(1.5 - 0.5*y*z*z); /* 2**-28 */
p.i[HIGH_HALF] -= n;
t = x*p.d;
return p.d*(1.5 - 0.5*p.d*t);
}