1 // FILE MODDATA.CC: Implementation of member functions for class moddata
2 //////////////////////////////////////////////////////////////////////////
3 //
4 // Copyright 1990-2012 John Cremona
5 //
6 // This file is part of the eclib package.
7 //
8 // eclib is free software; you can redistribute it and/or modify it
9 // under the terms of the GNU General Public License as published by the
10 // Free Software Foundation; either version 2 of the License, or (at your
11 // option) any later version.
12 //
13 // eclib is distributed in the hope that it will be useful, but WITHOUT
14 // ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
15 // FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
16 // for more details.
17 //
18 // You should have received a copy of the GNU General Public License
19 // along with eclib; if not, write to the Free Software Foundation,
20 // Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA
21 //
22 //////////////////////////////////////////////////////////////////////////
23
24 #include <eclib/moddata.h>
25 #include <eclib/arith.h>
26
level(long n,long neigs)27 level::level(long n, long neigs)
28 {
29 //cout<<"Creating a class level with n = " << n << endl;
30 modulus=n;
31 plist=pdivs(n); npdivs=plist.size();
32 dlist=posdivs(n); ndivs=dlist.size();
33 nap=neigs;
34 primelist=plist;
35 primevar pr; long p; p0=0;
36 while(primelist.size()<(unsigned)nap)
37 {
38 p=pr;
39 if (ndivides(p,modulus))
40 {
41 if(p0==0) p0=p;
42 primelist.push_back(p);
43 }
44 pr++;
45 }
46 sqfac=1;
47 for(long ip=0; ip<npdivs; ip++)
48 {
49 p = plist[ip];
50 if(::divides(p*p,n)) sqfac*=p;
51 }
52 long rootn=(long)(sqrt((double)n)+0.1); // rounded down
53 squarelevel=(n==rootn*rootn);
54 }
55
bezout_x(long aa,long bb,long & xx)56 long bezout_x(long aa, long bb, long& xx)
57 {long a,b,c,x,oldx,newx,q;
58 oldx = 1; x = 0; a = aa; b = bb;
59 while (b!=0)
60 { q = a/b;
61 c = a - q*b; a = b; b = c;
62 newx = oldx - q*x; oldx = x; x = newx;
63 }
64 if (a<0) {xx=-oldx; return(-a);}
65 else {xx= oldx; return( a);}
66 }
67
moddata(long n)68 moddata::moddata(long n) :level(n)
69 {
70 // cout << "In constructor moddata::moddata.\n";
71 long i,p,x,d,nd,nnoninv;
72 phi=psi=modulus;
73 for(i=0; i<npdivs; i++)
74 { p = plist[i];
75 phi -= phi/p;
76 psi += psi/p;
77 }
78 nsymb = psi;
79 nsymb1 = 2*modulus-phi;
80 nsymb2 = nsymb-nsymb1;
81 invlist.resize(modulus); //codes
82 noninvlist.resize(modulus-phi); //list of non-units
83 noninvdlist.resize(modulus-phi); //list of divisors for each nonunit
84 gcdtable.resize(modulus); //list of gcds mod N
85 unitdivlist.resize(modulus); //list of units s.t. u*res | N
86 nnoninv=0;
87 for (i=0; i<modulus; i++) //set up codes
88 { d = bezout_x(i,modulus,x);
89 gcdtable[i]=d;
90 if (d==1) {unitdivlist[i] = invlist[i] = reduce(x); }
91 else
92 {invlist[i]=-nnoninv;
93 noninvlist[nnoninv]=i;
94 noninvdlist[nnoninv]=-1;
95 if (d<modulus)
96 {
97 for (nd=0; (nd<ndivs)&&(dlist[nd]!=d); nd++) ;
98 noninvdlist[nnoninv]=nd;
99 }
100 nnoninv++;
101 if(::gcd(x,modulus)!=1) // adjust x so coprime to N
102 {
103 long m=modulus/d, mm, mpower, mmold,u,v;
104 mpower=mm=m; mmold=1;
105 while(mm!=mmold)
106 {
107 mpower=xmodmul(m,mpower,modulus);
108 mmold=mm;
109 mm=::gcd(mpower,modulus);
110 }
111
112 bezout(mm,modulus/mm,u,v);
113 // Must be careful of overflow!
114 x = (x*v)%mm;
115 x = (x*(modulus/mm))%modulus;
116 x = (x+(u*mm))%modulus;
117 }
118 unitdivlist[i]=x;
119 }
120 }
121 if (ndivs>0) {dstarts.resize(ndivs);}
122 }
123
display() const124 void moddata::display() const
125 {
126 cout << "Level = " << modulus << "\n";
127 cout << "Number of symbols = " << nsymb << "\n";
128 cout << ndivs << " non-trivial divisors: " << dlist << endl;
129 cout << npdivs << " prime divisors: " << plist << endl;
130 cout << "invlist: " << invlist << endl;
131 cout << "noninvlist: " << noninvlist << endl;
132 cout << "noninvdlist: " << noninvdlist << endl;
133 cout << "gcdtable: " << gcdtable << endl;
134 cout << "unitdivlist: " << unitdivlist << endl;
135 }
136
of_filename(long n,char c)137 string of_filename(long n, char c)
138 {
139 stringstream s;
140 s << getenv_with_default("OF_DIR","./newforms");
141 s << "/" << c << n;
142 return s.str();
143 }
144
nf_filename(long n,char c)145 string nf_filename(long n, char c)
146 {
147 stringstream s;
148 s << getenv_with_default("NF_DIR","./newforms");
149 s << "/" << c << n;
150 return s.str();
151 }
152
small_nf_filename(long n,char c)153 string small_nf_filename(long n, char c)
154 {
155 stringstream s;
156 s << getenv_with_default("SNF_DIR","./smallnf");
157 s << "/" << c << n;
158 return s.str();
159 }
160