[GiNaC-list] feature request

[Marko Riedel]

Greetings to all.

Continuing with my GiNaC project I have noticed that there was some 
functionality missing that other software offers. In particular, I 
require **factorization of integers** and a function to compute the 
**set of divisors** of an integer. I have implemented these using basic 
math such as might be used in a freshman course. It passes all my tests, 
but it really needs an implementation from a professional mathematician. 
Maybe you can file this with other feature requests you have received. 
An excerpt from my code follows below. I hope I have understood lists 
correctly that they are doubly linked and a list append of one element 
runs in constant time.

Best regards,

Marko

lst factorizeposint(numeric n)
{
   lst result = {}, ent;

   int pot = 0;
   while(irem(n, numeric(2)).is_zero()){
     pot++; n = iquo(n, numeric(2));
   }
   if(pot > 0){
     ent = {numeric(2), pot};
     result.append(ent);
   }

   numeric q = 3;

   while(q*q <= n){
     pot = 0;
     while(irem(n, q).is_zero()){
       pot++; n = iquo(n, q);
     }
     if(pot > 0){
       ent = {q, pot};
       result.append(ent);
     }
     q += 2;
   }
   if(n > 1){
     ent = {n, numeric(1)};
     result.append(ent);
   }

   return result;
}

numeric totient(numeric n)
{
   if(n==1) return 1;

   lst pfact = factorizeposint(n); numeric result = n;

   for(lst::const_iterator it = pfact.begin();
       it != pfact.end(); ++it){
     lst ent = ex_to<lst>(*it);
     result *= numeric(1)
       -numeric(1)/ex_to<numeric>(ent.op(0));
   }

   return result;
}


lst divisors(numeric n)
{
   if (n==1) return {1};

   lst pfact = factorizeposint(n), ent;

   numeric tau(1);

   lst::const_iterator it;
   for(it = pfact.begin(); it != pfact.end(); ++it){
     ent = ex_to<lst>(*it);
     tau *= numeric(1)+ex_to<numeric>(ent.op(1));
   }

   exvector res {};
   res.resize(tau.to_int(), numeric(1));

   int sofar = 1; int pos;
   for(it = pfact.begin(); it != pfact.end(); ++it){
     pos = sofar;
     ent = ex_to<lst>(*it);
     numeric p = ex_to<numeric>(ent.op(0));
     numeric v = ex_to<numeric>(ent.op(1));
     for(numeric vv(1); vv <= v; vv++){
       for(int q = 0; q < sofar; q++){
         res[pos++] = res[q]*pow(p,vv);
       }
     }
     sofar = pos;
   }

   lst result = {};
   for(pos=0; pos<sofar; pos++)
     result.append(res[pos]);

   return result;
}