June 2020 - GiNaC-list

Release 1.7.10
by Richard B. Kreckel 28 Jun '20

28 Jun '20

Hi, GiNaC 1.7.10 has just been released. The changes are: * Fix collect_common_factors() for hidden zero arguments. * Fix build on modern systems (C++20 and Python3). Take care! -richy. -- Richard B. Kreckel <https://in.terlu.de/~kreckel/>

1 0

C binding for GiNaC
by Xianwen Chen (陈贤文) 22 Jun '20

22 Jun '20

Dear all, I wonder if it is possible to use GiNaC in a C program? Is there a binding for that? I guess there may be one, because GiNaC is used in other software, for example, Octave and Sage Math. But I couldn't find relevant information. Yours sincerely, Xianwen

2 2

problem about collect_common_factors
by Feng Feng 22 Jun '20

22 Jun '20

Hi all, I got a strange problem about collect_common_factors, the problem can be reproduced by the following example: ------ symbol a("a"),b("b"),c("c"); symtab table; table["a"] = a; table["b"] = b; table["c"] = c; parser reader(table); ex res = reader("b*((c+c*(1+4*a+2*b+b^2)-c*(1+2*b+b^2)-(1+4*a)*c))-a*(c+c*(1+4*a+2*b+b^2)-c*(1+2*b+b^2)-(1+4*a)*c)"); cout << res << endl; cout << collect_common_factors(res) << endl; ------ the collect_common_factors will throw an error: division by zero. so is this a bug for collect_common_factors? Best regards! Feng

3 2

Issue with derivative of "abs"
by Pierangelo Masarati 16 Jun '20

16 Jun '20

Dear GiNaC users, we have been using GiNaC for years now, to provide our multibody dynamics solver, MBDyn <http://www.mbdyn.org/>, some symbolic differentiation capabilities (perhaps not the most appropriate use, but some users seem to like it). I'll make the story as short as possible. In one of the most interesting applications, the user needs to: - define the variables that represent strains (and strain rates), - the expression for the configuration-dependent generalized forces in short, a constitutive law. The expression's partial derivatives with respect to the strains and strain rates are then automatically computed by GiNaC, and, when needed, the numerical values are substituted and the expressions are evaluated. The related code snippet looks like > GiNaC::symbol gEps("v"); // parameter symbols > GiNaC::ex gExpr; // expressions > GiNaC::ex gExprDEps; // derivatives > > GiNaC::lst l; > > l.append(gEps); > gExpr = GiNaC::ex("v*abs(v)", l); > > gExprDEps = gExpr.diff(gEps); One user recently tried something like this: var = v expr = v*abs(v) # a turbulent viscous law the result was dexpr/dvar = abs(v)+1/2*v*(v+conjugate(v))*abs(v)^(-1) which fails when v = 0 because of the division by zero. I tried with "realsymbol" (since our values can only be real); things, however, didn't change much: dexpr/dvar = v^2*abs(v)^(-1)+abs(v) I realized that the issue is in the derivative of the abs() function; indeed, the proposed formula is rather general (despite failing when numerically evaluated with 0 as the argument). When the argument is real, a more robust representation would be "d(abs(v))/dv = sign(v)". In the case at hand, "d(v*abs(v))/dv = 2*abs(v)" or, as a second choice, "d(v*abs(v))/dv = abs(v) + v*sign(v)". First, I considered the possibility to implement the "sign()" function in GiNaC, but I am not familiar enough with its internals. If it makes sense, I would definitely recommend its introduction. However, I noticed that the "step()" function is implemented (and step(0) correctly evaluates to 1/2); the "sign()" function can then be obtained as "sign(x) = 2*step(x) - 1" when "x" is real. To this end, I modified GiNaC's code as reported here: > diff --git a/ginac/inifcns.cpp b/ginac/inifcns.cpp > index d68afbb5..ec74a435 100644 > --- a/ginac/inifcns.cpp > +++ b/ginac/inifcns.cpp > @@ -321,6 +321,8 @@ static ex abs_expand(const ex & arg, unsigned options) > static ex abs_expl_derivative(const ex & arg, const symbol & s) > { > ex diff_arg = arg.diff(s); > + if (arg.info(info_flags::real)) > + return diff_arg*(2*step(arg) - 1); > return > (diff_arg*arg.conjugate()+arg*diff_arg.conjugate())/2/abs(arg); > } which indeed gives the expected result: dexpr/dvar = v*(-1+2*step(v))+abs(v) I wonder whether my fix can be considered correct and, in case, whether it could be included in the distribution. Sincerely, p. -- Pierangelo Masarati Professore Ordinario di Costruzioni e Strutture Aerospaziali Dipartimento di Scienze e Tecnologie Aerospaziali Politecnico di Milano

2 1

Re: [GiNaC-list] Issue with derivative of "abs"
by Vladimir V. Kisil 16 Jun '20

