1 /**
2  *  \file complex.h
3  *  Class for Complex built on top of Number class
4  *
5  **/
6 #ifndef SYMENGINE_COMPLEX_H
7 #define SYMENGINE_COMPLEX_H
8 
9 #include <symengine/rational.h>
10 #include <symengine/symengine_exception.h>
11 
12 namespace SymEngine
13 {
14 //! ComplexBase Class for deriving all complex classes
15 class ComplexBase : public Number
16 {
17 public:
18     virtual RCP<const Number> real_part() const = 0;
19     virtual RCP<const Number> imaginary_part() const = 0;
20     bool is_re_zero() const;
21 };
22 
23 //! \return true if 'b' is any of the subclasses of ComplexBase
is_a_Complex(const Basic & b)24 inline bool is_a_Complex(const Basic &b)
25 {
26     return (b.get_type_code() == SYMENGINE_COMPLEX
27             || b.get_type_code() == SYMENGINE_COMPLEX_MPC
28             || b.get_type_code() == SYMENGINE_COMPLEX_DOUBLE);
29 }
30 
31 //! Complex Class
32 class Complex : public ComplexBase
33 {
34 public:
35     //! `real_` : Real part of the complex Number
36     //! `imaginary_` : Imaginary part of the complex Number
37     // Complex Number is of the form `real + i(imaginary)`
38     rational_class real_;
39     rational_class imaginary_;
40 
41 public:
42     IMPLEMENT_TYPEID(SYMENGINE_COMPLEX)
43     //! Constructor of Complex class
44     Complex(rational_class real, rational_class imaginary);
45     /*! Creates an instance of Complex if imaginary part is non-zero
46      * \param `re` must already be in rational_class canonical form
47      * \param `im` must already be in rational_class canonical form
48      * \return Complex or Rational depending on imaginary part.
49      * */
50     static RCP<const Number> from_mpq(const rational_class re,
51                                       const rational_class im);
52     //! \return true if canonical
53     bool is_canonical(const rational_class &real,
54                       const rational_class &imaginary) const;
55     //! \return size of the hash
56     virtual hash_t __hash__() const;
57     /*! Equality comparator
58      * \param o - Object to be compared with
59      * \return whether the 2 objects are equal
60      * */
61     virtual bool __eq__(const Basic &o) const;
62     virtual int compare(const Basic &o) const;
63     //! Get the real part of the complex number
64     virtual RCP<const Number> real_part() const;
65     //! Get the imaginary part of the complex number
66     virtual RCP<const Number> imaginary_part() const;
67     //! Get the conjugate of the complex number
68     virtual RCP<const Basic> conjugate() const;
69     //! \returns `false`
70     // False is returned because complex cannot be compared with zero
is_positive()71     inline virtual bool is_positive() const
72     {
73         return false;
74     }
75     //! \returns `false`
76     // False is returned because complex cannot be compared with zero
is_negative()77     inline virtual bool is_negative() const
78     {
79         return false;
80     }
81     //! \returns `true`
is_complex()82     inline virtual bool is_complex() const
83     {
84         return true;
85     }
86 
87     /*! Constructs Complex from re, im. If im is 0
88      * it will return a Rational instead.
89      * */
90     static RCP<const Number> from_two_rats(const Rational &re,
91                                            const Rational &im);
92 
93     /*! Constructs Complex from re, im. If im is 0
94      * it will return a Rational instead.
95      * */
96     static RCP<const Number> from_two_nums(const Number &re, const Number &im);
97 
98     //! \return `false` since `imaginary_` cannot be zero
is_zero()99     virtual bool is_zero() const
100     {
101         return false;
102     }
103     //! \return `false` since `imaginary_` cannot be zero
is_one()104     virtual bool is_one() const
105     {
106         return false;
107     }
108     //! \return `false` since `imaginary_` cannot be zero
is_minus_one()109     virtual bool is_minus_one() const
110     {
111         return false;
112     }
113 
114     /*! Add Complex
115      * \param other of type Complex
116      * */
addcomp(const Complex & other)117     inline RCP<const Number> addcomp(const Complex &other) const
118     {
119         return from_mpq(this->real_ + other.real_,
120                         this->imaginary_ + other.imaginary_);
121     }
122     /*! Add Complex
123      * \param other of type Rational
124      * */
addcomp(const Rational & other)125     inline RCP<const Number> addcomp(const Rational &other) const
126     {
127         return from_mpq(this->real_ + other.as_rational_class(),
128                         this->imaginary_);
129     }
130     /*! Add Complex
131      * \param other of type Integer
132      * */
addcomp(const Integer & other)133     inline RCP<const Number> addcomp(const Integer &other) const
134     {
135         return from_mpq(this->real_ + other.as_integer_class(),
136                         this->imaginary_);
137     }
138 
139     /*! Subtract Complex
140      * \param other of type Complex
141      * */
subcomp(const Complex & other)142     inline RCP<const Number> subcomp(const Complex &other) const
143     {
144         return from_mpq(this->real_ - other.real_,
145                         this->imaginary_ - other.imaginary_);
146     }
147     /*! Subtract Complex
148      * \param other of type Rational
149      * */
subcomp(const Rational & other)150     inline RCP<const Number> subcomp(const Rational &other) const
151     {
152         return from_mpq(this->real_ - other.as_rational_class(),
153                         this->imaginary_);
154     }
155     /*! Subtract Complex
156      * \param other of type Integer
157      * */
subcomp(const Integer & other)158     inline RCP<const Number> subcomp(const Integer &other) const
159     {
160         return from_mpq(this->real_ - other.as_integer_class(),
161                         this->imaginary_);
162     }
163     /*! Subtract Complex from other
164      * \param other of type Complex
165      * */
rsubcomp(const Complex & other)166     inline RCP<const Number> rsubcomp(const Complex &other) const
167     {
168         return from_mpq(other.real_ - this->real_,
169                         other.imaginary_ - this->imaginary_);
170     }
171     /*! Subtract Complex from other
172      * \param other of type Rational
173      * */
rsubcomp(const Rational & other)174     inline RCP<const Number> rsubcomp(const Rational &other) const
175     {
176         return from_mpq(other.as_rational_class() - this->real_,
177                         -this->imaginary_);
178     }
179     /*! Subtract Complex from other
180      * \param other of type Integer
181      * */
rsubcomp(const Integer & other)182     inline RCP<const Number> rsubcomp(const Integer &other) const
183     {
184         return from_mpq(other.as_integer_class() - this->real_,
185                         -this->imaginary_);
186     }
187 
188     /*! Multiply Complex
189      * \param other of type Complex
190      * */
mulcomp(const Complex & other)191     inline RCP<const Number> mulcomp(const Complex &other) const
192     {
193         return from_mpq(
194             this->real_ * other.real_ - this->imaginary_ * other.imaginary_,
195             this->real_ * other.imaginary_ + this->imaginary_ * other.real_);
196     }
197     /*! Multiply Complex
198      * \param other of type Rational
199      * */
mulcomp(const Rational & other)200     inline RCP<const Number> mulcomp(const Rational &other) const
201     {
202         return from_mpq(this->real_ * other.as_rational_class(),
203                         this->imaginary_ * other.as_rational_class());
204     }
205     /*! Multiply Complex
206      * \param other of type Integer
207      * */
mulcomp(const Integer & other)208     inline RCP<const Number> mulcomp(const Integer &other) const
209     {
210         return from_mpq(this->real_ * other.as_integer_class(),
211                         this->imaginary_ * other.as_integer_class());
212     }
213 
214     /*! Divide Complex
215      * \param other of type Complex
216      * */
divcomp(const Complex & other)217     inline RCP<const Number> divcomp(const Complex &other) const
218     {
219         rational_class modulus_sq_other
220             = other.real_ * other.real_ + other.imaginary_ * other.imaginary_;
221 
222         if (get_num(modulus_sq_other) == 0) {
223             rational_class modulus_sq_this
224                 = this->real_ * this->real_
225                   + this->imaginary_ * this->imaginary_;
226             if (get_num(modulus_sq_this) == 0) {
227                 return Nan;
228             } else {
229                 return ComplexInf;
230             }
231         } else {
232             return from_mpq((this->real_ * other.real_
233                              + this->imaginary_ * other.imaginary_)
234                                 / modulus_sq_other,
235                             (-this->real_ * other.imaginary_
236                              + this->imaginary_ * other.real_)
237                                 / modulus_sq_other);
238         }
239     }
240     /*! Divide Complex
241      * \param other of type Rational
242      * */
divcomp(const Rational & other)243     inline RCP<const Number> divcomp(const Rational &other) const
244     {
245         if (other.is_zero()) {
246             rational_class modulus_sq_this
247                 = this->real_ * this->real_
248                   + this->imaginary_ * this->imaginary_;
249 
250             if (get_num(modulus_sq_this) == 0) {
251                 return Nan;
252             } else {
253                 return ComplexInf;
254             }
255         } else {
256             return from_mpq(this->real_ / other.as_rational_class(),
257                             this->imaginary_ / other.as_rational_class());
258         }
259     }
260     /*! Divide Complex
261      * \param other of type Integer
262      * */
divcomp(const Integer & other)263     inline RCP<const Number> divcomp(const Integer &other) const
264     {
265         if (other.is_zero()) {
266             rational_class modulus_sq_this
267                 = this->real_ * this->real_
268                   + this->imaginary_ * this->imaginary_;
269 
270             if (get_num(modulus_sq_this) == 0) {
271                 return Nan;
272             } else {
273                 return ComplexInf;
274             }
275         } else {
276             return from_mpq(this->real_ / other.as_integer_class(),
277                             this->imaginary_ / other.as_integer_class());
278         }
279     }
280     /*! Divide other by the Complex
281      * \param other of type Integer
282      * */
rdivcomp(const Integer & other)283     inline RCP<const Number> rdivcomp(const Integer &other) const
284     {
285         rational_class modulus_sq_this
286             = this->real_ * this->real_ + this->imaginary_ * this->imaginary_;
287 
288         if (get_num(modulus_sq_this) == 0) {
289             if (other.is_zero()) {
290                 return Nan;
291             } else {
292                 return ComplexInf;
293             }
294         } else {
295             return from_mpq((this->real_ * other.as_integer_class())
296                                 / modulus_sq_this,
297                             (this->imaginary_ * (-other.as_integer_class()))
298                                 / modulus_sq_this);
299         }
300     }
301     /*! Pow Complex
302      * \param other of type Integer
303      * */
304     RCP<const Number> powcomp(const Integer &other) const;
305 
306     //! Converts the param `other` appropriately and then calls `addcomp`
add(const Number & other)307     virtual RCP<const Number> add(const Number &other) const
308     {
309         if (is_a<Rational>(other)) {
310             return addcomp(down_cast<const Rational &>(other));
311         } else if (is_a<Integer>(other)) {
312             return addcomp(down_cast<const Integer &>(other));
313         } else if (is_a<Complex>(other)) {
314             return addcomp(down_cast<const Complex &>(other));
315         } else {
316             return other.add(*this);
317         }
318     };
319     //! Converts the param `other` appropriately and then calls `subcomp`
sub(const Number & other)320     virtual RCP<const Number> sub(const Number &other) const
321     {
322         if (is_a<Rational>(other)) {
323             return subcomp(down_cast<const Rational &>(other));
324         } else if (is_a<Integer>(other)) {
325             return subcomp(down_cast<const Integer &>(other));
326         } else if (is_a<Complex>(other)) {
327             return subcomp(down_cast<const Complex &>(other));
328         } else {
329             return other.rsub(*this);
330         }
331     };
332     //! Converts the param `other` appropriately and then calls `rsubcomp`
rsub(const Number & other)333     virtual RCP<const Number> rsub(const Number &other) const
334     {
335         if (is_a<Rational>(other)) {
336             return rsubcomp(down_cast<const Rational &>(other));
337         } else if (is_a<Integer>(other)) {
338             return rsubcomp(down_cast<const Integer &>(other));
339         } else {
340             throw NotImplementedError("Not Implemented");
341         }
342     };
343     //! Converts the param `other` appropriately and then calls `mulcomp`
mul(const Number & other)344     virtual RCP<const Number> mul(const Number &other) const
345     {
346         if (is_a<Rational>(other)) {
347             return mulcomp(down_cast<const Rational &>(other));
348         } else if (is_a<Integer>(other)) {
349             return mulcomp(down_cast<const Integer &>(other));
350         } else if (is_a<Complex>(other)) {
351             return mulcomp(down_cast<const Complex &>(other));
352         } else {
353             return other.mul(*this);
354         }
355     };
356     //! Converts the param `other` appropriately and then calls `divcomp`
div(const Number & other)357     virtual RCP<const Number> div(const Number &other) const
358     {
359         if (is_a<Rational>(other)) {
360             return divcomp(down_cast<const Rational &>(other));
361         } else if (is_a<Integer>(other)) {
362             return divcomp(down_cast<const Integer &>(other));
363         } else if (is_a<Complex>(other)) {
364             return divcomp(down_cast<const Complex &>(other));
365         } else {
366             return other.rdiv(*this);
367         }
368     };
369     //! Converts the param `other` appropriately and then calls `rdivcomp`
rdiv(const Number & other)370     virtual RCP<const Number> rdiv(const Number &other) const
371     {
372         if (is_a<Integer>(other)) {
373             return rdivcomp(down_cast<const Integer &>(other));
374         } else {
375             throw NotImplementedError("Not Implemented");
376         }
377     };
378     //! Converts the param `other` appropriately and then calls `powcomp`
pow(const Number & other)379     virtual RCP<const Number> pow(const Number &other) const
380     {
381         if (is_a<Integer>(other)) {
382             return powcomp(down_cast<const Integer &>(other));
383         } else {
384             return other.rpow(*this);
385         }
386     };
387 
rpow(const Number & other)388     virtual RCP<const Number> rpow(const Number &other) const
389     {
390         throw NotImplementedError("Not Implemented");
391     };
392 };
393 
394 } // namespace SymEngine
395 
396 #endif
397