[BACK]Return to Fgfs.c CVS log [TXT][DIR] Up to [local] / OpenXM_contrib2 / asir2000 / engine

Diff for /OpenXM_contrib2/asir2000/engine/Fgfs.c between version 1.11 and 1.12

version 1.11, 2002/11/22 08:44:57 version 1.12, 2002/11/26 07:09:45
Line 1 
Line 1 
 /* $OpenXM: OpenXM_contrib2/asir2000/engine/Fgfs.c,v 1.10 2002/11/22 07:32:10 noro Exp $ */  /* $OpenXM: OpenXM_contrib2/asir2000/engine/Fgfs.c,v 1.11 2002/11/22 08:44:57 noro Exp $ */
   
 #include "ca.h"  #include "ca.h"
   
Line 12  void partial_sqfrsf(VL vl,V v,P f,P *r,DCP *dcp);
Line 12  void partial_sqfrsf(VL vl,V v,P f,P *r,DCP *dcp);
 void gcdsf(VL vl,P *pa,int k,P *r);  void gcdsf(VL vl,P *pa,int k,P *r);
 void mfctrsfmain(VL vl, P f, DCP *dcp);  void mfctrsfmain(VL vl, P f, DCP *dcp);
 void next_evaluation_point(int *mev,int n);  void next_evaluation_point(int *mev,int n);
 void estimatelc_sf(VL vl,VL rvl,P c,DCP dc,P *lcp);  void estimatelc_sf(VL vl,VL rvl,P c,DCP dc,int *mev,P *lcp);
 void mfctrsf_hensel(VL vl,VL rvl,P f,P pp0,P u0,P v0,P lcu,P lcv,P *up);  void mfctrsf_hensel(VL vl,VL rvl,P f,P pp0,P u0,P v0,P lcu,P lcv,int *mev,P *up);
 void substvp_sf(VL vl,VL rvl,P f,int *mev,P *r);  void substvp_sf(VL vl,VL rvl,P f,int *mev,P *r);
 void shift_sf(VL vl, VL rvl, P f, int *mev, int sgn, P *r);  void shift_sf(VL vl, VL rvl, P f, int *mev, int sgn, P *r);
 void adjust_coef_sf(VL vl,VL rvl,P lcu,P u0,P *r);  void adjust_coef_sf(VL vl,VL rvl,P lcu,P u0,int *mev,P *r);
 void extended_gcd_modyk(P u0,P v0,V x,V y,int dy,P *cu,P *cv);  void extended_gcd_modyk(P u0,P v0,V x,V y,int dy,P *cu,P *cv);
 void poly_to_gfsn_poly(VL vl,P f,V v,P *r);  void poly_to_gfsn_poly(VL vl,P f,V v,P *r);
 void gfsn_poly_to_poly(VL vl,P f,V v,P *r);  void gfsn_poly_to_poly(VL vl,P f,V v,P *r);
Line 515  void mfctrsfmain(VL vl, P f, DCP *dcp)
Line 515  void mfctrsfmain(VL vl, P f, DCP *dcp)
         P gcd,g,df,dfmin;          P gcd,g,df,dfmin;
         P pa[2];          P pa[2];
         P f0,pp0,spp0,c,c0,x,y,u,v,lcf,lcu,lcv,u0,v0,t,s;          P f0,pp0,spp0,c,c0,x,y,u,v,lcf,lcu,lcv,u0,v0,t,s;
         P ype,yme;          P ype,yme,fin;
         GFS ev,evy;          GFS ev,evy;
         P *fp0;          P *fp0;
         int *mev,*win;          int *mev,*win;
Line 578  void mfctrsfmain(VL vl, P f, DCP *dcp)
Line 578  void mfctrsfmain(VL vl, P f, DCP *dcp)
         }          }
         NEXT(tvl) = 0;          NEXT(tvl) = 0;
   
           fin = f;
         reorderp(nvl,vl,f,&g); f = g;          reorderp(nvl,vl,f,&g); f = g;
         vx = nvl->v;          vx = nvl->v;
         vy = NEXT(nvl)->v;          vy = NEXT(nvl)->v;
