_mode=POLY$
_cpu=1;
_language=en$
_modenum=NUMRATMP$

macro userdet[a,n]{

  private cr,ct,cs,cu,cv$
  private xx,yy,z,s0,s1$

  if(n==1) then{
     return a[1,1]$
  }$

  dim xx[1:n+1]$
  xx[1]=1$
  xx[2]=-a[1,1]$

  dim z[1:n+1,1:n]$

  z[1,1]=1$
  for ct=2 to n {
    z[1,ct]=0$
  }$

  z[2,2]=1$
  for ct=3 to n {
    z[2,ct]=0$
  }$

  dim s0[1:n-1]$
  
  for cr=2 to n-1 {
    z[2,1]=-a[cr,cr]$

    for cs=3 to cr+1 {
    
      z[cs,2:cr]=z[cs-1,1:cr-1]$

      if (cs==3) then{
        z[cs,1]=-(sum cv=1 to cr-1 {
	  return (a[cr,cv]*a[cv,cr])$
	})$
	if (cs<cr+1) then {
          s0[1:cr-1]=a[1:cr-1,1:cr-1]&*a[1:cr-1,cr:cr]$
        }$
      }
      else{
        z[cs,1]=-(sum cv=1 to cr-1 {
	  return (a[cr,cv]*s0[cv])$
	})$
	if (cs<cr+1) then {
          s0[1:cr-1]=a[1:cr-1,1:cr-1]&*s0[1:cr-1]$
	}$
      }$
      
    }$

    xx[1:cr+1]=z[1:cr+1,1:cr]&*xx[1:cr]$
  }$

  z[n+1,n]=-a[n,n]$
  for cs=3 to n+1 {

    if (cs==3) then{
      z[n+1,n-1]=-(sum cv=1 to n-1 {
        return (a[n,cv]*a[cv,n])$
      })$
      if (cs<n+1) then {
        s0[1:n-1]=a[1:n-1,1:n-1]&*a[1:n-1,n:n]$
      }$
    }
    else{
      z[n+1,n+2-cs]=-(sum cv=1 to n-1 {
        return (a[n,cv]*s0[cv])$
      })$
      if (cs<n+1) then {
        s0[1:n-1]=a[1:n-1,1:n-1]&*s0[1:n-1]$
      }$
    }$
      
  }$

  xx[n+1]=sum cv=1 to n {
    return (z[n+1,cv]*xx[cv])$
  }$
 
  if(mod(n,2)==1) then{
    return(-xx[n+1])$
  } else {
    return(xx[n+1])$
  }$

}$

macro sekimat[p1,m,p2,n]{

  private a,b,l,s,sekimat0,ii,j,q,tmp,f,g,h$

  if (m>=n) then{
    a=p1$
    b=p2$ 
    l=m$
    s=n$
  }
  else{
    a=p2$
    b=p1$
    l=n$
    s=m$
  }$

  dim sekimat0[1:l,1:l]$

  for ii=0 to l-s-1 {
    for j=0 to s {
      sekimat0[ii+1,j+ii+1]=coef_ext(b,(x,j))$
    }$
  }$

  b=x^(l-s)*b$

  div(a,x^l,q,f)$
  div(b,x^l,q,g)$
  h=coef_ext(a,(x,l))*g-coef_ext(b,(x,l))*f$
  for j=0 to l-1 {
    sekimat0[2*l-s-l+1,j+1]=coef_ext(h,(x,j))$
  }$

  for ii=l-1 to l-s+1 step(-1) {
    div(h,x^(l-1),q,tmp)$
    h=x*tmp+coef_ext(a,(x,ii))*g-coef_ext(b,(x,ii))*f$
    for j=0 to l-1 {
      sekimat0[2*l-s-ii+1,j+1]=coef_ext(h,(x,j))$
    }$
  }$ 
  
  return(sekimat0)$

}$

p=x^9+x^8+a7*x^7+a6*x^6+a5*x^5+a4*x^4+a3*x^3+a2*x^2+a1*x+a0$
q=deriv(p,x)$

time_s;
m9=%sekimat[p,9,q,8]$
w=%userdet[m9,9]$
size(w);
time_l;