comtool.c

Go to the documentation of this file.

 
/* #[ License : */
/*
 *   Copyright (C) 1984-2026 J.A.M. Vermaseren
 *   When using this file you are requested to refer to the publication
 *   J.A.M.Vermaseren "New features of FORM" math-ph/0010025
 *   This is considered a matter of courtesy as the development was paid
 *   for by FOM the Dutch physics granting agency and we would like to
 *   be able to track its scientific use to convince FOM of its value
 *   for the community.
 *
 *   This file is part of FORM.
 *
 *   FORM is free software: you can redistribute it and/or modify it under the
 *   terms of the GNU General Public License as published by the Free Software
 *   Foundation, either version 3 of the License, or (at your option) any later
 *   version.
 *
 *   FORM is distributed in the hope that it will be useful, but WITHOUT ANY
 *   WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
 *   FOR A PARTICULAR PURPOSE.  See the GNU General Public License for more
 *   details.
 *
 *   You should have received a copy of the GNU General Public License along
 *   with FORM.  If not, see <http://www.gnu.org/licenses/>.
 */
/* #] License : */ 
/*
    #[ Includes :
*/
 
#include "form3.h"
 
/*
    #] Includes : 
    #[ inicbufs :
*/
 
int inicbufs(void)
{
    int i, num = AC.cbufList.num;
    CBUF *C = cbuf;
    for ( i = 0; i < num; i++, C++ ) {
        if ( C->Buffer == 0 ) break;
    }
    if ( i >= num ) C = (CBUF *)FromList(&AC.cbufList);
    else num = i;
    C->BufferSize = 2000;
    C->Buffer = (WORD *)Malloc1(C->BufferSize*sizeof(WORD),"compiler buffer-1");
    C->Pointer = C->Buffer;
    C->Top = C->Buffer + C->BufferSize;
    C->maxlhs = 10;
    C->lhs = (WORD **)Malloc1(C->maxlhs*sizeof(WORD *),"compiler buffer-2");
    C->numlhs = 0;
    C->mnumlhs = 0;
    C->maxrhs = 25;
    C->rhs = (WORD **)Malloc1(C->maxrhs*(sizeof(WORD *)+2*sizeof(LONG)+2*sizeof(WORD)),"compiler buffer-3");
    C->CanCommu = (LONG *)(C->rhs+C->maxrhs);
    C->NumTerms = C->CanCommu+C->maxrhs;
    C->numdum = (WORD *)(C->NumTerms+C->maxrhs);
    C->dimension = C->numdum + C->maxrhs;
    C->numrhs = 0;
    C->mnumrhs = 0;
    C->rhs[0] = C->rhs[1] = C->Pointer;
    C->boomlijst = 0;
    RedoTree(C,C->maxrhs);
    ClearTree(num);
    return(num);
}
 
/*
    #] inicbufs : 
    #[ finishcbuf :
*/
 
void finishcbuf(WORD num)
{
    CBUF *C = cbuf+num;
    if ( C->Buffer ) M_free(C->Buffer,"compiler buffer-1");
    if ( C->rhs ) M_free(C->rhs,"compiler buffer-3");
    if ( C->lhs ) M_free(C->lhs,"compiler buffer-2");
    if ( C->boomlijst ) M_free(C->boomlijst,"boomlijst");
    C->Top = C->Pointer = C->Buffer = 0;
    C->rhs = C->lhs = 0;
    C->CanCommu = 0;
    C->NumTerms = 0;
    C->BufferSize = 0;
    C->boomlijst = 0;
    C->numlhs = C->numrhs = C->maxlhs = C->maxrhs = C->mnumlhs =
    C->mnumrhs = C->numtree = C->rootnum = C->MaxTreeSize = 0;
}
 
/*
    #] finishcbuf : 
    #[ clearcbuf :
*/
 