16 Jun '20

Dear Pierangelo, Thanks for the email. It is nice to see another usage of GiNaC! Your patch combines the best from both worlds as was discussed: https://www.ginac.de/pipermail/ginac-devel/2014-April/002105.html So I am voting in favour of this patch to be applied! Meanwhile, it seems the faster solution for your would be to clone abs function with the alteration towards the desired derivative. It is not difficult: 1. Cut & Paste everything for definition abs from GiNaC's inifcns.h and inifcns.cpp to your code. 2. Replace the derivative function by abs_expl_derivative() from your patch. 3. Search & replace abs -> myabs in this portion of code. Then the function myabs() shall be ready to be used. Best wishes, Vladimir -- Vladimir V. Kisil http://www.maths.leeds.ac.uk/~kisilv/ Book: Geometry of Mobius Transformations http://goo.gl/EaG2Vu Software: Geometry of cycles http://moebinv.sourceforge.net/ Jupyter: https://github.com/vvkisil/MoebInv-notebooks >>>>> On Tue, 16 Jun 2020 19:24:32 +0200, Pierangelo Masarati <pierangelo.masar\ ati(a)polimi.it> said: PM> Dear GiNaC users, PM> we have been using GiNaC for years now, to provide our multibody PM> dynamics solver, MBDyn <http://www.mbdyn.org/>, some symbolic PM> differentiation capabilities (perhaps not the most appropriate PM> use, but some users seem to like it). PM> I'll make the story as short as possible. In one of the most PM> interesting applications, the user needs to: PM> - define the variables that represent strains (and strain PM> rates), PM> - the expression for the configuration-dependent generalized PM> forces PM> in short, a constitutive law. The expression's partial PM> derivatives with respect to the strains and strain rates are PM> then automatically computed by GiNaC, and, when needed, the PM> numerical values are substituted and the expressions are PM> evaluated. PM> The related code snippet looks like >> GiNaC::symbol gEps("v"); // parameter symbols >> GiNaC::ex gExpr; // expressions GiNaC::ex gExprDEps; // >> derivatives >> >> GiNaC::lst l; >> >> l.append(gEps); gExpr = GiNaC::ex("v*abs(v)", l); >> >> gExprDEps = gExpr.diff(gEps); PM> One user recently tried something like this: PM> var = v PM> expr = v*abs(v) # a turbulent viscous law PM> the result was PM> dexpr/dvar = abs(v)+1/2*v*(v+conjugate(v))*abs(v)^(-1) PM> which fails when v = 0 because of the division by zero. I tried PM> with "realsymbol" (since our values can only be real); things, PM> however, didn't change much: PM> dexpr/dvar = v^2*abs(v)^(-1)+abs(v) PM> I realized that the issue is in the derivative of the abs() PM> function; indeed, the proposed formula is rather general PM> (despite failing when numerically evaluated with 0 as the PM> argument). When the argument is real, a more robust PM> representation would be "d(abs(v))/dv = sign(v)". In the case PM> at hand, "d(v*abs(v))/dv = 2*abs(v)" or, as a second choice, PM> "d(v*abs(v))/dv = abs(v) + v*sign(v)". PM> First, I considered the possibility to implement the "sign()" PM> function in GiNaC, but I am not familiar enough with its PM> internals. If it makes sense, I would definitely recommend its PM> introduction. PM> However, I noticed that the "step()" function is implemented PM> (and step(0) correctly evaluates to 1/2); the "sign()" function PM> can then be obtained as "sign(x) = 2*step(x) - 1" when "x" is PM> real. To this end, I modified GiNaC's code as reported here: >> diff --git a/ginac/inifcns.cpp b/ginac/inifcns.cpp index >> d68afbb5..ec74a435 100644 --- a/ginac/inifcns.cpp +++ >> b/ginac/inifcns.cpp @@ -321,6 +321,8 @@ static ex >> abs_expand(const ex & arg, unsigned options) static ex >> abs_expl_derivative(const ex & arg, const symbol & s) { >> ex diff_arg = arg.diff(s); + if >> (arg.info(info_flags::real)) + return >> diff_arg*(2*step(arg) - 1); return >> (diff_arg*arg.conjugate()+arg*diff_arg.conjugate())/2/abs(arg); >> } PM> which indeed gives the expected result: PM> dexpr/dvar = v*(-1+2*step(v))+abs(v) PM> I wonder whether my fix can be considered correct and, in case, PM> whether it could be included in the distribution. PM> Sincerely, p. PM> -- Pierangelo Masarati Professore Ordinario di Costruzioni e PM> Strutture Aerospaziali Dipartimento di Scienze e Tecnologie PM> Aerospaziali Politecnico di Milano PM> _______________________________________________ GiNaC-list PM> mailing list GiNaC-list(a)ginac.de PM> https://www.ginac.de/mailman/listinfo/ginac-list

1 0