C-XSC - A C++ Class Library for Extended Scientific Computing  2.5.4
lx_complex.cpp
1 /*
2 ** CXSC is a C++ library for eXtended Scientific Computing (V 2.5.4)
3 **
4 ** Copyright (C) 1990-2000 Institut fuer Angewandte Mathematik,
5 ** Universitaet Karlsruhe, Germany
6 ** (C) 2000-2014 Wiss. Rechnen/Softwaretechnologie
7 ** Universitaet Wuppertal, Germany
8 **
9 ** This library is free software; you can redistribute it and/or
10 ** modify it under the terms of the GNU Library General Public
11 ** License as published by the Free Software Foundation; either
12 ** version 2 of the License, or (at your option) any later version.
13 **
14 ** This library is distributed in the hope that it will be useful,
15 ** but WITHOUT ANY WARRANTY; without even the implied warranty of
16 ** MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
17 ** Library General Public License for more details.
18 **
19 ** You should have received a copy of the GNU Library General Public
20 ** License along with this library; if not, write to the Free
21 ** Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
22 */
23 
24 /* CVS $Id: lx_complex.cpp,v 1.9 2014/01/30 17:23:47 cxsc Exp $ */
25 
26 /*
27 ** F. Blomquist, University of Wuppertal, 19.09.2007;
28 */
29 
30 #include "lx_complex.hpp"
31 #include "lx_cinterval.hpp"
32 
33 namespace cxsc {
34 // ----------------------------------------------------------------------------
35 // ----- Functions and operators related to type lx_complex -------------------
36 // ----------------------------------------------------------------------------
37 
38  l_complex & l_complex::operator = (const lx_complex& a) throw()
39 {
40  l_real u,v;
41  u = Re(a); v = Im(a);
42  return *this = l_complex(u,v);
43 }
44 
45  complex & complex::operator = (const lx_complex& a) throw()
46 {
47  l_complex x;
48  complex y;
49 
50  x = a;
51  y = x;
52  return *this = y;
53 }
54 
55  std::string & operator >> (std::string& s, lx_complex& a) throw()
56 // Writes string s to variable a of type lx_complex;
57 // and returns an empty string s;
58 // Example: s = "({2,3} , {4,5})" delivers a value a
59 // with: a ~ 10^2 * 3.0 + i*10^4 * 5.0;
60 {
61  string su;
62  int i;
63  std::cout << "Halo 1" << std::endl;
64  s = skipwhitespacessinglechar (s, '(');
65  std::cout << "s = " << s << std::endl;
66  i = s.find("}");
67  std::cout << "i = " << i << std::endl;
68  su = s.substr(0,i+1);
69  std::cout << "su = " << su << std::endl;
70  su >> a.re;
71 
72  s.erase(0,i+1);
73  s = skipwhitespacessinglechar (s, ',');
74  std::cout << "s = " << s << std::endl;
75  s >> a.im;
76 
77  s = "";
78  return s;
79 }
80 
81  void operator >> (const std::string &s, lx_complex &a) throw()
82 {
83  // Writes string s to variable a of type lx_real;
84  std::string r(s);
85  r >> a;
86 }
87 
88  void operator >> (const char *s, lx_complex& a) throw()
89 {
90  std::string r(s);
91  r >> a;
92 }
93 
94  lx_real abs (const lx_complex& a) throw()
95 { return sqrtx2y2(a.re,a.im); }
96  lx_real abs2 (const lx_complex& a) throw()
97  { return a.re*a.re + a.im*a.im; }
98 
99 // -----------------------------------------------------------------------------
100 // --------- Elementary functions related to type lx_complex -------------------
101 // -----------------------------------------------------------------------------
102 
103 lx_complex sqr(const lx_complex& z) throw() { return (z*z); }
104 lx_complex sqrt(const lx_complex& z) throw()
105 { return mid(sqrt(lx_cinterval(z))); }
106 
107 lx_complex sqrt(const lx_complex& z,int n) throw()
108 { return mid( sqrt(lx_cinterval(z),n) ); }
109 
110 lx_complex exp(const lx_complex& z) throw()
111 { return mid(exp(lx_cinterval(z))); }
112 
113 lx_complex exp2(const lx_complex& z) throw()
114 { return mid(exp2(lx_cinterval(z))); }
115 
116 lx_complex exp10(const lx_complex& z) throw()
117 { return mid(exp10(lx_cinterval(z))); }
118 
119 lx_complex sin(const lx_complex& z) throw()
120 { return mid(sin(lx_cinterval(z))); }
121 
122 lx_complex cos(const lx_complex& z) throw()
123 { return mid(cos(lx_cinterval(z))); }
124 
125 lx_complex tan(const lx_complex& z) throw()
126 { return mid(tan(lx_cinterval(z))); }
127 
128 lx_complex cot(const lx_complex& z) throw()
129 { return mid(cot(lx_cinterval(z))); }
130 
131 lx_complex asin(const lx_complex& z) throw()
132 { return mid(asin(lx_cinterval(z))); }
133 
134 lx_complex acos(const lx_complex& z) throw()
135 { return mid(acos(lx_cinterval(z))); }
136 
137 lx_complex atan(const lx_complex& z) throw()
138 { return mid(atan(lx_cinterval(z))); }
139 
140 lx_complex acot(const lx_complex& z) throw()
141 { return mid(acot(lx_cinterval(z))); }
142 
143 lx_complex sinh(const lx_complex& z) throw()
144 { return mid(sinh(lx_cinterval(z))); }
145 
146 lx_complex cosh(const lx_complex& z) throw()
147 { return mid(cosh(lx_cinterval(z))); }
148 
149 lx_complex tanh(const lx_complex& z) throw()
150 { return mid(tanh(lx_cinterval(z))); }
151 
152 lx_complex coth(const lx_complex& z) throw()
153 { return mid(coth(lx_cinterval(z))); }
154 
155 lx_complex asinh(const lx_complex& z) throw()
156 { return mid(asinh(lx_cinterval(z))); }
157 
158 lx_complex acosh(const lx_complex& z) throw()
159 { return mid(acosh(lx_cinterval(z))); }
160 
161 lx_complex atanh(const lx_complex& z) throw()
162 { return mid(atanh(lx_cinterval(z))); }
163 
164 lx_complex acoth(const lx_complex& z) throw()
165 { return mid(acoth(lx_cinterval(z))); }
166 
167 // sqrt_all(c) computes a list of 2 values for all square roots of c
168 std::list<lx_complex> sqrt_all( const lx_complex& c )
169 {
170  lx_complex lc;
171  lc = sqrt(c);
172 
173  std::list<lx_complex> res;
174  res.push_back( lc );
175  res.push_back( -lc );
176 
177  return res;
178 } // end sqrt_all
179 
180 lx_real arg(const lx_complex& z) throw()
181 { return mid(arg(lx_cinterval(z))); }
182 
183 lx_real Arg(const lx_complex& z) throw()
184 { return mid(Arg(lx_cinterval(z))); }
185 
186 std::list<lx_complex> sqrt_all( const lx_complex& z, int n )
187  //
188 // sqrt_all(z,n) computes a list of n values approximating all n-th roots of z
189  //
190 // For n >=3, computing the optimal approximations is very expensive
191 // and thus not considered cost-effective.
192  //
193 // Hence, the polar form is used to calculate sqrt_all(z,n)
194  //
195 {
196  std::list<lx_complex> res;
197 
198  if( n == 0 )
199  {
200  res.push_back( lx_complex(l_real(1),l_real(0)) );
201  return res;
202  }
203  else if( n == 1 )
204  {
205  res.push_back(z);
206  return res;
207  }
208  else if( n == 2 ) return sqrt_all( z );
209  else
210  {
211  lx_real
212  arg_z = arg( z ), root_abs_z = sqrt( abs( z ), n );
213 
214  for(int k = 0; k < n; k++)
215  {
216  lx_real arg_k = ( arg_z + 2 * k * Pi_lx_real() ) / n;
217  res.push_back( lx_complex( root_abs_z * cos( arg_k ),
218  root_abs_z * sin( arg_k ) ) );
219  }
220  return res;
221  }
222 }
223  //
224 //-- end sqrt_all -------------------------------------------------------------
225 
226 lx_complex ln(const lx_complex& z) throw()
227 { return mid(ln(lx_cinterval(z))); }
228 
229 lx_complex log2(const lx_complex& z) throw()
230 { return mid(log2(lx_cinterval(z))); }
231 lx_complex log10(const lx_complex& z) throw()
232 { return mid(log10(lx_cinterval(z))); }
233 
234 lx_complex power_fast(const lx_complex& z, const real& n) throw()
235 {
236  if( n == 0 )
237  return lx_complex(l_real(1),l_real(0));
238  else
239  if( n == 1 ) return z;
240  else
241  if( n == -1 ) return 1 / z;
242  else
243  if( n == 2 ) return sqr(z);
244  else
245  {
246  lx_real abs_z = abs(z);
247  if( ((n < 0) && (abs_z == 0.0)) || !(Is_Integer(n)))
248  // z contains 0
249  cxscthrow (STD_FKT_OUT_OF_DEF(
250  "lx_complex power_fast(const lx_complex& z, const real& n ); z = 0 or n is not integer."));
251  if( abs_z == 0.0 )
252  return lx_complex(l_real(0),l_real(0));
253  else
254  {
255  lx_real arg_z = arg(z);
256  lx_real abs_z_n = exp( n * ln( abs_z ) );
257 
258  return lx_complex( abs_z_n * cos( n * arg_z ),
259  abs_z_n * sin( n * arg_z ) );
260  }
261  }
262 } // End: power_fast
263 
264 lx_complex power(const lx_complex& x, const real& n) throw()
265 {
266  if( !(Is_Integer(n)) )
267  // n is not an integer
268  cxscthrow(STD_FKT_OUT_OF_DEF(
269  "lx_complex power(const lx_complex& z, const real& n); n is not integer."));
270 
271  real zhi(2.0), N(n), r;
272  lx_complex one = lx_complex(l_real(1),l_real(0));
273  double dbl;
274  lx_complex y, neu, X(x);
275 
276  if (x == one) y = x;
277  else
278  if (N == 0.0) y = one;
279  else
280  {
281  if (N == 1) y = x;
282  else
283  if (N == 2) y = sqr(x);
284  else
285  {
286  if (N < 0)
287  {
288  X = 1.0/X;
289  N = -N;
290  }
291  // Initialisierung
292  if ( !Is_Integer(N/2) )
293  y = X;
294  else y = one;
295  neu = sqr(X); // neu = X*X;
296  do
297  {
298  dbl = _double(N/zhi);
299  dbl = floor(dbl);
300  r = (real) dbl;
301  if ( !Is_Integer( r/2 ) )
302  y *= neu;
303  zhi += zhi;
304  if (zhi <= N)
305  neu = sqr(neu); // neu = neu * neu;
306  } while (zhi <= N);
307  }
308  }
309  return y;
310 } // end power(z,n)
311 
312 lx_complex pow(const lx_complex& z, const lx_real& p) throw()
313 { return mid( pow( lx_cinterval(z) , lx_interval(p) ) ); }
314 
315 lx_complex pow(const lx_complex& z, const lx_complex& p) throw()
316 { return mid( pow( lx_cinterval(z) , lx_cinterval(p) ) ); }
317 
318 lx_complex sqrt1px2(const lx_complex& z) throw()
319 { return mid(sqrt1px2(lx_cinterval(z))); }
320 
321 lx_complex sqrt1mx2(const lx_complex& z) throw()
322 { return mid(sqrt1mx2(lx_cinterval(z))); }
323 
324 lx_complex sqrtx2m1(const lx_complex& z) throw()
325 { return mid(sqrtx2m1(lx_cinterval(z))); }
326 
327 lx_complex sqrtp1m1(const lx_complex& z) throw()
328 { return mid(sqrtp1m1(lx_cinterval(z))); }
329 
330 lx_complex expm1(const lx_complex& z) throw()
331 { return mid(expm1(lx_cinterval(z))); }
332 
333 lx_complex lnp1(const lx_complex& z) throw()
334 { return mid(lnp1(lx_cinterval(z))); }
335 
336 } // namespace cxsc
cxsc::Pi_lx_real
lx_real Pi_lx_real()
lx_real approximation for
Definition: lx_interval.cpp:8014
cxsc::ln
cinterval ln(const cinterval &z)
Calculates .
Definition: cimath.cpp:851
cxsc::power_fast
cinterval power_fast(const cinterval &z, int n)
Calculates .
Definition: cimath.cpp:1520
cxsc::lnp1
cinterval lnp1(const cinterval &z)
Calculates .
Definition: cimath.cpp:867
cxsc::atan
cinterval atan(const cinterval &z)
Calculates .
Definition: cimath.cpp:2938
cxsc::exp10
cinterval exp10(const cinterval &z)
Calculates .
Definition: cimath.cpp:172
cxsc::sqrt_all
std::list< cinterval > sqrt_all(const cinterval &z)
Calculates and returns all possible solutions.
Definition: cimath.cpp:1176
cxsc::exp2
cinterval exp2(const cinterval &z)
Calculates .
Definition: cimath.cpp:167
cxsc::tanh
cinterval tanh(const cinterval &z)
Calculates .
Definition: cimath.cpp:565
cxsc::Is_Integer
bool Is_Integer(const real &x)
Returns 1 if x is an integer value and if .
Definition: lx_real.inl:63
cxsc::asin
cinterval asin(const cinterval &z)
Calculates .
Definition: cimath.cpp:2311
cxsc::mid
cvector mid(const cimatrix_subv &mv)
Returns the middle of the matrix.
Definition: cimatrix.inl:739
cxsc::sqrtx2m1
cinterval sqrtx2m1(const cinterval &z)
Calculates .
Definition: cimath.cpp:1109
cxsc::cosh
cinterval cosh(const cinterval &z)
Calculates .
Definition: cimath.cpp:223
cxsc::sqrt1px2
cinterval sqrt1px2(const cinterval &z)
Calculates .
Definition: cimath.cpp:1071
cxsc::tan
cinterval tan(const cinterval &z)
Calculates .
Definition: cimath.cpp:393
cxsc::acoth
cinterval acoth(const cinterval &z)
Calculates .
Definition: cimath.cpp:3330
cxsc::sinh
cinterval sinh(const cinterval &z)
Calculates .
Definition: cimath.cpp:231
cxsc::coth
cinterval coth(const cinterval &z)
Calculates .
Definition: cimath.cpp:578
cxsc::sqrtp1m1
cinterval sqrtp1m1(const cinterval &z)
Calculates .
Definition: cimath.cpp:1054
cxsc::log2
cinterval log2(const cinterval &z)
Calculates .
Definition: cimath.cpp:898
cxsc::l_complex::operator=
l_complex & operator=(const l_real &lr)
Implementation of standard assigning operator.
Definition: l_complex.hpp:63
cxsc::sqrt
cinterval sqrt(const cinterval &z)
Calculates .
Definition: cimath.cpp:1007
cxsc::sqrt1mx2
cinterval sqrt1mx2(const cinterval &z)
Calculates .
Definition: cimath.cpp:1140
cxsc::pow
cinterval pow(const cinterval &z, const interval &p)
Calculates .
Definition: cimath.cpp:2074
cxsc::abs
ivector abs(const cimatrix_subv &mv)
Returns the absolute value of the matrix.
Definition: cimatrix.inl:737
cxsc::l_complex
The Multiple-Precision Data Type l_complex.
Definition: l_complex.hpp:46
cxsc::acosh
cinterval acosh(const cinterval &z)
Calculates .
Definition: cimath.cpp:2732
cxsc::Arg
interval Arg(const cinterval &z)
Calculates .
Definition: cimath.cpp:654
cxsc::sqrtx2y2
interval sqrtx2y2(const interval &x, const interval &y)
Calculates .
Definition: imath.cpp:80
cxsc::power
cinterval power(const cinterval &z, int n)
Calculates .
Definition: cimath.cpp:1941
cxsc::cos
cinterval cos(const cinterval &z)
Calculates .
Definition: cimath.cpp:207
cxsc::acot
cinterval acot(const cinterval &z)
Calculates .
Definition: cimath.cpp:3130
cxsc
The namespace cxsc, providing all functionality of the class library C-XSC.
Definition: cdot.cpp:29
cxsc::l_real
The Multiple-Precision Data Type l_real.
Definition: l_real.hpp:78
cxsc::sqr
cinterval sqr(const cinterval &z)
Calculates .
Definition: cimath.cpp:3342
cxsc::sin
cinterval sin(const cinterval &z)
Calculates .
Definition: cimath.cpp:215
cxsc::cot
cinterval cot(const cinterval &z)
Calculates .
Definition: cimath.cpp:538
cxsc::arg
interval arg(const cinterval &z)
Calculates .
Definition: cimath.cpp:741
cxsc::expm1
cinterval expm1(const cinterval &z)
Calculates .
Definition: cimath.cpp:177
cxsc::acos
cinterval acos(const cinterval &z)
Calculates .
Definition: cimath.cpp:2553
cxsc::exp
cinterval exp(const cinterval &z)
Calculates .
Definition: cimath.cpp:159
cxsc::atanh
cinterval atanh(const cinterval &z)
Calculates .
Definition: cimath.cpp:3317
cxsc::complex
The Scalar Type complex.
Definition: complex.hpp:50
cxsc::log10
cinterval log10(const cinterval &z)
Calculates .
Definition: cimath.cpp:903
cxsc::complex::operator=
complex & operator=(const real &r)
Implementation of standard assigning operator.
Definition: complex.inl:31
cxsc::asinh
cinterval asinh(const cinterval &z)
Calculates .
Definition: cimath.cpp:2718
cxsc::real
The Scalar Type real.
Definition: real.hpp:114