The first program calculates all intersections and outputs. The output format is that each row is
A, B, C, D, E, F, G, H [I, j, S, T]
The first eight numbers indicate that the intersection is
A + bz + cz ^ 2 + DZ ^ 3 + EZ ^ 4 + FZ ^ 5 + Gz ^ 6 + Hz ^ 7
The four numbers below indicate that the point is the intersection of the Z ^ I Z ^ J and the Z ^ s Z ^ t of the straight line.
# Include <stdio. h>
# Include <process. h>
_ Int64 factor (_ int64 A ,__ int64 B ){
If (B = 0) return;
Else Return Factor (B, A % B );
}
Class q {
_ Int64 _ up;
_ Int64 _ down;
Void normalize (){
_ Int64 FAC = factor (_ up, _ down );
_ Down/= FAC;
_ Up/= FAC;
If (_ down <0 ){
_ Down =-_ down;
_ Up =-_ up;
}
}
Public:
Q (): _ up (0), _ down (1 ){}
Q (_ int64 up ,__ int64 down): _ up (up), _ down (down ){}
Q (_ int64 X): _ up (x), _ down (1 ){}
Q & operator ++ = (const Q & X );
Q & operator-= (const Q & X );
Q & operator * = (const Q & X );
Q & operator/= (const Q & X );
Q operator-() const {q TMP = * This; TMP. _ up =-_ up; return TMP ;}
Q operator + (const Q & X) const {q TMP = * This; TMP + = x; return TMP ;}
Q operator-(const Q & X) const {q TMP = * This; TMP-= x; return TMP ;}
Q operator * (const Q & X) const {q TMP = * This; TMP * = x; return TMP ;}
Q operator/(const Q & X) const {q TMP = * This; tmp/= x; return TMP ;}
Bool iszero () const {return _ up = 0 ;}
Bool operator = (const Q & X) const {return _ up = x. _ up & _ down = x. _ down ;}
Void output () const {If (_ down = 1) printf ("% 4i64d", _ up); else printf ("% 2i64d/% i64d", _ up, _ down );}
};
Q & Q: Operator + = (const Q & X ){
_ Up = _ up * X. _ down + _ down * X. _ up;
_ Down * = x. _ down;
Normalize ();
Return * this;
}
Q & Q: Operator-= (const Q & X ){
_ Up = _ up * X. _ Down-_ down * X. _ up;
_ Down * = x. _ down;
Normalize ();
Return * this;
}
Q & Q: Operator * = (const Q & X ){
_ Up * = x. _ up;
_ Down * = x. _ down;
Normalize ();
Return * this;
}
Q & Q: Operator/= (const Q & X ){
_ Up * = x. _ down;
_ Down * = x. _ up;
Normalize ();
Return * this;
}
Q mod [9] = {q (1), Q (1), Q (0), Q (-1), Q (-1), Q (-1 ), q (0), Q (1), Q (1 )};
Void Mod (q a [15]) {
Int I, J;
For (I = 14; I> = 8; I --){
Q divider = A [I];
For (j = 0; j <= 8; j ++ ){
A [I-j]-= divider * mod [8-j];
}
}
}
Void dump (const q m [15] [15], const Q B [15]) {
Int I, J;
Printf ("// N // n ");
For (I = 0; I <15; I ++ ){
For (j = 0; j <14; j ++ ){
M [J] [I]. Output ();
Printf (",");
}
M [J] [I]. Output ();
Printf ("= ");
B [I]. Output ();
Printf ("// n ");
}
}
Void solveequation (q m [8] [8], Q B [15]) {
Int I, J, K;
For (I = 0; I <8; I ++ ){
// Dump (M, B );
For (j = I; j <8; j ++)
If (! M [I] [J]. iszero () break;
If (J! = I ){
If (j> = 8 ){
Printf ("There/'s no inverse for the matrix // n ");
Exit (-1 );
}
For (k = 0; k <8; k ++ ){
Q tmp = m [k] [I];
M [k] [I] = m [k] [J];
M [k] [J] = TMP;
}
Q tmp = B [I];
B [I] = B [J];
B [J] = TMP;
}
Q divider = m [I] [I];
For (j = I; j <8; j ++) m [J] [I]/= divider;
B [I]/= Divider;
For (j = 0; j <8; j ++ ){
If (j = I) continue;
Q mult = m [I] [J];
For (k = I; k <8; k ++) m [k] [J]-= mult * m [k] [I];
B [J]-= mult * B [I];
}
}
}
Void inverse (const q a [15], Q B [15]) {
Q m [8] [8];
Q tmp [15];
Int I, J;
For (I = 0; I <8; I ++ ){
For (j = 0; j <15; J ++) TMP [J] = Q (0 );
For (j = 0; j <8; j ++) TMP [I + J] = A [J];
MoD (TMP );
For (j = 0; j <8; j ++ ){
M [I] [J] = TMP [J];
}
}
B [0] = Q (1 );
For (I = 1; I <15; I ++) B [I] = Q (0 );
Solveequation (M, B );
}
Void multiple (const q a [15], const Q B [15], q c [15]) {
Int I, J;
For (I = 0; I <15; I ++) c= Q (0 );
For (I = 0; I <15; I ++) for (j = 0; j <15; J ++ ){
If (I + j <15 ){
C [I + J] + = A [I] * B [J];
} Else {
C [I + j-15]-= A [I] * B [J];
}
}
MoD (C );
}
Void output (const q r [15], int I, Int J, int S, int K ){
Int U;
For (u = 0; U <7; U ++ ){
R [u]. Output ();
Printf (",");
}
R [u]. Output ();
Printf ("[% d, % d] // n", I, j, S, k );
}
# Define Inc (z, x) do {//
Int TMP = (x) % 30; If (TMP <0) TMP + = 30 ;//
If (TMP <15) (z) [TMP] + = Q (1 );//
Else (z) [tmp-15]-= Q (1 );//
} While (0)
# Define Dec (z, x) do {//
Int TMP = (x) % 30; If (TMP <0) TMP + = 30 ;//
If (TMP <15) (z) [TMP]-= Q (1 );//
Else (z) [tmp-15] + = Q (1 );//
} While (0)
Void find (int I, Int J, int S, int t ){
Int K;
// Printf ("% d // n", I, j, S, T );
Q Delta [15], invdelta [15], R [15];
For (k = 0; k <15; k ++ ){
Delta [k] = 0;
}
INC (delta, J-T); Dec (Delta, I-t); Dec (delta, J-S); Inc (delta, I-s );
Dec (delta, T-j); Inc (delta, t-I); Inc (delta, S-j); Dec (delta, s-I );
MoD (DELTA );
Inverse (delta, invdelta );
For (k = 0; k <15; k ++) Delta [k] = 0;
INC (Delta, I-J + S); Dec (delta, I-j + T); Dec (delta,-I + J + S); Inc (delta, -I + J + t );
INC (delta, J + S-t); Dec (Delta, I + S-t); Dec (delta, J-S + T); Inc (delta, i-S + t );
MoD (DELTA );
Multiple (invdelta, Delta, R );
Output (R, I, j, S, T );
}
Int main (){
Int I, j, S, T;
For (I = 0; I <30; I ++) for (S = I + 1; S <30; S ++) for (j = S + 1; j <30; j ++) for (t = J + 1; t <30; t ++ ){
Find (I, j, S, T );
}
}
