00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023 #include <ppl-config.h>
00024
00025 #include "Grid.defs.hh"
00026 #include <cstddef>
00027
00028 namespace Parma_Polyhedra_Library {
00029
00030
00031
00032
00033
00034
00035
00036 bool
00037 Grid::lower_triangular(const Congruence_System& sys,
00038 const Dimension_Kinds& dim_kinds) {
00039 const dimension_type num_columns = sys.num_columns() - 1;
00040
00041
00042 if (sys.num_rows() > num_columns)
00043 return false;
00044
00045
00046
00047 dimension_type row = 0;
00048 for (dimension_type dim = num_columns; dim-- > 0; ) {
00049 if (dim_kinds[dim] == CON_VIRTUAL)
00050 continue;
00051 const Congruence& cg = sys[row];
00052 ++row;
00053
00054 if (cg[dim] <= 0)
00055 return false;
00056
00057 dimension_type col = dim;
00058 while (++col < num_columns)
00059 if (cg[col] != 0)
00060 return false;
00061 }
00062
00063
00064 return row == sys.num_rows();
00065 }
00066
00067
00068
00069
00070
00071
00072
00073 bool
00074 Grid::upper_triangular(const Grid_Generator_System& sys,
00075 const Dimension_Kinds& dim_kinds) {
00076 dimension_type num_columns = sys.space_dimension() + 1;
00077 dimension_type row = sys.num_rows();
00078
00079
00080 if (row > num_columns)
00081 return false;
00082
00083
00084 while (num_columns > 0) {
00085 --num_columns;
00086 if (dim_kinds[num_columns] == GEN_VIRTUAL)
00087 continue;
00088 const Grid_Generator& gen = sys[--row];
00089
00090 if (gen[num_columns] <= 0)
00091 return false;
00092
00093 dimension_type col = num_columns;
00094 while (col-- > 0)
00095 if (gen[col] != 0)
00096 return false;
00097 }
00098
00099
00100 return num_columns == row;
00101 }
00102
00103 void
00104 Grid::multiply_grid(const Coefficient& multiplier, Grid_Generator& gen,
00105 Grid_Generator_System& dest, const dimension_type num_rows,
00106 const dimension_type num_dims) {
00107 if (multiplier == 1)
00108 return;
00109
00110 if (gen.is_line())
00111
00112 for (dimension_type column = num_dims; column-- > 0; )
00113 gen[column] *= multiplier;
00114 else {
00115 assert(gen.is_parameter_or_point());
00116
00117 for (dimension_type index = num_rows; index-- > 0; ) {
00118 Grid_Generator& generator = dest[index];
00119 if (generator.is_parameter_or_point())
00120 for (dimension_type column = num_dims; column-- > 0; )
00121 generator[column] *= multiplier;
00122 }
00123 }
00124 }
00125
00126 void
00127 Grid::multiply_grid(const Coefficient& multiplier, Congruence& cg,
00128 Congruence_System& dest, const dimension_type num_rows,
00129 const dimension_type num_dims) {
00130 if (multiplier == 1)
00131 return;
00132
00133 if (cg.is_proper_congruence())
00134
00135 for (dimension_type index = num_rows; index-- > 0; ) {
00136 Congruence& congruence = dest[index];
00137 if (congruence.is_proper_congruence())
00138 for (dimension_type column = num_dims; column-- > 0; )
00139 congruence[column] *= multiplier;
00140 }
00141 else {
00142 assert(cg.is_equality());
00143
00144 for (dimension_type column = num_dims; column-- > 0; )
00145 cg[column] *= multiplier;
00146 }
00147 }
00148
00149
00150
00151
00152
00153 void
00154 Grid::conversion(Grid_Generator_System& source, Congruence_System& dest,
00155 Dimension_Kinds& dim_kinds) {
00156
00157
00158
00159 assert(upper_triangular(source, dim_kinds));
00160
00161
00162
00163
00164
00165
00166
00167 dimension_type dest_num_rows = 0;
00168 TEMP_INTEGER(diagonal_lcm);
00169 diagonal_lcm = 1;
00170 const dimension_type dims = source.space_dimension() + 1;
00171 dimension_type source_index = source.num_rows();
00172 for (dimension_type dim = dims; dim-- > 0; )
00173 if (dim_kinds[dim] == GEN_VIRTUAL)
00174
00175 ++dest_num_rows;
00176 else {
00177 --source_index;
00178 if (dim_kinds[dim] == PARAMETER) {
00179
00180
00181
00182 lcm_assign(diagonal_lcm, diagonal_lcm, source[source_index][dim]);
00183
00184 ++dest_num_rows;
00185 }
00186
00187 }
00188 assert(source_index == 0);
00189
00190
00191 if (diagonal_lcm == 0)
00192 throw std::runtime_error("PPL internal error: Grid::conversion:"
00193 " source matrix is singular.");
00194
00195 dest.resize_no_copy(dest_num_rows, dims + 1 );
00196
00197
00198
00199
00200 dimension_type dest_index = 0;
00201 source_index = source.num_rows();
00202 for (dimension_type dim = dims; dim-- > 0; ) {
00203 if (dim_kinds[dim] == LINE)
00204 --source_index;
00205 else {
00206 Congruence& cg = dest[dest_index];
00207 for (dimension_type j = dim; j-- > 0; )
00208 cg[j] = 0;
00209 for (dimension_type j = dim + 1; j < dims; ++j)
00210 cg[j] = 0;
00211
00212 if (dim_kinds[dim] == GEN_VIRTUAL) {
00213 cg[dims] = 0;
00214 cg[dim] = 1;
00215 }
00216 else {
00217 assert(dim_kinds[dim] == PARAMETER);
00218 --source_index;
00219 cg[dims] = 1;
00220 exact_div_assign(cg[dim], diagonal_lcm, source[source_index][dim]);
00221 }
00222 ++dest_index;
00223 }
00224 }
00225
00226 assert(source_index == 0);
00227 assert(dest_index == dest_num_rows);
00228 assert(lower_triangular(dest, dim_kinds));
00229
00230
00231
00232
00233
00234
00235
00236
00237 source_index = source.num_rows();
00238 dest_index = 0;
00239 TEMP_INTEGER(multiplier);
00240
00241 for (dimension_type dim = dims; dim-- > 0; ) {
00242 if (dim_kinds[dim] != GEN_VIRTUAL) {
00243 --source_index;
00244 const Coefficient& source_dim = source[source_index][dim];
00245
00246
00247
00248 for (dimension_type row = dest_index; row-- > 0; ) {
00249 Congruence& cg = dest[row];
00250
00251
00252
00253
00254 gcd_assign(multiplier, cg[dim], source_dim);
00255 exact_div_assign(multiplier, source_dim, multiplier);
00256 multiply_grid(multiplier, cg, dest, dest_num_rows, dims);
00257
00258 Coefficient& cg_dim = cg[dim];
00259 exact_div_assign(cg_dim, cg_dim, source_dim);
00260 }
00261 }
00262
00263
00264
00265
00266
00267
00268
00269 dimension_type tmp_source_index = source_index;
00270 if (dim_kinds[dim] != LINE)
00271 ++dest_index;
00272 for (dimension_type dim_prec = dim; dim_prec-- > 0; ) {
00273 if (dim_kinds[dim_prec] != GEN_VIRTUAL) {
00274 --tmp_source_index;
00275 const Coefficient& source_dim = source[tmp_source_index][dim];
00276
00277
00278
00279
00280
00281
00282
00283
00284 for (dimension_type row = dest_index; row-- > 0; ) {
00285 assert(row < dest_num_rows);
00286 Congruence& cg = dest[row];
00287 sub_mul_assign(cg[dim_prec], source_dim, cg[dim]);
00288 }
00289 }
00290 }
00291 }
00292
00293 const Coefficient& modulus = dest[dest_num_rows - 1][0];
00294 for (dimension_type row = dest_num_rows; row-- > 0; ) {
00295 Congruence& cg = dest[row];
00296 if (cg[dims] > 0)
00297
00298 cg[dims] = modulus;
00299 }
00300
00301 assert(lower_triangular(dest, dim_kinds));
00302
00303 #ifdef STRONG_REDUCTION
00304 for (dimension_type dim = dims, i = 0; dim-- > 0; )
00305 if (dim_kinds[dim] != CON_VIRTUAL)
00306
00307 reduce_reduced<Congruence_System, Congruence>
00308 (dest, dim, i++, 0, dim, dim_kinds, false);
00309 #endif
00310 }
00311
00312 void
00313 Grid::conversion(Congruence_System& source, Grid_Generator_System& dest,
00314 Dimension_Kinds& dim_kinds) {
00315
00316
00317
00318 assert(lower_triangular(source, dim_kinds));
00319
00320
00321
00322 dimension_type source_num_rows = 0, dest_num_rows = 0;
00323 TEMP_INTEGER(diagonal_lcm);
00324 diagonal_lcm = 1;
00325 dimension_type dims = source.num_columns() - 1;
00326 for (dimension_type dim = dims; dim-- > 0; )
00327 if (dim_kinds[dim] == CON_VIRTUAL)
00328
00329 ++dest_num_rows;
00330 else {
00331 if (dim_kinds[dim] == PROPER_CONGRUENCE) {
00332
00333
00334
00335 lcm_assign(diagonal_lcm, diagonal_lcm, source[source_num_rows][dim]);
00336
00337 ++dest_num_rows;
00338 }
00339
00340 ++source_num_rows;
00341 }
00342
00343
00344 if (diagonal_lcm == 0)
00345 throw std::runtime_error("PPL internal error: Grid::conversion:"
00346 " source matrix is singular.");
00347
00348 dest.set_index_first_pending_row(dest_num_rows);
00349 dest.resize_no_copy(dest_num_rows, dims + 1 );
00350
00351
00352
00353
00354
00355
00356
00357 dimension_type source_index = 0;
00358
00359 dimension_type dest_index = dest_num_rows - 1;
00360 for (dimension_type dim = dims; dim-- > 0; ) {
00361 if (dim_kinds[dim] == EQUALITY) {
00362 ++source_index;
00363 }
00364 else {
00365 Grid_Generator& g = dest[dest_index];
00366 for (dimension_type j = dim; j-- > 0; )
00367 g[j] = 0;
00368 for (dimension_type j = dim + 1; j < dims; ++j)
00369 g[j] = 0;
00370
00371 if (dim_kinds[dim] == CON_VIRTUAL) {
00372 g.set_is_line();
00373 g[dim] = 1;
00374 }
00375 else {
00376 assert(dim_kinds[dim] == PROPER_CONGRUENCE);
00377 g.set_is_parameter_or_point();
00378 exact_div_assign(g[dim], diagonal_lcm, source[source_index][dim]);
00379 ++source_index;
00380 }
00381 --dest_index;
00382 }
00383 }
00384
00385 assert(upper_triangular(dest, dim_kinds));
00386
00387
00388
00389
00390
00391
00392
00393
00394 source_index = source_num_rows;
00395 dest_index = 0;
00396 TEMP_INTEGER(reduced_source_dim);
00397
00398 for (dimension_type dim = 0; dim < dims; ++dim) {
00399 if (dim_kinds[dim] != CON_VIRTUAL) {
00400 --source_index;
00401 const Coefficient& source_dim = source[source_index][dim];
00402
00403
00404
00405 for (dimension_type row = dest_index; row-- > 0; ) {
00406 Grid_Generator& g = dest[row];
00407
00408
00409
00410
00411 gcd_assign(reduced_source_dim, g[dim], source_dim);
00412 exact_div_assign(reduced_source_dim, source_dim, reduced_source_dim);
00413 multiply_grid(reduced_source_dim, g, dest, dest_num_rows,
00414 dims + 1 );
00415
00416 Coefficient& g_dim = g[dim];
00417 exact_div_assign(g_dim, g_dim, source_dim);
00418 }
00419 }
00420
00421
00422
00423
00424
00425
00426
00427 dimension_type tmp_source_index = source_index;
00428 if (dim_kinds[dim] != EQUALITY)
00429 ++dest_index;
00430 for (dimension_type dim_fol = dim + 1; dim_fol < dims; ++dim_fol) {
00431 if (dim_kinds[dim_fol] != CON_VIRTUAL) {
00432 --tmp_source_index;
00433 const Coefficient& source_dim = source[tmp_source_index][dim];
00434
00435
00436
00437
00438
00439
00440
00441
00442 for (dimension_type row = dest_index; row-- > 0; ) {
00443 assert(row < dest_num_rows);
00444 Grid_Generator& g = dest[row];
00445 sub_mul_assign(g[dim_fol], source_dim, g[dim]);
00446 }
00447 }
00448 }
00449 }
00450
00451 assert(upper_triangular(dest, dim_kinds));
00452
00453 #ifdef STRONG_REDUCTION
00454 for (dimension_type dim = 0, i = 0; dim < dims; ++dim)
00455 if (dim_kinds[dim] != GEN_VIRTUAL)
00456
00457 reduce_reduced<Grid_Generator_System, Grid_Generator>
00458 (dest, dim, i++, dim, dims - 1, dim_kinds);
00459 #endif
00460
00461
00462
00463 const Coefficient& system_divisor = dest[0][0];
00464 for (dimension_type row = dest.num_rows() - 1, dim = dims;
00465 dim-- > 1; )
00466 switch (dim_kinds[dim]) {
00467 case PARAMETER:
00468 dest[row].set_divisor(system_divisor);
00469 case LINE:
00470 --row;
00471 case GEN_VIRTUAL:
00472 break;
00473 }
00474 }
00475
00476 }