void clearcbuf(WORD num)
{
    CBUF *C = cbuf+num;
    if ( C->boomlijst ) M_free(C->boomlijst,"boomlijst");
    C->Pointer = C->Buffer;
    C->numrhs = C->numlhs = 0;
    C->mnumlhs = 0;
    C->boomlijst = 0;
    C->mnumrhs = 0;
    C->rhs[0] = C->rhs[1] = C->Pointer;
    C->numtree = C->rootnum = C->MaxTreeSize = 0;
    RedoTree(C,C->maxrhs);
    ClearTree(num);
}
 
/*
    #] clearcbuf : 
    #[ DoubleCbuffer :
*/
 
WORD *DoubleCbuffer(int num, WORD *w,int par)
{
    CBUF *C = cbuf + num;
    LONG newsize = C->BufferSize*2;
    WORD *newbuffer = (WORD *)Malloc1(newsize*sizeof(WORD),"compiler buffer-4");
    WORD *w1, *w2;
    LONG offset, j, i;
    DUMMYUSE(par)
 
    w1 = C->Buffer; w2 = newbuffer;
    i = w - w1;
    j = i & 7;
    while ( --j >= 0 ) *w2++ = *w1++;
    i >>= 3;
    while ( --i >= 0 ) {
        *w2++ = *w1++; *w2++ = *w1++; *w2++ = *w1++; *w2++ = *w1++;
        *w2++ = *w1++; *w2++ = *w1++; *w2++ = *w1++; *w2++ = *w1++;
    }
    offset = newbuffer - C->Buffer;
    for ( i = 0; i <= C->numlhs; i++ ) C->lhs[i] += offset;
    for ( i = 1; i <= C->numrhs; i++ ) C->rhs[i] += offset;
    w1 = C->Buffer;
    C->Pointer += offset;
    C->Top = newbuffer + newsize;
    C->BufferSize = newsize;
    C->Buffer = newbuffer;
    M_free(w1,"DoubleCbuffer");
    return(w2);
}
 
/*
    #] DoubleCbuffer : 
    #[ AddLHS :
*/
 
WORD *AddLHS(int num)
{
    CBUF *C = cbuf + num;
    C->numlhs++;
    if ( C->numlhs >= (C->maxlhs-2) ) {
        WORD ***ppp = &(C->lhs);    /* to avoid compiler warning */
        if ( DoubleList((void ***)ppp,&(C->maxlhs),sizeof(WORD *),
        "statement lists") ) Terminate(-1);
    }
    C->lhs[C->numlhs] = C->Pointer;
    C->lhs[C->numlhs+1] = 0;
    return(C->Pointer);
}
 
/*
    #] AddLHS : 
    #[ AddRHS :
*/
 
