Richard Kreckel [Wed, 31 Jan 2018 09:04:25 +0000 (10:04 +0100)]
Improve gcd(a, b) where one argument is a power of a symbol.
The already implemented recursion
gcd(x^n, x*p(x)) -> x*gcd(x^(n-1), p(x))
is not ambitious enough: If p(x) has a factor of x, it just goes through
the same step again, and if p(x) has a factor which is a power of x, this
is reapeted many times.
Example:
gcd(x^n, expand(x^n*(1+x)))
used to go recurse through the gcd routine n times, which could
easily lead to a stack overflow for n about 10^5.
To improve the situation, introduce a special case for gcd where one
argument is a power of a symbol and just cancel the common factor.
This turned out to be the root cause of segfaults in matrix elimination
algorithms, as reported by Patrick Schulz and Vitaly Magerya:
Cf. <https://www.ginac.de/pipermail/ginac-list/2018-January/thread.html>
Richard Kreckel [Sun, 28 Jan 2018 18:59:54 +0000 (19:59 +0100)]
[CHECK] Add some more factorization exams.
Add the 15 multivariate polynomials from P. S. Wang's paper "An
Improved Multivariate Polynomial Factoring Algorithm" as exams.
Don't add them as benchmarks since timings vary considerably
depending on internal choice of variables. In any case, this should
never take hours. (But before 630db9a0b0 it did, occasionally.)
Richard Kreckel [Sun, 28 Jan 2018 18:23:50 +0000 (19:23 +0100)]
Add optional algorithm selection to matrix::inverse().
Following a proposal by Vitaly Magerya, a signature with a solve_algo
is added to class matrix. (Not using a default argument in order to
not break the ABI.)
Richard Kreckel [Wed, 24 Jan 2018 21:37:24 +0000 (22:37 +0100)]
Speed up special cases of square-free factorization.
Square-free factorization of polynomials containing a factor which is
a high power P of a symbol x did scale like O(P) in space and time.
This patch introduces a shortcut in Yun's algorithm, such that the
computation is only O(1) in space and time.
This makes it possible to compute, say sqrfree(x^P + x^(P+1)) =>
(1+x)*x^P with P=123456789.
Found this to be a bottleneck while debugging one of Vitaly Magerya's
examples.
Richard Kreckel [Wed, 15 Feb 2017 11:58:48 +0000 (12:58 +0100)]
Clean up combinatorial helpers.
* Convert all signatures from int to unsiged.
* Renamed partition_generator to partition_with_zero_parts_generator...
* ...and add class partition_generator which does not include zero parts.
* Rename ::current() to ::get() and make them compute vectors on first use.
Richard Kreckel [Sun, 6 Nov 2016 18:35:46 +0000 (19:35 +0100)]
Clean up some header files and fix compilation with MSC.
On MSC there is no <unistd.h>. It turns out that we need to #include <io.h>
in order to use close(3). This fixes the build failure introduced in 8305ec38.
Richard Kreckel [Sat, 1 Oct 2016 22:45:11 +0000 (00:45 +0200)]
Suppress compiler warnings about inconsistent override.
The introduction of the 'override' keyword in GiNaC in d5b86dd10 was
deliberately incomplete: Duplication of code in macros was avoided.
It turns out that some compilers (e.g. g++ -Wsuggest-override or clang++
-Wall) complain if the 'override' keyword is used inconsistently. Let's
make this consistent now, even if it leads to the introduction of two
more macros.
Vladimir Kisil [Thu, 11 Aug 2016 18:47:08 +0000 (20:47 +0200)]
Add proper functions to make clifford_bar() and clifford_star().
Previously operations clifford_bar() and clifford_star() called
conjugate() method. This results in reversion of all non-commutative
entries (not only Clifford units) and produced complex conjugation
of all non-real items. The new routine operates on Clifford units only.
Richard Kreckel [Thu, 21 Jul 2016 06:55:44 +0000 (08:55 +0200)]
[bugfix] Fix crash in basic::subs().
A regression in 1.7.0 was introduced with 1b8bcb06 in function
basic::subs_one_level(): implicitly casting *this to an ex first for
finding *this in m and later in the function's return statement caused
a crash in the second cast because *this was deleted in the first one.
After all, *this was dynamically allocated in basic::subs().
This bug was reported and hunted down by Mario Prausa.
Richard Kreckel [Mon, 4 Apr 2016 07:05:50 +0000 (09:05 +0200)]
Install ginac-excompiler in $LIBEXECDIR, not in $BINDIR.
...and make the GiNaC library aware of where it is installed.
The ginac-excompiler script is only invoked by GiNaC::compile_ex(...)
and serves no purpose on its own. In compliance with the FHS, it should
be installed in $LIBEXECDIR, not in $BINDIR. This also disburdens
distribution packagers from having to provide a manpage (which may be
required for all programs in $BINDIR).
The location for $LIBEXECDIR defaults to ${prefix}/libexec/. It may be
overwritten at configuration time. (Many distributions want to set it to
${prefix}/lib/ginac/.)
Richard Kreckel [Fri, 5 Feb 2016 23:47:08 +0000 (00:47 +0100)]
[bugfix] fix elusive bug in quo, rem,...
The power of two rational numeric objects needs not be rational. As a
result, some care is required when transforming (b^e)*l -> (b*l^(1/e))^e
for some rational e and l. This is a common transformation in order to
convert a polynomial over Q into a polynomial over Z when computing
their quotient, remainder, etc. Failure to be careful can potentially
introduce spurious non-rational numbers into rational polynomials and
make those operations fail. This patch avoids this transformation when
l^(1/e) is not a rational number.
Richard Kreckel [Fri, 5 Feb 2016 23:47:08 +0000 (00:47 +0100)]
[bugfix] fix elusive bug in quo, rem,...
The power of two rational numeric objects needs not be rational. As a
result, some care is required when transforming (b^e)*l -> (b*l^(1/e))^e
for some rational e and l. This is a common transformation in order to
convert a polynomial over Q into a polynomial over Z when computing
their quotient, remainder, etc. Failure to be careful can potentially
introduce spurious non-rational numbers into rational polynomials and
make those operations fail. This patch avoids this transformation when
l^(1/e) is not a rational number.
Richard Kreckel [Thu, 28 Jan 2016 21:45:56 +0000 (22:45 +0100)]
Remove 'level' argument of normal().
The 'level' argument was modeled after that of the eval() methods
(removed in 6c946d4c). It has never been very useful except for
confusing developers and it hasn't been documented in the tutorial.
Moreover, I have no indication that it has ever been used.
Richard Kreckel [Thu, 28 Jan 2016 21:11:46 +0000 (22:11 +0100)]
Remove 'level' argument of evalf().
The 'level' argument was modeled after that of the eval() methods
(removed in 6c946d4c). It has never been very useful except for
confusing developers. Moreover, I have no indication that it has
ever been used.
Richard Kreckel [Wed, 6 Jan 2016 21:34:23 +0000 (22:34 +0100)]
Avoid x^0 and Order(x^0) terms together in series expansion.
At the branch cut, the series expansions of log(), atan(), and atanh()
assembled pseries objects which contained both x^0 and Order(x^0) terms.
This patch removes the extra x^0 term if the order is 0. It also adds
a GINAC_ASSERT for these kinds of invariants to the pseries ctor and
simplifies the loops in pseries' degree(), ldegree(), eval() and evalf()
member functions.
Clarification on symmetries of metric of clifford object.
The metric of a clifford object may be non-symmetric. Even if
a metric is defined by a symmetric tensor, clifford object may
not be aware of the symmetry and it needs to be explicitly declared.
Signed-off-by: Vladimir V. Kisil <kisilv@maths.leeds.ac.uk>
Richard Kreckel [Thu, 17 Dec 2015 21:45:50 +0000 (22:45 +0100)]
Remove useless code in add::eval().
With commit ae6c004b, we have actually a more rigorous solution for the
bug fixed first (and wrongly) on September 23 2010 with commit 89d5356b.
This patch removes the now useless code and adds the new regression check
from master, just in case.
Richard Kreckel [Wed, 16 Dec 2015 20:22:36 +0000 (21:22 +0100)]
Make .eval() evaluate top-level only.
With the previous patch, some old workarounds have become unnecessary:
If all expressions are evaluated, any object which is an aggregate of
expressions will only ever have to evaluate the top level. As such, it
has become pointless to evaluate child objects of an expression prior
to doing its own term-rewriting. This patch removes the evaluation of
children from all GiNaC objects. It also removes the now superfluous
parameter 'level' of the eval methods.
Richard Kreckel [Wed, 16 Dec 2015 12:00:30 +0000 (13:00 +0100)]
Make add::eval(), mul::eval() work without compromise.
If a user asks to evaluate an object, the expectation is that the
library returns an evaluated, canonical expression. Everything else is
a bug. (It doesn't matter whether the user asks explicitly or by
assigning to an ex.) It turns out that it was possible to construct
objects which didn't fully evaluate. This patch fixes this by making
eval() a little bit more careful.
This obsoletes the need to go through combine_ex_with_coeff_to_pair()
when constructing an epvector that is then used to construct an add
or mul. (Alas, this was omitted frequently enough...)
Richard Kreckel [Mon, 30 Nov 2015 22:06:53 +0000 (23:06 +0100)]
Improve method of setting status_flags::dynallocated.
There seems to be no way to obsolete the need to mark an object derived
from basic and handled by ex as being 'on the heap', at least none that
doesn't have significant performance impact. Having said that, this
patch aims at making this process simpler and more intuitive.
Where, before, one would return from a function returning an ex with
return (new mul(a, b))->setflag(status_flags::dynallocated);
this patch lets us return with
return dynallocate<mul>(a, b);
which should be much clearer. In any case, it involves less typing.
The two points where the status_flags::dynallocated are set are now
* the dynallocate<B>(args...) template function and
* the virtual duplicate() member functions.
This patch rolls out the new functionality throughout the library.
Richard Kreckel [Sun, 29 Nov 2015 11:33:59 +0000 (12:33 +0100)]
Rename array of precomputed data in test suite.
Reason: C++17 may introduce a std::data<> template. Right now, the GCC 6.0
prerelease bails out at this code, when compiler with -std=c++17.
<http://en.cppreference.com/w/cpp/iterator/data>
Richard Kreckel [Sat, 28 Nov 2015 14:43:46 +0000 (15:43 +0100)]
Use neseted initializer lists to construct matrix objects.
Add constructor of initializer_list<initializer_list<ex>> to matrix.
Use this syntax where, previously, ctor from comma-separated list of
elements was used. Deprecate the ctor from comma-separated list.
Note: The output format '[[a,b],[c,d]]' and ginsh syntax are
unchanged because lists are printed '{a,b,c}' and a matrix is not a
list of lists.