Permission is granted to copy, distribute and/or modify this document under the terms of the GNU Free Documentation License, Version 1.2 or any later version published by the Free Software Foundation; with no Invariant Sections, no Front-Cover Texts, and no Back-Cover Texts. A copy of the license is included in the section entitled "GNU Free Documentation License".
The High Precision Arithmetic (HPA) library implements a high precision floating point arithmetic together with a comprehensive set of support functions. The general areas covered by these functions include:
The math library support includes evaluation of trigonometric, inverse trigonometric, hyperbolic, logarithm, and exponential functions at the same precision as the floating point math itself. The HPA library also supports high precision complex arithmetic and includes an Extended Precision Complex Math Library.
The HPA library comes from a branch of the source code of the CCMath library, which is a work by Daniel A. Atkinson. Daniel A. Atkinson was very kind to release the code of the CCMath Library under GNU Lesser General Public License. This made possible that Ivano Primi could modify, complete and redistribute this source code under the same terms.
So, the HPA (abbreviation of High Precision Arithmetic) Library is copyrighted by Ivano Primi <ivprimi (a) libero it> and Daniel A. Atkinson (the author of the original code). As it is for the source code of the CCMath Library, the source code of the HPA library is released 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.
The HPA library 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 library; if not, write to the Free Software Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA.
You can contact me by paper mail writing to the next address
Ivano Primi via Colle Cannetacce 50/A C.A.P. 00038 - Valmontone (ROMA) Italy .
If you prefer the electronic mail you can write to the address
<ivprimi (a) libero it> .
The functions forming the HPA library are all implemented in a portable fashion in the C language. The IEEE 754 standard for floating point hardware and software is assumed in the PC/Unix version of this library. The normal configuration of the library employs a floating point mantissa of 112 bits, or approximately 32 decimal digit precision. However, even higher precision is available as an option. An extended floating point number is represented as a combination of the following elements:
sign bit(s): 0 -> positive, 1 -> negative ;
exponent(e): 15-bits biased integer (bias=16383) ;
mantissa(m): 7 words of 16 bits length with the leading 1 explicitly represented .
Thus f = (-1)^s*2^[e-16383] *m , with 1 <= m < 2 . This format supports a dynamic range of:
2^16384 > f > 2^[-16383] or
1.19*10^4932 > f > 1.68*10^-[4932].
Special values of the exponent are:
all ones -> infinity (floating point overflow)
all zeros -> number = zero.
Underflow in operations is handled by a flush to zero. Thus, a number with the exponent zero and nonzero mantissa is invalid (not-a-number). From the point view of the HPA library, a complex number is simply a structure formed by two extended floating point numbers, representing respectively the real and the imaginary part of the complex number.
The HPA library is composed by two modules. The first one is formed
by the functions for the real arithmetic, that is to say by the functions
operating on real arguments. The second one is formed by all the functions
which manipulate complex arguments.
The lists of the functions which compose the HPA library
are in the header files xpre.h
and cxpre.h
.
xpre.h
is the header file for real arithmetic, it contains a definition
of the basic structure of an extended precision real number (struct xpr
),
and all the declarations of the functions provided by the library
to manipulate real arguments.
The numeric type struct xpr
can be used to declare and define
real variables, just as in
struct xpr q;
The size of a variable of struct xpr
type is given by
(2 x XDIM + 2) bytes, where XDIM is a constant defined in the file
xpre.h
(actually, in the file hpaconf.h
which is included by
xpre.h
).
cxpre.h
is the header file for complex arithmetic, it contains a definition
of the basic structure of an extended precision complex number
(struct cxpr
), and all the declarations of the functions provided
by the library to manipulate complex arguments.
The numeric type struct cxpr
can be used to declare and define
complex variables, just as in
struct cxpr q;
The size of a variable of struct cxpr
type is given by
(4 x XDIM + 4) bytes, where XDIM is the same constant as above.
Naturally, before declaring or defining variables of struct xpr
type and before using anyone of the functions declared in the header file
xpre.h
, you have to put the line
#include <xpre.h>
in your source code file.
Analogously, before declaring or defining variables of struct cxpr
type and before using anyone of the functions declared in the header file
cxpre.h
, you have to add the line
#include <cxpre.h>
to your source code file.
After including in your source code the header file xpre.h
or,
if you also need functions handling complex arguments, the header file
cxpre.h
(which automatically includes xpre.h
), you can start
to play with the HPA library, by defining all the variables and recalling
all the functions which are needed to do your computations.
In fact, the HPA library DOES NOT REQUIRE that a special initialization
routine must be called before any other function of the library.
Moreover, variables of struct xpr
or struct cxpr
type
DO NOT NEED to be initialized before they can be used.
So, with respect to these issues, the HPA library is different
from some other libraries for arbitrary precision computation,
like GNU MP(c) or MAPM(c). Actually, the HPA library is not for
arbitrary precision arithmetic, but only for high precision arithmetic,
if you know the difference.
When I wrote the HPA library, I tried to create a sort of namespace for all the identifiers used by the library. This has been achieved by sticking to the following rules:
xisNaN()
. Their names start
by x
(if they are defined in xpre.h
) or by cx
(if they are defined
in cxpre.h
) with the only exception of a few functions,
which however have a name ending by tox
or tocx
:
strtox(), strtocx(), atox(), atocx(),
dbltox(), dctocx(), flttox(), fctocx(),
inttox(), ictocx(), uinttox(), uctocx();
X
or CX
;
x
(real constants) or by cx
(complex constants) and the letter
which immediately follows this prefix is always uppercase, just as in
xZero
, xPi
, cxOne
, cxIU
(IU
stays for imaginary
unit);
int
type, whose name is xErrNo
.
A trivial program sampling the use of the HPA library, is given by:
#include <stdio.h> #include <xpre.h> int main (void) { struct xpr s; int i, n; do { printf ("Give me a number, n = ? \t"); scanf ("%d", &n); s = xZero; for (i = 0; i <= n; i++) s = xadd (s, xpr2(xOne, i), 0); printf ("The sum 2^0 + 2^1 + ... + 2^n is equal to\n"); xprxpr (s, 30); putchar ('\n'); } while (n > 0); return 0; }
This program takes in input from the user an integer value n
and prints on the screen the sum of the first n
powers of 2
.
In the program we use the functions xpr2()
and xprxpr()
.
xpr2(x, n)
, where n
is an integer, returns x* 2^n
, while
xprxpr(x, m)
, where m
is an integer, prints on the screen
the number x
with m
decimal digits after the dot (.) .
A last note: the HPA library is NOT thread safe. Some of the HPA internal data could get corrupted if multiple HPA functions are active at the same time. This is due to the fact that some functions of the HPA library use static variables to store information. So, the user should guarantee that only one thread is performing HPA functions. This can usually be achieved by a call to the operating system to obtain a semaphore, mutex, or critical code section so the operating system will guarantee that only one HPA thread will be active at a time.
During the use of the HPA library it could happen to pass
a function an illegal argument, that is to say an argument whose
value is against the mathematical definition of the function.
For instance, this occurs when a negative value is passed to
the function xsqrt()
. This function computes and returns
the square root of its argument, but the square root of a number
is defined only for non-negative numbers.
So, if x
is less than zero, then xsqrt(x)
can not
be computed and a mathematical error occurs (a so called domain error).
Another type of mathematical error occurs when the second argument
of the division function (xdiv()
) is zero: because it is impossible
to divide a number by zero, a division by zero error occurs.
What does it happen when a mathematical error is detected during the
execution of a function of the HPA library ?
Well, it depends upon the way the HPA library was compiled
when it was installed on the system where you are working.
If, during the installation process, the default setting was left
unchanged, then whenever a runtime error occurs within a function of the HPA
library, this function will set an external error indicator to
a suitable value. This value can be looked up later
to know what exactly went wrong.
The name of the variable of int
type used as error indicator
is xErrNo
. Before any function of the HPA library is
executed, the value of xErrNo
is 0
.
Then, when the first HPA function is called, if a mathematical error
occurs, then xErrNo
is set to a suitable positive value,
which indicates the exact type of the occurred error.
After, xErrNo
is modified if and only if, during the execution of
an HPA function, another mathematical error occurs. xErrNo
is never
reset to 0
by the HPA library, so one has to
zero xErrNo
before calling a function of the HPA library
in order to detect possible errors. A sample of it is given by:
#include <stdio.h> #include <xpre.h> extern int xErrNo; int main (void) { int n; struct xpr sr; do { printf ("Give me a number, n = ? \t"); scanf ("%d", &n); xErrNo = 0; sr = xsqrt (inttox (n)); if (xErrNo == 0) { printf ("The square root of %d is\n", n); xprxpr (sr, 30); putchar ('\n'); } else fprintf (stderr, "*** Error: Out of domain\n"); } while (n != 0); return 0; }
In this sample xErrNo
is reset to zero at each execution of the
do {...} while();
loop before the call to the xsqrt()
function.
However, the HPA library could be compiled to deal differently with
runtime errors.
For instance, in case of error a suitable message could be
printed onto stderr and the library could also cause the
termination of the calling program via EXIT(1).
At last, the library could also be compiled to ignore
at all any mathematical error (sigh !). To know how the routines of
the HPA library deal with errors is sufficient to examine the file
hpaconf.h
(which is automatically included by xpre.h
and
cxpre.h
). This file defines the macro:
xErrNo
is suitably set;
When the macro XERR_DFL is defined, the header file xpre.h
also defines the macros XENONE, XEDIV, XEDOM, XEBADEXP, XFPOFLOW and
XNERR:
#define XENONE 0 /* No error */ #define XEDIV 1 /* Division by zero */ #define XEDOM 2 /* Out of domain */ #define XEBADEXP 3 /* Bad exponent */ #define XFPOFLOW 4 /* Floating point overflow */ #define XNERR 4 /* Number of the non-null error codes */
that can be used, together with xErrNo
, to recover the exact type of
the error occurred during a call to a routine of the HPA library.
Now I will only concentrate on how compiling and building a program
which makes use of the HPA library.
Together with the HPA library is installed a little and simple
program called hpaconf
.
You can use it to quickly compile and build your programs.
If PREFIX is the root directory chosen to install the HPA library
(the default value for PREFIX is /usr/local), then hpaconf
should be installed inside PREFIX/bin. You can know where hpaconf
is installed by launching the command
which hpaconf
on your console or terminal.
In the following I will suppose that the directory PREFIX/bin
is included in your PATH environment variable (This is surely
true if the command which
was able to find hpaconf
).
hpaconf
recognizes four options:
-v
to return the version of HPA installed on your system
-c
to return the flags needed to compile using HPA
-l
to return the flags to link against HPA
-n
to print a newline at the end of the output.
The option -v
cannot be used together with the options -c
and -l
,
but it may be used together with -n
:
hpaconf -v
prints onto the standard output (console) the version of HPA installed on your system
hpaconf -v -n
or hpaconf -n -v
(order does not matter) behaves exactly the same but also prints a newline to force the following output to be written on the next line.
The options -c
and -l
cannot be used together with -v
, but they
can be used both at the same time and they can be accompanied by the option
-n
. Of course, order does not matter.
hpaconf -c
prints onto the standard output the flags needed to compile using HPA
hpaconf -l
prints onto the standard output the flags needed to link against HPA
hpaconf -c -l
or hpaconf -l -c
prints both the flags to compile using HPA and the flags to link against HPA.
If the -n
option is added, then the information printed is followed
by a newline. Why hpaconf
is so useful ?
I will give you an example. To compile the source
file example.c
you should tell the compiler where looking for
the header files of HPA and for the library itself; to do this
it is sufficient to specify the related paths with the options
-I and -L (at least if you are using GCC/G++ as C/C++ compiler).
However, in this way you are constrained to remember the path
where HPA was installed and this is quite uncomfortable.
With hpaconf
you can simply use the command
cc -c $(hpaconf -c) example.c
or
cc -c `hpaconf -c` example.c
to compile the file example.c
and obtain the object file example.o
.
Actually, the previous one is the right form of the command
for a shell sh compatible, like ash, bash or ksh.
If you are using another shell, probably the
right form to obtain the expansion of the command hpaconf -c
will be another one (see the manual of your preferred shell for this).
On GNU/Linux bash is the default shell for all users. If this is
not true for your machine, then ask your system administrator :)
Once you have obtained the object file example.o
, you may do the linkage
by using the command (for a shell sh-compatible):
cc example.o $(hpaconf -l) -o example
or
cc example.o `hpaconf -l` -o example
If you want, you may also compile and build at the same time by using (Do i need to repeat the next one is the right form only for a shell sh-compatible ?)
cc example.c $(hpaconf -c -l) -o example
or
cc example.c `hpaconf -c -l` -o example
which will compile example.c
and build the program example
.
Naturally, compiling and building at the same time is only
possible when the source code of your program is all contained in one
file.
The hpaconf
program may also be used to retrieve information
about the options used to compile the HPA library for the system
where you are working.
This information is displayed, together with hints about usage, when
hpaconf
is called with no options. This is the output obtained
on my personal machine:
ivano@darkstar[~]$ hpaconf *** Usage: hpaconf [-v] [-n] or *** hpaconf [-c] [-l] [-n] *** Meaning of the options: -v return the current version of the HPA library -c return flags to compile using HPA library -l return flags to link against HPA library -n insert a newline at the end ----- Features of the HPA library (including build options) ----- Size of an extended precision floating point value (in bytes): 16 Number of bits available for the sign: 1 Number of bits available for the exponent: 15 Number of bits available for the mantissa: 112 Decimal digits of accuracy: ~33 Dynamic range supported: 2^16384 > x > 2^(-16383) i.e. 1.19*10^4932 > x > 1.68*10^-(4932) In case of floating point error the global (extern) variable 'xErrNo' is suitably set
The first value shown after the header Features of the HPA library
is the size of a variable of struct xpr
type. When I installed the
HPA library on my machine, i chose to compile it by setting XDIM
to 7. So, a variable of type struct xpr
turns out to have a size of
16 = 2 * 7 + 2 bytes.
Since XDIM could have been set to another value on the system
where you are working (actually, XDIM could also have the values
11, 15, 19, 23, 27 and 31), the first value shown by hpaconf
could differ on your machine.
The next 2 values (bits available for sign and exponent) are
the same for all the installations.
The decimal digits of accuracy depend on the value of XDIM,
namely they increase together with XDIM (till to a maximum of
149 when XDIM is 31). As you can see, the actual value of XDIM
determines the accuracy provided by the mathematical functions of
the HPA library. Even if a larger value for XDIM implies
a greater accuracy, together with XDIM increases the memory
(and the time !) requested by the routines of the HPA library
to perform their computations.
The dynamic range supported by the HPA library is always the same
(or almost).
At end, the HPA library can be compiled to deal differently with
the error conditions (see previous section).
In the last line of the output of hpaconf
, you can find
information about treatment of the errors.
This information can also be retrieved by looking up the file
hpaconf.h
, as explained in the previous section.
This last way will turn out you very useful while writing your programs.
The file hpaconf.h
also defines the macro HPA_VERSION as a string
containing the version number of the HPA library currently in use.
The first module of the HPA library is made of functions for
Extended Precision Floating Point Arithmetic, functions of the
Extended Precision Math Library and applications of the Extended
Precision Arithmetic. They are all declared in the file xpre.h
together with some macros and numerical constants.
The header file xpre.h
also defines the structure xoutflags
:
struct xoutflags { short fmt, notat, sf, mfwd, lim; signed char padding, ldel, rdel; };
A structure of such kind is used by the output functions
xfout()
, xout()
and xsout()
to know how they have to print the
numbers.
The field notat
refers to the notation: it can
be equal to XOUT_SCIENTIFIC (scientific notation) or to XOUT_FIXED
(floating point notation). Both XOUT_SCIENTIFIC and XOUT_FIXED are
macros defined inside xpre.h
.
The field sf
refers to the sign: when sf
is not zero every non-negative
number is printed with a plus sign (+
) ahead.
The field mfwd
indicates the minimum field width to use in
printing numbers. When mfwd
is zero no minimum field width is used.
When mfwd
is negative, then the actual minimum field width is given
by -mfwd
and the printed number is left adjusted on the
field boundary. (The default is right justification).
lim
has a different meaning depending on the actual notation in use.
Together with the scientific notation, lim
gives the number of
decimal digits to the right of the decimal point
(lim+1
= total digits displayed). Otherwise, lim + 1
is the
number of significant digits displayed. When lim
is negative,
the default value (6
) is used.
At last, padding
defines the padding character to use together
with a non-zero minimum field width.
If padding
is negative, then the default padding char (i.e. the
blank character) is used.
The fields fmt
, ldel
and rdel
are ignored by the functions
xfout()
, xout()
and xsout()
. They are only used by the
functions cxfout()
, cxout()
and cxsout()
in order to
format and print complex numbers.
fmt
specifies the format to use in printing complex numbers.
The possible values for fmt
are XFMT_STD, XFMT_RAW and XFMT_ALT
(these macros are declared inside cxpre.h
!).
If fmt == XFMT_STD
, then the complex number (a, b)
is printed
using the notation a+bi
or a-bi
(depending on the sign of b
).
Naturally, a
and b
are printed under the rules above
exposed.
If fmt == XFMT_RAW
, then (a, b)
is printed in the form
a<two blank spaces>b
just as in
1.0 2.5
under the hypothesis that a = 1.0
and b = 2.5
.
At last, if fmt == XFMT_ALT
, then (a,b)
is printed
as
<left delimiter>a, b<right delimiter>
where <left_delimiter> and <right_delimiter> are the characters
given by the fields ldel
and rdel
respectively.
If ldel < 0
(rdel < 0
), then (
(or )
) is used as
default <left_delimiter> (<right_delimiter>).
Be careful ! None of the functions xfout()
, xout()
, xsout()
,
cxfout()
, cxout()
or cxsout()
adds a newline at the end
of the printed number !
The header file xpre.h
defines several constants. Between the
constants defined in xpre.h
there are those ones
corresponding to particular mathematical values:
extern const struct xpr xZero, xOne, xTwo, xTen; extern const struct xpr xPinf, xMinf, xNaN; extern const struct xpr xPi, xPi2, xPi4, xEe, xSqrt2; extern const struct xpr xLn2, xLn10, xLog2_e, xLog2_10, xLog10_e;
xZero (= 0), xOne (= 1), xTwo (= 2)and xTen (= 10) do not need a comment. xPi, xPi2, xPi4, xEe, xSqrt2, xLn2, xLn10, xLog2_e, xLog2_10, xLog10_e represent respectively the values PI (= 3.14159...), PI/2, PI/4, e (= 2.7182818...), square root of 2 (= 1.4142135...), natural logarithm of 2 and 10, base-2 logarithm of e and 10, 10-base logarithm of e.
xPinf, xMinf and xNan are special values: xPinf represents the value
+oo
(plus infty), xMinf the value -oo
(minus infty) and xNaN
is used to mean an invalid number (NaN stays for Not a Number).
xPinf and xMinf are usually returned by a function to signal a
floating point overflow (positive and negative respectively), while
xNaN is returned by the functions converting ASCII strings to floating
point numbers to indicate that the string given to them as argument
did not contain any valid number.
xPinf, xMinf and xNaN should never be used as arguments for functions,
since a such use has unpredictable results.
The arithmetic functions support the basic computation and input/output
operations needed for extended precision floating point mathematics. Some
of the operations supply capabilities designed to enhance the computational
efficiency of this arithmetic (e.g., xpwr
).
What follows is their complete list including the synopsis for each of them.
xadd
Add (subtract) two extended precision numbers.
struct xpr xadd(struct xpr s,struct xpr t,int f)
s
= structure containing first number;
t
= structure containing second number;
f
= control flag: if 0, then s
and t
are added,
else they are subtracted (s-t
).
The value returned by xadd()
is a structure containing the result
of the addition (subtraction). xadd()
can return xPinf or xMinf to
signal a floating point overflow.
xmul
Multiply two extended precision numbers.
struct xpr xmul(struct xpr s,struct xpr t)
s
= structure containing first number;
t
= structure containing second number.
The value returned by xmul()
is a structure containing the
product s*t
. It can be xPinf or xMinf in case of overflow.
xdiv
Divide one extended precision number by a second.
struct xpr xdiv(struct xpr s,struct xpr t)
s
= structure containing numerator;
t
= structure containing denominator.
The value returned by xdiv()
is a structure containing
the quotient s/t
.
xneg
Change sign (unary minus).
struct xpr xneg(struct xpr s)
s
= structure containing input number.
The value returned by xneg()
is a structure containing
its argument with the changed sign.
xabs
Compute absolute value.
struct xpr xabs(struct xpr s)
s
= structure containing input number.
The value returned by xabs()
is a structure containing
the absolute value of its argument.
x_exp
Extract the binary exponent.
int x_exp(const struct xpr *p)
p
= pointer to an extended precision number.
The value returned by x_exp()
is the binary exponent
(power of 2) of the number pointed to by p
.
x_neg
Test the sign of an extended precision number.
int x_neg(const struct xpr *p)
p
= pointer to an extended precision number.
The value returned by x_neg()
is a sign flag, with
0 meaning positive input, 1 negative input (the input
is, of course, the number pointed to by p
).
Note that neither x_exp() nor x_neg() alter the input number !
xpwr
Raise to integer powers.
struct xpr xpwr(struct xpr s,int n)
s
= structure containing input number;
n
= power desired.
The return value is a structure containing the n
th power
of the first argument.
xpr2
Multiplication by a power of 2.
struct xpr xpr2(struct xpr s,int m)
s
= structure containing input number;
m
= power of two desired.
The return value is a structure containing the product of
the first argument by the m
th power of two.
xpr2()
returns xZero in case of underflow, xPinf or
xMinf in case of overflow.
xpow
Power function.
struct xpr xpow (struct xpr x, struct xpr y)
x
= base;
y
= exponent.
The return value is a structure containing the power of the first argument raised to the second one. Note that the first argument must be positive, i.e. greater than zero.
xprcmp
Compare two extended precision numbers.
int xprcmp (const struct xpr *p, const struct xpr *q)
p
= pointer to first number;
q
= pointer to second number.
The value returned by xprcmp()
is a comparison flag, with
1 meaning *p
greater than *q
, 0 meaning *p
equal
to *q
and -1 meaning *p
less than *q
.
Note that the input numbers are not altered by xprcmp() !
xeq
Check if two numbers are or are not equal.
int xeq (struct xpr x1, struct xpr x2)
x1
= first number;
x2
= second number.
The return value is 0
if x1
and x2
are different, else a
non-null value.
xneq
Check if two numbers are or are not equal.
int xneq (struct xpr x1, struct xpr x2)
x1
= first number;
x2
= second number.
The return value is 0
if x1
and x2
are equal, else
a non-null value.
xgt
Check if a number is greater than another one.
int xgt (struct xpr x1, struct xpr x2)
x1
= first number;
x2
= second number.
The return value is 0
if x1
is less or equal to x2
,
else a non-null value.
xge
Check if a number is greater or equal to another one.
int xge (struct xpr x1, struct xpr x2)
x1
= first number;
x2
= second number.
The return value is 0
if x1
is less than x2
,
else a non-null value.
xlt
Check if a number is less than another one.
int xlt (struct xpr x1, struct xpr x2)
x1
= first number;
x2
= second number.
The return value is 0
if x1
is greater or equal to x2
,
else a non-null value.
xle
Check if a number is less or equal to another one.
int xle (struct xpr x1, struct xpr x2)
x1
= first number;
x2
= second number.
The return value is 0
if x1
is greater than x2
,
else a non-null value.
xisNaN
Check if a number is valid or not.
int xisNaN (const struct xpr *u)
u
= pointer to a structure containing a number.
The return value is 0
if *u
is a valid number,
else a non-null value.
Remark:
A number is considered invalid (not-a-number) when its exponent is zero but its mantissa isn't it.
xisPinf
Check if a number is equal to +oo
.
int xisPinf (const struct xpr *u)
u
= pointer to a structure containing a number.
The return value is 1
if *u
is equal to
xPinf (+oo
), 0
otherwise.
xisMinf
Check if a number is equal to -oo
.
int xisMinf (const struct xpr *u)
u
= pointer to a structure containing a number.
The return value is 1
if *u
is equal to
xMinf (-oo
), 0
otherwise.
xisordnumb
Check if a given number is an ordinary number.
int xisordnumb (const struct xpr *u)
u
= pointer to a structure containing a number.
The return value is 1
if *u
is a valid number
and is neither xPinf (+oo
) nor xMinf (-oo
),
else 0
.
xis0
Compare a number with zero.
int xis0 (const struct xpr *u)
u
= pointer to a structure containing a number.
The return value is 0
if *u
is not zero, else a non-zero value.
xnot0
Compare a number with zero.
int xnot0 (const struct xpr *u)
u
= pointer to a structure containing a number.
The return value is 0
if *u
is zero, else a non-zero value.
xsgn
Obtain the sign of a number.
int xsgn (const struct xpr *u)
u
= pointer to a structure containing a number.
The return value is 0
when *u
is zero or
is an invalid number (not-a-number), 1
if *u
is positive,
-1
if *u
is negative.
Remark:
xPinf is considered a positive value, xMinf a negative value.
xtodbl
Cast extended precision numbers to double precision ones.
double xtodbl(struct xpr s)
s
= structure containing extended precision input.
The return value is the double precision float corresponding
to s
.
dbltox
Convert double precision numbers to extended precision ones.
struct xpr dbltox(double y)
y
= double precision floating point input.
The return value is a structure containing extended equivalent
of y
.
xtoflt
float xtoflt(struct xpr s)
s
= structure containing extended precision input.
The return value is the single precision float corresponding
to s
.
flttox
Convert single precision numbers to extended precision ones.
struct xpr flttox(float y)
y
= single precision floating point input.
The return value is a structure containing extended equivalent
of y
.
inttox
Convert signed integers to extended precision numbers.
struct xpr inttox(long n)
n
= integer input.
The return value is a structure containing extended equivalent
of n
.
uinttox
Convert unsigned integers to extended precision numbers.
struct xpr uinttox(unsigned long n)
n
= integer input.
The return value is a structure containing extended equivalent
of n
.
strtox
Convert a floating point number, expressed as a decimal ASCII string in a form consistent with C, into the extended precision format.
struct xpr strtox (const char* s, char** endptr)
s
= pointer to null terminated ASCII string expressing a
decimal number;
endptr
= NULL or address of a pointer to char defined
outside strtox()
.
The value returned by strtox()
is a structure containing
the input number in the extended precision format.
Remarks:
The strtox()
function converts the initial portion of the string pointed to
by s
to its extended precision representation.
The expected form of the (initial portion of the) string is optional
leading white space as recognized by the standard library function
isspace()
, an optional plus (+
) or minus sign (-
) and
then a decimal number.
A decimal number consists of a nonempty sequence of decimal digits possibly
containing a radix character (decimal point, i.e. '.
'), optionally
followed by a decimal exponent. A decimal exponent consists of an
E
or e
, followed by an optional plus or minus sign, followed by
a non-empty sequence of decimal digits, and indicates multiplication by
a power of 10.
This function returns the converted value, if any. If the correct value would cause overflow, then xPinf or xMinf is returned, according to the sign of the value. If the correct value would cause underflow, xZero is returned. If no conversion is performed, xNaN is returned.
If endptr
is not NULL, a pointer to the character after the last character
used in the conversion is stored in the location referenced by
endptr
.
If no conversion is performed, the value of s
is
stored in the location referenced by endptr
.
atox
Convert a floating point number, expressed as a decimal ASCII string in a form consistent with C, into the extended precision format.
struct xpr atox(const char *s)
s
= pointer to null terminated ASCII string expressing a
decimal number.
The return value is a structure containing the input number in the extended precision format.
Remark:
The call atox(s)
is equivalent to strtox(s, NULL)
.
xfmod
This function is the extended precision analog of the standard library fmod function.
struct xpr xfmod(struct xpr s,struct xpr t,struct xpr* q)
s
= structure containing operand of fmod;
t
= structure containing base number (t!=0);
q
= pointer to store for output integer m.
The return value is the extended number with same sign as s
and absolute value less than that of t
, satisfying
s = m*t + x
if s*t>0
, or
s = -m*t + x
if s*t<0
.
xfrexp
This function is the extended precision analog of the standard library frexp function.
struct xpr xfrexp(struct xpr s,int *p)
s
= structure containing operand;
p
= pointer to store for output exponent e
.
The return value is the extended number satisfying
x = s*2^(-e)
with (-1 < x < +1)
.
xfrac
This function returns the fractional part of the input number.
struct xpr xfrac (struct xpr s)
s
= structure containing operand.
The return value is the fractional part of the number s
, with
the same sign as s
.
Remark:
The fractional part of the number s
is 0
if s
is an
integer number, else it is given by (-)0.xyz...
, where
xyz...
are the digits of s
following the radix character
(.
) in the decimal representation. xfrac(s)
has always
the same sign as s
.
xtrunc
This function returns the integer part of the input number.
struct xpr xtrunc (struct xpr s)
s
= structure containing operand.
The return value is the integer part of the number s
, with
the same sign as s
.
Remark:
The integer part of the number s
is given by (-)xyz...
, where
xyz...
are the digits of s
before the radix character (.
) in
its decimal representation.
xfix
Obtain the integer part of a number (2nd method).
struct xpr xfix (struct xpr s)
s
= structure containing operand.
Remark:
xfix()
is provided as an alternative to xtrunc()
. xfix()
tries
to take into account possible rounding errors and cancel their
effects. The use of xfix()
is strongly suggested whenever the argument
is presumed to be an integer number, but it is reasonable to expect
that some rounding errors make its actual value a bit different from
that one it should be.
For instance, when sizeof(struct xpr) == 64
, on my machine I
obtain
xtrunc (xdiv (inttox(100), inttox(100))) == xZero
while one could expect xOne
as result. This happens since
xdiv (inttox(100), inttox(100))
returns a number a bit lower
than 1
. On the other hand
xfix (xdiv (inttox(100), inttox(100))) == xOne .
Care that xfix()
introduces another type of rounding error
to give the correct answer in the cases similar to the previous
one. So, it is not always the right choice. :)
xround
Round an extended precision number to its nearest integer value (halfway cases are rounded away from zero).
struct xpr xround (struct xpr s)
s
= structure containing operand.
The return value is the integer value nearest to s
.
xceil
Round an extended precision number to the smallest integral value not less than it.
struct xpr xceil (struct xpr s)
s
= structure containing operand.
The return value is the smallest integral value not less than it.
xfloor
Round an extended precision number to the largest integral value not greater than it.
struct xpr xfloor (struct xpr s)
s
= structure containing operand.
The return value is the largest integral value not greater than it.
xpr_print
Print an extended precision number in scientific or floating point format onto a given file.
void xpr_print (FILE* stream, struct xpr u, int sc_not, int sign, int lim)
stream
= file where the number must be printed;
u
= structure containing number to be printed;
sc_not
= zero to mean floating point format, not zero to mean scientific format;
sign
= not zero to put a plus sign (+
) before the number when
it is non-negative (in case of negative number a minus sign (-
)
is always printed even if the parameter sign
is zero);
lim
= number of decimal digits to the right of the
decimal point (lim+1
= total digits displayed) in
case of scientific format, else number of significant digits - 1
(lim+1
= total of significant digits).
xpr_asprint
Convert an extended precision number to a string.
char* xpr_asprint (struct xpr u, int sc_not, int sign, int lim)
u
= structure containing number to be printed;
sc_not
= zero to mean floating point format, not zero to mean scientific format;
sign
= not zero to put a plus sign (+
) before the number when
it is non-negative (in case of negative number a minus sign (-
)
is always printed even if the parameter sign
is zero);
lim
= number of decimal digits to the right of the
decimal point (lim+1
= total digits displayed) in
case of scientific format, else number of significant digits - 1
(lim+1
= total of significant digits).
xpr_asprint()
returns the string with the converted number.
The memory for this string is calloc'ed inside the function.
xprxpr
Print an extended precision number in scientific format.
void xprxpr(struct xpr u,int lim)
u
= structure containing number to be printed;
lim
= number of decimal digits to the right of the
decimal point (lim+1
= total digits displayed).
Remark:
xprxpr(u, lim)
is equivalent to
xpr_print(stdout, u, 1, 0, lim)
xtoa
This function converts an extended precision number to a string. Scientific format is always used.
char* xtoa (struct xpr u,int lim)
u
= structure containing number to be printed;
lim
= number of decimal digits to the right of the
decimal point (lim+1
= total digits displayed).
Remark:
xtoa(u, lim)
is equivalent to
xpr_asprint(u, 1, 0, lim)
xbprint
Print an extended precision number in binary format.
void xbprint (FILE* stream, struct xpr u)
stream
= file where the number must be printed;
u
= structure containing number to be printed.
The xbprint()
function supports a bit oriented analysis of
rounding error effects. It always prints a newline (\n
)
at the end of the binary string.
xprint
Print an extended precision number as a string of hexadecimal numbers.
void xprint(FILE* stream, struct xpr u)
stream
= file where the number must be printed;
u
= structure containing number to be printed.
The xprint()
function supports a bit oriented analysis of
rounding error effects. It always prints a newline (\n
)
at the end.
xfout
Print an extended precision number on file according to a given set of I/O flags.
int xfout (FILE * stream, struct xoutflags ofs, struct xpr x)
stream
= file where the number must be printed;
ofs
= structure containing all the I/O flags;
x
= structure containing number to be printed.
The return value is 0
in case of success,
-1
to mean a failure.
Remark:
For the definition of struct xoutflags
and the meaning
of its fields see section "Real Arithmetic".
xfout()
does not add any newline
at the end of the printed number.
xout
Print an extended precision number on stdout
according to
a given set of I/O flags.
int xout (struct xoutflags ofs, struct xpr x)
ofs
= structure containing all the I/O flags;
x
= structure containing number to be printed.
The return value is 0
in case of success, -1
to mean a failure.
Remark:
For the definition of struct xoutflags
and the meaning
of its fields see section "Real Arithmetic".
xout()
does not add any newline
at the end of the printed number.
xsout
Write an extended precision number on a string according to a given set of I/O flags.
int xsout (char* s, unsigned long n, struct xoutflags ofs, struct xpr x)
s
= pointer to a buffer of characters (char);
n
= size of the buffer;
ofs
= structure containing all the I/O flags;
x
= structure containing number to be printed.
The return value is the number of the non-null characters written
to the buffer or, if it is greater or equal than n
, which would have been
written to the buffer if enough space had been available.
Remark:
For the definition of struct xoutflags
and the meaning of its fields
see section "Real Arithmetic".
xsout()
always adds a null character ('\0') at the end of
the written number.
xsout()
does not write more than n
bytes (including
the trailing '\0'). Thus, a return value of n
or more
means that the output was truncated. In this case, the contents of the
buffer pointed to by the first argument of xsout()
are
COMPLETELY unreliable.
The Extended Precision Math Library provides the elementary functions normally supported in a C math library. They are designed to provide full precision accuracy.
xsqrt
Compute the square root of an extended precision number.
struct xpr xsqrt(struct xpr z)
z
= structure containing the input number.
The return value is a structure containing the square root of the input number. Input of a negative argument results in a domain error.
xexp
Compute the exponential function.
struct xpr xexp(struct xpr z)
z
= structure containing function argument.
The return value is e
(the base of natural logarithms)
raised to the power of z
.
xexp2
Compute the base-2 exponential function.
struct xpr xexp2(struct xpr z)
z
= structure containing function argument.
The return value is 2
raised to the power of z
.
xexp10
Compute the base-10 exponential function.
struct xpr xexp10(struct xpr z)
z
= structure containing function argument.
The return value is 10
raised to the power of z
.
xlog
Compute natural (base e
) logarithms.
struct xpr xlog(struct xpr z)
z
= structure containing function argument.
This function returns the natural logarithm of its argument. A null or negative argument results in a domain error.
xlog2
Compute base-2 logarithms.
struct xpr xlog2(struct xpr z)
z
= structure containing function argument.
This function returns the base-2 logarithm of its argument. A null or negative argument results in a domain error.
xlog10
Compute base-10 logarithms.
struct xpr xlog10(struct xpr z)
z
= structure containing function argument.
This function returns the base-10 logarithm of its argument. A null or negative argument results in a domain error.
xtan
Tangent function.
struct xpr xtan(struct xpr z)
z
= structure containing function argument.
The return value is the tangent of z
, where z
is
given in radians.
xtan(z)
returns xPinf
when z
is equal to xPi2
(up to an integer multiple of xPi
),
xMinf
when z
is equal to -xPi2
(up to an
integer multiple of xPi
). In both cases a
domain error is produced.
xcos
Cosine function.
struct xpr xcos(struct xpr z)
z
= structure containing function argument.
The return value is the cosine of z
, where z
is
given in radians.
xsin
Sine function.
struct xpr xsin(struct xpr z)
z
= structure containing function argument.
The return value is the sine of z
, where z
is
given in radians.
xatan
Arc tangent function.
struct xpr xatan(struct xpr z)
z
= structure containing function argument.
This function returns the arc tangent of z
in radians and
the value is mathematically defined to be between -Pi/2
and Pi/2 (inclusive).
xasin
Arc sine function.
struct xpr xasin(struct xpr z)
z
= structure containing function argument.
This function returns the arc sine of z
in radians and
the value is mathematically defined to be between -Pi/2
and Pi/2 (inclusive). If z
falls outside the range
-1
to 1
, a domain error is produced.
xacos
Arc cosine function.
struct xpr xacos(struct xpr z)
z
= structure containing function argument.
This function returns the arc cosine of z
in radians and
the value is mathematically defined to be between zero
and Pi (inclusive). If z
falls outside the range
-1
to 1
, a domain error is produced.
xtanh
Hyperbolic tangent function.
struct xpr xtanh(struct xpr z)
z
= structure containing function argument.
The return value is the hyperbolic tangent of z
.
xcosh
Hyperbolic cosine function.
struct xpr xcosh(struct xpr z)
z
= structure containing function argument.
The return value is the hyperbolic cosine of z
.
xsinh
Hyperbolic sine function.
struct xpr xsinh(struct xpr z)
z
= structure containing function argument.
The return value is the hyperbolic sine of z
.
xatanh
Hyperbolic arc tangent function.
struct xpr xatanh(struct xpr z)
z
= structure containing function argument.
This function returns the hyperbolic arc tangent of z
.
If the absolute value of z
is greater than 1
,
then a domain error is produced.
xasinh
Hyperbolic arc sine function.
struct xpr xasinh(struct xpr z)
z
= structure containing function argument.
This function returns the hyperbolic arc sine of z
.
xacosh
Hyperbolic arc cosine function.
struct xpr xacosh(struct xpr z)
z
= structure containing function argument.
This function returns the hyperbolic arc cosine of z
.
If z
is less than 1
, a domain error is produced.
The Tchebycheff expansion supplied with the library can be used to
compute the Tchebycheff expansion coefficients of a function to an accuracy
of 32 digits at least. This ability is useful in developing high accuracy function
approximations, since the effect of rounding error on coefficients used in
double precision can effectively be eliminated with these inputs.
The functions provided to this purpose are xchcof()
and
xevtch()
.
xchcof
Compute the Tchebycheff expansion coefficients of a specified function f(x)
.
struct xpr* xchcof(int m,struct xpr (*xfunc)(struct xpr))
m
= index of the last coefficient (the computed coefficients will
be m+1
, indexed from 0
to m
);
xfunc
= pointer to user defined function returning extended
precision values of the function f()
;
The return value is the array of the computed coefficients. One has that
f(x) = c[0]/2 + Sum(k=1 to m) c[k]*Tk(x)
where Tk(x)
is the kth Tchebycheff polynomial.
Remarks:
The memory needed by the returned array is malloc'ed
inside the function xchcof()
.
To avoid memory leaks it should be explicitly freed through a call to free().
In case of insufficient memory xchcof()
will return NULL.
If m <= XMAX_DEGREE (= 50)
, then the array returned by xchcof()
will have exactly m+1
elements, indexed from 0
to m
.
But if m > XMAX_DEGREE
, then xchcof()
will behave as if
m
was equal to XMAX_DEGREE
. So, if m > XMAX_DEGREE
,
then the array returned by xchcof()
will have only
XMAX_DEGREE + 1
elements, indexed from 0
to XMAX_DEGREE
.
In other words, a value of m
greater than XMAX_DEGREE
is
ignored and replaced by XMAX_DEGREE
!
Remembering this fact is really important in order to avoid
buffer overruns on the array returned by xchcof()
!
XMAX_DEGREE
is a macro declared inside the header file
xpre.h
.
xevtch
Evaluate an extended precision Tchebycheff expansion.
struct xpr xevtch(struct xpr z,struct xpr *a,int m)
z
= structure containing function argument;
a
= structure array containing expansion coefficients;
m
= maximum index of coefficient array (dimension=m+1
).
The return value is the number given by the formula
f(z) = Sum(k=0 to m) a[k]*Tk(z)
where Tk(x)
is the kth Tchebycheff polynomial.
The second module of the HPA library is formed by
functions for Extended Precision Complex Arithmetic and functions
of the Extended Precision Complex Math Library.
They are all declared in the file cxpre.h
together with some
macros and numerical constants.
The header file cxpre.h
defines the constants
cxZero (= 0), cxOne (= 1) and cxIU (= 1i):
extern const struct cxpr cxZero; extern const struct cxpr cxOne; extern const struct cxpr cxIU;
which do not need a comment.
The functions for complex arithmetic support the basic computation
and input/output operations needed for extended precision complex mathematics.
Some of the operations supply capabilities designed to enhance the
computational efficiency of this arithmetic (e.g., cxpwr
and cxpow
).
What follows is their complete list including the synopsis for each of them.
cxreset
Make a new complex number from its real and imaginary parts.
struct cxpr cxreset (struct xpr re, struct xpr im)
re
= structure containing real part;
im
= structure containing imaginary part.
The value returned by cxreset()
is the complex
number having re
as its real part, im
as its imaginary part.
Remark:
cxreset()
is also available in the form of a macro:
#define CXRESET(re, im) (struct cxpr){re, im}
cxconv
Convert a real number into a complex one.
struct cxpr cxconv (struct xpr x)
re
= structure containing real part.
The value returned by cxconv()
is the complex
number having x
as its real part, zero
as its imaginary part.
Remark:
cxconv()
is also available in the form of a macro:
#define CXCONV(x) (struct cxpr){x, xZero}
cxre
Obtain the real part of a complex number.
struct xpr cxre (struct cxpr z)
z
= structure containing a complex number.
The value returned by cxre(z)
is the real part of z
.
Remark:
cxre()
is also available in the form of a macro:
#define CXRE(z) (z).re
cxim
Obtain the imaginary part of a complex number.
struct xpr cxim (struct cxpr z)
z
= structure containing a complex number.
The value returned by cxim(z)
is the imaginary part of z
.
Remark:
cxim()
is also available in the form of a macro:
#define CXIM(z) (z).im
cxswap
Swapping the real and the imaginary part of a complex number.
struct cxpr cxswap (struct cxpr z)
z
= structure containing a complex number.
The value returned by cxswap(z)
is the complex
number {cxim(z), cxre(z)}
.
Remark:
cxswap()
is also available in the form of a macro:
#define CXSWAP(z) (struct cxpr){(z).im, (z).re}
cxconj
Calculate the complex conjugate of a complex number.
struct cxpr cxconj (struct cxpr z)
z
= structure containing a complex number.
The cxconj()
function returns the complex conjugate value of
its argument, that is the value obtained by changing the sign of the
imaginary part.
cxneg
struct cxpr cxneg (struct cxpr z)
z
= structure containing a complex number.
cxneg(z)
returns the complex number -z
: if
z = a+ib
, then cxneg(z)
returns -a-ib
.
cxinv
Obtain the reciprocal of a complex number.
struct cxpr cxinv (struct cxpr z)
z
= structure containing a complex number.
The value returned by cxinv()
is the
reciprocal of its argument. If z
is zero,
then a division-by-zero error is produced.
cxabs
Calculate the absolute value of a complex number.
struct xpr cxabs (struct cxpr z)
z
= structure containing a complex number.
cxabs(z)
returns the absolute value of the complex number z
.
The result is a real number.
cxarg
Calculate the argument of a complex number.
struct xpr cxarg (struct cxpr z)
z
= structure containing a complex number.
cxarg(z)
returns the argument or phase angle
of the complex number z
. The result is a real number
in the range [-xPi,xPi).
cxadd
Add (subtract) two extended precision complex numbers.
struct cxpr cxadd (struct cxpr z1, struct cxpr z2, int k)
z1
= structure containing first number;
z2
= structure containing second number;
k
= control flag: if 0, then z1
and z2
are added,
else they are subtracted (z1-z2
).
The value returned by cxadd()
is a structure containing the result
of the addition (subtraction).
cxsum
Add two extended precision complex numbers.
struct cxpr cxsum (struct cxpr z1, struct cxpr z2)
z1
= structure containing first number;
z2
= structure containing second number.
The value returned by cxsum()
is a structure containing the result
of the addition z1 + z2
.
cxsub
Subtract two extended precision complex numbers.
struct cxpr cxsub (struct cxpr z1, struct cxpr z2)
z1
= structure containing first number;
z2
= structure containing second number.
The value returned by cxsub()
is a structure containing the result
of the subtraction z1 - z2
.
cxmul
Multiply two extended precision complex numbers.
struct cxpr cxmul (struct cxpr z1, struct cxpr z2)
z1
= structure containing first number;
z2
= structure containing second number.
The value returned by cxmul()
is a structure containing the
product z1 * z2
.
cxrmul
Multiply a complex number by a real one.
struct cxpr cxrmul (struct xpr c, struct cxpr z)
c
= structure containing a real number;
z
= structure containing a complex number.
The value returned by cxrmul()
is a structure containing the
product c * z
.
cxdrot
Multiply a complex number by the imaginary unit.
struct cxpr cxdrot (struct cxpr z)
z
= structure containing a complex number.
The value returned by cxdrot(z)
is the product of
z
by the imaginary unit.
cxrrot
Multiply a complex number by -1i
, where
1i
is the imaginary unit.
struct cxpr cxrrot (struct cxpr z)
z
= structure containing a complex number.
The value returned by cxrrot(z)
is the product of
z
by the negative imaginary unit.
cxdiv
Divide two extended precision complex numbers.
struct cxpr cxdiv (struct cxpr z1, struct cxpr z2)
z1
= structure containing first number;
z2
= structure containing second number.
The value returned by cxdiv(z1, z2)
is a structure containing
the quotient z1 / z2
. If z2
is zero, then a
division-by-zero error is produced.
cxgdiv
Gaussian division between complex numbers having real and imaginary part both integer.
struct cxpr cxgdiv (struct cxpr z1, struct cxpr z2)
z1
= structure containing first number;
z2
= structure containing second number.
After eventually rounding z1
and z2
by recalling
cxround()
(see below) on them, cxgdiv(z1, z2)
returns
a structure containing the quotient of the gaussian
division of z1
by z2
. If z2
is zero, then a
division-by-zero error is produced.
If you do not know what gaussian division means, probably you will never need this function :)
cxgmod
Remainder of the Gaussian division.
struct cxpr cxgmod (struct cxpr z1, struct cxpr z2)
z1
= structure containing first number;
z2
= structure containing second number.
After eventually rounding z1
and z2
by recalling
cxround()
(see below) on them, cxgmod(z1, z2)
returns
a structure containing the remainder of the gaussian
division of z1
by z2
. If z2
is zero, then a
division-by-zero error is produced.
If you do not know what gaussian division means, probably you will never need this function :)
cxidiv
Integer division.
struct cxpr cxidiv (struct cxpr z1, struct cxpr z2)
z1
= structure containing first number;
z2
= structure containing second number.
After eventually rounding z1
and z2
by recalling
cxround()
(see below) on them, cxidiv(z1, z2)
returns
a structure containing the quotient of the integer
division of z1
by z2
. If z2
is zero, then a
division-by-zero error is produced.
cxidiv()
is a smooth extension of the integer division
between real numbers.
cxmod
Remainder of the integer division .
struct cxpr cxmod (struct cxpr z1, struct cxpr z2)
After eventually rounding z1
and z2
by recalling
cxround()
(see below) on them, cxmod(z1, z2)
returns
a structure containing the remainder of the integer
division of z1
by z2
. If z2
is zero, then a
division-by-zero error is produced.
cxpwr
Raise to integer powers.
struct cxpr cxpwr (struct cxpr z, int n)
z
= structure containing input number;
n
= power desired.
The return value is a structure containing the n
th power
of the first argument.
cxpow
Power function.
struct cxpr cxpow (struct cxpr z1, struct cxpr z2)
z1
= base;
z2
= exponent.
The return value is a structure containing the
power of the first argument raised to the second one.
Note that the modulus of the first argument must be positive,
i.e. greater than zero, if the real part of z2
is
less or equal than zero, otherwise a bad-exponent error
will be produced.
cxsqr
Square of a number.
struct cxpr cxsqr (struct cxpr z)
This function returns the square of its argument.
cxroot
Nth root of a complex number.
struct cxpr cxroot (struct cxpr z, int i, int n)
cxroot(z,i,n)
returns the i
th branch of
the n
th root of z
.
If n
is zero or n
is negative and the
modulus of z
is zero, then a bad-exponent error
is produced.
cxsqrt
Principal branch of the square root of a complex number.
struct cxpr cxsqrt (struct cxpr z)
This function returns the principal branch of the square root of its argument.
cxprcmp
Compare two extended precision complex numbers.
struct cxprcmp_res cxprcmp (const struct cxpr* z1, const struct cxpr* z2)
z1
= pointer to first number;
z2
= pointer to second number.
The value returned by cxprcmp()
is a structure formed by two
comparison flags:
struct cxprcmp_res { int re, im; };
If the .re
field of the returned structure is:
+1, then z1->re
is greater than z2->re
,
0, then z1->re
is equal to z2->re
,
-1, then z1->re
is less than z2->re
.
The meaning of the .im
field is the same but with
respect to z1->im
and z2->im
respectively.
Note that the input numbers are not altered by cxprcmp() !
cxis0
Compare a complex number with zero.
int cxis0 (const struct cxpr* z)
z
= pointer to an extended precision complex number.
The return value is 0
if *z
is not zero, else a non-zero value.
cxnot0
Compare a complex number with zero.
int cxnot0 (const struct cxpr* z)
z
= pointer to an extended precision complex number.
The return value is 0
if *z
is zero, else a non-zero value.
cxeq
Check if two complex numbers are or are not equal.
int cxeq (struct cxpr z1, struct cxpr z2)
z1
= first number;
z2
= second number.
The return value is 0
if z1
and z2
are different, else a
non-null value.
cxneq
Check if two complex numbers are or are not equal.
int cxneq (struct cxpr z1, struct cxpr z2)
z1
= first number;
z2
= second number.
The return value is 0
if z1
and z2
are equal, else
a non-null value.
cxgt
Check if a complex number is greater than another one.
int cxgt (struct cxpr z1, struct cxpr z2)
z1
= first number;
z2
= second number.
The return value is 0
if z1
is not greater than z2
,
else a non-null value.
cxge
Check if a complex number is greater or equal to another one.
int cxge (struct cxpr z1, struct cxpr z2)
z1
= first number;
z2
= second number.
The return value is 0
if z1
is not greater or equal to z2
,
else a non-null value.
cxlt
Check if a complex number is less than another one.
int cxlt (struct cxpr z1, struct cxpr z2)
z1
= first number;
z2
= second number.
The return value is 0
if z1
is not less than z2
,
else a non-null value.
cxle
Check if a complex number is less or equal to another one.
int cxle (struct cxpr z1, struct cxpr z2)
z1
= first number;
z2
= second number.
The return value is 0
if z1
is not less or equal to z2
,
else a non-null value.
dctocx
Convert a double precision complex number to an extended precision number.
struct cxpr dctocx (double re, double im)
re
= real part of the double precision complex number;
im
= imaginary part of the double precision complex number.
The returned value of dctocx(re,im)
is the
extended precision equivalent of the complex number
(re, im)
.
cxtodc
Convert an extended precision complex number to a double precision complex number.
void cxtodc (const struct cxpr *z, double *re, double *im)
z
= pointer to an extended precision complex number;
re
= pointer to a double precision number;
im
= pointer to a double precision number.
cxtodc()
stores in its second and last argument respectively
the real and the imaginary part of the number
pointed to by its first argument.
fctocx
Convert a single precision complex number to an extended precision number.
struct cxpr fctocx (float re, float im)
re
= real part of the single precision complex number;
im
= imaginary part of the single precision complex number.
The returned value of fctocx(re,im)
is the
extended precision equivalent of the complex number
(re, im)
.
cxtofc
Convert an extended precision complex number to a single precision complex number.
void cxtofc (const struct cxpr *z, float *re, float *im)
z
= pointer to an extended precision complex number;
re
= pointer to a single precision number;
im
= pointer to a single precision number.
cxtofc()
stores in its second and last argument respectively
the real and the imaginary part of the number
pointed to by its first argument.
ictocx
Convert an integer complex number to an extended precision number.
struct cxpr ictocx (long re, long im)
re
= real part of the integer complex number;
im
= imaginary part of the integer complex number.
The returned value of ictocx(re,im)
is the
extended precision equivalent of the complex number
(re, im)
.
uctocx
Convert an integer complex number to an extended precision number.
struct cxpr uctocx (unsigned long re, unsigned long im)
re
= real part of the integer complex number;
im
= imaginary part of the integer complex number.
The returned value of uctocx(re,im)
is the
extended precision equivalent of the complex number
(re, im)
.
strtocx
Convert a floating point complex number, expressed as a decimal ASCII string in a form consistent with C, into the extended precision format.
struct cxpr strtocx (const char *s, char **endptr)
s
= pointer to null terminated ASCII string expressing a
complex number;
endptr
= NULL or address of a pointer to char defined
outside strtocx()
.
The value returned by strtocx()
is a structure containing
the input number in the extended precision format.
Remarks:
The strtocx()
function converts the initial portion of
the string pointed to
by s
to its extended precision representation.
The expected form of the (initial portion of the) string is optional
leading white space as recognized by the standard library function
isspace()
, an optional plus (+
) or minus sign (-
) and
then a decimal number.
A decimal number consists of a nonempty sequence of decimal digits possibly
containing a radix character (decimal point, i.e. '.
'), optionally
followed by a decimal exponent. A decimal exponent consists of an
E
or e
, followed by an optional plus or minus sign, followed by
a non-empty sequence of decimal digits, and indicates multiplication by
a power of 10.
After this decimal number there can be an i
character or,
alternatively, some optional white spaces,
an optional plus (+
) or minus sign (-
) and
then another decimal number followed by an i
character.
Examples of valid representations of complex numbers are:
"12"
,"34.56"
,".7895i"
,"-34.56-7.23i"
,
"-45.7 +23.4i"
.
This function returns the converted value, if any. If the correct value for the real or/and the imaginary part would cause overflow, then xPinf or xMinf is returned in the corresponding field, according to the sign of the value. If the correct value would cause underflow, xZero is returned. If no conversion is performed, xNaN is returned.
If endptr
is not NULL, a pointer to the character after the last character
used in the conversion is stored in the location referenced by
endptr
.
If no conversion is performed, the value of s
is
stored in the location referenced by endptr
.
atocx
Convert a floating point complex number, expressed as a decimal ASCII string in a form consistent with C, into the extended precision format.
struct cxpr atocx (const char *s)
s
= pointer to null terminated ASCII string expressing a
complex number.
The return value is a structure containing the input number in the extended precision format.
Remark:
The call atocx(s)
is equivalent to strtocx(s, NULL)
.
cxpr_asprint
Convert an extended precision complex number to a string.
char *cxpr_asprint (struct cxpr z, int sc_not, int sign, int lim)
z
= structure containing number to be printed;
sc_not
= zero to mean floating point notation, not zero to mean scientific notation;
sign
= not zero to put a plus sign (+
) before the real part
of the number when it is non-negative (in case of negative real part
a minus sign (-
) is always printed even if the parameter sign
is zero);
lim
= number of decimal digits to the right of the
decimal point (lim+1
= total digits displayed) in
case of scientific notation, else number of significant digits - 1
(lim+1
= total of significant digits).
cxpr_asprint()
returns the string with the converted number.
The memory for this string is calloc'ed inside the function.
cxpr_asprint()
uses always the format "a+bi"
in
the conversion.
cxtoa
This function converts an extended precision complex number to a string. Scientific notation is always used for both real and imaginary part.
char *cxtoa (struct cxpr z, int lim)
z
= structure containing number to be printed;
lim
= number of decimal digits to the right of the
decimal point (lim+1
= total digits displayed).
Remark:
cxtoa(z, lim)
is equivalent to
cxpr_asprint(z, 1, 0, lim)
cxfrac
Fractional part of both real and imaginary part.
struct cxpr cxfrac (struct cxpr z)
cxfrac(z)
returns {xfrac(z.re), xfrac(z.im)}
.
cxtrunc
Integer part of both real and imaginary part.
struct cxpr cxtrunc (struct cxpr z)
cxtrunc(z)
returns {xtrunc(z.re), xtrunc(z.im)}
.
cxfix
Integer part of both real and imaginary part (2nd method).
struct cxpr cxfix (struct cxpr z)
cxfix(z)
returns {xfix(z.re), xfix(z.im)}
.
cxround
Rounding real and imaginary part to the nearest integer values (halfway cases are rounded away from zero).
struct cxpr cxround (struct cxpr z)
cxround(z)
returns {xround(z.re), xround(z.im)}
.
cxfloor
Rounding real and imaginary part to the largest integral values not greater than them.
struct cxpr cxfloor (struct cxpr z)
cxfloor(z)
returns {xfloor(z.re), xfloor(z.im)}
.
cxceil
Rounding real and imaginary part to the smallest integral values not less than them.
struct cxpr cxceil (struct cxpr z)
cxceil(z)
returns {xceil(z.re), xceil(z.im)}
.
cxpr_print
Print an extended precision complex number in scientific or floating point notation onto a given file.
void cxpr_print (FILE * stream, struct cxpr z, int sc_not, int sign, int lim)
stream
= file where the number must be printed;
z
= structure containing number to be printed;
sc_not
= zero to mean floating point notation, not zero to mean scientific notation;
sign
= not zero to put a plus sign (+
) before the real part
of the number when it is non-negative (in case of negative real part
a minus sign (-
) is always printed even if the parameter sign
is zero);
lim
= number of decimal digits to the right of the
decimal point (lim+1
= total digits displayed) in
case of scientific notation, else number of significant digits - 1
(lim+1
= total of significant digits).
cxprcxpr
Print an extended precision complex number in scientific notation.
void cxprcxpr (struct cxpr z, int m)
z
= structure containing number to be printed;
lim
= number of decimal digits to the right of the
decimal point (lim+1
= total digits displayed).
Remark:
cxprcxpr(z, lim)
is equivalent to
cxpr_print(stdout, z, 1, 0, lim)
cxprint
Print an extended precision complex number as a couple of strings of hexadecimal numbers.
void cxprint (FILE * stream, struct cxpr z)
stream
= file where the number must be printed;
z
= structure containing number to be printed.
The cxprint()
function supports a bit oriented analysis of
rounding error effects. It always prints a newline (\n
)
at the end.
cxfout
Print an extended precision complex number on file according to a given set of I/O flags.
int cxfout (FILE * stream, struct xoutflags ofs, struct cxpr z)
stream
= file where the number must be printed;
ofs
= structure containing all the I/O flags;
z
= structure containing number to be printed.
The return value is 0
in case of success,
-1
to mean a failure.
Remark:
For the definition of struct xoutflags
and the meaning
of its fields see section "Real Arithmetic".
cxfout()
does not add any newline
at the end of the printed number.
cxout
Print an extended precision complex number on stdout
according to
a given set of I/O flags.
int cxout (struct xoutflags ofs, struct cxpr z)
ofs
= structure containing all the I/O flags;
z
= structure containing number to be printed.
The return value is 0
in case of success, -1
to mean a failure.
Remark:
For the definition of struct xoutflags
and the meaning
of its fields see section "Real Arithmetic".
cxout()
does not add any newline
at the end of the printed number.
cxsout
Write an extended precision complex number on a string according to a given set of I/O flags.
unsigned long cxsout (char *s, unsigned long n, struct xoutflags ofs, struct cxpr z)
s
= pointer to a buffer of characters (char);
n
= size of the buffer;
ofs
= structure containing all the I/O flags;
z
= structure containing number to be printed.
The return value is the number of the non-null characters written
to the buffer or, if it is greater or equal than n
, which would have been
written to the buffer if enough space had been available.
Remark:
For the definition of struct xoutflags
and the meaning of its fields
see section "Real Arithmetic".
cxsout()
always adds a null character ('\0') at the end of
the written number.
cxsout()
does not write more than n
bytes (including
the trailing '\0'). Thus, a return value of n
or more
means that the output was truncated. In this case, the contents of the
buffer pointed to by the first argument of cxsout()
are
COMPLETELY unreliable.
These functions provide the elementary function evaluations normally supported in a complex math library. They are designed to provide full precision accuracy.
cxexp
Compute the complex exponential function.
struct cxpr cxexp (struct cxpr z)
z
= structure containing function argument.
The return value is e
(the base of natural logarithms)
raised to the power of z
.
cxexp2
Compute the base-2 complex exponential function.
struct cxpr cxexp2 (struct cxpr z)
z
= structure containing function argument.
The return value is 2
raised to the power of z
.
cxexp10
Compute the base-10 complex exponential function.
struct cxpr cxexp10 (struct cxpr z)
z
= structure containing function argument.
The return value is 10
raised to the power of z
.
cxlog
Compute natural (base e
) logarithm of a complex number.
struct cxpr cxlog (struct cxpr z)
z
= structure containing function argument.
This function returns the natural logarithm of its argument. A null argument results in a domain error. The imaginary part of the result is chosen in the interval [-xPi,xPi).
cxlog2
Compute base-2 logarithm of a complex number.
struct cxpr cxlog2 (struct cxpr z)
z
= structure containing function argument.
This function returns the base-2 logarithm of its argument. A null argument results in a domain error.
cxlog10
Compute base-10 logarithm of a complex number.
struct cxpr cxlog10 (struct cxpr z)
z
= structure containing function argument.
This function returns the base-10 logarithm of its argument. A null argument results in a domain error.
cxlog_sqrt
Compute natural (base e
) logarithm of the
principal branch of the square root of a complex number.
struct cxpr cxlog_sqrt (struct cxpr z)
z
= structure containing function argument.
This function returns the natural logarithm of the principal branch of the square root of its argument. A null argument results in a domain error. The imaginary part of the result is chosen in the interval [-xPi2,xPi2).
cxtan
Complex tangent function.
struct cxpr cxtan (struct cxpr z)
z
= structure containing function argument.
The return value is the tangent of the complex number
stored in z
.
cxtan(z)
yields a domain error when the imaginary
part of z
is null and the real part of it is equal
to xPi2
up to an integer multiple of xPi
.
cxcos
Complex cosine function.
struct cxpr cxcos (struct cxpr z)
z
= structure containing function argument.
The return value is the cosine of the complex number
stored in z
.
cxsin
Complex sine function.
struct cxpr cxsin (struct cxpr z)
z
= structure containing function argument.
The return value is the sine of the complex number
stored in z
.
cxatan
Complex arc tangent function.
struct cxpr cxatan (struct cxpr z)
z
= structure containing function argument.
This function returns the arc tangent of the complex
number z
, i.e. a number w
such that
z = tan(w)
.
If the real part of z
is null and the
imaginary part is equal to +1
or -1
,
then a domain error is produced.
cxasin
Complex arc sine function.
struct cxpr cxasin (struct cxpr z)
z
= structure containing function argument.
This function returns the arc sine of the complex
number z
, i.e. a number w
such that
z = sin(w)
.
cxacos
Complex arc cosine function.
struct cxpr cxacos (struct cxpr z)
z
= structure containing function argument.
This function returns the arc cosine of the complex
number z
, i.e. a number w
such that
z = cos(w)
.
cxtanh
Complex hyperbolic tangent function.
struct cxpr cxtanh (struct cxpr z)
z
= structure containing function argument.
The return value is the hyperbolic tangent
of the complex number stored in z
.
cxtanh(z)
yields a domain error when the real
part of z
is null and the imaginary part of it is equal
to xPi2
up to an integer multiple of xPi
.
cxcosh
Complex hyperbolic cosine function.
struct cxpr cxcosh (struct cxpr z)
z
= structure containing function argument.
The return value is the hyperbolic cosine of
the complex number stored in z
.
cxsinh
Complex hyperbolic sine function.
struct cxpr cxsinh (struct cxpr z)
z
= structure containing function argument.
The return value is the hyperbolic sine of
the complex number stored in z
.
cxatanh
Complex hyperbolic arc tangent function.
struct cxpr cxatanh (struct cxpr z)
z
= structure containing function argument.
This function returns the hyperbolic arc tangent of
the complex number z
, i.e. a number w
such that
z = tanh(w)
.
If the imaginary part of z
is null and the
real part of it is equal to +1
or -1
,
then a domain error is produced.
cxasinh
Complex hyperbolic arc sine function.
struct cxpr cxasinh (struct cxpr z)
z
= structure containing function argument.
This function returns the hyperbolic arc sine of
the complex number z
, i.e. a number w
such that
z = sinh(w)
.
cxacosh
Complex hyperbolic arc cosine function.
struct cxpr cxacosh (struct cxpr z)
z
= structure containing function argument.
This function returns the hyperbolic arc cosine of
the complex number z
, i.e. a number w
such that
z = cosh(w)
.
FINAL REMARK:
The header file cxpre.h
also defines the macros
cxconvert
, cxdiff
, cxprod
and cxipow
:
#define cxconvert cxconv #define cxdiff cxsub #define cxprod cxmul #define cxipow cxpwr
allowing to use cxconvert
as synonym of cxconv
,
cxdiff
as synonym of cxsub
, cxprod
in place
of cxmul
and cxipow
for cxpwr
.
The HPA library supplies a C++ wrapper allowing to perform high precision computations with the same syntax of the normal code. Using the C++ wrapper allows you to do things like:
// Compute the factorial of the integer n xreal factorial(xreal n) { xreal i; xreal product = 1; for (i=2; i <= n; i++) product *= i; return product; }
Technically, the C++ wrapper is contained in a separate
library. However, this module is distributed together
with the HPA library as an extension to it.
Therefore, I decided to add a section to this manual
in order to describe the C++ interface of this
extension.
This one is formed by two classes, called xreal
and xcomplex
respectively.
The type xreal
can be used to declare or define
real variables. Moreover, it defines the mathematical
operations and functions which can be used with them.
The same is for the type xcomplex
with respect to
complex variables and to their manipulation.
Having variables of xreal
and xcomplex
type
you can use the same syntax as if you were just writing normal code,
but all the computations will be performed with the
high precision math library. Of course, this possibility
greatly simplifies the writing of code.
Note the difference between
struct xpr x, y; char buffer[256]; while ( (fgets (buffer, 256, stdin)) ) { x = atox (buffer); printf ("The square root of 2 * %s + 1 is ", buffer); y = xadd (xpr2(x, 1), xOne, 0); xprxpr (xsqrt (y), 30); putchar ('\n'); }
and its C++ version:
xreal x; while ( x.getfrom(cin) > 0 ) cout << "The square root of 2 * " << x << " + 1 is " << sqrt(2 * x + 1) << endl;
which is much more compact and easier to understand.
The use of the C++ wrapper requires however a recent
and ANSI-compliant C++ compiler (for instance g++ 3.x).
Moreover, before declaring or defining variables of xreal
type and before using anyone of the functions or operators
declared in the header file xreal.h
, you have to put the line
#include <xreal.h>
in your source code file.
Analogously, before declaring or defining variables of xcomplex
type and before using anyone of the functions or operators
declared in the header file xcomplex.h
, you have to add the line
#include <xcomplex.h>
to your source code file.
After that, it is recommendable to add the directive
using namespace HPA;
since all the objects of the C++ wrapper are declared
or defined within the namespace HPA.
Alternatively, you should always use the prefix
HPA::
in front of any identifiers coming
from the files xreal.h
or xcomplex.h
.
This is in order to tell the compiler where looking for
the classes xreal
, xcomplex
and all the related
stuff.
This code comes from a real world program:
#include<iostream> #include<xreal.h> using namespace HPA; int main (void) { xreal x; while ( x.getfrom(cin) > 0 ) cout << "The square root of 2 * " << x << " + 1 is " << sqrt(2 * x + 1) << endl; return 0; }
It is remarkable that the header files xreal.h
and
xcomplex.h
make also directly available the standard classes
ostream
, istream
and string
(see sections "The xreal class"
and "The xcomplex class").
When you have to compile and build a program making use
of the C++ wrapper to the HPA library, you can do it by following
the same instructions given in the section Compiling and linking,
by only taking care to use hpaxxconf
in place of
hpaconf
, just as in:
c++ -c $(hpaxxconf -c) example.cc
c++ -c `hpaxxconf -c` example.cc
to compile the file example.cc
and obtain the object file example.o
,
or in
c++ example.o $(hpaxxconf -l) -o example
c++ example.o `hpaxxconf -l` -o example
to do the linkage. If you want, you may also compile and build at the same time by using
c++ example.cc $(hpaxxconf -c -l) -o example
c++ example.cc `hpaxxconf -c -l` -o example
Naturally, all these samples assume that you are working
with bash or with another shell sh-compatible.
In any case, the synopsis of hpaxxconf
is the same as
for hpaconf
.
Warning:
The command hpaxxconf
is available only if
the C++ wrapper was built together with the
HPA library when this one was installed on
the system where you are working. Actually
it is possible to install the HPA library without
its C++ wrapper. In this case, if you try to recall
the command hpaxxconf
, the shell will print
the error message "Command not found"
or
something like that.
The interface of the xreal
class is
contained in the header file xreal.h
,
whose contents you can find here suitably commented:
#ifndef _XREAL_H_ #define _XREAL_H_ #include <xpre.h> #include <iostream> #include <string> using std::ostream; using std::istream; using std::string; namespace HPA { class xreal { // << and >> are used respectively for the output and the // input of extended precision numbers. // The input operator >> actually reads a double precision // number and then converts it to an extended precision // number. This can have undesirable rounding effects. // To avoid them, use the input function // xreal::getfrom() (see below). friend ostream& operator<< (ostream& os, const xreal& x); friend istream& operator>> (istream& is, xreal& x); // +, -, *, / are the usual arithmetic operators friend xreal operator+ (const xreal& x1, const xreal& x2); friend xreal operator- (const xreal& x1, const xreal& x2); friend xreal operator* (const xreal& x1, const xreal& x2); friend xreal operator/ (const xreal& x1, const xreal& x2); // x % n is equal to x * pow (2,n) friend xreal operator% (const xreal& x1, int n); // ==, !=, <=, >=, <, > are the usual comparison operators friend int operator== (const xreal& x1, const xreal& x2); friend int operator!= (const xreal& x1, const xreal& x2); friend int operator<= (const xreal& x1, const xreal& x2); friend int operator>= (const xreal& x1, const xreal& x2); friend int operator< (const xreal& x1, const xreal& x2); friend int operator> (const xreal& x1, const xreal& x2); // sget (s, n, x) tries to read an extended precision // number from the string 's' starting from the position // 'n'. The retrieved number is converted and stored in // 'x'. The return value is the number of characters // composing the decimal representation of this number // as read from 's'. For instance, if s == "12.34dog" and // n == 0, then 'x' is set to 12.34 and the return value // is 5. // If the portion of 's' starting from the position 'n' // can not be converted to a number, then 'x' is set to // xNAN and 0 is returned. // If the exactly converted value would cause overflow, // then xINF or x_INF is returned, according to the sign // of the value. // If 'n' is greater or equal to the length of 's', then 0 // is returned while 'x' is set to xZERO. friend unsigned long sget (string s, unsigned long startptr, xreal& x); // bget (buff, x) tries to read an extended precision // number from the buffer pointed to by 'buff'. // The retrieved number is converted and stored in 'x'. // The return value is a pointer to the character after // the last character used in the conversion. // For instance, if 'buff' is a pointer to the buffer // "12.34dog", then 'x' is set to 12.34 and the return // value is a pointer to "dog" (that is to say, it points // to the character 'd'). // If the initial portion of the string pointed to by 'buff' // can not be converted to a number, then 'x' is set to xNAN // and 'buff' is returned. // If the exactly converted value would cause overflow, // then xINF or x_INF is returned, according to the sign // of the value. // If 'buff' is NULL (0), then an error message is printed // on 'cerr' (standard error device). friend const char* bget (const char* buff, xreal& x); // compare (x1, x2) returns //+1 to mean x1 > x2 // 0 to mean x1 == x2 //-1 to mean x1 < x2 friend int compare (const xreal& x1, const xreal& x2); //isNaN (x) returns 1 when x == xNAN, else 0 friend int isNaN (const xreal& x); // The following functions do not need a particular comment: // each of them is defined as the corresponding function // of the standard math library, that is to say the function // from <cmath> having the same name. // However qfmod(), sfmod(), frac() and fix() do not have // counterparts in the standard math library. // With respect to fmod(), qfmod() requires one more // argument, where it stores the quotient of the division of // the first argument by the second one. // sfmod (x,&p) stores in the integer variable pointed to by // 'p' the integer part of 'x' and, at the same time, it // returns the fractional part of 'x'. // The usage of sfmod() is strongly discouraged. // frac() returns the fractional part of its argument. // At last, fix() is a frontend to the xfix() // function (see section "Extended Precision Floating // Point Arithmetic"). friend xreal abs (const xreal& s); friend xreal frexp (const xreal& s, int *p); friend xreal qfmod (const xreal& s, const xreal& t, xreal& q); friend xreal fmod (const xreal& s, const xreal& t); friend xreal sfmod (const xreal& s, int *p); friend xreal frac (const xreal& x); friend xreal trunc (const xreal& x); friend xreal round (const xreal& x); friend xreal ceil (const xreal& x); friend xreal floor (const xreal& x); friend xreal fix (const xreal& x); friend xreal tan (const xreal& x); friend xreal sin (const xreal& x); friend xreal cos (const xreal& x); friend xreal atan (const xreal& a); friend xreal asin (const xreal& a); friend xreal acos (const xreal& a); friend xreal sqrt (const xreal& u); friend xreal exp (const xreal& u); friend xreal exp2 (const xreal& u); friend xreal exp10 (const xreal& u); friend xreal log (const xreal& u); friend xreal log2 (const xreal& u); friend xreal log10 (const xreal& u); friend xreal tanh (const xreal& v); friend xreal sinh (const xreal& v); friend xreal cosh (const xreal& v); friend xreal atanh (const xreal& v); friend xreal asinh (const xreal& v); friend xreal acosh (const xreal& v); friend xreal pow (const xreal& x, const xreal& y); public: // Various constructors. They allow to define // an extended precision number in several ways. // Moreover, they allow for conversions from other // numeric types. xreal (const struct xpr* px = &xZero); xreal (struct xpr x); xreal (double x); xreal (float x); xreal (int n); xreal (long n); xreal (unsigned int u); xreal (unsigned long u); // This constructor requires a special comment. When // only the first argument is present, the initial portion // of the string pointed to by this argument is converted // into an extended precision number, if a conversion is // possible. If no conversion is possible, then the // returned number is xNAN. When the second argument is // present and is not null, it must be the address of a // valid pointer to 'char'. // Before returning, the constructor will set this pointer // so that it points to the character after the last // character used in the conversion. xreal (const char* str, char** endptr = 0); xreal (string str); xreal (const xreal& x); // Assignment operators. They do not require // any explanation with the only exception of '%=', // which combines a '%' operation with an assignment. // So, x %= n is equivalent to x *= pow(2,n) . xreal& operator= (const xreal& x); xreal& operator+= (const xreal& x); xreal& operator-= (const xreal& x); xreal& operator*= (const xreal& x); xreal& operator/= (const xreal& x); xreal& operator%= (int n); // Increment and decrement operators. Both prefixed // and postfixed versions are defined. xreal& operator++ (); xreal& operator-- (); xreal& operator++ (int dummy); xreal& operator-- (int dummy); // Destructor. You will never have to recall it // explicitly in your code. ~xreal (void); // Integer exponent power. For any extended precision // number 'x', x(n) is equal to 'x' raised to the power // of 'n'. xreal operator() (int n) const; // This is the usual unary minus. xreal operator-() const; // For any extended precision number 'x', !x evaluates to 1 // when 'x' is null, else it evaluates to 0. int operator!() const; // x.isneg() returns 1 if 'x' is negative, else it // returns 0. int isneg() const; // x.exp() returns the exponent part // of the binary representation of 'x'. int exp() const; // Functions for conversion. x._2double(), x._2float(), // x._2xpr() and x._2string() convert the extended precision // number 'x' in a double precision number, in a single // precision number, in a structure of // type 'xpr' and in a string respectively. double _2double () const; float _2float() const; struct xpr _2xpr() const; string _2string() const; // The member function xreal::getfrom() can be used to // recover an extended precision number from an input // stream. The input stream is passed as argument to the // function. // The return value is 0 in case of input error (in case // of End-Of-File for instance). // When it starts to process its input, this function drops // all the eventual leading white spaces. // After reading the first non space character, it // continues to read from the input stream until it finds // a white space or reaches the End-Of-File. // Then it tries to convert into an extended // precision number the (initial portion of the) string just // read. // If no conversion can be performed, then x.getfrom(is) // sets 'x' to the value xNAN. // If the exactly converted value would cause overflow, // then 'x' is set to xINF or x_INF, according to the sign // of the correct value. int getfrom (istream& is); // The member function xreal::print() can be used to print // an extended precision number onto an output stream. // The output stream is passed to the function as first // argument. The next three arguments have the same meanings // of the fields 'notat', 'sf' and 'lim' of the // structure 'xoutflags' (see section "Real Arithmetic"). int print (ostream& os, int sc_not, int sign, int lim) const; // The function call x.asprint(sc_not, sign, lim) returns // a buffer of characters with the representation, in form // of a decimal ASCII string, of the extended precision // number 'x'. The arguments 'sc_not', 'sign' and 'lim' are // used to format the string. // They have the same meanings of the fields 'notat', 'sf' // and 'lim' of the structure 'xoutflags' (see section // "Real Arithmetic"). // The buffer returned by this function is malloc'ed inside // the function. In case of insufficient memory, the null // pointer is returned. char* asprint (int sc_not, int sign, int lim) const; // The following static functions are used to set // or get the values of the fields of the structure // 'xreal::ioflags'. This structure is a static member // variable of the class 'xreal' and it is used by the // output operator << to know how format its second // argument. The meaning of the fields of the structure // 'xreal::ioflags' is explained in the section // "Real arithmetic". // xreal::set_notation (which) sets to 'which' the value // of 'xreal::ioflags.notat' . static void set_notation (short notat); // xreal::set_signflag (which) sets to 'which' the value // of 'xreal::ioflags.sf' . static void set_signflag (short onoff); // xreal::set_mfwd (which) sets to 'which' the value // of 'xreal::ioflags.mfwd' . static void set_mfwd (short wd); // xreal::set_lim (which) sets to 'which' the value // of 'xreal::ioflags.lim' . static void set_lim (short lim); // xreal::set_padding (which) sets to 'which' the value // of 'xreal::ioflags.padding' . static void set_padding (signed char ch); // xreal::get_notation () returns the current value // of 'xreal::ioflags.notat' . static short get_notation (void); // xreal::get_signflag () returns the current value // of 'xreal::ioflags.sf' . static short get_signflag (void); // xreal::get_mfwd () returns the current value // of 'xreal::ioflags.mfwd' . static short get_mfwd (void); // xreal::get_lim () returns the current value // of 'xreal::ioflags.lim' . static short get_lim (void); // xreal::get_padding () returns the current value // of 'xreal::ioflags.padding' . static signed char get_padding (void); private: struct xpr br; /* binary representation */ static struct xoutflags ioflags; /* output flags */ }; // xmatherrcode() returns the current value of the global // variable 'xErrNo' (see section // "Dealing with runtime errors") if // this variable is defined. // Otherwise xmatherrcode() returns -1. int xmatherrcode (); // clear_xmatherr() resets to 0 the value of the global // variable 'xErrNo' (see section "Dealing with runtime // errors") if this variable is defined. // Otherwise, clear_xmatherr() prints a suitable warning // on 'cerr' (standard error device). void clear_xmatherr (); // Some useful constants: // xZERO == 0 // xONE == 1 // xTWO == 2 // xTEN == 10 // xINF == +INF // x_INF == -INF // xNAN == Not-A-Number // xPI == Pi Greek // xPI2 == Pi / 2 // xPI4 == Pi / 4 // xEE == e (base of natural logarithms) // xSQRT2 == square root of 2 // xLN2 == natural logarithm of 2 // xLN10 == natural logarithm of 10 // xLOG2_E == base-2 logarithm of e // xLOG2_10 == base-2 logarithm of 10 // xLOG10_E == base-10 logarithm of e extern const xreal xZERO, xONE, xTWO, xTEN; extern const xreal xINF, x_INF, xNAN; extern const xreal xPI, xPI2, xPI4, xEE, xSQRT2; extern const xreal xLN2, xLN10, xLOG2_E, xLOG2_10, xLOG10_E; } /* End namespace HPA */ #endif /* _XREAL_H_ */
The interface of the xcomplex
class is
contained in the header file xcomplex.h
,
whose contents you can find here suitably commented:
#ifndef _XCOMPLEX_H_ #define _XCOMPLEX_H_ #include <cxpre.h> #include <iostream> #include <string> #include "xreal.h" using std::istream; using std::ostream; using std::string; namespace HPA { struct double_complex { double re, im; }; struct float_complex { float re, im; }; class xcomplex { // << and >> are used respectively for the output and the // input of extended precision complex numbers. // The input operator >> actually reads a couple of double // precision numbers and then converts it into // an extended precision complex number. This can have // undesirable rounding effects. To avoid them, use the // input function xcomplex::getfrom() (see below). friend ostream& operator<< (ostream& os, const xcomplex& x); friend istream& operator>> (istream& is, xcomplex& x); // +, -, *, / are the usual arithmetic operators friend xcomplex operator+ (const xcomplex& x1, const xcomplex& x2); friend xcomplex operator- (const xcomplex& x1, const xcomplex& x2); friend xcomplex operator* (const xcomplex& x1, const xcomplex& x2); friend xcomplex operator/ (const xcomplex& x1, const xcomplex& x2); // z % n is equal to z * pow (2,n) friend xcomplex operator% (const xcomplex& x1, int n); // ==, !=, <=, >=, <, > are the usual comparison operators friend int operator== (const xcomplex& x1, const xcomplex& x2); friend int operator!= (const xcomplex& x1, const xcomplex& x2); friend int operator<= (const xcomplex& x1, const xcomplex& x2); friend int operator>= (const xcomplex& x1, const xcomplex& x2); friend int operator< (const xcomplex& x1, const xcomplex& x2); friend int operator> (const xcomplex& x1, const xcomplex& x2); // sget (s, n, x) tries to read an extended precision // complex number from the string 's' starting from the // position 'n'. The retrieved number is converted and // stored in 'x'. The return value is the number of // characters composing the decimal representation of this // number as read from 's'. // For instance, if s == "12.34+6.7idog" and n == 0, // then 'x' is set to 12.34+6.7i and the return value // is 10. // If the portion of 's' starting from the position 'n' can // not be converted to a number, then 'x' is set to // xNAN + xNANi and 0 is returned. // If the exactly converted value would cause overflow in // the real or/and imaginary part, then the real or/and the // imaginary part of 'x' is set to xINF or x_INF, according // to the signs of the correct value. // If 'n' is greater or equal to the length of 's', then 0 // is returned while 'x' is set to cxZERO. friend unsigned long sget (string s, unsigned long startptr, xcomplex& x); // bget (buff, x) tries to read an extended precision // complex number from the buffer pointed to by 'buff'. // The retrieved number is converted and stored in 'x'. // The return value is a pointer to the character after // the last character used in the conversion. // For instance, if 'buff' is a pointer to the buffer // "12.34+6.7idog", then 'x' is set to 12.34+6.7i and // the return value is a pointer to "dog" (that is to say, // it points to the character 'd'). // If the initial portion of the string pointed to by 'buff' // can not be converted to a number, then 'x' is set to // xNAN + xNANi and 'buff' is returned. // If the exactly converted value would cause overflow // in the real or/and imaginary part, then the real or/and // the imaginary part of 'x' is set to xINF or x_INF, // according to the signs of the correct value. // If 'buff' is NULL (0), then an error message is printed // on 'cerr' (standard error device). friend const char* bget (const char* buff, xcomplex& x); // rmul (x,z) (here 'x' is a real number) returns the // product x * z. // It is faster than the * operator. friend xcomplex rmul (const xreal& x, const xcomplex& z); // After eventually rounding 'z1' and 'z2' by recalling // round() (see below) on them, gdiv(z1, z2) returns // the quotient of the gaussian division of 'z1' by 'z2'. // If you do not know what gaussian division means, probably // you will never need this function :) friend xcomplex gdiv (const xcomplex& x1, const xcomplex& x2); // After eventually rounding 'z1' and 'z2' by recalling // round() (see below) on them, gmod(z1, z2) returns // the remainder of the gaussian division of 'z1' by 'z2'. // If you do not know what gaussian division means, probably // you will never need this function :) friend xcomplex gmod (const xcomplex& x1, const xcomplex& x2); // idiv() is a wrapper to cxidiv() (see section // "Extended Precision Complex Arithmetic"). friend xcomplex idiv (const xcomplex& x1, const xcomplex& x2); // mod() is a wrapper to cxmod() (see section // "Extended Precision Complex Arithmetic"). friend xcomplex mod (const xcomplex& x1, const xcomplex& x2); // conj() returns the complex conjugate of its argument. friend xcomplex conj (const xcomplex& z); // inv() returns the complex reciprocal of its argument. // So inv(z) == 1/z . friend xcomplex inv (const xcomplex& z); // swap(z) returns the complex number {z.im, z.re}. friend xcomplex swap (const xcomplex& z); // Multiplication by 1i (imaginary unit). friend xcomplex drot (const xcomplex& z); // Multiplication by -1i friend xcomplex rrot (const xcomplex& z); // abs() returns the absolute value (or modulus) of its // argument. // The return value of abs() is then an 'xreal' number. friend xreal abs (const xcomplex& s); // arg(z) returns the phase angle (or argument) // of the complex number 'z'. // The return value of arg() is an 'xreal' number // in the range [-xPI, xPI). // If 'z' is null, then a domain-error is produced. friend xreal arg (const xcomplex& s); // The next six functions have the same // meanings of the corresponding real functions, // but they act upon both the real // and the imaginary part of their argument. friend xcomplex frac (const xcomplex& x); friend xcomplex trunc (const xcomplex& x); friend xcomplex round (const xcomplex& x); friend xcomplex ceil (const xcomplex& x); friend xcomplex floor (const xcomplex& x); friend xcomplex fix (const xcomplex& x); // sqr() returns the square of its argument. friend xcomplex sqr (const xcomplex& u); // sqrt() returns the principal branch of the square root // of its argument. friend xcomplex sqrt (const xcomplex& u); // root (u,i,n) returns the 'i'th branch of the 'n'th root // of 'u'. If 'n' is zero or 'n' is negative and 'u' is // null, then a bad-exponent error is produced. friend xcomplex root (const xcomplex& u, int i, int n); // These functions do not require any comment, except that // tan() and tanh() yield a domain-error in the same cases // as cxtan() and cxtanh() respectively. friend xcomplex exp (const xcomplex& u); friend xcomplex exp2 (const xcomplex& u); friend xcomplex exp10 (const xcomplex& u); friend xcomplex tan (const xcomplex& x); friend xcomplex sin (const xcomplex& x); friend xcomplex cos (const xcomplex& x); friend xcomplex tanh (const xcomplex& v); friend xcomplex sinh (const xcomplex& v); friend xcomplex cosh (const xcomplex& v); // Natural, base-2 and base-10 logarithm of a complex // number. // A null argument results in a domain-error. // The imaginary part of the return value of log() is always // in the interval [-xPI,xPI). friend xcomplex log (const xcomplex& u); friend xcomplex log2 (const xcomplex& u); friend xcomplex log10 (const xcomplex& u); // log_sqrt(u) returns the natural logarithm of the // principal branch of the square root of 'u'. // A null argument results in a domain-error. // The imaginary part of the return value of log_sqrt() // is always in the interval [-xPI2,xPI2). friend xcomplex log_sqrt (const xcomplex& u); // These functions are self-explanatory. atan(z) // yields a domain-error if the real part of 'z' is null and // the imaginary part is equal to '+1' or '-1'. // Analogously, atanh(z) yields a domain-error if the // imaginary part of 'z' is null and the // real part is equal to '+1' or '-1'. friend xcomplex atan (const xcomplex& a); friend xcomplex asin (const xcomplex& a); friend xcomplex acos (const xcomplex& a); friend xcomplex atanh (const xcomplex& v); friend xcomplex asinh (const xcomplex& v); friend xcomplex acosh (const xcomplex& v); // The return value of pow() is the power of the first // argument raised to the second one. // Note that the modulus of the first argument must be // positive, i.e. greater than zero, if the real part of // the second argument is less or equal than zero, otherwise // a bad-exponent error will be produced. friend xcomplex pow (const xcomplex& x, const xcomplex& y); public: // Various constructors. They allow to define // an extended precision complex number in several ways. // Moreover, they allow for conversions from other // numeric types. xcomplex (const struct cxpr* px = &cxZero); xcomplex (struct cxpr x); xcomplex (struct xpr x, struct xpr y = xZero); xcomplex (xreal x, xreal y = xZERO); xcomplex (double x, double y = 0.0); xcomplex (float x, float y = 0.0); xcomplex (int m, int n = 0); xcomplex (long m, long n = 0); xcomplex (unsigned int u, unsigned int v = 0U); xcomplex (unsigned long u, unsigned long v = 0U); // This constructor requires a special comment. When // only the first argument is present, the initial portion // of the string pointed to by this argument is converted // into an extended precision complex number, if a // conversion is possible. If no conversion is possible, // then the returned number is xNAN + xNANi. // When the second argument is present and is not null, // it must be the address of a valid pointer to 'char'. // Before returning, the constructor will set this pointer // so that it points to the character after the last // character used in the conversion. xcomplex (const char* str, char** endptr = 0); xcomplex (string str); xcomplex (const xcomplex& x); // Assignment operators. They do not require // any explanation with the only exception of '%=', // which combines a '%' operation with an assignment. // So, x %= n is equivalent to x *= pow(2,n) . xcomplex& operator= (const xcomplex& x); xcomplex& operator+= (const xcomplex& x); xcomplex& operator-= (const xcomplex& x); xcomplex& operator*= (const xcomplex& x); xcomplex& operator*= (const xreal& x); xcomplex& operator/= (const xcomplex& x); xcomplex& operator%= (int n); // Increment and decrement operators. Both prefixed // and postfixed versions are defined. These operators // only act on the real part of their argument. xcomplex& operator++ (); xcomplex& operator-- (); xcomplex& operator++ (int dummy); xcomplex& operator-- (int dummy); // Destructor. You will never have to recall it // explicitly in your code. ~xcomplex (void); // Integer exponent power. For any extended precision // complex number 'x', x(n) is equal to 'x' raised to // the power of 'n'. xcomplex operator() (int n) const; // This is the usual unary minus. xcomplex operator-() const; // For any extended precision complex number 'x', // !x evaluates to 1 when // 'x' is null, else it evaluates to 0. int operator!() const; // Functions for conversion. x._2dcomplex(), x._2fcomplex(), // x._2cxpr() and x._2string() convert the extended // precision complex number 'x' in a double precision // complex number, in a single precision complex // number, in a structure of type 'cxpr' and in a // string respectively. double_complex _2dcomplex () const; float_complex _2fcomplex() const; struct cxpr _2cxpr() const; string _2string() const; // For any extended precision complex number 'z', // z.real() and z.imag() return, respectively, the // real and the imaginary part of 'z' in the form // of an extended precision number. xreal real () const; xreal imag () const; // For any extended precision complex number 'z', // z._real() and z._imag() return, respectively, the // real and the imaginary part of 'z' in the form // of a structure of 'xpr' type. struct xpr _real () const; struct xpr _imag () const; // For any extended precision complex number 'z', // z.dreal() and z.dimag() return, respectively, the // real and the imaginary part of 'z' in the form // of a double precision number. double dreal () const; double dimag () const; // For any extended precision complex number 'z', // z.freal() and z.fimag() return, respectively, the // real and the imaginary part of 'z' in the form // of a single precision number. double freal () const; double fimag () const; // For any extended precision complex number 'z', // z.sreal() and z.simag() return, respectively, the // real and the imaginary part of 'z' in the form // of a string. string sreal () const; string simag () const; // The next functions allow to set (or reset) // the real and the imaginary part of a complex number. void real (const xreal& x); void imag (const xreal& x); void real (struct xpr x); void imag (struct xpr x); void real (const struct xpr* px); void imag (const struct xpr* px); void real (double x); void imag (double x); void real (float x); void imag (float x); void real (int x); void imag (int x); void real (long x); void imag (long x); void real (unsigned int x); void imag (unsigned int x); void real (unsigned long x); void imag (unsigned long x); void real (const char* str, char** endptr = 0); void imag (const char* str, char** endptr = 0); void real (string str); void imag (string str); // The member function xcomplex::getfrom() can be used // to recover an extended precision complex number from an // input stream. The input stream is passed as argument to // the function. // The return value is 0 in case of input error (in case of // End-Of-File for instance). // When it starts to process its input, this function drops // all the eventual leading white spaces. // After reading the first non space character, it continues // to read from the input stream until it finds a white // space or reaches the End-Of-File. // Then it tries to convert into an extended // precision complex number the (initial portion of the) // string just read. // If no conversion can be performed, then x.getfrom(is) // sets 'x' to the value xNAN + xNANi. // If the exactly converted value would cause overflow in // the real or/and in the imaginary part, then the real part // or/and the imaginary part of 'x' is set to xINF or x_INF, // according to the signs of the correct value. int getfrom (istream& is); // The member function xcomplex::print() can be used to // print an extended precision complex number onto an output // stream. The output stream is passed to the function as // first argument. The next three arguments have the same // meanings of the fields 'notat', 'sf' // and 'lim' of the structure 'xoutflags' (see section // "Real Arithmetic"). int print (ostream& os, int sc_not, int sign, int lim) const; // The function call x.asprint(sc_not, sign, lim) returns // a buffer of characters with the representation, // in form of a decimal ASCII string, // of the extended precision complex number 'x'. // The arguments 'sc_not', 'sign' and 'lim' are used // to format the string. // They have the same meanings of the fields 'notat', 'sf' // and 'lim' of the structure 'xoutflags' (see section // "Real Arithmetic"). // The buffer returned by this function is malloc'ed inside // the function. In case of insufficient memory, the null // pointer is returned. char* asprint (int sc_not, int sign, int lim) const; // The following static functions are used to set // or get the values of the fields of the structure // 'xcomplex::ioflags'. This structure is a static member // variable of the class 'xcomplex' and it is used by // the output operator << to know how format its second // argument. The meaning of the // fields of the structure 'xcomplex::ioflags' is explained // in the section "Real arithmetic". // xcomplex::set_fmt (which) sets to 'which' the value // of 'xcomplex::ioflags.fmt' . static void set_fmt (short format); // xcomplex::set_notation (which) sets to 'which' the value // of 'xcomplex::ioflags.notat' . static void set_notation (short notat); // xcomplex::set_signflag (which) sets to 'which' the value // of 'xcomplex::ioflags.sf' . static void set_signflag (short onoff); // xcomplex::set_mfwd (which) sets to 'which' the value // of 'xcomplex::ioflags.mfwd' . static void set_mfwd (short wd); // xcomplex::set_lim (which) sets to 'which' the value // of 'xcomplex::ioflags.lim' . static void set_lim (short lim); // xcomplex::set_padding (which) sets to 'which' the value // of 'xcomplex::ioflags.padding' . static void set_padding (signed char ch); // xcomplex::set_ldelim (ch) sets to 'ch' the value // of 'xcomplex::ioflags.ldel' . static void set_ldelim (signed char ch); // xcomplex::set_rdelim (ch) sets to 'ch' the value // of 'xcomplex::ioflags.rdel' . static void set_rdelim (signed char ch); // xcomplex::get_fmt () returns the current value // of 'xcomplex::ioflags.fmt' . static short get_fmt (void); // xcomplex::get_notation () returns the current value // of 'xcomplex::ioflags.notat' . static short get_notation (void); // xcomplex::get_signflag () returns the current value // of 'xcomplex::ioflags.sf' . static short get_signflag (void); // xcomplex::get_mfwd () returns the current value // of 'xcomplex::ioflags.mfwd' . static short get_mfwd (void); // xcomplex::get_lim () returns the current value // of 'xcomplex::ioflags.lim' . static short get_lim (void); // xcomplex::get_padding () returns the current value // of 'xcomplex::ioflags.padding' . static signed char get_padding (void); // xcomplex::get_ldelim () returns the current value // of 'xcomplex::ioflags.ldel' . static signed char get_ldelim (void); // xcomplex::get_rdelim () returns the current value // of 'xcomplex::ioflags.rdel' . static signed char get_rdelim (void); private: struct cxpr br; /* binary representation */ static struct xoutflags ioflags; /* output flags */ }; // Some useful constants: // cxZERO == 0 // cxONE == 1 // cxI == 1i extern const xcomplex cxZERO, cxONE, cxI; } /* End namespace HPA */ #define xi cxI #define xj cxI #define _i cxI #define _j cxI #endif /* _XCOMPLEX_H_ */
Firstly I feel like to thank Daniel A. Atkinson, since without its work the HPA library would not exist. Next I have to thank Aurelio Marinho Jargas (<verde (a) aurelio net>), author of txt2tags (http://txt2tags.sf.net), a wonderful text formatting and conversion tool, which I extensively used in writing this document. Last but not least, I want to thank all the people till now involved in the Free Software community, starting from those ones directly involved in the GNU project (http://www.gnu.org). Without their great work, this little one would never be done.
GNU Free Documentation License Version 1.2, November 2002 Copyright (C) 2000,2001,2002 Free Software Foundation, Inc. 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA Everyone is permitted to copy and distribute verbatim copies of this license document, but changing it is not allowed. 0. PREAMBLE The purpose of this License is to make a manual, textbook, or other functional and useful document "free" in the sense of freedom: to assure everyone the effective freedom to copy and redistribute it, with or without modifying it, either commercially or noncommercially. Secondarily, this License preserves for the author and publisher a way to get credit for their work, while not being considered responsible for modifications made by others. This License is a kind of "copyleft", which means that derivative works of the document must themselves be free in the same sense. It complements the GNU General Public License, which is a copyleft license designed for free software. We have designed this License in order to use it for manuals for free software, because free software needs free documentation: a free program should come with manuals providing the same freedoms that the software does. But this License is not limited to software manuals; it can be used for any textual work, regardless of subject matter or whether it is published as a printed book. We recommend this License principally for works whose purpose is instruction or reference. 1. APPLICABILITY AND DEFINITIONS This License applies to any manual or other work, in any medium, that contains a notice placed by the copyright holder saying it can be distributed under the terms of this License. Such a notice grants a world-wide, royalty-free license, unlimited in duration, to use that work under the conditions stated herein. The "Document", below, refers to any such manual or work. Any member of the public is a licensee, and is addressed as "you". You accept the license if you copy, modify or distribute the work in a way requiring permission under copyright law. A "Modified Version" of the Document means any work containing the Document or a portion of it, either copied verbatim, or with modifications and/or translated into another language. A "Secondary Section" is a named appendix or a front-matter section of the Document that deals exclusively with the relationship of the publishers or authors of the Document to the Document's overall subject (or to related matters) and contains nothing that could fall directly within that overall subject. (Thus, if the Document is in part a textbook of mathematics, a Secondary Section may not explain any mathematics.) The relationship could be a matter of historical connection with the subject or with related matters, or of legal, commercial, philosophical, ethical or political position regarding them. The "Invariant Sections" are certain Secondary Sections whose titles are designated, as being those of Invariant Sections, in the notice that says that the Document is released under this License. If a section does not fit the above definition of Secondary then it is not allowed to be designated as Invariant. The Document may contain zero Invariant Sections. If the Document does not identify any Invariant Sections then there are none. The "Cover Texts" are certain short passages of text that are listed, as Front-Cover Texts or Back-Cover Texts, in the notice that says that the Document is released under this License. A Front-Cover Text may be at most 5 words, and a Back-Cover Text may be at most 25 words. A "Transparent" copy of the Document means a machine-readable copy, represented in a format whose specification is available to the general public, that is suitable for revising the document straightforwardly with generic text editors or (for images composed of pixels) generic paint programs or (for drawings) some widely available drawing editor, and that is suitable for input to text formatters or for automatic translation to a variety of formats suitable for input to text formatters. A copy made in an otherwise Transparent file format whose markup, or absence of markup, has been arranged to thwart or discourage subsequent modification by readers is not Transparent. An image format is not Transparent if used for any substantial amount of text. A copy that is not "Transparent" is called "Opaque". Examples of suitable formats for Transparent copies include plain ASCII without markup, Texinfo input format, LaTeX input format, SGML or XML using a publicly available DTD, and standard-conforming simple HTML, PostScript or PDF designed for human modification. Examples of transparent image formats include PNG, XCF and JPG. Opaque formats include proprietary formats that can be read and edited only by proprietary word processors, SGML or XML for which the DTD and/or processing tools are not generally available, and the machine-generated HTML, PostScript or PDF produced by some word processors for output purposes only. The "Title Page" means, for a printed book, the title page itself, plus such following pages as are needed to hold, legibly, the material this License requires to appear in the title page. For works in formats which do not have any title page as such, "Title Page" means the text near the most prominent appearance of the work's title, preceding the beginning of the body of the text. A section "Entitled XYZ" means a named subunit of the Document whose title either is precisely XYZ or contains XYZ in parentheses following text that translates XYZ in another language. (Here XYZ stands for a specific section name mentioned below, such as "Acknowledgements", "Dedications", "Endorsements", or "History".) To "Preserve the Title" of such a section when you modify the Document means that it remains a section "Entitled XYZ" according to this definition. The Document may include Warranty Disclaimers next to the notice which states that this License applies to the Document. These Warranty Disclaimers are considered to be included by reference in this License, but only as regards disclaiming warranties: any other implication that these Warranty Disclaimers may have is void and has no effect on the meaning of this License. 2. VERBATIM COPYING You may copy and distribute the Document in any medium, either commercially or noncommercially, provided that this License, the copyright notices, and the license notice saying this License applies to the Document are reproduced in all copies, and that you add no other conditions whatsoever to those of this License. You may not use technical measures to obstruct or control the reading or further copying of the copies you make or distribute. However, you may accept compensation in exchange for copies. If you distribute a large enough number of copies you must also follow the conditions in section 3. You may also lend copies, under the same conditions stated above, and you may publicly display copies. 3. COPYING IN QUANTITY If you publish printed copies (or copies in media that commonly have printed covers) of the Document, numbering more than 100, and the Document's license notice requires Cover Texts, you must enclose the copies in covers that carry, clearly and legibly, all these Cover Texts: Front-Cover Texts on the front cover, and Back-Cover Texts on the back cover. Both covers must also clearly and legibly identify you as the publisher of these copies. The front cover must present the full title with all words of the title equally prominent and visible. You may add other material on the covers in addition. Copying with changes limited to the covers, as long as they preserve the title of the Document and satisfy these conditions, can be treated as verbatim copying in other respects. If the required texts for either cover are too voluminous to fit legibly, you should put the first ones listed (as many as fit reasonably) on the actual cover, and continue the rest onto adjacent pages. If you publish or distribute Opaque copies of the Document numbering more than 100, you must either include a machine-readable Transparent copy along with each Opaque copy, or state in or with each Opaque copy a computer-network location from which the general network-using public has access to download using public-standard network protocols a complete Transparent copy of the Document, free of added material. If you use the latter option, you must take reasonably prudent steps, when you begin distribution of Opaque copies in quantity, to ensure that this Transparent copy will remain thus accessible at the stated location until at least one year after the last time you distribute an Opaque copy (directly or through your agents or retailers) of that edition to the public. It is requested, but not required, that you contact the authors of the Document well before redistributing any large number of copies, to give them a chance to provide you with an updated version of the Document. 4. MODIFICATIONS You may copy and distribute a Modified Version of the Document under the conditions of sections 2 and 3 above, provided that you release the Modified Version under precisely this License, with the Modified Version filling the role of the Document, thus licensing distribution and modification of the Modified Version to whoever possesses a copy of it. In addition, you must do these things in the Modified Version: A. Use in the Title Page (and on the covers, if any) a title distinct from that of the Document, and from those of previous versions (which should, if there were any, be listed in the History section of the Document). You may use the same title as a previous version if the original publisher of that version gives permission. B. List on the Title Page, as authors, one or more persons or entities responsible for authorship of the modifications in the Modified Version, together with at least five of the principal authors of the Document (all of its principal authors, if it has fewer than five), unless they release you from this requirement. C. State on the Title page the name of the publisher of the Modified Version, as the publisher. D. Preserve all the copyright notices of the Document. E. Add an appropriate copyright notice for your modifications adjacent to the other copyright notices. F. Include, immediately after the copyright notices, a license notice giving the public permission to use the Modified Version under the terms of this License, in the form shown in the Addendum below. G. Preserve in that license notice the full lists of Invariant Sections and required Cover Texts given in the Document's license notice. H. Include an unaltered copy of this License. I. Preserve the section Entitled "History", Preserve its Title, and add to it an item stating at least the title, year, new authors, and publisher of the Modified Version as given on the Title Page. If there is no section Entitled "History" in the Document, create one stating the title, year, authors, and publisher of the Document as given on its Title Page, then add an item describing the Modified Version as stated in the previous sentence. J. Preserve the network location, if any, given in the Document for public access to a Transparent copy of the Document, and likewise the network locations given in the Document for previous versions it was based on. These may be placed in the "History" section. You may omit a network location for a work that was published at least four years before the Document itself, or if the original publisher of the version it refers to gives permission. K. For any section Entitled "Acknowledgements" or "Dedications", Preserve the Title of the section, and preserve in the section all the substance and tone of each of the contributor acknowledgements and/or dedications given therein. L. Preserve all the Invariant Sections of the Document, unaltered in their text and in their titles. Section numbers or the equivalent are not considered part of the section titles. M. Delete any section Entitled "Endorsements". Such a section may not be included in the Modified Version. N. Do not retitle any existing section to be Entitled "Endorsements" or to conflict in title with any Invariant Section. O. Preserve any Warranty Disclaimers. If the Modified Version includes new front-matter sections or appendices that qualify as Secondary Sections and contain no material copied from the Document, you may at your option designate some or all of these sections as invariant. To do this, add their titles to the list of Invariant Sections in the Modified Version's license notice. These titles must be distinct from any other section titles. You may add a section Entitled "Endorsements", provided it contains nothing but endorsements of your Modified Version by various parties--for example, statements of peer review or that the text has been approved by an organization as the authoritative definition of a standard. You may add a passage of up to five words as a Front-Cover Text, and a passage of up to 25 words as a Back-Cover Text, to the end of the list of Cover Texts in the Modified Version. Only one passage of Front-Cover Text and one of Back-Cover Text may be added by (or through arrangements made by) any one entity. If the Document already includes a cover text for the same cover, previously added by you or by arrangement made by the same entity you are acting on behalf of, you may not add another; but you may replace the old one, on explicit permission from the previous publisher that added the old one. The author(s) and publisher(s) of the Document do not by this License give permission to use their names for publicity for or to assert or imply endorsement of any Modified Version. 5. COMBINING DOCUMENTS You may combine the Document with other documents released under this License, under the terms defined in section 4 above for modified versions, provided that you include in the combination all of the Invariant Sections of all of the original documents, unmodified, and list them all as Invariant Sections of your combined work in its license notice, and that you preserve all their Warranty Disclaimers. The combined work need only contain one copy of this License, and multiple identical Invariant Sections may be replaced with a single copy. If there are multiple Invariant Sections with the same name but different contents, make the title of each such section unique by adding at the end of it, in parentheses, the name of the original author or publisher of that section if known, or else a unique number. Make the same adjustment to the section titles in the list of Invariant Sections in the license notice of the combined work. In the combination, you must combine any sections Entitled "History" in the various original documents, forming one section Entitled "History"; likewise combine any sections Entitled "Acknowledgements", and any sections Entitled "Dedications". You must delete all sections Entitled "Endorsements". 6. COLLECTIONS OF DOCUMENTS You may make a collection consisting of the Document and other documents released under this License, and replace the individual copies of this License in the various documents with a single copy that is included in the collection, provided that you follow the rules of this License for verbatim copying of each of the documents in all other respects. You may extract a single document from such a collection, and distribute it individually under this License, provided you insert a copy of this License into the extracted document, and follow this License in all other respects regarding verbatim copying of that document. 7. AGGREGATION WITH INDEPENDENT WORKS A compilation of the Document or its derivatives with other separate and independent documents or works, in or on a volume of a storage or distribution medium, is called an "aggregate" if the copyright resulting from the compilation is not used to limit the legal rights of the compilation's users beyond what the individual works permit. When the Document is included in an aggregate, this License does not apply to the other works in the aggregate which are not themselves derivative works of the Document. If the Cover Text requirement of section 3 is applicable to these copies of the Document, then if the Document is less than one half of the entire aggregate, the Document's Cover Texts may be placed on covers that bracket the Document within the aggregate, or the electronic equivalent of covers if the Document is in electronic form. Otherwise they must appear on printed covers that bracket the whole aggregate. 8. TRANSLATION Translation is considered a kind of modification, so you may distribute translations of the Document under the terms of section 4. Replacing Invariant Sections with translations requires special permission from their copyright holders, but you may include translations of some or all Invariant Sections in addition to the original versions of these Invariant Sections. You may include a translation of this License, and all the license notices in the Document, and any Warranty Disclaimers, provided that you also include the original English version of this License and the original versions of those notices and disclaimers. In case of a disagreement between the translation and the original version of this License or a notice or disclaimer, the original version will prevail. If a section in the Document is Entitled "Acknowledgements", "Dedications", or "History", the requirement (section 4) to Preserve its Title (section 1) will typically require changing the actual title. 9. TERMINATION You may not copy, modify, sublicense, or distribute the Document except as expressly provided for under this License. Any other attempt to copy, modify, sublicense or distribute the Document is void, and will automatically terminate your rights under this License. However, parties who have received copies, or rights, from you under this License will not have their licenses terminated so long as such parties remain in full compliance. 10. FUTURE REVISIONS OF THIS LICENSE The Free Software Foundation may publish new, revised versions of the GNU Free Documentation License from time to time. Such new versions will be similar in spirit to the present version, but may differ in detail to address new problems or concerns. See http://www.gnu.org/copyleft/. Each version of the License is given a distinguishing version number. If the Document specifies that a particular numbered version of this License "or any later version" applies to it, you have the option of following the terms and conditions either of that specified version or of any later version that has been published (not as a draft) by the Free Software Foundation. If the Document does not specify a version number of this License, you may choose any version ever published (not as a draft) by the Free Software Foundation. ADDENDUM: How to use this License for your documents To use this License in a document you have written, include a copy of the License in the document and put the following copyright and license notices just after the title page: Copyright (c) YEAR YOUR NAME. Permission is granted to copy, distribute and/or modify this document under the terms of the GNU Free Documentation License, Version 1.2 or any later version published by the Free Software Foundation; with no Invariant Sections, no Front-Cover Texts, and no Back-Cover Texts. A copy of the license is included in the section entitled "GNU Free Documentation License". If you have Invariant Sections, Front-Cover Texts and Back-Cover Texts, replace the "with...Texts." line with this: with the Invariant Sections being LIST THEIR TITLES, with the Front-Cover Texts being LIST, and with the Back-Cover Texts being LIST. If you have Invariant Sections without Cover Texts, or some other combination of the three, merge those two alternatives to suit the situation. If your document contains nontrivial examples of program code, we recommend releasing these examples in parallel under your choice of free software license, such as the GNU General Public License, to permit their use in free software.