WORD *AddRHS(int num, int type)
{
    LONG fullsize, *lold, newsize;
    int i;
    WORD **old, *wold;
    CBUF *C;
restart:;
    C = cbuf + num;
    if ( C->numrhs >= (C->maxrhs-2) ) {
        if ( C->maxrhs == 0 ) newsize = 100;
        else                  newsize = C->maxrhs * 2;
        if ( newsize > MAXCOMBUFRHS ) newsize = MAXCOMBUFRHS;
        if ( newsize == C->maxrhs ) {
            if ( AC.tablefilling ) {
                TABLES T = functions[AC.tablefilling].tabl;
/*
                We add a compiler buffer, change a few settings and continue.
*/
                if ( T->buffersfill >= T->bufferssize ) {
                    int new1 = 2*T->bufferssize;
                    WORD *nbufs = (WORD *)Malloc1(new1*sizeof(WORD),"Table compile buffers");
                    for ( i = 0; i < T->buffersfill; i++ )
                            nbufs[i] = T->buffers[i];
                    for ( ; i < new1; i++ ) nbufs[i] = 0;
                    M_free(T->buffers,"Table compile buffers");
                    T->buffers = nbufs;
                    T->bufferssize = new1;
                }
                T->buffers[T->buffersfill++] = T->bufnum = inicbufs();
                AC.cbufnum = num = T->bufnum;
                goto restart;
            }
            else {
                MesPrint("@Compiler buffer overflow. Try to make modules smaller");
                Terminate(-1);
            }
        }
        old = C->rhs;
        fullsize = newsize * (sizeof(WORD *) + 2*sizeof(LONG) + 2*sizeof(WORD));
        C->rhs = (WORD **)Malloc1(fullsize,"subexpression lists");
        for ( i = 0; i < C->maxrhs; i++ ) C->rhs[i] = old[i];
        lold = C->CanCommu; C->CanCommu = (LONG *)(C->rhs+newsize);
        for ( i = 0; i < C->maxrhs; i++ ) C->CanCommu[i] = lold[i];
        lold = C->NumTerms; C->NumTerms = (LONG *)(C->rhs+2*newsize);
        for ( i = 0; i < C->maxrhs; i++ ) C->NumTerms[i] = lold[i];
        wold = C->numdum; C->numdum = (WORD *)(C->NumTerms+newsize);
        for ( i = 0; i < C->maxrhs; i++ ) C->numdum[i] = wold[i];
        wold = C->dimension; C->dimension = (WORD *)(C->numdum+newsize);
        for ( i = 0; i < C->maxrhs; i++ ) C->dimension[i] = wold[i];
        if ( old ) M_free(old,"subexpression lists");
        C->maxrhs = newsize;
        if ( type == 0 ) RedoTree(C,C->maxrhs);
    }
    C->numrhs++;
    C->CanCommu[C->numrhs] = 0;
    C->NumTerms[C->numrhs] = 0;
    C->numdum[C->numrhs] = 0;
    C->dimension[C->numrhs] = 0;
    C->rhs[C->numrhs] = C->Pointer;
    return(C->Pointer);
}
 
/*
    #] AddRHS : 
    #[ AddNtoL :
*/
 
int AddNtoL(int n, WORD *array)
{
    int i;
    CBUF *C = cbuf+AC.cbufnum;
#ifdef COMPBUFDEBUG
    MesPrint("LH: %a",n,array);
#endif
    AddLHS(AC.cbufnum);
    while ( C->Pointer+n >= C->Top ) DoubleCbuffer(AC.cbufnum,C->Pointer,1);
    for ( i = 0; i < n; i++ ) *(C->Pointer)++ = *array++;
    return(0);
}
 
/*
    #] AddNtoL : 
    #[ AddNtoC :
 
    Commentary: added the bufnum on 14-sep-2010 to make the whole a bit
    more flexible (JV). Still to do with AddNtoL.
*/
 
int AddNtoC(int bufnum, int n, WORD *array,int par)
{
    int i;
    WORD *w;
    CBUF *C = cbuf+bufnum;
#ifdef COMPBUFDEBUG
    MesPrint("RH: %a",n,array);
#endif
    while ( C->Pointer+n+1 >= C->Top ) DoubleCbuffer(bufnum,C->Pointer,50+par);
    w = C->Pointer;
    for ( i = 0; i < n; i++ ) *w++ = *array++;
    C->Pointer = w;
    return(0);
}
 
/*
    #] AddNtoC : 
    #[ InsTree :
 
    Routines for balanced tree searching and insertion.
    Compared to Knuth we have a parent link. This minimizes the
    number of compares. That is better for anything that is more
    complicated than just single numbers.
    There are no provisions for removing elements from the tree.
    The routines are:
    void RedoTree(size) Re-allocates the tree space. There will
                        be MaxTreeSize = size elements.
    void ClearTree()    Prunes the tree down to the root element.
    int InsTree(int,int)Searches for the requested element. If not found it
                        will allocate a new element, balance the tree if
                        necessary and return the called number.
                        If it was in the tree, it returns the tree 'value'.
 
    Commentary: added the bufnum on 14-sep-2010 to make the whole a bit
    more flexible (JV).
*/
static COMPTREE comptreezero = {0,0,0,0,0,0};
 
