bash-2.05$ gp
Reading GPRC: /home/his/.gprc ...Done.

                  GP/PARI CALCULATOR Version 2.1.3 (released)
                       i386 running netbsd 32-bit version
                (readline v1.0 enabled, extended help available)

                       Copyright (C) 2000 The PARI Group

PARI/GP is free software, covered by the GNU General Public License, and 
comes WITHOUT ANY WARRANTY WHATSOEVER.

Type ? for help, \q to quit.
Type ?12 for how to get moral (and possibly technical) support.

   realprecision = 28 significant digits
   seriesprecision = 16 significant terms
   format = g0.28

parisize = 4000000, primelimit = 500000
(20:56) gp> read("modular.gp")
(20:56) gp> f(n,x,y)=modularpoly(n)
(20:57) gp> h(n,x)=f(n,x,x)
(20:57) gp> h(2,x)
%1 = -x^4 + 2978*x^3 + 40449375*x^2 + 17496000000*x - 157464000000000
(20:57) gp> factor(%1)
%2 = 
[x - 8000 1]

[x - 1728 1]

[x + 3375 2]

(20:58) gp> h(3,x)
%3 = -x^6 + 4464*x^5 + 2585778176*x^4 + 17800519680000*x^3 - 769939996672000000*x^2 + 3710851743744000000000*x
(20:58) gp> factor(%3)
%4 = 
[x - 54000 1]

[x - 8000 2]

[x 1]

[x + 32768 2]

(20:58) gp> h(5,x)
%5 = -x^10 + 7440*x^9 + 1665990262720*x^8 + 215757860427776000*x^7 - 440440798293848579637248*x^6 + 53797234800359280738891202560*x^5 + 4726025910884027749483397649530880*x^4 + 73669962723556137647021587795909017600*x^3 - 250688456991364600842741491417948646014976*x^2 + 106548661606848850900840320546713018302464000*x + 141359947154721358697753474691071362751004672000
(20:58) gp> factor(%5)
%6 = 
[x - 287496 2]

[x - 1728 2]

[x + 32768 2]

[x + 884736 2]

[x^2 - 1264000*x - 681472000 1]

(20:58) gp> read("modularj.gp")
(21:09) gp> jtau(q,20)
%12 = q^-1 + 744 + 196884*q + 21493760*q^2 + 864299970*q^3 + 20245856256*q^4 + 333202640600*q^5 + 4252023300096*q^6 + 44656994071935*q^7 + 401490886656000*q^8 + 3176440229784420*q^9 + 22567393309593600*q^10 + 146211911499519294*q^11 + 874313719685775360*q^12 + 4872010111798142520*q^13 + 25497827389410525184*q^14 + 126142916465781843075*q^15 + 593121772421445058560*q^16 + 2662842413150775245160*q^17 + 11459912788444786513920*q^18 + 47438786801234168813250*q^19 + 189449976248893390028800*q^20 + O(q^21)
(21:17) gp> jj(q)=q^-1 + 744 + 196884*q + 21493760*q^2 + 864299970*q^3 + 20245856256*q^4 +333202640600*q^5 + 4252023300096*q^6 + 44656994071935*q^7 + 401490886656000*q^8 + 3176440229784420*q^9 + 22567393309593600*q^10 + 146211911499519294*q^11 + 874313719685775360*q^12 + 4872010111798142520*q^13 + 25497827389410525184*q^14 + 126142916465781843075*q^15 + 593121772421445058560*q^16 + 2662842413150775245160*q^17 + 11459912788444786513920*q^18 + 47438786801234168813250*q^19 + 189449976248893390028800*q^20;
(21:19) gp> jj(exp(2*Pi*I*(1+sqrt(7)*I)/2))
%13 = -3374.9999999999999999999999999999999999999999999999999998856624613044205965884625288270274004090445953212181667679775862020009363088534758991298325198485262471379951400536062555295511486472284515051346 - 1.643695324883795048 E-198*I
(21:19) gp> jj(exp(2*Pi*I*(1+sqrt(11)*I)/2))
%14 = -32767.999999999999999999999999999999999999999999999999999999999999999999999993133468617050631512551101796957437365426617294534005695910741194363444230540165227972584459136599541226509772299005855526867 - 4.103013251097894468 E-197*I
(21:19) gp> jj(exp(2*Pi*I*(1+sqrt(19)*I)/2))
%15 = -884735.99999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999058878433097618412945226933946506836566969433846979449275685440161531010816124450586376531012 - 1.084507834830101640202014833 E-195*I
(21:20) gp> quit;
Good bye!
bash-2.05$