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
33namespace cxsc {
34// ----------------------------------------------------------------------------
35// ----- Functions and operators related to type lx_complex -------------------
36// ----------------------------------------------------------------------------
37
38 l_complex & l_complex::operator = (const lx_complex& a) noexcept
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) noexcept
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) noexcept
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) noexcept
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) noexcept
89{
90 std::string r(s);
91 r >> a;
92}
93
94 lx_real abs (const lx_complex& a) noexcept
95{ return sqrtx2y2(a.re,a.im); }
96 lx_real abs2 (const lx_complex& a) noexcept
97 { return a.re*a.re + a.im*a.im; }
98
99// -----------------------------------------------------------------------------
100// --------- Elementary functions related to type lx_complex -------------------
101// -----------------------------------------------------------------------------
102
103lx_complex sqr(const lx_complex& z) noexcept { return (z*z); }
104lx_complex sqrt(const lx_complex& z) noexcept
105{ return mid(sqrt(lx_cinterval(z))); }
106
107lx_complex sqrt(const lx_complex& z,int n) noexcept
108{ return mid( sqrt(lx_cinterval(z),n) ); }
109
110lx_complex exp(const lx_complex& z) noexcept
111{ return mid(exp(lx_cinterval(z))); }
112
113lx_complex exp2(const lx_complex& z) noexcept
114{ return mid(exp2(lx_cinterval(z))); }
115
116lx_complex exp10(const lx_complex& z) noexcept
117{ return mid(exp10(lx_cinterval(z))); }
118
119lx_complex sin(const lx_complex& z) noexcept
120{ return mid(sin(lx_cinterval(z))); }
121
122lx_complex cos(const lx_complex& z) noexcept
123{ return mid(cos(lx_cinterval(z))); }
124
125lx_complex tan(const lx_complex& z) noexcept
126{ return mid(tan(lx_cinterval(z))); }
127
128lx_complex cot(const lx_complex& z) noexcept
129{ return mid(cot(lx_cinterval(z))); }
130
131lx_complex asin(const lx_complex& z) noexcept
132{ return mid(asin(lx_cinterval(z))); }
133
134lx_complex acos(const lx_complex& z) noexcept
135{ return mid(acos(lx_cinterval(z))); }
136
137lx_complex atan(const lx_complex& z) noexcept
138{ return mid(atan(lx_cinterval(z))); }
139
140lx_complex acot(const lx_complex& z) noexcept
141{ return mid(acot(lx_cinterval(z))); }
142
143lx_complex sinh(const lx_complex& z) noexcept
144{ return mid(sinh(lx_cinterval(z))); }
145
146lx_complex cosh(const lx_complex& z) noexcept
147{ return mid(cosh(lx_cinterval(z))); }
148
149lx_complex tanh(const lx_complex& z) noexcept
150{ return mid(tanh(lx_cinterval(z))); }
151
152lx_complex coth(const lx_complex& z) noexcept
153{ return mid(coth(lx_cinterval(z))); }
154
155lx_complex asinh(const lx_complex& z) noexcept
156{ return mid(asinh(lx_cinterval(z))); }
157
158lx_complex acosh(const lx_complex& z) noexcept
159{ return mid(acosh(lx_cinterval(z))); }
160
161lx_complex atanh(const lx_complex& z) noexcept
162{ return mid(atanh(lx_cinterval(z))); }
163
164lx_complex acoth(const lx_complex& z) noexcept
165{ return mid(acoth(lx_cinterval(z))); }
166
167// sqrt_all(c) computes a list of 2 values for all square roots of c
168std::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
180lx_real arg(const lx_complex& z) noexcept
181{ return mid(arg(lx_cinterval(z))); }
182
183lx_real Arg(const lx_complex& z) noexcept
184{ return mid(Arg(lx_cinterval(z))); }
185
186std::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
226lx_complex ln(const lx_complex& z) noexcept
227{ return mid(ln(lx_cinterval(z))); }
228
229lx_complex log2(const lx_complex& z) noexcept
230{ return mid(log2(lx_cinterval(z))); }
231lx_complex log10(const lx_complex& z) noexcept
232{ return mid(log10(lx_cinterval(z))); }
233
234lx_complex power_fast(const lx_complex& z, const real& n) noexcept
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
264lx_complex power(const lx_complex& x, const real& n) noexcept
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
312lx_complex pow(const lx_complex& z, const lx_real& p) noexcept
313{ return mid( pow( lx_cinterval(z) , lx_interval(p) ) ); }
314
315lx_complex pow(const lx_complex& z, const lx_complex& p) noexcept
316{ return mid( pow( lx_cinterval(z) , lx_cinterval(p) ) ); }
317
318lx_complex sqrt1px2(const lx_complex& z) noexcept
319{ return mid(sqrt1px2(lx_cinterval(z))); }
320
321lx_complex sqrt1mx2(const lx_complex& z) noexcept
322{ return mid(sqrt1mx2(lx_cinterval(z))); }
323
324lx_complex sqrtx2m1(const lx_complex& z) noexcept
325{ return mid(sqrtx2m1(lx_cinterval(z))); }
326
327lx_complex sqrtp1m1(const lx_complex& z) noexcept
328{ return mid(sqrtp1m1(lx_cinterval(z))); }
329
330lx_complex expm1(const lx_complex& z) noexcept
331{ return mid(expm1(lx_cinterval(z))); }
332
333lx_complex lnp1(const lx_complex& z) noexcept
334{ return mid(lnp1(lx_cinterval(z))); }
335
336} // namespace cxsc
The Scalar Type complex.
Definition complex.hpp:50
complex & operator=(const real &r) noexcept
Implementation of standard assigning operator.
Definition complex.inl:31
The Multiple-Precision Data Type l_complex.
Definition l_complex.hpp:46
l_complex & operator=(const l_real &lr) noexcept
Implementation of standard assigning operator.
Definition l_complex.hpp:63
The Multiple-Precision Data Type l_real.
Definition l_real.hpp:78
The Scalar Type real.
Definition real.hpp:114
The namespace cxsc, providing all functionality of the class library C-XSC.
Definition cdot.cpp:29
cinterval sqrtp1m1(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:1054
cinterval exp2(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:167
cinterval sqrt1mx2(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:1140
cinterval asinh(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:2718
cinterval coth(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:578
cinterval log2(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:898
cinterval power(const cinterval &z, int n) noexcept
Calculates .
Definition cimath.cpp:1941
cinterval log10(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:903
cinterval ln(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:851
bool Is_Integer(const real &x)
Returns 1 if x is an integer value and if .
Definition lx_real.inl:63
cinterval pow(const cinterval &z, const interval &p) noexcept
Calculates .
Definition cimath.cpp:2074
cinterval sinh(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:231
cinterval asin(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:2311
cinterval tan(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:393
cinterval exp10(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:172
interval arg(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:741
std::list< cinterval > sqrt_all(const cinterval &z)
Calculates and returns all possible solutions.
Definition cimath.cpp:1176
cinterval acos(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:2553
cinterval sqrtx2m1(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:1109
cinterval acosh(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:2732
cinterval cosh(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:223
cinterval cos(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:207
cinterval sqrt1px2(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:1071
cinterval exp(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:159
cinterval tanh(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:565
cinterval expm1(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:177
cinterval cot(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:538
ivector abs(const cimatrix_subv &mv) noexcept
Returns the absolute value of the matrix.
Definition cimatrix.inl:737
cinterval sqrt(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:1007
cinterval power_fast(const cinterval &z, int n) noexcept
Calculates .
Definition cimath.cpp:1520
cinterval acot(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:3130
cinterval sqr(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:3342
cinterval lnp1(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:867
cvector mid(const cimatrix_subv &mv) noexcept
Returns the middle of the matrix.
Definition cimatrix.inl:739
cinterval atan(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:2938
cinterval atanh(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:3317
interval Arg(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:654
cinterval acoth(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:3330
interval sqrtx2y2(const interval &x, const interval &y) noexcept
Calculates .
Definition imath.cpp:80
lx_real Pi_lx_real() noexcept
lx_real approximation for
cinterval sin(const cinterval &z) noexcept
Calculates .
Definition cimath.cpp:215