int InsTree(int bufnum, int h)
{
    CBUF *C = cbuf + bufnum;
    COMPTREE *boomlijst = C->boomlijst, *q = boomlijst + C->rootnum, *p, *s;
    WORD *v1, *v2, *v3;
    int ip, iq, is;
 
    if ( C->numtree + 1 >= C->MaxTreeSize ) {
        if ( C->MaxTreeSize == 0 ) {
            COMPTREE *root;
            C->MaxTreeSize = 125;
            C->boomlijst = (COMPTREE *)Malloc1((C->MaxTreeSize+1)*sizeof(COMPTREE),
                "ClearInsTree");
            root = C->boomlijst;
            C->numtree = 0;
            C->rootnum = 0;
            root->left = -1;        
            root->right = -1;
            root->parent = -1;
            root->blnce = 0;
            root->value = -1;
            root->usage = 0;
            for ( ip = 1; ip < C->MaxTreeSize; ip++ ) { C->boomlijst[ip] = comptreezero; }
        }
        else {
            is = C->MaxTreeSize * 2;
            s  = (COMPTREE *)Malloc1((is+1)*sizeof(COMPTREE),"InsTree");
            for ( ip = 0; ip < C->MaxTreeSize; ip++ ) { s[ip] = C->boomlijst[ip]; }
            for ( ip = C->MaxTreeSize; ip <= is; ip++ ) { s[ip] = comptreezero; }
            if ( C->boomlijst ) M_free(C->boomlijst,"InsTree");
            C->boomlijst = s;
            C->MaxTreeSize = is;
        }
        boomlijst = C->boomlijst;
        q = boomlijst + C->rootnum;
    }
 
    if ( q->right == -1 ) { /* First element */
        C->numtree++;
        s = boomlijst+C->numtree;
        q->right = C->numtree;
        s->parent = C->rootnum;
        s->left = s->right = -1;
        s->blnce = 0;
        s->value = h;
        s->usage = 1;
        return(h);
    }
    ip = q->right;
    while ( ip >= 0 ) {
        p = boomlijst + ip;
        v1 = C->rhs[p->value]; v2 = v3 = C->rhs[h];
        while ( *v3 ) v3 += *v3;  /* find the 0 that indicates end-of-expr */
        while ( *v1 == *v2 && v2 < v3 ) { v1++; v2++; }
        if ( *v1 > *v2 ) {
            iq = p->right;
            if ( iq >= 0 ) { ip = iq; }
            else {
                C->numtree++;
                is = C->numtree; 
                p->right = is;
                s = boomlijst + is;
                s->parent = ip; s->left = s->right = -1;
                s->blnce = 0;   s->value = h; s->usage = 1;
                p->blnce++;
                if ( p->blnce == 0 ) return(h);
                goto balance;
            }
        }
        else if ( *v1 < *v2 ) {
            iq = p->left;
            if ( iq >= 0 ) { ip = iq; }
            else {
                C->numtree++;
                is = C->numtree;
                s = boomlijst+is;
                p->left = is;
                s->parent = ip; s->left = s->right = -1;
                s->blnce = 0;   s->value = h; s->usage = 1;
                p->blnce--;
                if ( p->blnce == 0 ) return(h);
                goto balance;
            }
        }
        else {
            p->usage++;
            return(p->value);
        }
    }
/* INTERNAL_ERROR_EXCL_START */
    MesPrint("!>We vallen uit de boom!");
    Terminate(-1);
/* INTERNAL_ERROR_EXCL_STOP */
    return(h);
balance:;
    for (;;) {
        p = boomlijst + ip;
        iq = p->parent;
        if ( iq == C->rootnum ) break;
        q = boomlijst + iq;
        if ( ip == q->left ) q->blnce--;
        else                 q->blnce++;
        if ( q->blnce == 0 ) break;
        if ( q->blnce == -2 ) {
            if ( p->blnce == -1 ) { /* single rotation */
                q->left = p->right;
                p->right = iq;
                p->parent = q->parent;
                q->parent = ip;
                if ( boomlijst[p->parent].left == iq ) boomlijst[p->parent].left = ip;
                else                                   boomlijst[p->parent].right = ip;
                if ( q->left >= 0 ) boomlijst[q->left].parent = iq;
                q->blnce = p->blnce = 0;
            }
            else {  /* double rotation */
                s = boomlijst + is;
                q->left = s->right;
                p->right = s->left;
                s->right = iq;
                s->left = ip;
                if ( p->right >= 0 ) boomlijst[p->right].parent = ip;
                if ( q->left >= 0 ) boomlijst[q->left].parent = iq;
                s->parent = q->parent;
                q->parent = is;
                p->parent = is;
                if ( boomlijst[s->parent].left == iq )
                       boomlijst[s->parent].left = is;
                else   boomlijst[s->parent].right = is;
                if ( s->blnce > 0 ) { q->blnce = s->blnce = 0; p->blnce = -1; }
                else if ( s->blnce < 0 ) { p->blnce = s->blnce = 0; q->blnce = 1; }
                else { p->blnce = s->blnce = q->blnce = 0; }
            }
            break;
        }
        else if ( q->blnce == 2 ) {
            if ( p->blnce == 1 ) {  /* single rotation */
                q->right = p->left;
                p->left = iq;
                p->parent = q->parent;
                q->parent = ip;
                if ( boomlijst[p->parent].left == iq ) boomlijst[p->parent].left = ip;
                else                                   boomlijst[p->parent].right = ip;
                if ( q->right >= 0 ) boomlijst[q->right].parent = iq;
                q->blnce = p->blnce = 0;
            }
            else {  /* double rotation */
                s = boomlijst + is;
                q->right = s->left;
                p->left = s->right;
                s->left = iq;
                s->right = ip;
                if ( p->left >= 0 ) boomlijst[p->left].parent = ip;
                if ( q->right >= 0 ) boomlijst[q->right].parent = iq;
                s->parent = q->parent;
                q->parent = is;
                p->parent = is;
                if ( boomlijst[s->parent].left == iq ) boomlijst[s->parent].left = is;
                else                                   boomlijst[s->parent].right = is;
                if ( s->blnce < 0 ) { q->blnce = s->blnce = 0; p->blnce = 1; }
                else if ( s->blnce > 0 ) { p->blnce = s->blnce = 0; q->blnce = -1; }
                else { p->blnce = s->blnce = q->blnce = 0; }
            }
            break;
        }
        is = ip; ip = iq;
    }
    return(h);
}
 