Line 626  void mfctrsfmain(VL vl, P f, DCP *dcp)
Line 627  void mfctrsfmain(VL vl, P f, DCP *dcp)
   
         if ( !NEXT(dc) ) {          if ( !NEXT(dc) ) {
                 /* f is irreducible */                  /* f is irreducible */
                 NEWDC(dc); DEG(dc) = ONE; COEF(dc) = f; NEXT(dc) = 0;                  NEWDC(dc); DEG(dc) = ONE; COEF(dc) = fin; NEXT(dc) = 0;
                 *dcp = dc;                  *dcp = dc;
                 return;                  return;
         }          }
Line 645  void mfctrsfmain(VL vl, P f, DCP *dcp)
Line 646  void mfctrsfmain(VL vl, P f, DCP *dcp)
         lcf = COEF(DC(f));          lcf = COEF(DC(f));
         mfctrsf(nvl,lcf,&dct);          mfctrsf(nvl,lcf,&dct);
         /* skip the first element (= a number) */          /* skip the first element (= a number) */
         dct = NEXT(dct);          lcfdc = NEXT(dct);
   
         /* shift lcfdc; c <- c(X+mev) */  
         for ( lcfdc = 0; dct; dct = NEXT(dct) ) {  
                 NEXTDC(lcfdc,dcs);  
                 DEG(dcs) = DEG(dct);  
                 shift_sf(nvl,rvl,COEF(dct),mev,1,&COEF(dcs));  
         }  
         NEXT(dcs) = 0;  
   
         /* np = number of bivariate factors */          /* np = number of bivariate factors */
         for ( np = 0, dct = dc; dct; dct = NEXT(dct), np++ );          for ( np = 0, dct = dc; dct; dct = NEXT(dct), np++ );
         fp0 = (P *)ALLOCA((np+1)*sizeof(P));          fp0 = (P *)ALLOCA((np+1)*sizeof(P));
