share/affine/macbasic.mac

   1 /*-*- Mode: MACSYMA; Package: MACSYMA -*-*/
   2 1$
   3 dotexptsimp:false$
   4 dotscrules:true$
   5 matrix_element_mult:"."$
   6 algebraic:true$
   7 fast_solve:true$
   8 aaaa:create_list (concat(aa,i),i,1,105)$
   9 bbbb:create_list (concat(bb,i),i,1,35)$
  10 cccc:create_list (concat(cc,i),i,1,35)$
  11 parr:create_list (concat(par,i),i,0,20)$
  12 dddd:create_list (concat(dd,i),i,1,35)$
  13 eeee:create_list (concat(ee,i),i,1,35)$
  14 declare_scalar_list(dddd)$
  15 declare_scalar_list(eeee)$
  16 tellrat(%zeta3^2=-1-%zeta3,%zeta8^4=-1,
  17         %sqrt2^2=2,
  18         %zeta4^2+1,
  19         %zeta5^4=-(+%zeta5^3+%zeta5^2+%zeta5+1),
  20         %zeta6^2-%zeta6+1,
  21         %zeta7^6+%zeta7^5+%zeta7^4+%zeta7^3+%zeta7^2+%zeta7+1,
  22         %zeta9^6+%zeta9^3+1,
  23         %zeta24^8=%zeta24^4-1)$
  24 declare_scalar_list([a,%eps,%zeta3,%zeta4,
  25         %sqrt2,%zeta6,%zeta7,
  26         %zeta8,%zeta9,%zeta5,
  27         %i,b,c,d,e,%alpha,%alpha4,%alpha2,%alpha3,%beta,%zeta24])$
  28 declare_scalar_list(parr)$
  29 declare_scalar_list(bbbb)$
  30 declare_scalar_list(cccc)$
  31 declare_scalar_list(aaaa)$
  32 declare_weights(s2,4,t,3,w,1,v,1,x2,2,s,1,q1,5,q2,5,r,4,z5,5,z3,1,z2,1,z4,4,x3,3,x5,2,x1,2,x0,3,x10,3,u,1,x13,3,x03,3,s4,4,x04,4,y,1,s6,6)$
  33 declare_order_weights(y,1)$
  34 deg2_mat:matrix([a00*x.x+a01*x.y+a10*y.x+a11*y.y,sc00*x.x+sc01*x.y+sc10*y.x+sc11*y.y],
  35        [b00*x.x+b01*x.y+b10*y.x+b11*y.y,sd00*x.x+sd01*x.y+sd10*y.x+sd11*y.y])$
  36
  37 deg2_mat:matrix(
  38        [a0000*x.x+a0010*x.y+a0100*y.x+a0110*y.y,a0001*x.x+a0011*x.y+a0101*y.x+a0111*y.y],
  39        [a1000*x.x+a1010*x.y+a1100*y.x+a1110*y.y,a1001*x.x+a1011*x.y+a1101*y.x+a1111*y.y])$
  40 declare_scalar_list(delete(x,(delete(y,listofvars(deg2_mat)))))$
  41 rtx:matrix([x],[y])$
  42 lftx:matrix([x,y])$
  43 declare_scalar_list(delete(x,(delete(y,list_of_variables(deg2_mat)))))$
  44 fast_solve:true$
  45 central_elements(variables,deg):=block([expop:100],
  46         temp:[],
  47         monoms:mono(variables,deg),
  48         monoms_higher:mono(variables,1+deg),
  49         number_of_monoms:length(monoms),
  50         if not (listp(aaaa) or length(monoms) > length(aaaa)) then
  51           aaaa: create_list(concat(aa,i),i,1,number_of_monoms),
  52         declare_scalar_list(aaaa),
  53         f:0,i:0,
  54         for u in monoms do (i:i+1,
  55           f: u*?nth(i,aaaa)+f),
  56         for v in variables do (
  57           temp : append (temp,[com(f,v)])),
  58         eqns:Extract_linear_equations (temp,monoms_higher),
  59         if fast_solve then ans : fast_linsolve (eqns,firstn(length(monoms),aaaa)) else
  60         ans:linsolve(eqns,firstn (length( monoms), aaaa)),
  61         return(sublis(ans,f)))$
  62
  63
  64 skew_centralizer(variables,elements,zet,deg):=block([expop:100,monoms_higher:[]],
  65         temp:[],
  66         monoms:mono(variables,deg),
  67         for j in map(nc_degree,elements)
  68         do( monoms_higher:append(monoms_higher,mono(variables,j+deg))),
  69         number_of_monoms:length(monoms),
  70         if not( listp(aaaa) or length(monoms) > length(aaaa)) then
  71           aaaa: create_list(concat(aa,i),i,1,number_of_monoms),
  72         declare_scalar_list(aaaa),
  73         f:0,i:0,
  74         for u in monoms do (i:i+1,
  75           f: u*?nth(i,aaaa)+f),
  76         for v in elements do (
  77            temp : append (temp,[com(f,v,zet)])),
  78
  79         print(temp,monoms_higher),
  80         eqns:Extract_linear_equations (temp,monoms_higher),
  81         equations:eqns,
  82         if fast_solve then ans : fast_linsolve (eqns,firstn(length(monoms),aaaa)) else
  83         ans:linsolve(eqns,firstn (length( monoms), aaaa)),
  84         return(sublis(ans,f)))$
  85
  86 dot_simplifications:false$
  87 central_elements_and_check_dot_simplifications(variables,deg):=block([expop:100],
  88         check_overlapping_dot_simplifications(1+deg,true),
  89         central_elements(variables,deg))$
  90
  91 free_ideal_generated_by(a_list,variables):=block([temp],temp:[],
  92         for v in variables do
  93         (temp:append(temp,create_list(v.u,u,a_list),
  94         create_List(u.v,u,a_list))),temp)$
  95
  96  /*  %beta^2=-%beta*%zeta8-%zeta8^2 ,%alpha2^2=-2*%alpha2-2,
  97         %alpha3^3=-(-5*%alpha3^2+9*%alpha3-1),%alpha^3-23*%alpha^2+210*%alpha+539)*/
  98 rtx3:
  99                           matrix([x],[y],[z])$
 100  p3:
 101            matrix([p11,p12,p13],[p21,p22,p23],[p31,p32,p33])$
 102 declare_scalar_list(list_matrix_entries(p3))$
 103 new_fast_dotsimp:true$
 104
 105 tchirn(f,x):=ratsimp(sublis([x = x-ratcoef(f,x,2)/(ratcoef(f,x,3)*3)],f))$
 106
 107                               /* algorithm for j of ternary cubic
 108 eliminate the y^3 term by z=z+eps*y where AA1*EPS^3+3*BB2*EPS^2+3*CC1*EPS+DD0
 109 If have y^2*L+y*L'*L''=h(x,z)
 110 then can make replacement like Z=L,X=x,Y=y to make it y^2*z+y*L'*l''
 111  then have y^2*z* to deal with?*/
 112       /* for a binary quartic put one of the roots at infinity to get cubic*/
 113 j(g2,g3):=ratsimp(g2^3/(g2^3-27*g3^2))$
 114 /* note that this should really be multiplied by 1728 */
 115 jcubic(g,x):=block([ fac,f],f:tchirn(g,x),fac:1/4*ratcoef(f,x,3),
 116         if fac=0 then (print("not degree 3"),infinity)
 117         else if poly_discriminant(ratsimp(numerator(g)),x)=0
 118          then (print("discriminant is zero", infinity))
 119          else
 120           ratsimp(j(-ratcoef(f/fac,x),-sublis([x = 0],f/fac))))$
 121
 122 /*j-invariant of  y^2=x^4+bb*x^2+cc (with 1728 factor) */
 123
 124  y^2=x^4+bb*x^2+cc;
 125
 126 jinv(bb,cc):=16*(12*c^2+bb^2)^3/(c^2*(2*c-bb)^2*(2*c+bb)^2)$
 127
 128
 129
 130                               /* cases  April 20 1985*/
 131
 132 %alpha-1$
 133  m300:matrix([x,0,0],[0,0,0],[0,0,0])$
 134
 135 %beta-1$
 136  m030:matrix([0,0,0],[0,y,0],[0,0,0])$
 137
 138 %gamma-1$
 139  m003: matrix([0,0,0],[0,0,0],[0,0,z])$
 140
 141
 142
 143 %alpha^2*%beta-1$
 144 m210:matrix([y,x/%alpha,0],[%alpha*x,0,0],[0,0,0])$
 145
 146 (%alpha*%beta^2-1)$
 147 m120:
 148 matrix([0,%beta*y,0],[y/%beta,x,0],[0,0,0])$
 149
 150 %alpha^2*%gamma-1$
 151 m201:
 152 matrix([z,0,x/%alpha],[0,0,0],[%alpha*x,0,0])$
 153
 154
 155 %alpha*%gamma^2-1$
 156 m102:matrix([0,0,%gamma*z],[0,0,0],[z/%gamma,0,x])$
 157
 158
 159 %gamma*%beta^2-1$
 160 m021:matrix([0,0,0],[0,z,y/%beta],[0,%beta*y,0])$
 161 hhh(x):= x^2+3$
 162 %beta*%gamma^2-1$
 163 m012:matrix([0,0,0],[0,0,%gamma*z],[0,z/%gamma,y])$
 164
 165
 166
 167 %alpha*%beta*%gamma-1$        /* agree with manuscript Mar 85*/
 168 m111p:matrix([0,%beta*z,0],[0,0,x],[%alpha*%beta*y,0,0])$
 169 m111: matrix([0,0,y],[%alpha*z,0,0],[0,%alpha*%beta*x,0])$
 170 type3a:sublis([%beta=1,%alpha=1],m300+m030+m003+a*m111+b*m111p);
 171 type3b:sublis([%beta=1,%alpha=1,%gamma=-1],m210+m120+m102+a*m012);
 172
 173
 174