Karim BELABAS on Wed, 2 Jul 2003 23:36:46 +0200 (MEST)

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

Re: using zbrent in C mode

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

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

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! ]


  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...

Karim Belabas                     Tel: (+33) (0)1 69 15 57 48
Dép. de Mathématiques, Bât. 425   Fax: (+33) (0)1 69 15 60 19
Université Paris-Sud              http://www.math.u-psud.fr/~belabas/
F-91405 Orsay (France)            http://www.parigp-home.de/  [PARI/GP]