Annotation of OpenXM_contrib2/asir2000/lib/solve, Revision 1.3
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
1.3 ! noro 26: * e-mail at risa-admin@sec.flab.fujitsu.co.jp of the detailed specification
1.2 noro 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: *
1.3 ! noro 48: * $OpenXM: OpenXM_contrib2/asir2000/lib/solve,v 1.2 2000/08/21 08:31:43 noro Exp $
1.2 noro 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>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to solve CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib2 / asir2000 / lib