Is GiNaC really that dirty?
Chris Dams
chrisd at sci.kun.nl
Tue Jul 17 12:28:13 CEST 2001
Hello everybody,
I have been playing around with ginac a bit, but the impression that I get
from the tutorial is that one has to write really dirty code to get
something done. Since it is not my intention to start a flame war based on
vague sentences like the previous one, I will give an example of a
function I want to write. The answer I am hoping for is something like:
``No, you stupid, you don't need to write things like that. It can be done
much more beautifull, namely like this: ......''
To be more concrete, I want to write the function
ex checkdeg(ex f,ex var,int n)
which expects that *f* is a polynomial in expanded form in a variable in
*var* and perhaps other variables. The return value is the same polynomial
except that all terms in which *var* occurs in a power higher than *n* are
discarded. Maybe you want to compare the function I have written to
accomplish this with the FORM statement
id var^m?{>n}=0;
which does the same thing. ``Comparing'' especially with respect to its
length and elegance.
I want to make it a *map_function*. Now I have a problem since
*map_function*s only have one argument. The solution I came up with looks
like
struct Checkdeg:public map_function
{ ex var;
int n;
ex operator()(const ex&f,ex v,int nint)
{ var=v;
n=nint;
return (*this)(f);
}
ex operator()(const ex&f)
{
... Here comes the code that actually does something ...
}
} checkdeg;
which looks somewhat queer. I really would like to be able to say
something like
return f.map(*this)(var,n)
or maybe like this
return f.map((*this)(var,n))
Well, you get the idea. I want to map a function with more than just one
argument.
Now after having ``solved'' the problem of writing the *map_function*, now
let us consider the the function
ex operator()(const ex&f)
This function first checks whether *f* is an addition and then uses the
*map* function to tell ginac that it should treat each term seperately. No
problem here anymore. Then we have the case were *f* is not an addition. I
really would want to write someting like
lst repls;
return f.match(wild(0)*pow(var,wild(1)),repls)
if(repls.op(1)>n)
f.subs(var==0);
but then I discover that *pow(var,m)* is not mathched since it is not a
product. I really cannot imagine any reason whatsoever why not to
automatically consider something that is not a product as a product with
a single factor in patern matching against a product. It would make code
so much more elegant, and cannot be really difficult or inefficient. The
same thing with the power. The manual mentions something like ten times
that one may not use the pattern *pow(var,wild(1))* to match the
expression *var*. Why on earth not? The cynical remark would be that
if something needs to be said ten times, it can hardly be something you
want to hear.
A solution to this problem would seem to be using the *.has* method, but
then again, this one does not return the thing that matches the wildcard.
So I have to use the *.find* method. Then I have to interpret the list it
returns. This I have to do by writing an *if* like
if(repls.op(0).op(1)>n)
to check the first element in the list, which contains an expression like
*var^m* and get the *m* out of this. This is may be only slightly
inelegant, but I think it is inelegant. Why does the *.find* method return
the entire match expression, why not simply the thing that is equal to the
wildcard. It would change my *if* into
if(repls.op(0)>n)
which is somewhat shorter (but in complicated expressions could be very
much sorter). In my opinion the user can generally be expected to know the
expression againts which he matches, so why does he have to dissemine that
one again when he reads the result from the match? It simply
does not give any new information. Furthermore it is inconsistent since
the match method indeed only returns the values of the wildcards. The only
objection I could imagine is that the *find* function could return more
than one match. This is true, but in this case one could return
{{value of $0,value of $1},{value of $0,value of $1}}
Indeed, now I would have to write down exactly the same test, but the
advantage is that I don't need to know that the exponent of a power is
stored in *.op(1)*, which I really don't want to and shouldn't need to
know as an ordinary user.
After these considerations I felt I had to write something like
ex operator()(const ex&f)
{ if(is_a<add>(f))
return f.map((*this)(var,n));
switch(n)
{ case 0: return f.subs(var==0);
case 1: return f.subs(pow(var,wild())==0);
default:
lst repls;
if(f.find(pow(var,wild()),repls))
{ if(repls.op(0).op(1)>n)
return f.subs(var==0);
else
return f;
}
else
return f;
}
}
for my function. Now I guess the whole *struct Checkdeg* is at least
thrice as long as the thing I saw before me in my dreams. Does it really
have to be that dirty?
Greetings,
Chris Dams
More information about the GiNaC-list
mailing list