Return to tk_nm_dn.rr CVS log

Up to [local] / OpenXM / src / asir-contrib / packages / src

generate a sample problem: L=gtt_ekn.prob1(3,5,10 | factor=1, factor_row=3);
get a list of degrees of rational functions:  tk_nm_dn([x/(1+x),(x^3+1)/x,0],x);

/* $OpenXM: OpenXM/src/asir-contrib/packages/src/tk_nm_dn.rr,v 1.3 2019/02/15 05:22:44 takayama Exp $ */
import("names.rr")$

#define USE_MODULE 1

#ifdef USE_MODULE
module tk_nm_dn$
static Tk_nm_dn_assert$
#else
extern Tk_nm_dn_assert$
Tk_nm_dn_assert=1$
#endif

#ifdef USE_MODULE
localf number_nm_dn_0$
localf mul_list$
localf number_nm_dn_list$
localf number_nm_dn$
localf poly_nm_dn_list$
localf poly_nm_dn$
localf set_assert$
localf nm_dn$
localf rdeg$
localf nm_dn_assert$
#else
#endif

/*
  F=G/D;
  map(number_nm_dn_0,M=map(nm,newmat(2,2,mm(a,-10,5,13/51)[0])));
*/
def number_nm_dn_0(F) {
  if (F==0) return [0,1];
  if (type(F)<2) {
    return [nm(F),dn(F)];
  }
  G=ptozp(F |factor=0);  // F=G[0]*G[1];
  D=dn(G[1]); G=F*D; 
  // G=F*D;
  return [G,D];
}
def mul_list(V,S) {  
  if (type(V) < 4) return(red(S*V));
  if (type(V) == 4) {
    return map(mul_list,V,S);
  }
  return map(red,S*V);
}
/*
  number_nm_dn_list([[x/2,(x+y)/3],(x+1)/7]);
 */
def number_nm_dn_list(F) {
  if (type(F) < 4) {
    return(number_nm_dn_0(F));
  }
  D=1;
  for (I=0; I<length(F);I++)  {
    D1=number_nm_dn(F[I])[1];
    D=ilcm(D,D1);
  }
  G=[];
  for (I=0; I<length(F); I++) {
    H=number_nm_dn(mul_list(F[I],D));
    if (type(H[1])>1) {print("Error: number_nm_dn_list"); debug();}
    G=cons(H[0],G);
  }
  return([reverse(G),D]);
}

/*
nm_dn([a-1163/56,-551/560]);
*/
def number_nm_dn(F) {
  if (type(F) < 4) return number_nm_dn_0(F);
  if (type(F)==4) return number_nm_dn_list(F);
  if ((type(F)==5) || (type(F)==6)) {
    G=number_nm_dn_list(matrix_matrix_to_list(F));
    return [matrix_list_to_matrix(G[0]),G[1]];
  }
  print("Error: number_nm_dn");
  debug();
}

def poly_nm_dn_list(F) {
  if (type(F) < 4) {
    return(poly_nm_dn(F));
  }
  D=1;
  for (I=0; I<length(F);I++)  {
    D1=poly_nm_dn(F[I])[1];
    D=lcm(D,D1);
  }
  G=[];
  for (I=0; I<length(F); I++) {
    H=poly_nm_dn(mul_list(F[I],D));
    if (type(H[1])>1) {print("Error: poly_nm_dn_list");debug();}
    G=cons(H[0],G);
  }
  return([reverse(G),D]);
}

/*
H=poly_nm_dn(newmat(2,2,[[x/(y+1),y],[y/x,y/(y+1/2)]]));
H2=poly_nm_dn(newmat(2,2,[[x/(y+1),y],[y/x,0]]));
*/
def poly_nm_dn(F) {
  if (type(F)<2) return [F,1];
  else if (type(F)<4) return [nm(F),dn(F)];
  else if (type(F)==4) return poly_nm_dn_list(F);
  else if ((type(F)==5) || (type(F)==6)) {
    G=poly_nm_dn_list(matrix_matrix_to_list(F));
    return [matrix_list_to_matrix(G[0]),G[1]];
  }
  print("Error: poly_nm_dn");
  debug();
}

def set_assert(L) {
  Tk_nm_dn_assert=L;
}
/*
 set_assert(1);
 nm_dn([1/2*x-1,0,5]);
 nm_dn([1/2*x-1,0,1/3]);
*/
def nm_dn(F) 
"nm_dn(F) returns [numerator (list, vector) M,denominator D]. F=M/D. Example: nm_dn([1/2*x-1,0,1/3]);"
{
  G=poly_nm_dn(F);
  H=number_nm_dn(G[0]);
  Ans=[H[0],H[1]*G[1]];
  if (Tk_nm_dn_assert) {
     nm_dn_assert(F,Ans);
  }
  return(Ans);
}

def nm_dn_assert(F,Ans) {
  G=Ans[0]; D=Ans[1];
  if (type(G)<4) {
    if (red(G/D - F) != 0) {
      print("Assert fails:"); debug();
      return(-1);
    }
    return(0);
  }
  if (type(F) == 4) {
    for (I=0; I<length(F); I++) {
      return(nm_dn_assert(F[I],[Ans[0][I],D]));
    }
  }
  G=matrix_list_to_matrix(G);
  E=map(red,G/D-F);
  if (matrix_is_zero(E)) return(0);
  else {
    print("Assert fails:"); debug();
    return(-1);
  }
}

def rdeg(F,V) {
  if (type(F) <3) return(deg(F,V));
  if (type(F)==3) return([deg(nm(F),V),deg(dn(F),V)]);
  if (type(F)>=4) return(map(rdeg,F,V));
}

#ifdef USE_MODULE
endmodule$
#endif


end$