Return to mvol.cpp CVS log

Up to [local] / OpenXM / src / ox_cdd

Annotation of OpenXM/src/ox_cdd/mvol.cpp, Revision 1.1

1.1     ! noro        1: /*     $Id: mvol.cpp,v 1.4 2004/11/12 04:37:15 fujiwara Exp $  */
        !             2:
        !             3: /*
        !             4: * ##  origianl  #################################################
        !             5: MixedVol: Polynomial system mixed volume calculation program.
        !             6: * Version 1.0 Copyright (C) 2003  Tangan Gao, T. Y. Li, Xing Li, and Mengnien Wu.
        !             7: * ###############################################################
        !             8: */
        !             9:
        !            10: #include <cstdlib>
        !            11: #include <string>
        !            12: #include <ctime>
        !            13: #include <fstream>
        !            14: #include <iostream>
        !            15:
        !            16: #include "GetOption/getopt.h"
        !            17: #include "PolyReader/PolynomialSystemReader.h"
        !            18: #include "PolyReader/PolynomialException.h"
        !            19: #include "PreProcess/Pre4MV.h"
        !            20: #include "MixedVol/MixedVol.h"
        !            21:
        !            22: extern int debug_print;
        !            23:
        !            24: extern "C" int mvol(int nVar,int *Spt1stIdx,int **Spt)
        !            25: {
        !            26:        int i, j=0, k;
        !            27:        int nPts = 0;
        !            28:        int nSpt;
        !            29:
        !            30:        nPts = Spt1stIdx[nVar];
        !            31:
        !            32:        for (i=nVar-1; i>=0; i--){
        !            33:                Spt1stIdx[i] = Spt1stIdx[i+1] - Spt1stIdx[i];
        !            34:        }
        !            35:
        !            36:        // Quick return if 1-variable system or less than two terms
        !            37:        k = -1;
        !            38:        for (i=0; i<nVar; i++){
        !            39:                if ( Spt1stIdx[i+1]-Spt1stIdx[i] < 2 )  {
        !            40:                        k = i;
        !            41:                        break;
        !            42:                }
        !            43:        }
        !            44:        if ( k >= 0 ){
        !            45:                if( debug_print >= 2 ){
        !            46:                        cout << "The " << k+1 << "-th support has less than 2 points" << endl;
        !            47:                }
        !            48:                return 0;   // end of the case: too few terms
        !            49:        }
        !            50:
        !            51:        if ( nVar == 1 ){
        !            52:                int kmin, kmax;
        !            53:                kmin = INT_MAX/2;
        !            54:                kmax = -INT_MAX/2;
        !            55:                for (i=0; i<Spt1stIdx[1]; i++){
        !            56:                        kmax = max(kmax, Spt[i][0]);
        !            57:                        kmin = min(kmin, Spt[i][0]);
        !            58:                }
        !            59:                if( debug_print >= 2 ){
        !            60:                        cout << "A 1-dimensional support, its mixed volume is " << kmax-kmin << endl;
        !            61:                }
        !            62:                return (kmax-kmin);   // end of the case: 1-variable
        !            63:        }
        !            64:        for (i=0; i<Spt1stIdx[nVar]; i++){
        !            65:                int kmin, kmax;
        !            66:                kmin = INT_MAX/2;
        !            67:                kmax = -INT_MAX/2;
        !            68:                k=0;
        !            69:                for (j=0; j<nVar; j++){
        !            70:                        kmax = max(kmax, Spt[i][j]);
        !            71:                        kmin = min(kmin, Spt[i][j]);
        !            72:                        k += abs(Spt[i][j]);
        !            73:                }
        !            74:                if ( kmin != 0 || kmax > 1 || k > 1 ){
        !            75:                        j = -1;
        !            76:                        break;
        !            77:                }
        !            78:        }
        !            79:        if ( j != -1 ){
        !            80:                if( debug_print >= 2 ){
        !            81:                        cout << "A linear support, its mixed volume <= 1" << endl;
        !            82:                }
        !            83:                return 1;
        !            84:        }
        !            85:        // end of quick return
        !            86:
        !            87:        // To preprocess the support
        !            88:
        !            89:        int* SptType = new int [nVar];
        !            90:        int* SptVtx1stIdx = new int [nVar+1];
        !            91:        int** SptVtx = new int* [nPts];
        !            92:        SptVtx[0] = new int [nVar*nPts];
        !            93:        for (i=1; i<nPts; i++) SptVtx[i] = SptVtx[i-1] + nVar;
        !            94:        int* NuIdx2OldIdx = new int [nPts];
        !            95:        nSpt = nVar;
        !            96:
        !            97:        Pre4MV(nVar,nSpt,SptType,Spt,Spt1stIdx,SptVtx,SptVtx1stIdx,NuIdx2OldIdx);
        !            98:
        !            99:        // To calculate the mixed volume of the support
        !           100:
        !           101:        srand( unsigned(time(0)) );
        !           102:        double* lft = new double[SptVtx1stIdx[nSpt]];
        !           103:        for (j=0; j<SptVtx1stIdx[nSpt]; j++){
        !           104:                lft[j] = 1.0+(3*double(rand())-double(rand()))/RAND_MAX;
        !           105:        }
        !           106:
        !           107:        int CellSize = nSpt;
        !           108:        for (j=0; j<nSpt; j++ ){
        !           109:                CellSize += SptType[j];
        !           110:        }
        !           111:        CellStack  MCells( CellSize );
        !           112:        int MVol;
        !           113:
        !           114:        MixedVol(nVar,nSpt,SptType,SptVtx1stIdx,SptVtx,lft,MCells,MVol);
        !           115:
        !           116:        if( debug_print >= 2 ){
        !           117:                cout << "The mixed volume is " << MVol << endl;
        !           118:        }
        !           119:        // Clean the memory
        !           120:
        !           121:        while( ! MCells.IsEmpty() ) MCells.Pop();
        !           122:
        !           123:        delete [] SptType;
        !           124:        delete [] SptVtx1stIdx;
        !           125:        delete [] SptVtx[0];
        !           126:        delete [] SptVtx;
        !           127:        delete [] NuIdx2OldIdx;
        !           128:        delete [] lft;
        !           129:
        !           130:        return MVol;
        !           131: }