Igor Schein on Wed, 2 Jul 2003 23:30:16 -0400


[Date Prev] [Date Next] [Thread Prev] [Thread Next] [Date Index] [Thread Index]

Re: using zbrent in C mode


On Wed, Jul 02, 2003 at 11:36:46PM +0200, Karim BELABAS wrote:
> On Wed, 2 Jul 2003, Bill Allombert wrote:
> > On Wed, Jul 02, 2003 at 08:58:55PM +0200, Olivier Ramare wrote:
> >>   I'll be able to program in pari-C one day ...
> >> Just now, can anyone tell me how to program that
> >>
> >> x=solve(t=a,b,sin(t))
> >
> > Short answer: You cannot.
> >
> > Long answer: You can try to use lisexpr():
> 
> Intermediate answer: it was a 2 minutes job to convert the existing
> zbrent into
> 
> GEN
> zbrent(void *dat, GEN (*f)(GEN, void*), GEN a, GEN b, long prec)
> 
> where f(x) is computed by f(x, dat) [ to cater for functions requiring
> external parameters: dat being a pointer to an arbitrary data structure,
> typecast to void *, reinterpreted within f(). Possibly dat = NULL ].
> 
> This will be cleaner and much faster than using lisexpr() for simple
> functions.
> 
> It's not GP2C-friendly [ iterators would be better ], but at least it is
> usable in libpari, and requires tiny code changes  [ I wrote a few iterators
> when fiddling with intnum(), the resulting code would be a nightmare to
> maintain ]
> 
> It's not mandatory that f be stack-clean [ although it won't hurt if it is! ]
> 
> Example:
> 
>   GEN myfun(GEN x, void *dat) { return gsub(gsqr(x), gdeux); }
> 
>   GEN test() { return zbrent(NULL, &myfun, gzero, gdeux, DEFAULTPREC); }
> 
> Warning: all this is untested, but it should work. Provided you update from
> CVS, that is...

Originally posted on pari-users.

Here're some instant regression cases:

? solve(x=-1,0,sin(x)-x);
  ***   forbidden assignment t_SMALL --> t_REAL.
? solve(x=0,1,deriv(x));
  ***   forbidden  ***   unknown type 52.

Thanks

Igor