1 // -*- C++ -*-
2 //==============================================================================================
3 //
4 // This file is part of LiDIA --- a library for computational number theory
5 //
6 // Copyright (c) 1994--2001 the LiDIA Group. All rights reserved.
7 //
8 // See http://www.informatik.tu-darmstadt.de/TI/LiDIA/
9 //
10 //----------------------------------------------------------------------------------------------
11 //
12 // $Id$
13 //
14 // Author : Patrick Theobald (PT)
15 // Changes : See CVS log
16 //
17 //==============================================================================================
18
19
20 #ifndef LIDIA_SPARSE_FP_MATRIX_ALGORITHMS_H_GUARD_
21 #define LIDIA_SPARSE_FP_MATRIX_ALGORITHMS_H_GUARD_
22
23
24 #ifndef LIDIA_RANDOM_GENERATOR_H_GUARD_
25 # include "LiDIA/random_generator.h"
26 #endif
27 #ifndef LIDIA_POLYNOMIAL_H_GUARD_
28 # include "LiDIA/polynomial.h"
29 #endif
30 #ifndef LIDIA_ARITH_INL_GUARD_
31 # include "LiDIA/arith.inl"
32 #endif
33 #ifndef LIDIA_MODULAR_OPERATIONS_INL_GUARD_
34 # include "LiDIA/modular_operations.inl"
35 #endif
36 #ifndef LIDIA_MATH_VECTOR_H_GUARD_
37 # include "LiDIA/math_vector.h"
38 #endif
39 #ifndef LIDIA_BASE_POWER_PRODUCT_H_GUARD_
40 # include "LiDIA/base_power_product.h"
41 #endif
42 #ifndef LIDIA_FP_POLYNOMIAL_H_GUARD_
43 # include "LiDIA/Fp_polynomial.h"
44 #endif
45
46
47
48 #ifdef LIDIA_NAMESPACE
49 namespace LiDIA {
50 # define IN_NAMESPACE_LIDIA
51 #endif
52
53
54
55 template< class T, class MATRIX_TYPE >
56 class sparse_fp_matrix_algorithms
57 {
58 protected:
59
60 //
61 // internal functions
62 //
63
64 T * berlekamp_massay(T *s, lidia_size_t n, const T &mod) const;
65 const T multiply_special(const math_vector< T > &w,
66 const math_vector< T > &v,
67 const T &mod) const;
68 void multiply_special(math_vector< T > &w,
69 const base_power_product< ring_matrix < T >, lidia_size_t > &B,
70 const T &mod) const;
71 void multiply_special(math_vector< T > &w,
72 const MATRIX_TYPE &B,
73 const T &mod) const;
74 void multiply_special(math_vector< T > &v,
75 const math_vector< T > &w,
76 const MATRIX_TYPE &B,
77 const T &mod) const;
78
79 public:
80
81 //
82 // constructor
83 //
84
sparse_fp_matrix_algorithms()85 sparse_fp_matrix_algorithms() {}
86
87 //
88 // destructor
89 //
90
~sparse_fp_matrix_algorithms()91 ~sparse_fp_matrix_algorithms() {}
92
93 //
94 // column step form
95 //
96
97 int STF(MATRIX_TYPE &, const T &) const;
98
99 //
100 // extended column step form
101 //
102
STF_extended(MATRIX_TYPE &,const T &)103 const T * STF_extended(MATRIX_TYPE &, const T &) const
104 {
105 return NULL;
106 }
107
108 //
109 // rank
110 //
111
112 lidia_size_t rank(MATRIX_TYPE &, const T &) const;
113
114 //
115 // rank and linearly independent rows or columns
116 //
117
118 lidia_size_t *lininr(MATRIX_TYPE &, const T &) const;
119 lidia_size_t *lininc(MATRIX_TYPE &, const T &) const;
120
121 //
122 // adjoint matrix
123 //
124
125 void adj(MATRIX_TYPE &, const T &) const;
126
127 //
128 // determinant
129 //
130
131 const T det(MATRIX_TYPE &A, const T &mod) const;
132 const T det(const base_power_product< ring_matrix < T >, lidia_size_t > &A, const T &mod) const;
133
134 //
135 // Hessenberg form
136 //
137
138 void HBF(MATRIX_TYPE &, const T &) const;
139
140 //
141 // characteristic polynomial
142 //
143
144 T *charpoly(MATRIX_TYPE &, const T &) const;
145 T *charpoly(const base_power_product< ring_matrix < T >, lidia_size_t > &A, const T &mod) const;
146
147 //
148 // lanczos algorithm
149 //
150
151 bool lanczos(const MATRIX_TYPE &A, math_vector< T > &x, const math_vector< T > &b, T &mod) const;
152 T lanczos_ZmZ(const MATRIX_TYPE &A, math_vector< T > &x, const math_vector< T > &b, T &mod) const;
153
154 //
155 // wiedemann algorithm
156 //
157
158 bool wiedemann(const ring_matrix< T > &A,
159 math_vector< T > &x,
160 const math_vector< T > &b,
161 T &mod) const;
162 bool wiedemann(const base_power_product< ring_matrix < T >, lidia_size_t > &A,
163 math_vector< T > &x,
164 const math_vector< T > &b, T &mod) const;
165
166 //
167 // conjugate_gradient
168 //
169
170 bool conjugate_gradient(const ring_matrix< T > &A,
171 math_vector< T > &x,
172 const math_vector< T > &b,
173 T &mod) const;
174
175 };
176
177
178
179 template< class T >
Fp_polynomial_convert(const Fp_polynomial & res,const T &)180 inline T * Fp_polynomial_convert(const Fp_polynomial &res, const T&)
181 {
182 T *bigres = new T[res.degree() + 1];
183 for (lidia_size_t i = 0; i <= res.degree(); i++)
184 bigres[i] = res[i];
185 return bigres;
186 }
187
188
189
Fp_polynomial_convert(const Fp_polynomial & res,const long &)190 inline long *Fp_polynomial_convert(const Fp_polynomial &res, const long &)
191 {
192 long *bigres = new long[res.degree() + 1];
193 for (lidia_size_t i = 0; i <= res.degree(); i++)
194 res[i].longify(bigres[i]);
195 return bigres;
196 }
197
198
199
200 #ifdef LIDIA_NAMESPACE
201 } // end of namespace LiDIA
202 # undef IN_NAMESPACE_LIDIA
203 #endif
204
205
206
207 #endif // LIDIA_SPARSE_FP_MATRIX_ALGORITHMS_H_GUARD_
208