1  Solving from the operators of multiplication

1. The eigenvalues of the linear operator M_p (resp. M_p^t) are {p(z₁), ...,p(z_d)}.
2. The common eigenvectors of (M_xi^t) are (up to a scalar) 1_z1, ..., 1_zd.
3. When m=n, the common eigenvectors of (M_xi) are (up to a scalar) J(x)e₁, ..., J(x)e_d, where J(x) is the Jacobian of p₁, ..., p_n.

1. Compute the transposed of the matrix of multiplication M_a for a in A.
2. Compute its eigenvectors v_i= (v_i,1, v_i,x1, ..., v_i,xn, ...) for i=1,...,d.
3. For i=1,...,d, compute and output

   zi = (

   v   i,x1 vi,1

   , ...,

   v   i,xn vi,1

   ).

Implementation:

<B><PRE>
 E =Solve(M,n,Eigen()) 
</PRE></B>

2  Solving from matrix of univariate polynomials

• If the system of n equations f₁=···=f_n=0, wher f_i is considered as a polynomial in x₁, ..., x_n-1, with parameter x_n has a solution, then the ``resultant'' (to be precised) vanishes. So any multiple of the resultant must vanish. Let M(x_n) be a resultant matrix of the system (see X), where we consider the variable x_n as hidden in the coefficients field. By construction, M(x_n) is the matrix of multiples of (f₁,...,f_n) in a monomial basis w=(x^a)_{a in F}. Thus we have w   M(x_n) =0 at a root of our system (and usually by construction of M, w ¹ 0). If we solve det(M(x_n))=0, we obtain the values of the coordinnate x_n of the roots. If we compute the corresponding w, we obtain the other coordinates.
• Another approach consists in adding a new polynomial, eg f=u₀+u₁x₁+ ... + u_nx_n, where u₀ is an indeterminate. Then we compute a resultant matrix M(u₀) of these n+1 polynomials, and solve the generalized eigenproblem w  M(u₀)=0.

Implementation:

<B><PRE>
    E = Solve$<$ MatC$>$ (Mx, Compagnon()) 
    </PRE></B>

3  Computing the matrix of multiplication from resultant matrices

1. its rows are indexed by monomials (x^a)_{a in F},
2. the set of monomials B₀ = (x^a)_{a in E0} indexing the rows of the block (A B) is a basis of A=R/(f₁,...,f_m),
3. the columns of (

   A

   C

   ) represents the elements x^a f₀ for a in E₀, expressed as linear combinations of the monomials (x^b)_{b in F},
4. the columns of (

   B

   D

   ) represent some multiples of the polynomials f₁, ...,f_m, expressed as linear combinations of the monomials (x^b)_{b in F},
5. the block D is invertible.

Implementation:

   M= USchur(S);

4  An example based on resultant matrices

#include <list>
#include "mpoly.H"
#include "mpoly/resultant.H"
#include "lapack.H"
#include "linalg.H"

typedef Monom<double, dynamicexp<'x'> > Mon;
typedef MPoly<vector<Mon>, Dlex<Mon> > Pol;
typedef MatDense<lapack<double> > Matr;
typedef VectStd<array1d<double> > Vect;
typedef VectStd<array1d<complex<double> > > VectC;

int main (int argc, char **argv)
{

  int n;
  cin >>n;
  
  array1d<Pol> l(n);
  for(int i=0; i< n;i++) cin >> l[i];

  cout <<"Polynomes:"<<endl;
  for(int i=0; i< n;i++) cout<<l[i]<<";"<<endl;

  int d=1;
  for(int i=1; i< n;i++) d*=Degree(l[i]);
  cout <<"Bezout bound:"<<d<<endl;

  Matr S=Macaulay<lapack<double> >(l,'N');
  cout <<"Macaulay matrix of size "<< S.nbrow()<<endl;

  for(unsigned int i=d;  i<S.nbrow(); i++) 
    if(S(i,i)==0) 
      cout<<"Probleme d'ordre des polynomes "<<i<<" "<<i <<endl; 
  
  Matr D = S(Range(d,S.nbrow()-1),Range(d,S.nbcol()-1));
   
  cout <<"Lower diag. block of size "<< D.nbrow()<<endl;

  cout <<"Singular values of D: "<<endl;
  Vect Sv= Svd(D);
  cout << Sv<<endl;

  //  cout <<"Singular values of S: "<<endl;
  //  cout <<Eval<Vect>(Svd(S))<<endl;
  
  Matr U= USchur(S,d);
  cout <<"Eigenvalues of the Schur complement: "<<endl;
  cout << Eigenval(U)<<endl;
  
}