mingwrt/mingwex/complex/catan_generic.c

   1 /*
   2  * catan_generic.c
   3  *
   4  * $Id$
   5  *
   6  * Compute the complex arctan corresponding to a complex tangent value;
   7  * a generic implementation for catan(), catanf(), and catanl().
   8  *
   9  * Written by Keith Marshall <keithmarshall@users.sourceforge.net>
  10  * This is an adaptation of an original contribution by Danny Smith
  11  * Copyright (C) 2003-2005, 2015, MinGW.org Project
  12  *
  13  *
  14  * Permission is hereby granted, free of charge, to any person obtaining a
  15  * copy of this software and associated documentation files (the "Software"),
  16  * to deal in the Software without restriction, including without limitation
  17  * the rights to use, copy, modify, merge, publish, distribute, sublicense,
  18  * and/or sell copies of the Software, and to permit persons to whom the
  19  * Software is furnished to do so, subject to the following conditions:
  20  *
  21  * The above copyright notice, this permission notice, and the following
  22  * disclaimer shall be included in all copies or substantial portions of
  23  * the Software.
  24  *
  25  * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
  26  * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
  27  * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT.  IN NO EVENT SHALL
  28  * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
  29  * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
  30  * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OF OR OTHER
  31  * DEALINGS IN THE SOFTWARE.
  32  *
  33  *
  34  * This is a generic implementation for all of the catan(), catanl(),
  35  * and catanh() functions; each is to be compiled separately, i.e.
  36  *
  37  *   gcc -D FUNCTION=catan  -o catan.o  catan_generic.c
  38  *   gcc -D FUNCTION=catanl -o catanl.o catan_generic.c
  39  *   gcc -D FUNCTION=catanf -o catanf.o catan_generic.c
  40  *
  41  */
  42 #include <math.h>
  43 #include <complex.h>
  44
  45 #ifndef FUNCTION
  46 /* If user neglected to specify it, the default compilation is for
  47  * the catan() function.
  48  */
  49 # define FUNCTION catan
  50 #endif
  51
  52 #define argtype_catan  double
  53 #define argtype_catanl long double
  54 #define argtype_catanf float
  55
  56 #define PASTE(PREFIX,SUFFIX)   PREFIX##SUFFIX
  57 #define mapname(PREFIX,SUFFIX) PASTE(PREFIX,SUFFIX)
  58
  59 #define ARGTYPE(FUNCTION) PASTE(argtype_,FUNCTION)
  60 #define ARGCAST(VALUE)    mapname(argcast_,FUNCTION)(VALUE)
  61
  62 #define argcast_catan(VALUE)  VALUE
  63 #define argcast_catanl(VALUE) PASTE(VALUE,L)
  64 #define argcast_catanf(VALUE) PASTE(VALUE,F)
  65
  66 #define mapfunc(NAME) mapname(mapfunc_,FUNCTION)(NAME)
  67
  68 #define mapfunc_catan(NAME)  NAME
  69 #define mapfunc_catanl(NAME) PASTE(NAME,l)
  70 #define mapfunc_catanf(NAME) PASTE(NAME,f)
  71
  72 /* Define the generic function implementation.
  73  */
  74 ARGTYPE(FUNCTION) complex FUNCTION( ARGTYPE(FUNCTION) complex Z )
  75 {
  76   ARGTYPE(FUNCTION) x = __real__ Z;     /* real part of Z saved as x */
  77   ARGTYPE(FUNCTION) y = __imag__ Z;     /* imaginary part of Z saved as iy */
  78
  79   /* Having thus saved the original value of Z as separate real and
  80    * imaginary parts, we may reuse Z to accumulate the arctan value,
  81    * noting that...
  82    */
  83   if( isfinite( mapfunc(hypot)( x, y )) )
  84   {
  85     /* ...for any complex number with finite modulus, its argtan may
  86      * be computed as:
  87      *
  88      *   catan(Z) = (clog(1.0 + iZ) - clog(1.0 - iZ)) / 2.0i
  89      *
  90      * For computational convenience, we may introduce an additional
  91      * intermediate accumulator variable...
  92      */
  93     ARGTYPE(FUNCTION) complex tmp;
  94
  95     /* ...in which we initially compute the (1.0 - iZ) term...
  96      */
  97     __real__ tmp = ARGCAST(1.0) + y;
  98     __imag__ tmp = -x;
  99
 100     /* ...while computing (1.0 + iZ) directly in Z itself.
 101      */
 102     __real__ Z = ARGCAST(1.0) - y;
 103     __imag__ Z = x;
 104
 105     /* We may then compute and combine the complex logarithms
 106      * of this pair of sub-expressions, while we simultaneously
 107      * perform the REAL division by 2.0...
 108      */
 109     tmp = ARGCAST(0.5) * (mapfunc(clog)( Z ) - mapfunc(clog)( tmp ));
 110
 111     /* ...and finally, the division by i has the effect of
 112      * exchanging the real and imaginary parts of the result,
 113      * with a change of sign in the resultant imaginary part.
 114      */
 115     __real__ Z = __imag__ tmp;
 116     __imag__ Z = -(__real__ tmp);
 117   }
 118
 119   else
 120   { /* Conversely, if the modulus of Z is in any way undefined,
 121      * we return a real value, with modulus of PI / 2.0, and with
 122      * the same sign as the original real part of Z.
 123      */
 124     __real__ Z = mapfunc(copysign)( ARGCAST(M_PI_2), x );
 125     __imag__ Z = ARGCAST(0.0);
 126   }
 127
 128   /* In either case, Z now represents the arctan of its original
 129    * value, which we are required to return.
 130    */
 131   return Z;
 132 }
 133
 134 /* $RCSfile$: end of file */