1/*
2 * Copyright (c) 1999
3 * Silicon Graphics Computer Systems, Inc.
4 *
5 * Copyright (c) 1999
6 * Boris Fomitchev
7 *
8 * This material is provided "as is", with absolutely no warranty expressed
9 * or implied. Any use is at your own risk.
10 *
11 * Permission to use or copy this software for any purpose is hereby granted
12 * without fee, provided the above notices are retained on all copies.
13 * Permission to modify the code and to distribute modified code is granted,
14 * provided the above notices are retained, and a notice that the code was
15 * modified is included with the above copyright notice.
16 *
17 */
18#include "stlport_prefix.h"
19
20
21// Trigonometric and hyperbolic functions for complex<float>,
22// complex<double>, and complex<long double>
23#include <complex>
24#include <cfloat>
25#include <cmath>
26
27_STLP_BEGIN_NAMESPACE
28
29
30//----------------------------------------------------------------------
31// helpers
32#if defined (__sgi)
33  static const union { unsigned int i; float f; } float_ulimit = { 0x42b2d4fc };
34  static const float float_limit = float_ulimit.f;
35  static union {
36    struct { unsigned int h; unsigned int l; } w;
37    double d;
38  } double_ulimit = { 0x408633ce, 0x8fb9f87d };
39  static const double double_limit = double_ulimit.d;
40  static union {
41    struct { unsigned int h[2]; unsigned int l[2]; } w;
42    long double ld;
43  } ldouble_ulimit = {0x408633ce, 0x8fb9f87e, 0xbd23b659, 0x4e9bd8b1};
44#  if !defined (_STLP_NO_LONG_DOUBLE)
45  static const long double ldouble_limit = ldouble_ulimit.ld;
46#  endif
47#else
48#  if defined (M_LN2) && defined (FLT_MAX_EXP)
49  static const float float_limit = float(M_LN2 * FLT_MAX_EXP);
50  static const double double_limit = M_LN2 * DBL_MAX_EXP;
51#  else
52  static const float float_limit = ::log(FLT_MAX);
53  static const double double_limit = ::log(DBL_MAX);
54#  endif
55#  if !defined (_STLP_NO_LONG_DOUBLE)
56#    if defined (M_LN2l)
57  static const long double ldouble_limit = M_LN2l * LDBL_MAX_EXP;
58#    else
59  static const long double ldouble_limit = ::log(LDBL_MAX);
60#    endif
61#  endif
62#endif
63
64
65//----------------------------------------------------------------------
66// sin
67template <class _Tp>
68static complex<_Tp> sinT(const complex<_Tp>& z) {
69  return complex<_Tp>(::sin(z._M_re) * ::cosh(z._M_im),
70                      ::cos(z._M_re) * ::sinh(z._M_im));
71}
72
73_STLP_DECLSPEC complex<float> _STLP_CALL sin(const complex<float>& z)
74{ return sinT(z); }
75
76_STLP_DECLSPEC complex<double> _STLP_CALL sin(const complex<double>& z)
77{ return sinT(z); }
78
79#if !defined (_STLP_NO_LONG_DOUBLE)
80_STLP_DECLSPEC complex<long double> _STLP_CALL sin(const complex<long double>& z)
81{ return sinT(z); }
82#endif
83
84//----------------------------------------------------------------------
85// cos
86template <class _Tp>
87static complex<_Tp> cosT(const complex<_Tp>& z) {
88  return complex<_Tp>(::cos(z._M_re) * ::cosh(z._M_im),
89                     -::sin(z._M_re) * ::sinh(z._M_im));
90}
91
92_STLP_DECLSPEC complex<float> _STLP_CALL cos(const complex<float>& z)
93{ return cosT(z); }
94
95_STLP_DECLSPEC complex<double> _STLP_CALL cos(const complex<double>& z)
96{ return cosT(z); }
97
98#if !defined (_STLP_NO_LONG_DOUBLE)
99_STLP_DECLSPEC complex<long double> _STLP_CALL cos(const complex<long double>& z)
100{ return cosT(z); }
101#endif
102
103//----------------------------------------------------------------------
104// tan
105template <class _Tp>
106static complex<_Tp> tanT(const complex<_Tp>& z, const _Tp& Tp_limit) {
107  _Tp re2 = 2.f * z._M_re;
108  _Tp im2 = 2.f * z._M_im;
109
110  if (::abs(im2) > Tp_limit)
111    return complex<_Tp>(0.f, (im2 > 0 ? 1.f : -1.f));
112  else {
113    _Tp den = ::cos(re2) + ::cosh(im2);
114    return complex<_Tp>(::sin(re2) / den, ::sinh(im2) / den);
115  }
116}
117
118_STLP_DECLSPEC complex<float> _STLP_CALL tan(const complex<float>& z)
119{ return tanT(z, float_limit); }
120
121_STLP_DECLSPEC complex<double> _STLP_CALL tan(const complex<double>& z)
122{ return tanT(z, double_limit); }
123
124#if !defined (_STLP_NO_LONG_DOUBLE)
125_STLP_DECLSPEC complex<long double> _STLP_CALL tan(const complex<long double>& z)
126{ return tanT(z, ldouble_limit); }
127#endif
128
129//----------------------------------------------------------------------
130// sinh
131template <class _Tp>
132static complex<_Tp> sinhT(const complex<_Tp>& z) {
133  return complex<_Tp>(::sinh(z._M_re) * ::cos(z._M_im),
134                      ::cosh(z._M_re) * ::sin(z._M_im));
135}
136
137_STLP_DECLSPEC complex<float> _STLP_CALL sinh(const complex<float>& z)
138{ return sinhT(z); }
139
140_STLP_DECLSPEC complex<double> _STLP_CALL sinh(const complex<double>& z)
141{ return sinhT(z); }
142
143#if !defined (_STLP_NO_LONG_DOUBLE)
144_STLP_DECLSPEC complex<long double> _STLP_CALL sinh(const complex<long double>& z)
145{ return sinhT(z); }
146#endif
147
148//----------------------------------------------------------------------
149// cosh
150template <class _Tp>
151static complex<_Tp> coshT(const complex<_Tp>& z) {
152  return complex<_Tp>(::cosh(z._M_re) * ::cos(z._M_im),
153                      ::sinh(z._M_re) * ::sin(z._M_im));
154}
155
156_STLP_DECLSPEC complex<float> _STLP_CALL cosh(const complex<float>& z)
157{ return coshT(z); }
158
159_STLP_DECLSPEC complex<double> _STLP_CALL cosh(const complex<double>& z)
160{ return coshT(z); }
161
162#if !defined (_STLP_NO_LONG_DOUBLE)
163_STLP_DECLSPEC complex<long double> _STLP_CALL cosh(const complex<long double>& z)
164{ return coshT(z); }
165#endif
166
167//----------------------------------------------------------------------
168// tanh
169template <class _Tp>
170static complex<_Tp> tanhT(const complex<_Tp>& z, const _Tp& Tp_limit) {
171  _Tp re2 = 2.f * z._M_re;
172  _Tp im2 = 2.f * z._M_im;
173  if (::abs(re2) > Tp_limit)
174    return complex<_Tp>((re2 > 0 ? 1.f : -1.f), 0.f);
175  else {
176    _Tp den = ::cosh(re2) + ::cos(im2);
177    return complex<_Tp>(::sinh(re2) / den, ::sin(im2) / den);
178  }
179}
180
181_STLP_DECLSPEC complex<float> _STLP_CALL tanh(const complex<float>& z)
182{ return tanhT(z, float_limit); }
183
184_STLP_DECLSPEC complex<double> _STLP_CALL tanh(const complex<double>& z)
185{ return tanhT(z, double_limit); }
186
187#if !defined (_STLP_NO_LONG_DOUBLE)
188_STLP_DECLSPEC complex<long double> _STLP_CALL tanh(const complex<long double>& z)
189{ return tanhT(z, ldouble_limit); }
190#endif
191
192_STLP_END_NAMESPACE
193