/*
    #] InsTree : 
    #[ FindTree :
 
    Routines for balanced tree searching.
    Is like InsTree but without the insertions.
    Returns -1 if the element is not in the tree.
    The advantage of this routine over InsTree is that this routine
    can be run in parallel.
*/
 
int FindTree(int bufnum, WORD *subexpr)
{
    CBUF *C = cbuf + bufnum;
    COMPTREE *boomlijst = C->boomlijst, *q = boomlijst + C->rootnum, *p;
    WORD *v1, *v2, *v3;
    int ip, iq;
 
    ip = q->right;
    while ( ip >= 0 ) {
        p = boomlijst + ip;
        v1 = C->rhs[p->value]; v2 = v3 = subexpr;
        while ( *v3 ) v3 += *v3;  /* find the 0 that indicates end-of-expr */
        while ( *v1 == *v2 && v2 < v3 ) { v1++; v2++; }
        if ( *v1 > *v2 ) {
            iq = p->right;
            if ( iq >= 0 ) { ip = iq; }
            else { return(-1); }
        }
        else if ( *v1 < *v2 ) {
            iq = p->left;
            if ( iq >= 0 ) { ip = iq; }
            else { return(-1); }
        }
        else {
            p->usage++;
            return(p->value);
        }
    }
    return(-1);
}
 
/*
    #] FindTree : 
    #[ RedoTree :
*/
 
void RedoTree(CBUF *C, int size)
{
    COMPTREE *newboomlijst;
    int i;
    newboomlijst = (COMPTREE *)Malloc1((size+1)*sizeof(COMPTREE),"newboomlijst");
    if ( C->boomlijst ) {
        if ( C->MaxTreeSize > size ) C->MaxTreeSize = size;
        for ( i = 0; i < C->MaxTreeSize; i++ ) newboomlijst[i] = C->boomlijst[i];
        M_free(C->boomlijst,"boomlijst");
    }
    C->boomlijst = newboomlijst;
    C->MaxTreeSize = size;
}
 
/*
    #] RedoTree : 
    #[ ClearTree :
*/
 
void ClearTree(int i)
{
    CBUF *C = cbuf + i;
    COMPTREE *root = C->boomlijst;
    if ( root ) {
        C->numtree = 0;
        C->rootnum = 0;
        root->left = -1;        
        root->right = -1;
        root->parent = -1;
        root->blnce = 0;
        root->value = -1;
        root->usage = 0;
    }       
}
 
/*
    #] ClearTree : 
    #[ IniFbuffer :
*/
int IniFbuffer(WORD bufnum)
{
    CBUF *C = cbuf + bufnum;
    COMPTREE *root;
    int i;
    LONG fullsize;
    C->maxrhs = AM.fbuffersize;
    C->MaxTreeSize = AM.fbuffersize;
 
    /*
     * Note that bufnum is a return value of inicbufs(). So C has been already
     * initialized. (TU 20 Dec 2011)
     */
    if ( C->boomlijst ) M_free(C->boomlijst, "IniFbuffer-tree");
    if ( C->rhs ) M_free(C->rhs, "IniFbuffer-rhs");
 
    C->boomlijst = (COMPTREE *)Malloc1((C->MaxTreeSize+1)*sizeof(COMPTREE),"IniFbuffer-tree");
    root = C->boomlijst;
    C->numtree = 0;
    C->rootnum = 0;
    root->left = -1;        
    root->right = -1;
    root->parent = -1;
    root->blnce = 0;
    root->value = -1;
    root->usage = 0;
    for ( i = 1; i < C->MaxTreeSize; i++ ) { C->boomlijst[i] = comptreezero; }
 
    fullsize = (C->maxrhs+1) * (sizeof(WORD *) + 2*sizeof(LONG) + 2*sizeof(WORD));
    C->rhs = (WORD **)Malloc1(fullsize,"IniFbuffer-rhs");
    C->CanCommu = (LONG *)(C->rhs+C->maxrhs);
    C->NumTerms = (LONG *)(C->rhs+2*C->maxrhs);
    C->numdum = (WORD *)(C->NumTerms+C->maxrhs);
    C->dimension = (WORD *)(C->numdum+C->maxrhs);
 
    return(0);
}
 
/*
    #] IniFbuffer : 
    #[ numcommute :
 
    Returns the number of non-commuting terms in the expression
*/
 
LONG numcommute(WORD *terms, LONG *numterms)
{
    LONG num = 0;
    WORD *t, *m;
    *numterms = 0;
    while ( *terms ) {
        *numterms += 1;
        t = terms + 1;
        GETSTOP(terms,m);
        while ( t < m ) {
            if ( *t >= FUNCTION ) {
                if ( functions[*t-FUNCTION].commute ) { num++; break; }
            }
            t += t[1];
        }
        terms = terms + *terms;
    }
    return(num);
}
 
/*
    #] numcommute : 
*/