Line 664  void mfctrsfmain(VL vl, P f, DCP *dcp)
Line 657  void mfctrsfmain(VL vl, P f, DCP *dcp)
         l = tl = (P *)ALLOCA((np+1)*sizeof(P));          l = tl = (P *)ALLOCA((np+1)*sizeof(P));
         win = W_ALLOC(np+1);          win = W_ALLOC(np+1);
   
         /* f <- f(X+mev) */  
         shift_sf(nvl,rvl,f,mev,1,&s); f = s;  
   
         for ( k = 1, win[0] = 1, --np; ; ) {          for ( k = 1, win[0] = 1, --np; ; ) {
                 itogfs(1,&u0);                  itogfs(1,&u0);
                 /* u0 = product of selected factors */                  /* u0 = product of selected factors */
Line 675  void mfctrsfmain(VL vl, P f, DCP *dcp)
Line 665  void mfctrsfmain(VL vl, P f, DCP *dcp)
                 }                  }
                 /* we have to consider the content */                  /* we have to consider the content */
                 /* f0 = c0*u0*v0 */                  /* f0 = c0*u0*v0 */
                 mulp(nvl,LC(u0),c0,&c); estimatelc_sf(nvl,rvl,c,lcfdc,&lcu);                  mulp(nvl,LC(u0),c0,&c); estimatelc_sf(nvl,rvl,c,lcfdc,mev,&lcu);
                 divsp(nvl,pp0,u0,&v0);                  divsp(nvl,pp0,u0,&v0);
                 mulp(nvl,LC(v0),c0,&c); estimatelc_sf(nvl,rvl,c,lcfdc,&lcv);                  mulp(nvl,LC(v0),c0,&c); estimatelc_sf(nvl,rvl,c,lcfdc,mev,&lcv);
                 mfctrsf_hensel(nvl,rvl,f,pp0,u0,v0,lcu,lcv,&u);                  mfctrsf_hensel(nvl,rvl,f,pp0,u0,v0,lcu,lcv,mev,&u);
                 if ( u ) {                  if ( u ) {
                         /* save the factor */                          /* save the factor */
                         reorderp(vl,nvl,u,&t);                          reorderp(vl,nvl,u,&t);
                         /* x -> x-mev, y -> y-evy */                          /* y -> y-evy */
                         shift_sf(vl,rvl,t,mev,-1,&s); substp(vl,s,vy,yme,tl++);                          substp(vl,t,vy,yme,tl++);
   
                         /* update f,pp0 */                          /* update f,pp0 */
                         divsp(nvl,f,u,&t); f = t;                          divsp(nvl,f,u,&t); f = t;
Line 714  void mfctrsfmain(VL vl, P f, DCP *dcp)
Line 704  void mfctrsfmain(VL vl, P f, DCP *dcp)
                 }                  }
         }          }
         reorderp(vl,nvl,f,&t);          reorderp(vl,nvl,f,&t);
         /* x -> x-mev, y -> y-evy */          /* y -> y-evy */
         shift_sf(vl,rvl,t,mev,-1,&s); substp(vl,s,vy,yme,tl++);          substp(vl,t,vy,yme,tl++);
         *tl = 0;          *tl = 0;
         for ( dc0 = 0, i = 0; l[i]; i++ ) {          for ( dc0 = 0, i = 0; l[i]; i++ ) {
                 NEXTDC(dc0,dc); DEG(dc) = ONE; COEF(dc) = l[i];                  NEXTDC(dc0,dc); DEG(dc) = ONE; COEF(dc) = l[i];
Line 742  void next_evaluation_point(int *e,int n)
Line 732  void next_evaluation_point(int *e,int n)
   
 /*  /*
  * dc : f1^E1*...*fk^Ek   * dc : f1^E1*...*fk^Ek
  * find e1,...,ek s.t. fi(0)^ei | c   * find e1,...,ek s.t. fi(mev)^ei | c
  * and return f1^e1*...*fk^ek   * and return f1^e1*...*fk^ek
  * vl = (vx,vy,rvl)   * vl = (vx,vy,rvl)
  */   */
   
 void estimatelc_sf(VL vl,VL rvl,P c,DCP dc,P *lcp)  void estimatelc_sf(VL vl,VL rvl,P c,DCP dc,int *mev,P *lcp)
 {  {
         DCP dct;          DCP dct;
         P r,c1,c2,t,s,f;          P r,c1,c2,t,s,f;
Line 758  void estimatelc_sf(VL vl,VL rvl,P c,DCP dc,P *lcp)
Line 748  void estimatelc_sf(VL vl,VL rvl,P c,DCP dc,P *lcp)
                 if ( NUM(COEF(dct)) )                  if ( NUM(COEF(dct)) )
                         continue;                          continue;
                 /* constant part */                  /* constant part */
                 substvp_sf(vl,rvl,COEF(dct),0,&f);                  substvp_sf(vl,rvl,COEF(dct),mev,&f);
                 d = QTOS(DEG(dct));                  d = QTOS(DEG(dct));
                 for ( i = 0, c1 = c; i < d; i++ )                  for ( i = 0, c1 = c; i < d; i++ )
                         if ( !divtp(vl,c1,f,&c2) )                          if ( !divtp(vl,c1,f,&c2) )
Line 819  void shift_sf(VL vl, VL rvl, P f, int *mev, int sgn, P
Line 809  void shift_sf(VL vl, VL rvl, P f, int *mev, int sgn, P
  * pp(f(0)) = u0*v0   * pp(f(0)) = u0*v0
  */   */
   
 void mfctrsf_hensel(VL vl,VL rvl,P f,P pp0,P u0,P v0,P lcu,P lcv,P *up)  void mfctrsf_hensel(VL vl,VL rvl,P f,P pp0,P u0,P v0,P lcu,P lcv,int *mev,P *up)
 {  {
         VL tvl,onevl;          VL tvl,onevl;
         P t,s,w,u,v,ff,si,wu,wv,fj,cont;          P t,s,w,u,v,ff,si,wu,wv,fj,cont;
Line 828  void mfctrsf_hensel(VL vl,VL rvl,P f,P pp0,P u0,P v0,P
Line 818  void mfctrsf_hensel(VL vl,VL rvl,P f,P pp0,P u0,P v0,P
         int dy,n,i,dbd,nv,j;          int dy,n,i,dbd,nv,j;
         int *md;          int *md;
         P *uh,*vh;          P *uh,*vh;
         P x,du0,dv0,m,q,r;          P x,du0,dv0,m,q,r,fin;
         P *cu,*cv;          P *cu,*cv;
         GFSN inv;          GFSN inv;
   
         /* adjust coeffs */          /* adjust coeffs */
         /* u0 = am x^m+ ... -> lcu*x^m + a(m-1)*(lcu(0)/am)*x^(m-1)+... */          /* u0 = am x^m+ ... -> lcu*x^m + a(m-1)*(lcu(mev)/am)*x^(m-1)+... */
         /* v0 = bm x^l+ ... -> lcv*x^l + b(l-1)*(lcv(0)/bl)*x^(l-1)+... */          /* v0 = bm x^l+ ... -> lcv*x^l + b(l-1)*(lcv(mev)/bl)*x^(l-1)+... */
         /* f                -> lcu*lcv*x^(m+l)+... */          /* f                -> lcu*lcv*x^(m+l)+... */
         adjust_coef_sf(vl,rvl,lcu,u0,&u);          adjust_coef_sf(vl,rvl,lcu,u0,mev,&u);
         adjust_coef_sf(vl,rvl,lcv,v0,&v);          adjust_coef_sf(vl,rvl,lcv,v0,mev,&v);
         mulp(vl,lcu,lcv,&t); divsp(vl,t,LC(f),&m); mulp(vl,m,f,&t); f = t;          mulp(vl,lcu,lcv,&t); divsp(vl,t,LC(f),&m); mulp(vl,m,f,&t); f = t;
   
           /* f <- f(X+mev), u <- u(X+mev), v <- v(X+mev) */
           fin = f;
           shift_sf(vl,rvl,f,mev,1,&s); f = s;
           shift_sf(vl,rvl,u,mev,1,&s); u = s;
           shift_sf(vl,rvl,v,mev,1,&s); v = s;
   
         vx = vl->v; vy = NEXT(vl)->v;          vx = vl->v; vy = NEXT(vl)->v;
         n = getdeg(vx,f);          n = getdeg(vx,f);
         dy = getdeg(vy,f)+1;          dy = getdeg(vy,f)+1;
Line 930  void mfctrsf_hensel(VL vl,VL rvl,P f,P pp0,P u0,P v0,P
Line 926  void mfctrsf_hensel(VL vl,VL rvl,P f,P pp0,P u0,P v0,P
                 }                  }
                 if ( !wu ) {                  if ( !wu ) {
                         gfsn_poly_to_poly(vl,u,vy,&t);                          gfsn_poly_to_poly(vl,u,vy,&t);
                         if ( divtp(vl,f,t,&q) ) {                          shift_sf(vl,rvl,t,mev,-1,&s);
                                 cont_pp_mv_sf(vl,onevl,t,&cont,up);                          if ( divtp(vl,fin,s,&q) ) {
                                   cont_pp_mv_sf(vl,onevl,s,&cont,up);
                                 return;                                  return;
                         }                          }
                 }                  }
                 if ( !wv ) {                  if ( !wv ) {
                         gfsn_poly_to_poly(vl,v,vy,&t);                          gfsn_poly_to_poly(vl,v,vy,&t);
                         if ( divtp(vl,f,t,&q) ) {                          shift_sf(vl,rvl,t,mev,-1,&s);
                           if ( divtp(vl,fin,s,&q) ) {
                                 cont_pp_mv_sf(vl,onevl,q,&cont,up);                                  cont_pp_mv_sf(vl,onevl,q,&cont,up);
                                 return;                                  return;
                         }                          }
Line 949  void mfctrsf_hensel(VL vl,VL rvl,P f,P pp0,P u0,P v0,P
Line 947  void mfctrsf_hensel(VL vl,VL rvl,P f,P pp0,P u0,P v0,P
         }          }
 }  }
   
 void adjust_coef_sf(VL vl,VL rvl,P lcu,P u0,P *r)  void adjust_coef_sf(VL vl,VL rvl,P lcu,P u0,int *mev,P *r)
 {  {
         P lcu0,cu;          P lcu0,cu;
         DCP dc0,dcu,dc;          DCP dc0,dcu,dc;
   
         substvp_sf(vl,rvl,lcu,0,&lcu0);          substvp_sf(vl,rvl,lcu,mev,&lcu0);
         divsp(vl,lcu0,LC(u0),&cu);          divsp(vl,lcu0,LC(u0),&cu);
         for ( dc0 = 0, dcu = DC(u0); dcu; dcu = NEXT(dcu) ) {          for ( dc0 = 0, dcu = DC(u0); dcu; dcu = NEXT(dcu) ) {
                 if ( !dc0 ) {                  if ( !dc0 ) {
Legend:
Removed from v.1.11  
changed lines
  Added in v.1.12
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>