GiNaC 1.8.5
power.cpp
Go to the documentation of this file.
1
5/*
6 * GiNaC Copyright (C) 1999-2023 Johannes Gutenberg University Mainz, Germany
7 *
8 * This program is free software; you can redistribute it and/or modify
9 * it under the terms of the GNU General Public License as published by
10 * the Free Software Foundation; either version 2 of the License, or
11 * (at your option) any later version.
12 *
13 * This program is distributed in the hope that it will be useful,
14 * but WITHOUT ANY WARRANTY; without even the implied warranty of
15 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16 * GNU General Public License for more details.
17 *
18 * You should have received a copy of the GNU General Public License
19 * along with this program; if not, write to the Free Software
20 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
21 */
22
23#include "power.h"
24#include "expairseq.h"
25#include "add.h"
26#include "mul.h"
27#include "ncmul.h"
28#include "numeric.h"
29#include "constant.h"
30#include "operators.h"
31#include "inifcns.h" // for log() in power::derivative()
32#include "matrix.h"
33#include "indexed.h"
34#include "symbol.h"
35#include "lst.h"
36#include "archive.h"
37#include "utils.h"
38#include "relational.h"
39#include "compiler.h"
40
41#include <iostream>
42#include <limits>
43#include <stdexcept>
44#include <vector>
45#include <algorithm>
46
47namespace GiNaC {
48
51 print_func<print_latex>(&power::do_print_latex).
52 print_func<print_csrc>(&power::do_print_csrc).
53 print_func<print_python>(&power::do_print_python).
54 print_func<print_python_repr>(&power::do_print_python_repr).
55 print_func<print_csrc_cl_N>(&power::do_print_csrc_cl_N))
56
57
58// default constructor
60
61power::power() { }
62
64// other constructors
66
67// all inlined
68
70// archiving
72
73void power::read_archive(const archive_node &n, lst &sym_lst)
74{
75 inherited::read_archive(n, sym_lst);
76 n.find_ex("basis", basis, sym_lst);
77 n.find_ex("exponent", exponent, sym_lst);
78}
79
81{
82 inherited::archive(n);
83 n.add_ex("basis", basis);
84 n.add_ex("exponent", exponent);
85}
86
88// functions overriding virtual functions from base classes
90
91// public
92
93void power::print_power(const print_context & c, const char *powersymbol, const char *openbrace, const char *closebrace, unsigned level) const
94{
95 // Ordinary output of powers using '^' or '**'
96 if (precedence() <= level)
97 c.s << openbrace << '(';
99 c.s << powersymbol;
100 c.s << openbrace;
102 c.s << closebrace;
103 if (precedence() <= level)
104 c.s << ')' << closebrace;
105}
106
107void power::do_print_dflt(const print_dflt & c, unsigned level) const
108{
109 if (exponent.is_equal(_ex1_2)) {
110
111 // Square roots are printed in a special way
112 c.s << "sqrt(";
113 basis.print(c);
114 c.s << ')';
115
116 } else
117 print_power(c, "^", "", "", level);
118}
119
120void power::do_print_latex(const print_latex & c, unsigned level) const
121{
122 if (is_exactly_a<numeric>(exponent) && ex_to<numeric>(exponent).is_negative()) {
123
124 // Powers with negative numeric exponents are printed as fractions
125 c.s << "\\frac{1}{";
126 power(basis, -exponent).eval().print(c);
127 c.s << '}';
128
129 } else if (exponent.is_equal(_ex1_2)) {
130
131 // Square roots are printed in a special way
132 c.s << "\\sqrt{";
133 basis.print(c);
134 c.s << '}';
135
136 } else
137 print_power(c, "^", "{", "}", level);
138}
139
140static void print_sym_pow(const print_context & c, const symbol &x, int exp)
141{
142 // Optimal output of integer powers of symbols to aid compiler CSE.
143 // C.f. ISO/IEC 14882:2011, section 1.9 [intro.execution], paragraph 15
144 // to learn why such a parenthesation is really necessary.
145 if (exp == 1) {
146 x.print(c);
147 } else if (exp == 2) {
148 x.print(c);
149 c.s << "*";
150 x.print(c);
151 } else if (exp & 1) {
152 x.print(c);
153 c.s << "*";
154 print_sym_pow(c, x, exp-1);
155 } else {
156 c.s << "(";
157 print_sym_pow(c, x, exp >> 1);
158 c.s << ")*(";
159 print_sym_pow(c, x, exp >> 1);
160 c.s << ")";
161 }
162}
163
164void power::do_print_csrc_cl_N(const print_csrc_cl_N& c, unsigned level) const
165{
166 if (exponent.is_equal(_ex_1)) {
167 c.s << "recip(";
168 basis.print(c);
169 c.s << ')';
170 return;
171 }
172 c.s << "expt(";
173 basis.print(c);
174 c.s << ", ";
176 c.s << ')';
177}
178
179void power::do_print_csrc(const print_csrc & c, unsigned level) const
180{
181 // Integer powers of symbols are printed in a special, optimized way
183 (is_a<symbol>(basis) || is_a<constant>(basis))) {
184 int exp = ex_to<numeric>(exponent).to_int();
185 if (exp > 0)
186 c.s << '(';
187 else {
188 exp = -exp;
189 c.s << "1.0/(";
190 }
191 print_sym_pow(c, ex_to<symbol>(basis), exp);
192 c.s << ')';
193
194 // <expr>^-1 is printed as "1.0/<expr>" or with the recip() function of CLN
195 } else if (exponent.is_equal(_ex_1)) {
196 c.s << "1.0/(";
197 basis.print(c);
198 c.s << ')';
199
200 // Otherwise, use the pow() function
201 } else {
202 c.s << "pow(";
203 basis.print(c);
204 c.s << ',';
206 c.s << ')';
207 }
208}
209
210void power::do_print_python(const print_python & c, unsigned level) const
211{
212 print_power(c, "**", "", "", level);
213}
214
215void power::do_print_python_repr(const print_python_repr & c, unsigned level) const
216{
217 c.s << class_name() << '(';
218 basis.print(c);
219 c.s << ',';
221 c.s << ')';
222}
223
224bool power::info(unsigned inf) const
225{
226 switch (inf) {
235 case info_flags::real:
238 return (flags & status_flags::expanded);
246 return true;
248 return false;
252 return true;
253 } else {
256 return false;
257 }
258 }
259 }
260 return inherited::info(inf);
261}
262
263size_t power::nops() const
264{
265 return 2;
266}
267
268ex power::op(size_t i) const
269{
270 GINAC_ASSERT(i<2);
271
272 return i==0 ? basis : exponent;
273}
274
276{
277 const ex &mapped_basis = f(basis);
278 const ex &mapped_exponent = f(exponent);
279
280 if (!are_ex_trivially_equal(basis, mapped_basis)
281 || !are_ex_trivially_equal(exponent, mapped_exponent))
282 return dynallocate<power>(mapped_basis, mapped_exponent);
283 else
284 return *this;
285}
286
287bool power::is_polynomial(const ex & var) const
288{
289 if (basis.is_polynomial(var)) {
290 if (basis.has(var))
291 // basis is non-constant polynomial in var
293 else
294 // basis is constant in var
295 return !exponent.has(var);
296 }
297 // basis is a non-polynomial function of var
298 return false;
299}
300
301int power::degree(const ex & s) const
302{
303 if (is_equal(ex_to<basic>(s)))
304 return 1;
305 else if (is_exactly_a<numeric>(exponent) && ex_to<numeric>(exponent).is_integer()) {
306 if (basis.is_equal(s))
307 return ex_to<numeric>(exponent).to_int();
308 else
309 return basis.degree(s) * ex_to<numeric>(exponent).to_int();
310 } else if (basis.has(s))
311 throw(std::runtime_error("power::degree(): undefined degree because of non-integer exponent"));
312 else
313 return 0;
314}
315
316int power::ldegree(const ex & s) const
317{
318 if (is_equal(ex_to<basic>(s)))
319 return 1;
320 else if (is_exactly_a<numeric>(exponent) && ex_to<numeric>(exponent).is_integer()) {
321 if (basis.is_equal(s))
322 return ex_to<numeric>(exponent).to_int();
323 else
324 return basis.ldegree(s) * ex_to<numeric>(exponent).to_int();
325 } else if (basis.has(s))
326 throw(std::runtime_error("power::ldegree(): undefined degree because of non-integer exponent"));
327 else
328 return 0;
329}
330
331ex power::coeff(const ex & s, int n) const
332{
333 if (is_equal(ex_to<basic>(s)))
334 return n==1 ? _ex1 : _ex0;
335 else if (!basis.is_equal(s)) {
336 // basis not equal to s
337 if (n == 0)
338 return *this;
339 else
340 return _ex0;
341 } else {
342 // basis equal to s
343 if (is_exactly_a<numeric>(exponent) && ex_to<numeric>(exponent).is_integer()) {
344 // integer exponent
345 int int_exp = ex_to<numeric>(exponent).to_int();
346 if (n == int_exp)
347 return _ex1;
348 else
349 return _ex0;
350 } else {
351 // non-integer exponents are treated as zero
352 if (n == 0)
353 return *this;
354 else
355 return _ex0;
356 }
357 }
358}
359
375{
377 return *this;
378
379 const numeric *num_basis = nullptr;
380 const numeric *num_exponent = nullptr;
381
382 if (is_exactly_a<numeric>(basis)) {
383 num_basis = &ex_to<numeric>(basis);
384 }
385 if (is_exactly_a<numeric>(exponent)) {
386 num_exponent = &ex_to<numeric>(exponent);
387 }
388
389 // ^(x,0) -> 1 (0^0 also handled here)
390 if (exponent.is_zero()) {
391 if (basis.is_zero())
392 throw (std::domain_error("power::eval(): pow(0,0) is undefined"));
393 else
394 return _ex1;
395 }
396
397 // ^(x,1) -> x
399 return basis;
400
401 // ^(0,c1) -> 0 or exception (depending on real value of c1)
402 if (basis.is_zero() && num_exponent) {
403 if ((num_exponent->real()).is_zero())
404 throw (std::domain_error("power::eval(): pow(0,I) is undefined"));
405 else if ((num_exponent->real()).is_negative())
406 throw (pole_error("power::eval(): division by zero",1));
407 else
408 return _ex0;
409 }
410
411 // ^(1,x) -> 1
412 if (basis.is_equal(_ex1))
413 return _ex1;
414
415 // power of a function calculated by separate rules defined for this function
416 if (is_exactly_a<function>(basis))
417 return ex_to<function>(basis).power(exponent);
418
419 // Turn (x^c)^d into x^(c*d) in the case that x is positive and c is real.
420 if (is_exactly_a<power>(basis) && basis.op(0).info(info_flags::positive) && basis.op(1).info(info_flags::real))
421 return dynallocate<power>(basis.op(0), basis.op(1) * exponent);
422
423 if ( num_exponent ) {
424
425 // ^(c1,c2) -> c1^c2 (c1, c2 numeric(),
426 // except if c1,c2 are rational, but c1^c2 is not)
427 if ( num_basis ) {
428 const bool basis_is_crational = num_basis->is_crational();
429 const bool exponent_is_crational = num_exponent->is_crational();
430 if (!basis_is_crational || !exponent_is_crational) {
431 // return a plain float
432 return dynallocate<numeric>(num_basis->power(*num_exponent));
433 }
434
435 const numeric res = num_basis->power(*num_exponent);
436 if (res.is_crational()) {
437 return res;
438 }
439 GINAC_ASSERT(!num_exponent->is_integer()); // has been handled by now
440
441 // ^(c1,n/m) -> *(c1^q,c1^(n/m-q)), 0<(n/m-q)<1, q integer
442 if (basis_is_crational && exponent_is_crational
443 && num_exponent->is_real()
444 && !num_exponent->is_integer()) {
445 const numeric n = num_exponent->numer();
446 const numeric m = num_exponent->denom();
447 numeric r;
448 numeric q = iquo(n, m, r);
449 if (r.is_negative()) {
450 r += m;
451 --q;
452 }
453 if (q.is_zero()) { // the exponent was in the allowed range 0<(n/m)<1
454 if (num_basis->is_rational() && !num_basis->is_integer()) {
455 // try it for numerator and denominator separately, in order to
456 // partially simplify things like (5/8)^(1/3) -> 1/2*5^(1/3)
457 const numeric bnum = num_basis->numer();
458 const numeric bden = num_basis->denom();
459 const numeric res_bnum = bnum.power(*num_exponent);
460 const numeric res_bden = bden.power(*num_exponent);
461 if (res_bnum.is_integer())
462 return dynallocate<mul>(dynallocate<power>(bden,-*num_exponent),res_bnum).setflag(status_flags::evaluated);
463 if (res_bden.is_integer())
464 return dynallocate<mul>(dynallocate<power>(bnum,*num_exponent),res_bden.inverse()).setflag(status_flags::evaluated);
465 }
466 return this->hold();
467 } else {
468 // assemble resulting product, but allowing for a re-evaluation,
469 // because otherwise we'll end up with something like
470 // (7/8)^(4/3) -> 7/8*(1/2*7^(1/3))
471 // instead of 7/16*7^(1/3).
472 return pow(basis, r.div(m)) * pow(basis, q);
473 }
474 }
475 }
476
477 // ^(^(x,c1),c2) -> ^(x,c1*c2)
478 // (c1, c2 numeric(), c2 integer or -1 < c1 <= 1 or (c1=-1 and c2>0),
479 // case c1==1 should not happen, see below!)
480 if (is_exactly_a<power>(basis)) {
481 const power & sub_power = ex_to<power>(basis);
482 const ex & sub_basis = sub_power.basis;
483 const ex & sub_exponent = sub_power.exponent;
484 if (is_exactly_a<numeric>(sub_exponent)) {
485 const numeric & num_sub_exponent = ex_to<numeric>(sub_exponent);
486 GINAC_ASSERT(num_sub_exponent!=numeric(1));
487 if (num_exponent->is_integer() || (abs(num_sub_exponent) - (*_num1_p)).is_negative() ||
488 (num_sub_exponent == *_num_1_p && num_exponent->is_positive())) {
489 return dynallocate<power>(sub_basis, num_sub_exponent.mul(*num_exponent));
490 }
491 }
492 }
493
494 // ^(*(x,y,z),c1) -> *(x^c1,y^c1,z^c1) (c1 integer)
495 if (num_exponent->is_integer() && is_exactly_a<mul>(basis)) {
496 return expand_mul(ex_to<mul>(basis), *num_exponent, false);
497 }
498
499 // (2*x + 6*y)^(-4) -> 1/16*(x + 3*y)^(-4)
500 if (num_exponent->is_integer() && is_exactly_a<add>(basis)) {
501 numeric icont = basis.integer_content();
502 const numeric lead_coeff =
503 ex_to<numeric>(ex_to<add>(basis).seq.begin()->coeff).div(icont);
504
505 const bool canonicalizable = lead_coeff.is_integer();
506 const bool unit_normal = lead_coeff.is_pos_integer();
507 if (canonicalizable && (! unit_normal))
508 icont = icont.mul(*_num_1_p);
509
510 if (canonicalizable && (icont != *_num1_p)) {
511 const add& addref = ex_to<add>(basis);
512 add & addp = dynallocate<add>(addref);
514 addp.overall_coeff = ex_to<numeric>(addp.overall_coeff).div_dyn(icont);
515 for (auto & i : addp.seq)
516 i.coeff = ex_to<numeric>(i.coeff).div_dyn(icont);
517
518 const numeric c = icont.power(*num_exponent);
519 if (likely(c != *_num1_p))
520 return dynallocate<mul>(dynallocate<power>(addp, *num_exponent), c);
521 else
522 return dynallocate<power>(addp, *num_exponent);
523 }
524 }
525
526 // ^(*(...,x;c1),c2) -> *(^(*(...,x;1),c2),c1^c2) (c1, c2 numeric(), c1>0)
527 // ^(*(...,x;c1),c2) -> *(^(*(...,x;-1),c2),(-c1)^c2) (c1, c2 numeric(), c1<0)
528 if (is_exactly_a<mul>(basis)) {
529 GINAC_ASSERT(!num_exponent->is_integer()); // should have been handled above
530 const mul & mulref = ex_to<mul>(basis);
531 if (!mulref.overall_coeff.is_equal(_ex1)) {
532 const numeric & num_coeff = ex_to<numeric>(mulref.overall_coeff);
533 if (num_coeff.is_real()) {
534 if (num_coeff.is_positive()) {
535 mul & mulp = dynallocate<mul>(mulref);
536 mulp.overall_coeff = _ex1;
538 return dynallocate<mul>(dynallocate<power>(mulp, exponent),
539 dynallocate<power>(num_coeff, *num_exponent));
540 } else {
541 GINAC_ASSERT(num_coeff.compare(*_num0_p)<0);
542 if (!num_coeff.is_equal(*_num_1_p)) {
543 mul & mulp = dynallocate<mul>(mulref);
544 mulp.overall_coeff = _ex_1;
546 return dynallocate<mul>(dynallocate<power>(mulp, exponent),
547 dynallocate<power>(abs(num_coeff), *num_exponent));
548 }
549 }
550 }
551 }
552 }
553
554 // ^(nc,c1) -> ncmul(nc,nc,...) (c1 positive integer, unless nc is a matrix)
555 if (num_exponent->is_pos_integer() &&
557 !is_a<matrix>(basis)) {
558 return ncmul(exvector(num_exponent->to_int(), basis));
559 }
560 }
561
562 return this->hold();
563}
564
566{
567 ex ebasis = basis.evalf();
568 ex eexponent;
569
570 if (!is_exactly_a<numeric>(exponent))
571 eexponent = exponent.evalf();
572 else
573 eexponent = exponent;
574
575 return dynallocate<power>(ebasis, eexponent);
576}
577
579{
580 const ex ebasis = basis.evalm();
581 const ex eexponent = exponent.evalm();
582 if (is_a<matrix>(ebasis)) {
583 if (is_exactly_a<numeric>(eexponent)) {
584 return dynallocate<matrix>(ex_to<matrix>(ebasis).pow(eexponent));
585 }
586 }
587 return dynallocate<power>(ebasis, eexponent);
588}
589
590bool power::has(const ex & other, unsigned options) const
591{
593 return basic::has(other, options);
594 if (!is_a<power>(other))
595 return basic::has(other, options);
597 !other.op(1).info(info_flags::integer))
598 return basic::has(other, options);
600 other.op(1).info(info_flags::posint) &&
601 ex_to<numeric>(exponent) > ex_to<numeric>(other.op(1)) &&
602 basis.match(other.op(0)))
603 return true;
605 other.op(1).info(info_flags::negint) &&
606 ex_to<numeric>(exponent) < ex_to<numeric>(other.op(1)) &&
607 basis.match(other.op(0)))
608 return true;
609 return basic::has(other, options);
610}
611
612// from mul.cpp
613extern bool tryfactsubs(const ex &, const ex &, int &, exmap&);
614
615ex power::subs(const exmap & m, unsigned options) const
616{
617 const ex &subsed_basis = basis.subs(m, options);
618 const ex &subsed_exponent = exponent.subs(m, options);
619
620 if (!are_ex_trivially_equal(basis, subsed_basis)
621 || !are_ex_trivially_equal(exponent, subsed_exponent))
622 return dynallocate<power>(subsed_basis, subsed_exponent);
623
625 return subs_one_level(m, options);
626
627 for (auto & it : m) {
628 int nummatches = std::numeric_limits<int>::max();
629 exmap repls;
630 if (tryfactsubs(*this, it.first, nummatches, repls)) {
631 ex anum = it.second.subs(repls, subs_options::no_pattern);
632 ex aden = it.first.subs(repls, subs_options::no_pattern);
633 ex result = (*this) * pow(anum/aden, nummatches);
634 return (ex_to<basic>(result)).subs_one_level(m, options);
635 }
636 }
637
638 return subs_one_level(m, options);
639}
640
642{
643 return inherited::eval_ncmul(v);
644}
645
647{
648 // conjugate(pow(x,y))==pow(conjugate(x),conjugate(y)) unless on the
649 // branch cut which runs along the negative real axis.
651 ex newexponent = exponent.conjugate();
652 if (are_ex_trivially_equal(exponent, newexponent)) {
653 return *this;
654 }
655 return dynallocate<power>(basis, newexponent);
656 }
658 ex newbasis = basis.conjugate();
659 if (are_ex_trivially_equal(basis, newbasis)) {
660 return *this;
661 }
662 return dynallocate<power>(newbasis, exponent);
663 }
664 return conjugate_function(*this).hold();
665}
666
668{
669 // basis == a+I*b, exponent == c+I*d
670 const ex a = basis.real_part();
671 const ex c = exponent.real_part();
672 if (basis.is_equal(a) && exponent.is_equal(c) &&
674 // Re(a^c)
675 return *this;
676 }
677
678 const ex b = basis.imag_part();
680 // Re((a+I*b)^c) w/ c ∈ ℤ
681 long N = ex_to<numeric>(c).to_long();
682 // Use real terms in Binomial expansion to construct
683 // Re(expand(pow(a+I*b, N))).
684 long NN = N > 0 ? N : -N;
685 ex numer = N > 0 ? _ex1 : pow(pow(a,2) + pow(b,2), NN);
686 ex result = 0;
687 for (long n = 0; n <= NN; n += 2) {
688 ex term = binomial(NN, n) * pow(a, NN-n) * pow(b, n) / numer;
689 if (n % 4 == 0) {
690 result += term; // sign: I^n w/ n == 4*m
691 } else {
692 result -= term; // sign: I^n w/ n == 4*m+2
693 }
694 }
695 return result;
696 }
697
698 // Re((a+I*b)^(c+I*d))
699 const ex d = exponent.imag_part();
700 return pow(abs(basis),c) * exp(-d*atan2(b,a)) * cos(c*atan2(b,a)+d*log(abs(basis)));
701}
702
704{
705 // basis == a+I*b, exponent == c+I*d
706 const ex a = basis.real_part();
707 const ex c = exponent.real_part();
708 if (basis.is_equal(a) && exponent.is_equal(c) &&
710 // Im(a^c)
711 return 0;
712 }
713
714 const ex b = basis.imag_part();
716 // Im((a+I*b)^c) w/ c ∈ ℤ
717 long N = ex_to<numeric>(c).to_long();
718 // Use imaginary terms in Binomial expansion to construct
719 // Im(expand(pow(a+I*b, N))).
720 long p = N > 0 ? 1 : 3; // modulus for positive sign
721 long NN = N > 0 ? N : -N;
722 ex numer = N > 0 ? _ex1 : pow(pow(a,2) + pow(b,2), NN);
723 ex result = 0;
724 for (long n = 1; n <= NN; n += 2) {
725 ex term = binomial(NN, n) * pow(a, NN-n) * pow(b, n) / numer;
726 if (n % 4 == p) {
727 result += term; // sign: I^n w/ n == 4*m+p
728 } else {
729 result -= term; // sign: I^n w/ n == 4*m+2+p
730 }
731 }
732 return result;
733 }
734
735 // Im((a+I*b)^(c+I*d))
736 const ex d = exponent.imag_part();
737 return pow(abs(basis),c) * exp(-d*atan2(b,a)) * sin(c*atan2(b,a)+d*log(abs(basis)));
738}
739
740// protected
741
745{
746 if (is_a<numeric>(exponent)) {
747 // D(b^r) = r * b^(r-1) * D(b) (faster than the formula below)
748 const epvector newseq = {expair(basis, exponent - _ex1), expair(basis.diff(s), _ex1)};
749 return dynallocate<mul>(std::move(newseq), exponent);
750 } else {
751 // D(b^e) = b^e * (D(e)*ln(b) + e*D(b)/b)
752 return *this * (exponent.diff(s)*log(basis) + exponent*basis.diff(s)*pow(basis, _ex_1));
753 }
754}
755
756int power::compare_same_type(const basic & other) const
757{
758 GINAC_ASSERT(is_exactly_a<power>(other));
759 const power &o = static_cast<const power &>(other);
760
761 int cmpval = basis.compare(o.basis);
762 if (cmpval)
763 return cmpval;
764 else
765 return exponent.compare(o.exponent);
766}
767
768unsigned power::return_type() const
769{
770 return basis.return_type();
771}
772
774{
775 return basis.return_type_tinfo();
776}
777
778ex power::expand(unsigned options) const
779{
780 if (is_a<symbol>(basis) && exponent.info(info_flags::integer)) {
781 // A special case worth optimizing.
783 return *this;
784 }
785
786 // (x*p)^c -> x^c * p^c, if p>0
787 // makes sense before expanding the basis
788 if (is_exactly_a<mul>(basis) && !basis.info(info_flags::indefinite)) {
789 const mul &m = ex_to<mul>(basis);
790 exvector prodseq;
791 epvector powseq;
792 prodseq.reserve(m.seq.size() + 1);
793 powseq.reserve(m.seq.size() + 1);
794 bool possign = true;
795
796 // search for positive/negative factors
797 for (auto & cit : m.seq) {
798 ex e=m.recombine_pair_to_ex(cit);
800 prodseq.push_back(pow(e, exponent).expand(options));
801 else if (e.info(info_flags::negative)) {
802 prodseq.push_back(pow(-e, exponent).expand(options));
803 possign = !possign;
804 } else
805 powseq.push_back(cit);
806 }
807
808 // take care on the numeric coefficient
809 ex coeff=(possign? _ex1 : _ex_1);
810 if (m.overall_coeff.info(info_flags::positive) && m.overall_coeff != _ex1)
811 prodseq.push_back(pow(m.overall_coeff, exponent));
812 else if (m.overall_coeff.info(info_flags::negative) && m.overall_coeff != _ex_1) {
813 prodseq.push_back(pow(-m.overall_coeff, exponent));
814 coeff = -coeff;
815 } else
816 coeff *= m.overall_coeff;
817
818 // If positive/negative factors are found, then extract them.
819 // In either case we set a flag to avoid the second run on a part
820 // which does not have positive/negative terms.
821 if (prodseq.size() > 0) {
822 ex newbasis = dynallocate<mul>(std::move(powseq), coeff);
823 ex_to<basic>(newbasis).setflag(status_flags::purely_indefinite);
824 return dynallocate<mul>(std::move(prodseq)) * pow(newbasis, exponent);
825 } else
826 ex_to<basic>(basis).setflag(status_flags::purely_indefinite);
827 }
828
829 const ex expanded_basis = basis.expand(options);
830 const ex expanded_exponent = exponent.expand(options);
831
832 // x^(a+b) -> x^a * x^b
833 if (is_exactly_a<add>(expanded_exponent)) {
834 const add &a = ex_to<add>(expanded_exponent);
835 exvector distrseq;
836 distrseq.reserve(a.seq.size() + 1);
837 for (auto & cit : a.seq) {
838 distrseq.push_back(pow(expanded_basis, a.recombine_pair_to_ex(cit)));
839 }
840
841 // Make sure that e.g. (x+y)^(2+a) expands the (x+y)^2 factor
842 if (ex_to<numeric>(a.overall_coeff).is_integer()) {
843 const numeric &num_exponent = ex_to<numeric>(a.overall_coeff);
844 long int_exponent = num_exponent.to_int();
845 if (int_exponent > 0 && is_exactly_a<add>(expanded_basis))
846 distrseq.push_back(expand_add(ex_to<add>(expanded_basis), int_exponent, options));
847 else
848 distrseq.push_back(pow(expanded_basis, a.overall_coeff));
849 } else
850 distrseq.push_back(pow(expanded_basis, a.overall_coeff));
851
852 // Make sure that e.g. (x+y)^(1+a) -> x*(x+y)^a + y*(x+y)^a
853 ex r = dynallocate<mul>(distrseq);
854 return r.expand(options);
855 }
856
857 if (!is_exactly_a<numeric>(expanded_exponent) ||
858 !ex_to<numeric>(expanded_exponent).is_integer()) {
859 if (are_ex_trivially_equal(basis,expanded_basis) && are_ex_trivially_equal(exponent,expanded_exponent)) {
860 return this->hold();
861 } else {
862 return dynallocate<power>(expanded_basis, expanded_exponent).setflag(options == 0 ? status_flags::expanded : 0);
863 }
864 }
865
866 // integer numeric exponent
867 const numeric & num_exponent = ex_to<numeric>(expanded_exponent);
868 long int_exponent = num_exponent.to_long();
869
870 // (x+y)^n, n>0
871 if (int_exponent > 0 && is_exactly_a<add>(expanded_basis))
872 return expand_add(ex_to<add>(expanded_basis), int_exponent, options);
873
874 // (x*y)^n -> x^n * y^n
875 if (is_exactly_a<mul>(expanded_basis))
876 return expand_mul(ex_to<mul>(expanded_basis), num_exponent, options, true);
877
878 // cannot expand further
879 if (are_ex_trivially_equal(basis,expanded_basis) && are_ex_trivially_equal(exponent,expanded_exponent))
880 return this->hold();
881 else
882 return dynallocate<power>(expanded_basis, expanded_exponent).setflag(options == 0 ? status_flags::expanded : 0);
883}
884
886// new virtual functions which can be overridden by derived classes
888
889// none
890
892// non-virtual functions in this class
894
897ex power::expand_add(const add & a, long n, unsigned options)
898{
899 // The special case power(+(x,...y;x),2) can be optimized better.
900 if (n==2)
901 return expand_add_2(a, options);
902
903 // method:
904 //
905 // Consider base as the sum of all symbolic terms and the overall numeric
906 // coefficient and apply the binomial theorem:
907 // S = power(+(x,...,z;c),n)
908 // = power(+(+(x,...,z;0);c),n)
909 // = sum(binomial(n,k)*power(+(x,...,z;0),k)*c^(n-k), k=1..n) + c^n
910 // Then, apply the multinomial theorem to expand all power(+(x,...,z;0),k):
911 // The multinomial theorem is computed by an outer loop over all
912 // partitions of the exponent and an inner loop over all compositions of
913 // that partition. This method makes the expansion a combinatorial
914 // problem and allows us to directly construct the expanded sum and also
915 // to re-use the multinomial coefficients (since they depend only on the
916 // partition, not on the composition).
917 //
918 // multinomial power(+(x,y,z;0),3) example:
919 // partition : compositions : multinomial coefficient
920 // [0,0,3] : [3,0,0],[0,3,0],[0,0,3] : 3!/(3!*0!*0!) = 1
921 // [0,1,2] : [2,1,0],[1,2,0],[2,0,1],... : 3!/(2!*1!*0!) = 3
922 // [1,1,1] : [1,1,1] : 3!/(1!*1!*1!) = 6
923 // => (x + y + z)^3 =
924 // x^3 + y^3 + z^3
925 // + 3*x^2*y + 3*x*y^2 + 3*y^2*z + 3*y*z^2 + 3*x*z^2 + 3*x^2*z
926 // + 6*x*y*z
927 //
928 // multinomial power(+(x,y,z;0),4) example:
929 // partition : compositions : multinomial coefficient
930 // [0,0,4] : [4,0,0],[0,4,0],[0,0,4] : 4!/(4!*0!*0!) = 1
931 // [0,1,3] : [3,1,0],[1,3,0],[3,0,1],... : 4!/(3!*1!*0!) = 4
932 // [0,2,2] : [2,2,0],[2,0,2],[0,2,2] : 4!/(2!*2!*0!) = 6
933 // [1,1,2] : [2,1,1],[1,2,1],[1,1,2] : 4!/(2!*1!*1!) = 12
934 // (no [1,1,1,1] partition since it has too many parts)
935 // => (x + y + z)^4 =
936 // x^4 + y^4 + z^4
937 // + 4*x^3*y + 4*x*y^3 + 4*y^3*z + 4*y*z^3 + 4*x*z^3 + 4*x^3*z
938 // + 6*x^2*y^2 + 6*y^2*z^2 + 6*x^2*z^2
939 // + 12*x^2*y*z + 12*x*y^2*z + 12*x*y*z^2
940 //
941 // Summary:
942 // r = 0
943 // for k from 0 to n:
944 // f = c^(n-k)*binomial(n,k)
945 // for p in all partitions of n with m parts (including zero parts):
946 // h = f * multinomial coefficient of p
947 // for c in all compositions of p:
948 // t = 1
949 // for e in all elements of c:
950 // t = t * a[e]^e
951 // r = r + h*t
952 // return r
953
954 epvector result;
955 // The number of terms will be the number of combinatorial compositions,
956 // i.e. the number of unordered arrangements of m nonnegative integers
957 // which sum up to n. It is frequently written as C_n(m) and directly
958 // related with binomial coefficients: binomial(n+m-1,m-1).
959 size_t result_size = binomial(numeric(n+a.nops()-1), numeric(a.nops()-1)).to_long();
960 if (!a.overall_coeff.is_zero()) {
961 // the result's overall_coeff is one of the terms
962 --result_size;
963 }
964 result.reserve(result_size);
965
966 // Iterate over all terms in binomial expansion of
967 // S = power(+(x,...,z;c),n)
968 // = sum(binomial(n,k)*power(+(x,...,z;0),k)*c^(n-k), k=1..n) + c^n
969 for (int k = 1; k <= n; ++k) {
970 numeric binomial_coefficient; // binomial(n,k)*c^(n-k)
971 if (a.overall_coeff.is_zero()) {
972 // degenerate case with zero overall_coeff:
973 // apply multinomial theorem directly to power(+(x,...z;0),n)
974 binomial_coefficient = 1;
975 if (k < n) {
976 continue;
977 }
978 } else {
979 binomial_coefficient = binomial(numeric(n), numeric(k)) * pow(ex_to<numeric>(a.overall_coeff), numeric(n-k));
980 }
981
982 // Multinomial expansion of power(+(x,...,z;0),k)*c^(n-k):
983 // Iterate over all partitions of k with exactly as many parts as
984 // there are symbolic terms in the basis (including zero parts).
985 partition_with_zero_parts_generator partitions(k, a.seq.size());
986 do {
987 const std::vector<unsigned>& partition = partitions.get();
988 // All monomials of this partition have the same number of terms and the same coefficient.
989 const unsigned msize = std::count_if(partition.begin(), partition.end(), [](int i) { return i > 0; });
990 const numeric coeff = multinomial_coefficient(partition) * binomial_coefficient;
991
992 // Iterate over all compositions of the current partition.
993 composition_generator compositions(partition);
994 do {
995 const std::vector<unsigned>& exponent = compositions.get();
996 epvector monomial;
997 monomial.reserve(msize);
999 for (unsigned i = 0; i < exponent.size(); ++i) {
1000 const ex & r = a.seq[i].rest;
1001 GINAC_ASSERT(!is_exactly_a<add>(r));
1002 GINAC_ASSERT(!is_exactly_a<power>(r) ||
1003 !is_exactly_a<numeric>(ex_to<power>(r).exponent) ||
1004 !ex_to<numeric>(ex_to<power>(r).exponent).is_pos_integer() ||
1005 !is_exactly_a<add>(ex_to<power>(r).basis) ||
1006 !is_exactly_a<mul>(ex_to<power>(r).basis) ||
1007 !is_exactly_a<power>(ex_to<power>(r).basis));
1008 GINAC_ASSERT(is_exactly_a<numeric>(a.seq[i].coeff));
1009 const numeric & c = ex_to<numeric>(a.seq[i].coeff);
1010 if (exponent[i] == 0) {
1011 // optimize away
1012 } else if (exponent[i] == 1) {
1013 // optimized
1014 monomial.emplace_back(expair(r, _ex1));
1015 if (c != *_num1_p)
1016 factor = factor.mul(c);
1017 } else { // general case exponent[i] > 1
1018 monomial.emplace_back(expair(r, exponent[i]));
1019 if (c != *_num1_p)
1020 factor = factor.mul(c.power(exponent[i]));
1021 }
1022 }
1023 result.emplace_back(expair(mul(std::move(monomial)).expand(options), factor));
1024 } while (compositions.next());
1025 } while (partitions.next());
1026 }
1027
1028 GINAC_ASSERT(result.size() == result_size);
1029 if (a.overall_coeff.is_zero()) {
1030 return dynallocate<add>(std::move(result)).setflag(status_flags::expanded);
1031 } else {
1032 return dynallocate<add>(std::move(result), ex_to<numeric>(a.overall_coeff).power(n)).setflag(status_flags::expanded);
1033 }
1034}
1035
1036
1039ex power::expand_add_2(const add & a, unsigned options)
1040{
1041 epvector result;
1042 size_t result_size = (a.nops() * (a.nops()+1)) / 2;
1043 if (!a.overall_coeff.is_zero()) {
1044 // the result's overall_coeff is one of the terms
1045 --result_size;
1046 }
1047 result.reserve(result_size);
1048
1049 auto last = a.seq.end();
1050
1051 // power(+(x,...,z;c),2)=power(+(x,...,z;0),2)+2*c*+(x,...,z;0)+c*c
1052 // first part: ignore overall_coeff and expand other terms
1053 for (auto cit0=a.seq.begin(); cit0!=last; ++cit0) {
1054 const ex & r = cit0->rest;
1055 const ex & c = cit0->coeff;
1056
1057 GINAC_ASSERT(!is_exactly_a<add>(r));
1058 GINAC_ASSERT(!is_exactly_a<power>(r) ||
1059 !is_exactly_a<numeric>(ex_to<power>(r).exponent) ||
1060 !ex_to<numeric>(ex_to<power>(r).exponent).is_pos_integer() ||
1061 !is_exactly_a<add>(ex_to<power>(r).basis) ||
1062 !is_exactly_a<mul>(ex_to<power>(r).basis) ||
1063 !is_exactly_a<power>(ex_to<power>(r).basis));
1064
1065 if (c.is_equal(_ex1)) {
1066 if (is_exactly_a<mul>(r)) {
1067 result.emplace_back(expair(expand_mul(ex_to<mul>(r), *_num2_p, options, true),
1068 _ex1));
1069 } else {
1070 result.emplace_back(expair(dynallocate<power>(r, _ex2),
1071 _ex1));
1072 }
1073 } else {
1074 if (is_exactly_a<mul>(r)) {
1075 result.emplace_back(expair(expand_mul(ex_to<mul>(r), *_num2_p, options, true),
1076 ex_to<numeric>(c).power_dyn(*_num2_p)));
1077 } else {
1078 result.emplace_back(expair(dynallocate<power>(r, _ex2),
1079 ex_to<numeric>(c).power_dyn(*_num2_p)));
1080 }
1081 }
1082
1083 for (auto cit1=cit0+1; cit1!=last; ++cit1) {
1084 const ex & r1 = cit1->rest;
1085 const ex & c1 = cit1->coeff;
1086 result.emplace_back(expair(mul(r,r1).expand(options),
1087 _num2_p->mul(ex_to<numeric>(c)).mul_dyn(ex_to<numeric>(c1))));
1088 }
1089 }
1090
1091 // second part: add terms coming from overall_coeff (if != 0)
1092 if (!a.overall_coeff.is_zero()) {
1093 for (auto & i : a.seq)
1094 result.push_back(a.combine_pair_with_coeff_to_pair(i, ex_to<numeric>(a.overall_coeff).mul_dyn(*_num2_p)));
1095 }
1096
1097 GINAC_ASSERT(result.size() == result_size);
1098
1099 if (a.overall_coeff.is_zero()) {
1100 return dynallocate<add>(std::move(result)).setflag(status_flags::expanded);
1101 } else {
1102 return dynallocate<add>(std::move(result), ex_to<numeric>(a.overall_coeff).power(2)).setflag(status_flags::expanded);
1103 }
1104}
1105
1108ex power::expand_mul(const mul & m, const numeric & n, unsigned options, bool from_expand)
1109{
1110 GINAC_ASSERT(n.is_integer());
1111
1112 if (n.is_zero()) {
1113 return _ex1;
1114 }
1115
1116 // do not bother to rename indices if there are no any.
1120 // Leave it to multiplication since dummy indices have to be renamed
1122 (get_all_dummy_indices(m).size() > 0) && n.is_positive()) {
1123 ex result = m;
1125 sort(va.begin(), va.end(), ex_is_less());
1126
1127 for (int i=1; i < n.to_int(); i++)
1128 result *= rename_dummy_indices_uniquely(va, m);
1129 return result;
1130 }
1131
1132 epvector distrseq;
1133 distrseq.reserve(m.seq.size());
1134 bool need_reexpand = false;
1135
1136 for (auto & cit : m.seq) {
1137 expair p = m.combine_pair_with_coeff_to_pair(cit, n);
1138 if (from_expand && is_exactly_a<add>(cit.rest) && ex_to<numeric>(p.coeff).is_pos_integer()) {
1139 // this happens when e.g. (a+b)^(1/2) gets squared and
1140 // the resulting product needs to be reexpanded
1141 need_reexpand = true;
1142 }
1143 distrseq.push_back(p);
1144 }
1145
1146 const mul & result = dynallocate<mul>(std::move(distrseq), ex_to<numeric>(m.overall_coeff).power_dyn(n));
1147 if (need_reexpand)
1148 return ex(result).expand(options);
1149 if (from_expand)
1150 return result.setflag(status_flags::expanded);
1151 return result;
1152}
1153
1155
1156} // namespace GiNaC
Interface to GiNaC's sums of expressions.
Archiving of GiNaC expressions.
#define GINAC_ASSERT(X)
Assertion macro for checking invariances.
Definition: assertion.h:33
Sum of expressions.
Definition: add.h:32
ex recombine_pair_to_ex(const expair &p) const override
Form an ex out of an expair, using the corresponding semantics.
Definition: add.cpp:564
expair combine_pair_with_coeff_to_pair(const expair &p, const ex &c) const override
Definition: add.cpp:550
This class stores all properties needed to record/retrieve the state of one object of class basic (or...
Definition: archive.h:49
This class is the ABC (abstract base class) of GiNaC's class hierarchy.
Definition: basic.h:105
const basic & clearflag(unsigned f) const
Clear some status_flags.
Definition: basic.h:291
const basic & setflag(unsigned f) const
Set some status_flags.
Definition: basic.h:288
virtual bool has(const ex &other, unsigned options=0) const
Test for occurrence of a pattern.
Definition: basic.cpp:280
friend class ex
Definition: basic.h:108
unsigned flags
of type status_flags
Definition: basic.h:302
ex subs_one_level(const exmap &m, unsigned options) const
Helper function for subs().
Definition: basic.cpp:585
const basic & hold() const
Stop further evaluation.
Definition: basic.cpp:887
bool is_equal(const basic &other) const
Test for syntactic equality.
Definition: basic.cpp:863
virtual int compare_same_type(const basic &other) const
Returns order relation between two objects of same type.
Definition: basic.cpp:719
Generate all compositions of a partition of an integer n, starting with the compositions which has no...
Definition: utils.h:395
const std::vector< unsigned > & get() const
Definition: utils.h:460
Wrapper template for making GiNaC classes out of STL containers.
Definition: container.h:73
Lightweight wrapper for GiNaC's symbolic objects.
Definition: ex.h:72
numeric integer_content() const
Compute the integer content (= GCD of all numeric coefficients) of an expanded polynomial.
Definition: normal.cpp:318
bool match(const ex &pattern) const
Check whether expression matches a specified pattern.
Definition: ex.cpp:95
bool is_polynomial(const ex &vars) const
Check whether expression is a polynomial.
Definition: ex.cpp:241
ex diff(const symbol &s, unsigned nth=1) const
Compute partial derivative of an expression.
Definition: ex.cpp:86
ex expand(unsigned options=0) const
Definition: ex.cpp:73
bool is_equal(const ex &other) const
Definition: ex.h:345
int degree(const ex &s) const
Definition: ex.h:173
bool has(const ex &pattern, unsigned options=0) const
Definition: ex.h:151
ex evalf() const
Definition: ex.h:121
ex conjugate() const
Definition: ex.h:146
unsigned return_type() const
Definition: ex.h:230
return_type_t return_type_tinfo() const
Definition: ex.h:231
ex imag_part() const
Definition: ex.h:148
ex subs(const exmap &m, unsigned options=0) const
Definition: ex.h:841
bool info(unsigned inf) const
Definition: ex.h:132
int compare(const ex &other) const
Definition: ex.h:322
bool is_zero() const
Definition: ex.h:213
void print(const print_context &c, unsigned level=0) const
Print expression to stream.
Definition: ex.cpp:56
ex op(size_t i) const
Definition: ex.h:136
int ldegree(const ex &s) const
Definition: ex.h:174
ex real_part() const
Definition: ex.h:147
ex evalm() const
Definition: ex.h:122
ex coeff(const ex &s, int n=1) const
Definition: ex.h:175
A pair of expressions.
Definition: expair.h:38
ex coeff
second member of pair, must be numeric
Definition: expair.h:91
size_t nops() const override
Number of operands/members.
Definition: expairseq.cpp:202
epvector seq
Definition: expairseq.h:127
@ expand_rename_idx
used internally by mul::expand()
Definition: flags.h:34
@ algebraic
enable algebraic matching
Definition: flags.h:43
@ integer_polynomial
Definition: flags.h:256
@ cinteger_polynomial
Definition: flags.h:257
@ crational_polynomial
Definition: flags.h:259
@ rational_polynomial
Definition: flags.h:258
Product of expressions.
Definition: mul.h:32
Non-commutative product of expressions.
Definition: ncmul.h:33
This class is a wrapper around CLN-numbers within the GiNaC class hierarchy.
Definition: numeric.h:82
bool is_pos_integer() const
True if object is an exact integer greater than zero.
Definition: numeric.cpp:1161
const numeric & mul_dyn(const numeric &other) const
Numerical multiplication method.
Definition: numeric.cpp:957
bool is_crational() const
True if object is an exact rational number, may even be complex (denominator may be unity).
Definition: numeric.cpp:1243
bool is_rational() const
True if object is an exact rational number, may even be complex (denominator may be unity).
Definition: numeric.cpp:1201
const numeric real() const
Real part of a number.
Definition: numeric.cpp:1339
long to_long() const
Converts numeric types to machine's long.
Definition: numeric.cpp:1313
int compare(const numeric &other) const
This method establishes a canonical order on all numbers.
Definition: numeric.cpp:1104
bool is_positive() const
True if object is not complex and greater than zero.
Definition: numeric.cpp:1136
bool is_real() const
True if object is a real integer, rational or float (but not complex).
Definition: numeric.cpp:1208
const numeric numer() const
Numerator.
Definition: numeric.cpp:1356
bool is_integer() const
True if object is a non-complex integer.
Definition: numeric.cpp:1154
const numeric power(const numeric &other) const
Numerical exponentiation.
Definition: numeric.cpp:900
const numeric denom() const
Denominator.
Definition: numeric.cpp:1387
const numeric mul(const numeric &other) const
Numerical multiplication method.
Definition: numeric.cpp:880
bool is_equal(const numeric &other) const
Definition: numeric.cpp:1122
int to_int() const
Converts numeric types to machine's int.
Definition: numeric.cpp:1303
const numeric inverse() const
Inverse of a number.
Definition: numeric.cpp:1053
bool is_zero() const
True if object is zero.
Definition: numeric.cpp:1129
const numeric div(const numeric &other) const
Numerical division method.
Definition: numeric.cpp:890
Generate all bounded combinatorial partitions of an integer n with exactly m parts (including zero pa...
Definition: utils.h:327
const std::vector< unsigned > & get() const
Definition: utils.h:337
Exception class thrown when a singularity is encountered.
Definition: numeric.h:70
This class holds a two-component object, a basis and and exponent representing exponentiation.
Definition: power.h:39
static ex expand_mul(const mul &m, const numeric &n, unsigned options, bool from_expand=false)
Expand factors of m in m^n where m is a mul and n is an integer.
Definition: power.cpp:1108
void do_print_dflt(const print_dflt &c, unsigned level) const
Definition: power.cpp:107
void do_print_csrc(const print_csrc &c, unsigned level) const
Definition: power.cpp:179
static ex expand_add(const add &a, long n, unsigned options)
expand a^n where a is an add and n is a positive integer.
Definition: power.cpp:897
int degree(const ex &s) const override
Return degree of highest power in object s.
Definition: power.cpp:301
void read_archive(const archive_node &n, lst &syms) override
Read (a.k.a.
Definition: power.cpp:73
ex subs(const exmap &m, unsigned options=0) const override
Substitute a set of objects by arbitrary expressions.
Definition: power.cpp:615
int ldegree(const ex &s) const override
Return degree of lowest power in object s.
Definition: power.cpp:316
static ex expand_add_2(const add &a, unsigned options)
Special case of power::expand_add.
Definition: power.cpp:1039
ex real_part() const override
Definition: power.cpp:667
void archive(archive_node &n) const override
Save (a.k.a.
Definition: power.cpp:80
ex derivative(const symbol &s) const override
Implementation of ex::diff() for a power.
Definition: power.cpp:744
friend class mul
Definition: power.h:42
ex map(map_function &f) const override
Construct new expression by applying the specified function to all sub-expressions (one level only,...
Definition: power.cpp:275
ex eval() const override
Perform automatic term rewriting rules in this class.
Definition: power.cpp:374
ex basis
Definition: power.h:105
void do_print_python(const print_python &c, unsigned level) const
Definition: power.cpp:210
void print_power(const print_context &c, const char *powersymbol, const char *openbrace, const char *closebrace, unsigned level) const
Definition: power.cpp:93
void do_print_python_repr(const print_python_repr &c, unsigned level) const
Definition: power.cpp:215
ex exponent
Definition: power.h:106
unsigned precedence() const override
Return relative operator precedence (for parenthezing output).
Definition: power.h:53
ex op(size_t i) const override
Return operand/member at position i.
Definition: power.cpp:268
power(const ex &lh, const ex &rh)
Definition: power.h:48
void do_print_latex(const print_latex &c, unsigned level) const
Definition: power.cpp:120
ex conjugate() const override
Definition: power.cpp:646
ex coeff(const ex &s, int n=1) const override
Return coefficient of degree n in object s.
Definition: power.cpp:331
bool has(const ex &other, unsigned options=0) const override
Test for occurrence of a pattern.
Definition: power.cpp:590
ex evalm() const override
Evaluate sums, products and integer powers of matrices.
Definition: power.cpp:578
return_type_t return_type_tinfo() const override
Definition: power.cpp:773
size_t nops() const override
Number of operands/members.
Definition: power.cpp:263
ex expand(unsigned options=0) const override
Expand expression, i.e.
Definition: power.cpp:778
bool info(unsigned inf) const override
Information about the object.
Definition: power.cpp:224
ex eval_ncmul(const exvector &v) const override
Definition: power.cpp:641
void do_print_csrc_cl_N(const print_csrc_cl_N &c, unsigned level) const
Definition: power.cpp:164
unsigned return_type() const override
Definition: power.cpp:768
ex imag_part() const override
Definition: power.cpp:703
bool is_polynomial(const ex &var) const override
Check whether this is a polynomial in the given variables.
Definition: power.cpp:287
ex evalf() const override
Evaluate object numerically.
Definition: power.cpp:565
Base class for print_contexts.
Definition: print.h:103
Context for C source output using CLN numbers.
Definition: print.h:182
Base context for C source output.
Definition: print.h:158
Context for default (ginsh-parsable) output.
Definition: print.h:115
Context for latex-parsable output.
Definition: print.h:123
Context for python-parsable output.
Definition: print.h:139
Context for python pretty-print output.
Definition: print.h:131
@ expanded
.expand(0) has already done its job (other expand() options ignore this flag)
Definition: flags.h:204
@ evaluated
.eval() has already done its job
Definition: flags.h:203
@ hash_calculated
.calchash() has already done its job
Definition: flags.h:205
@ no_pattern
disable pattern matching
Definition: flags.h:51
@ algebraic
enable algebraic substitutions
Definition: flags.h:53
Basic CAS symbol.
Definition: symbol.h:39
Definition of optimizing macros.
#define likely(cond)
Definition: compiler.h:32
Interface to GiNaC's constant types and some special constants.
Interface to sequences of expression pairs.
unsigned options
Definition: factor.cpp:2475
vector< int > k
Definition: factor.cpp:1435
size_t n
Definition: factor.cpp:1432
size_t c
Definition: factor.cpp:757
ex x
Definition: factor.cpp:1610
size_t r
Definition: factor.cpp:757
mvec m
Definition: factor.cpp:758
size_t last
Definition: factor.cpp:1434
Interface to GiNaC's indexed expressions.
Interface to GiNaC's initially known functions.
Definition of GiNaC's lst.
Interface to symbolic matrices.
Interface to GiNaC's products of expressions.
Definition: add.cpp:38
const numeric * _num_1_p
Definition: utils.cpp:351
const ex _ex2
Definition: utils.cpp:389
const numeric pow(const numeric &x, const numeric &y)
Definition: numeric.h:251
const ex _ex1_2
Definition: utils.cpp:381
std::map< ex, ex, ex_is_less > exmap
Definition: basic.h:50
std::vector< expair > epvector
expair-vector
Definition: expairseq.h:33
const numeric abs(const numeric &x)
Absolute value.
Definition: numeric.cpp:2315
bool is_negative(const numeric &x)
Definition: numeric.h:269
const ex _ex1
Definition: utils.cpp:385
bool are_ex_trivially_equal(const ex &e1, const ex &e2)
Compare two objects of class quickly without doing a deep tree traversal.
Definition: ex.h:699
const numeric binomial(const numeric &n, const numeric &k)
The Binomial coefficients.
Definition: numeric.cpp:2143
print_func< print_dflt >(&diracone::do_print). print_func< print_latex >(&diracone
Definition: clifford.cpp:51
const numeric * _num2_p
Definition: utils.cpp:388
const numeric exp(const numeric &x)
Exponential function.
Definition: numeric.cpp:1439
const numeric cos(const numeric &x)
Numeric cosine (trigonometric function).
Definition: numeric.cpp:1470
const ex _ex_1
Definition: utils.cpp:352
GINAC_IMPLEMENT_REGISTERED_CLASS_OPT(add, expairseq, print_func< print_context >(&add::do_print). print_func< print_latex >(&add::do_print_latex). print_func< print_csrc >(&add::do_print_csrc). print_func< print_tree >(&add::do_print_tree). print_func< print_python_repr >(&add::do_print_python_repr)) add
Definition: add.cpp:40
const numeric iquo(const numeric &a, const numeric &b)
Numeric integer quotient.
Definition: numeric.cpp:2404
bool is_pos_integer(const numeric &x)
Definition: numeric.h:275
const numeric log(const numeric &x)
Natural logarithm.
Definition: numeric.cpp:1450
const numeric sin(const numeric &x)
Numeric sine (trigonometric function).
Definition: numeric.cpp:1461
const numeric * _num1_p
Definition: utils.cpp:384
ex numer(const ex &thisex)
Definition: ex.h:760
ex factor(const ex &poly, unsigned options)
Interface function to the outside world.
Definition: factor.cpp:2576
bool is_integer(const numeric &x)
Definition: numeric.h:272
static void print_sym_pow(const print_context &c, const symbol &x, int exp)
Definition: power.cpp:140
lst rename_dummy_indices_uniquely(const exvector &va, const exvector &vb)
Similar to above, where va and vb are the same and the return value is a list of two lists for substi...
Definition: indexed.cpp:1460
GINAC_IMPLEMENT_REGISTERED_CLASS_OPT_T(lst, basic, print_func< print_context >(&lst::do_print). print_func< print_tree >(&lst::do_print_tree)) template<> bool lst GINAC_BIND_UNARCHIVER(lst)
Specialization of container::info() for lst.
Definition: lst.cpp:42
const ex _ex0
Definition: utils.cpp:369
bool tryfactsubs(const ex &origfactor, const ex &patternfactor, int &nummatches, exmap &repls)
Definition: mul.cpp:672
std::vector< ex > exvector
Definition: basic.h:48
exvector get_all_dummy_indices(const ex &e)
Returns all dummy indices from the exvector.
Definition: indexed.cpp:1435
const numeric multinomial_coefficient(const std::vector< unsigned > &p)
Compute the multinomial coefficient n!/(p1!*p2!*...*pk!) where n = p1+p2+...+pk, i....
Definition: utils.cpp:60
const numeric * _num0_p
Definition: utils.cpp:367
Interface to GiNaC's non-commutative products of expressions.
Makes the interface to the underlying bignum package available.
Interface to GiNaC's overloaded operators.
Interface to GiNaC's symbolic exponentiation (basis^exponent).
Interface to relations between expressions.
Function object for map().
Definition: basic.h:85
To distinguish between different kinds of non-commutative objects.
Definition: registrar.h:44
Interface to GiNaC's symbolic objects.
Interface to several small and furry utilities needed within GiNaC but not of any interest to the use...

This page is part of the GiNaC developer's reference. It was generated automatically by doxygen. For an introduction, see the tutorial.