Return to toric.c CVS log

Up to [local] / OpenXM_contrib / TiGERS_0.9

Annotation of OpenXM_contrib/TiGERS_0.9/toric.c, Revision 1.1

1.1     ! maekawa     1: /*
        !             2: ** toric.c                                 Birk Huber, 4/99
        !             3: ** -- support for computing first grobner bases of toric ideals.
        !             4: ** -- includes functions gset_rlex_rgb() and gset_compute_colon()
        !             5: **    which once tested and working will move back to the gset.c file.
        !             6: **
        !             7: **
        !             8: ** TiGERS,  Toric Groebner Basis Enumeration by Reverse Search
        !             9: ** copyright (c) 1999  Birk Huber
        !            10: **
        !            11: */
        !            12: #include<stdio.h>
        !            13: #include<stdlib.h>
        !            14: #include "utils.h"
        !            15: #include "gset.h"
        !            16: #include "matrices.h"
        !            17: #include "Rsimp.h"
        !            18:
        !            19:
        !            20:
        !            21:
        !            22: /*
        !            23: ** void gset_compute_colon(gset g,int lv):
        !            24: **    Compute deg rev lex rgb (lv least) for colon of g with variable lv.
        !            25: **    - set ring_lv=lv.
        !            26: **    - convert g to rgb wrt deg_rev_lex as implemented in binomial.h.
        !            27: **    - divide each binomial by highest power of var lv possible.
        !            28: **    - autoreduce.
        !            29: **
        !            30: */
        !            31: #define min(m,n) (((m)<(n))?(m):(n))
        !            32: void gset_compute_colon(gset g,int lv){
        !            33:   binomial ptr;
        !            34:   int lold,jk;
        !            35:
        !            36:   /* progress report */
        !            37:   fprintf(stderr,"taking colon of J with %c:\n",'a'+lv);
        !            38:
        !            39:   /* set lex last variable to be lv */
        !            40:   lold=ring_lv;
        !            41:   ring_lv=lv;
        !            42:
        !            43:   /* compute grobner basis wrt to this deg rev lex term order */
        !            44:   gset_rgb(g,monomial_rlexcomp);
        !            45:
        !            46:   /* now devide each binomial by highest power of lvar which divides it */
        !            47:   for(ptr=gset_first(g); ptr!=0; ptr=binomial_next(ptr)){
        !            48:     jk=min(binomial_lead(ptr)[lv],binomial_trail(ptr)[lv]);
        !            49:     binomial_lead(ptr)[lv]-=jk;
        !            50:     binomial_trail(ptr)[lv]-=jk;
        !            51:   }
        !            52:
        !            53:   /* result is a grobner basis, reduce it to make it an rgb */
        !            54:   gset_autoreduce(g);
        !            55:
        !            56:   /* restore original value of lvar*/
        !            57:   ring_lv=lold;
        !            58: }
        !            59:
        !            60: /*
        !            61: ** gset gset_toric_ideal(int **M,int m, int n):
        !            62: **  Given an mxn integer matrix M -- compute an rgb for the toric ideal I_M.
        !            63: **  Uses repeated colon computations.
        !            64: */
        !            65: gset gset_toric_ideal(int **M,int m, int n){
        !            66:    int i,j,jk,crk;
        !            67:    binomial b;
        !            68:    int **V,**S,**new_imatrix();
        !            69:    gset g;
        !            70:
        !            71:    /* set ring dimension*/
        !            72:    if (ring_N!=n){
        !            73:     fprintf(stderr,"ERROR gset_toric_ideal():");
        !            74:     fprintf(stderr," matrix dimensions incompatible with ring\n");
        !            75:     return 0;
        !            76:    }
        !            77:    /* set wieghting of A */
        !            78:    for(i=0;i<n;i++){
        !            79:     ring_weight[i]=0;
        !            80:     for(j=0;j<m;j++) ring_weight[i]+=M[j][i];
        !            81:    }
        !            82:
        !            83:    /* reserve space for S and V */
        !            84:    S=new_imatrix(n,m);
        !            85:    V=new_imatrix(n,n);
        !            86:
        !            87:    /* copy transpose of M to S*/
        !            88:    for(i=0;i<n;i++)for(j=0;j<m;j++) IMref(S,i,j)=IMref(M,j,i);
        !            89:
        !            90:    /* compute hermite normal form of M,S,V and corank of M*/
        !            91:    crk=ihermite(S,V,m,n);
        !            92:
        !            93:    /* create new gset and read off equations corresponding to rows of V*/
        !            94:    g=gset_new();
        !            95:    for(i=crk;i<n;i++){
        !            96:      b=binomial_new();
        !            97:        for(j=0;j<n;j++){
        !            98:          jk=IMref(V,i,j);
        !            99:          if (jk>0){
        !           100:              binomial_lead(b)[j]=jk; binomial_trail(b)[j]=0;}
        !           101:          else {binomial_trail(b)[j]=-1*jk; binomial_lead(b)[j]=0;}
        !           102:         }
        !           103:      gset_insert(g,b);
        !           104:    }
        !           105:
        !           106:    /* compute successive colon ideals */
        !           107:    for(i=0;i<n;i++)gset_compute_colon(g,i);
        !           108:
        !           109:    /*unset weighting */
        !           110:    for(i=0;i<n;i++) ring_weight[i]=1;
        !           111:
        !           112:    /* free space used by matrices */
        !           113:    free_imatrix(M);
        !           114:    free_imatrix(S);
        !           115:    free_imatrix(V);
        !           116:
        !           117:    /* make sure result is the lex gbasis */
        !           118:    gset_rgb(g,monomial_lexcomp);
        !           119:
        !           120: return g;
        !           121: }
        !           122:
        !           123:
        !           124:
        !           125:
        !           126:
        !           127:
        !           128:
        !           129:
        !           130:
        !           131:
        !           132:
        !           133: