[BACK]Return to solve CVS log [TXT][DIR] Up to [local] / OpenXM_contrib2 / asir2000 / lib

Annotation of OpenXM_contrib2/asir2000/lib/solve, Revision 1.2

1.2     ! noro        1: /*
        !             2:  * Copyright (c) 1994-2000 FUJITSU LABORATORIES LIMITED
        !             3:  * All rights reserved.
        !             4:  *
        !             5:  * FUJITSU LABORATORIES LIMITED ("FLL") hereby grants you a limited,
        !             6:  * non-exclusive and royalty-free license to use, copy, modify and
        !             7:  * redistribute, solely for non-commercial and non-profit purposes, the
        !             8:  * computer program, "Risa/Asir" ("SOFTWARE"), subject to the terms and
        !             9:  * conditions of this Agreement. For the avoidance of doubt, you acquire
        !            10:  * only a limited right to use the SOFTWARE hereunder, and FLL or any
        !            11:  * third party developer retains all rights, including but not limited to
        !            12:  * copyrights, in and to the SOFTWARE.
        !            13:  *
        !            14:  * (1) FLL does not grant you a license in any way for commercial
        !            15:  * purposes. You may use the SOFTWARE only for non-commercial and
        !            16:  * non-profit purposes only, such as academic, research and internal
        !            17:  * business use.
        !            18:  * (2) The SOFTWARE is protected by the Copyright Law of Japan and
        !            19:  * international copyright treaties. If you make copies of the SOFTWARE,
        !            20:  * with or without modification, as permitted hereunder, you shall affix
        !            21:  * to all such copies of the SOFTWARE the above copyright notice.
        !            22:  * (3) An explicit reference to this SOFTWARE and its copyright owner
        !            23:  * shall be made on your publication or presentation in any form of the
        !            24:  * results obtained by use of the SOFTWARE.
        !            25:  * (4) In the event that you modify the SOFTWARE, you shall notify FLL by
        !            26:  * e-mail at risa-admin@flab.fujitsu.co.jp of the detailed specification
        !            27:  * for such modification or the source code of the modified part of the
        !            28:  * SOFTWARE.
        !            29:  *
        !            30:  * THE SOFTWARE IS PROVIDED AS IS WITHOUT ANY WARRANTY OF ANY KIND. FLL
        !            31:  * MAKES ABSOLUTELY NO WARRANTIES, EXPRESSED, IMPLIED OR STATUTORY, AND
        !            32:  * EXPRESSLY DISCLAIMS ANY IMPLIED WARRANTY OF MERCHANTABILITY, FITNESS
        !            33:  * FOR A PARTICULAR PURPOSE OR NONINFRINGEMENT OF THIRD PARTIES'
        !            34:  * RIGHTS. NO FLL DEALER, AGENT, EMPLOYEES IS AUTHORIZED TO MAKE ANY
        !            35:  * MODIFICATIONS, EXTENSIONS, OR ADDITIONS TO THIS WARRANTY.
        !            36:  * UNDER NO CIRCUMSTANCES AND UNDER NO LEGAL THEORY, TORT, CONTRACT,
        !            37:  * OR OTHERWISE, SHALL FLL BE LIABLE TO YOU OR ANY OTHER PERSON FOR ANY
        !            38:  * DIRECT, INDIRECT, SPECIAL, INCIDENTAL, PUNITIVE OR CONSEQUENTIAL
        !            39:  * DAMAGES OF ANY CHARACTER, INCLUDING, WITHOUT LIMITATION, DAMAGES
        !            40:  * ARISING OUT OF OR RELATING TO THE SOFTWARE OR THIS AGREEMENT, DAMAGES
        !            41:  * FOR LOSS OF GOODWILL, WORK STOPPAGE, OR LOSS OF DATA, OR FOR ANY
        !            42:  * DAMAGES, EVEN IF FLL SHALL HAVE BEEN INFORMED OF THE POSSIBILITY OF
        !            43:  * SUCH DAMAGES, OR FOR ANY CLAIM BY ANY OTHER PARTY. EVEN IF A PART
        !            44:  * OF THE SOFTWARE HAS BEEN DEVELOPED BY A THIRD PARTY, THE THIRD PARTY
        !            45:  * DEVELOPER SHALL HAVE NO LIABILITY IN CONNECTION WITH THE USE,
        !            46:  * PERFORMANCE OR NON-PERFORMANCE OF THE SOFTWARE.
        !            47:  *
        !            48:  * $OpenXM: OpenXM_contrib2/asir2000/lib/solve,v 1.1.1.1 1999/12/03 07:39:11 noro Exp $
        !            49: */
1.1       noro       50: def kenzan(El,Sl)
                     51: {
                     52:        for ( Tl = El; Tl != []; Tl = cdr(Tl) ) {
                     53:                if ( substv(car(Tl),Sl) ) {
                     54:                        print("kenzan : error");
                     55:                        return 0;
                     56:                }
                     57:        }
                     58:        print("kenzan : ok");
                     59:        return 1;
                     60: }
                     61:
                     62: def substv(P,Sl)
                     63: {
                     64:        for ( A = P; Sl != []; Sl = cdr(Sl) )
                     65:                A = subst(A,car(car(Sl)),car(cdr(car(Sl))));
                     66:        return A;
                     67: }
                     68:
                     69: def co(X,V,D)
                     70: {
                     71:        for ( I = 0; I < D; I++ )
                     72:                X = diff(X,V);
                     73:        return sdiv(subst(X,V,0),fac(D));
                     74: }
                     75:
                     76: def solve(El,Vl)
                     77: /*
                     78:  * El : list of linear forms
                     79:  * Vl : list of variable
                     80:  */
                     81: {
                     82:        N = length(El); M = length(Vl);
                     83:        Mat = newmat(N,M+1);
                     84:        W = newvect(M+1); Index = newvect(N); Vs = newvect(M);
                     85:        for ( I = 0, Tl = Vl; I < M; Tl = cdr(Tl), I++ )
                     86:                Vs[I] = car(Tl);
                     87:        for ( I = 0, Tl = El; I < N; Tl = cdr(Tl), I++ ) {
                     88:                ltov(car(Tl),Vl,W);
                     89:                for ( J = 0; J <= M; J++ )
                     90:                        Mat[I][J] = W[J];
                     91:        }
                     92:        Tl = solvemain(Mat,Index,N,M); L = car(Tl); D = car(cdr(Tl));
                     93:        if ( L < 0 )
                     94:                return [];
                     95:        for ( I = L - 1, S = []; I >= 0; I-- ) {
                     96:                for ( J = Index[I]+1, A = 0; J < M; J++ ) {
                     97:                        A += Mat[I][J]*Vs[J];
                     98:                }
                     99:                S = cons([Vs[Index[I]],-red((A+Mat[I][M])/D)],S);
                    100:        }
                    101:        if ( kenzan(El,S) )
                    102:                return S;
                    103:        else
                    104:                return [];
                    105:        return S;
                    106: }
                    107:
                    108: def solvemain(Mat,Index,N,M)
                    109: /*
                    110:  *     Mat : matrix of size Nx(M+1)
                    111:  *     Index : vector of length N
                    112:  */
                    113: {
                    114:        for ( J = 0, L = 0, D = 1; J < M; J++ ) {
                    115:                for ( I = L; I < N && !Mat[I][J]; I++ );
                    116:                if ( I == N )
                    117:                        continue;
                    118:                Index[L] = J;
                    119:                for ( K = 0; K <= M; K++ ) {
                    120:                        T = Mat[I][K]; Mat[I][K] = Mat[L][K]; Mat[L][K] = T;
                    121:                }
                    122:                for ( I = L + 1, V = Mat[L][J]; I < N; I++ )
                    123:                        for ( K = J, U = Mat[I][J]; K <= M; K++ )
                    124:                                Mat[I][K] = sdiv(Mat[I][K]*V-Mat[L][K]*U,D);
                    125:                D = V; L++;
                    126:        }
                    127:        for ( I = L; I < N; I++ )
                    128:                for ( J = 0; J <= M; J++ )
                    129:                        if ( Mat[I][J] )
                    130:                                return -1;
                    131:        for ( I = L - 2, W = newvect(M+1); I >= 0; I-- ) {
                    132:                for ( J = 0; J <= M; J++ )
                    133:                        W[J] = 0;
                    134:                for ( G = I + 1; G < L; G++ )
                    135:                        for ( H = Index[G], U = Mat[I][H]; H <= M; H++ )
                    136:                                W[H] += Mat[G][H]*U;
                    137:                for ( J = Index[I], U = Mat[I][J]; J <= M; J++ )
                    138:                        Mat[I][J] = sdiv(Mat[I][J]*D-W[J],U);
                    139:        }
                    140:        return [L,D];
                    141: }
                    142:
                    143: def length(L)
                    144: {
                    145:        for ( I = 0; L != []; L = cdr(L), I++ );
                    146:        return I;
                    147: }
                    148:
                    149: def ltov(P,VL,W)
                    150: {
                    151:        for ( I = 0, L = VL; L != []; L = cdr(L), I++ ) {
                    152:                W[I] = co(P,car(L),1); P -= W[I]*car(L);
                    153:        }
                    154:        W[I] = P;
                    155: }
                    156: end$

FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>