ADD file via upload
This commit is contained in:
parent
f15c452a9c
commit
fc0848f3e5
|
|
@ -0,0 +1,632 @@
|
|||
#include "conf.h"
|
||||
#include "monomial.h"
|
||||
#include "polynomial.h"
|
||||
#include "matrix.h"
|
||||
#include "polynomial_store.h"
|
||||
#include "Store.h"
|
||||
int grevlex_interval[D + 1][N] = { 0 };
|
||||
int grevlex_degree_interval[D + 1] = { 0 };
|
||||
int** mult_res;
|
||||
int part_varis[M][1 << D] = { 0 };
|
||||
int d; //当前处理多项式的次数
|
||||
double gauss_time = 0;
|
||||
double append_time = 0;
|
||||
double make_time = 0;
|
||||
double setReducer_time = 0;
|
||||
long long mult_count;// = 0;
|
||||
void set_part_varis(int mono_order)
|
||||
{
|
||||
Monomial m;
|
||||
m.from_order(mono_order);
|
||||
for (int i = 0; i < (1 << m.degree); i++)
|
||||
{
|
||||
Monomial m_part;
|
||||
char part_varis_1[D] = { 0 };
|
||||
int part_index = 0;
|
||||
int temp = i;
|
||||
for (int j = 0; j < m.degree; j++)
|
||||
{
|
||||
if (temp & 1)
|
||||
{
|
||||
part_varis_1[part_index++] = m.factor[j];
|
||||
}
|
||||
temp >>= 1;
|
||||
}
|
||||
m_part.from_factor(part_varis_1, part_index);
|
||||
part_varis[mono_order][i] = m_part.order;
|
||||
}
|
||||
}
|
||||
bool isfactor(int a, int b)
|
||||
{
|
||||
for (int i = 0; i < (1 << D); i++)
|
||||
{
|
||||
if (part_varis[b][i] == a)
|
||||
return true;
|
||||
}
|
||||
return false;
|
||||
}
|
||||
int get_part_varis(char varis1[], char part_varis1[], int get_binary)
|
||||
{
|
||||
int part_index = 0;
|
||||
int i = 0;
|
||||
while (get_binary > 0)
|
||||
{
|
||||
if (get_binary & 1)
|
||||
part_varis1[part_index++] = varis1[i];
|
||||
i++;
|
||||
get_binary >>= 1;
|
||||
}
|
||||
return part_index;
|
||||
}
|
||||
void Initial()
|
||||
{
|
||||
//初始化次数
|
||||
d = 0;
|
||||
//初始化grevlex_interval
|
||||
int now_interval = 0;
|
||||
for (int i = 0; i <= D; i++)
|
||||
{
|
||||
now_interval++;
|
||||
grevlex_degree_interval[i] = now_interval;
|
||||
grevlex_interval[i][0] = now_interval;
|
||||
for (int j = 1; j < N - i; j++)
|
||||
{
|
||||
now_interval += comb(N - j, i);
|
||||
grevlex_interval[i][j] = now_interval;
|
||||
}
|
||||
for (int j = N - i; j < N; j++)
|
||||
grevlex_interval[i][j] = 1 << 30; //置为足够大的数表示取不到
|
||||
}
|
||||
//初始化乘法表
|
||||
mult_res = new int* [grevlex_degree_interval[D - 1]];
|
||||
mult_res[0] = new int[grevlex_degree_interval[D]];
|
||||
for (int i = 0; i < grevlex_degree_interval[D]; i++)
|
||||
mult_res[0][i] = i;
|
||||
for (int i = 1; i < D; i++)
|
||||
{
|
||||
for (int j = grevlex_degree_interval[i - 1]; j < grevlex_degree_interval[i]; j++)
|
||||
{
|
||||
Monomial a;
|
||||
a.from_order(j);
|
||||
mult_res[j] = new int[grevlex_degree_interval[D - i]];
|
||||
for (int k = 0; k < grevlex_degree_interval[D - i]; k++)
|
||||
{
|
||||
Monomial b;
|
||||
b.from_order(k);
|
||||
Monomial ab = a * b;
|
||||
mult_res[j][k] = ab.order;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
for (int i = 0; i < M; i++)
|
||||
set_part_varis(i);
|
||||
}
|
||||
|
||||
Polynomial_Store ps(MAX_POLYS);
|
||||
Temp_Polynomial_Store ps_temp;
|
||||
Matrix matrix;
|
||||
Polynomial_Input ps_input;
|
||||
//Polynomial_Store ps_temp2(30000);
|
||||
MPTS induce_polys;
|
||||
MPTS cross_polys;
|
||||
Store system_store;
|
||||
bool cross_mono_disable[__PARALLEL][monomials(N, D - 1)];
|
||||
|
||||
//设置消元子
|
||||
void set_reducer(int degree, int mode = 0)
|
||||
{
|
||||
auto setR_start = high_resolution_clock::now();
|
||||
//遍历系统存储各索引，查看其中非稠密多项式会不会有因子
|
||||
Monomial temp_leader, temp_part;
|
||||
char part_varis[D] = { 0 };
|
||||
int part_degree;
|
||||
//模式0：先导入交叉多项式产生的新首项，搜索系统中各项看哪些项可以设置乘积多项式
|
||||
if (mode == 0)
|
||||
{
|
||||
//清除当前系统中乘积多项式
|
||||
//for (int i = 0; i < grevlex_degree_interval[degree]; i++)
|
||||
parallel_for(0, grevlex_degree_interval[degree], [&](int i) {
|
||||
if (system_store[i].type == 0) //统计系统中存储的低次多项式1的数量
|
||||
{
|
||||
return;
|
||||
}
|
||||
system_store[i].mul_mono = -1;
|
||||
system_store[i].mul_poly = -1;
|
||||
system_store[i].reductor = -1;
|
||||
system_store[i].type = -1;
|
||||
system_store[i].poly = NULL;
|
||||
});
|
||||
|
||||
//将产生新首项的交叉多项式记入系统
|
||||
for (int i = 0; i < cross_polys.polys; i++)
|
||||
{
|
||||
Monomial m1, m2, mult;
|
||||
m1.from_order(cross_polys[i].first);
|
||||
m2.from_order(cross_polys[i].second);
|
||||
mult = m1 * m2;
|
||||
if (system_store[mult.order].type == -1)
|
||||
{
|
||||
system_store[mult.order].type = 1;
|
||||
system_store[mult.order].reductor = cross_polys[i].first;
|
||||
system_store[mult.order].mul_poly = cross_polys[i].first;
|
||||
system_store[mult.order].mul_mono = cross_polys[i].second;
|
||||
cross_polys[i].first = -1;
|
||||
cross_polys[i].second = -1;
|
||||
}
|
||||
}
|
||||
|
||||
for (int i = 0; i < grevlex_degree_interval[degree]; i++)
|
||||
{
|
||||
if (system_store[i].type == 0 || system_store[i].type == 1)
|
||||
continue;
|
||||
temp_leader.from_order(i);
|
||||
for (int j = 0; j < (1 << temp_leader.degree) - 1; j++)
|
||||
{
|
||||
part_degree = get_part_varis(temp_leader.factor, part_varis, j);
|
||||
temp_part.from_factor(part_varis, part_degree);
|
||||
if (system_store[temp_part.order].type == 0)
|
||||
{
|
||||
if (system_store[i].mul_poly < 0)
|
||||
//if (temp_part.order < system_store[i].mul_poly || system_store[i].mul_poly < 0)
|
||||
{
|
||||
system_store[i].type = 1;
|
||||
system_store[i].reductor = temp_part.order;
|
||||
system_store[i].mul_poly = temp_part.order;
|
||||
system_store[i].mul_mono = (temp_leader / temp_part).order;
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
if (mode == 1)
|
||||
{
|
||||
for (int i = 0; i < grevlex_degree_interval[degree]; i++)
|
||||
{
|
||||
system_store[i].mul_poly = -1;
|
||||
system_store[i].mul_mono = -1;
|
||||
temp_leader.from_order(i);
|
||||
for (int j = 0; j < (1 << temp_leader.degree) - 1; j++)
|
||||
{
|
||||
part_degree = get_part_varis(temp_leader.factor, part_varis, j);
|
||||
temp_part.from_factor(part_varis, part_degree);
|
||||
if (system_store[temp_part.order].type == 0)
|
||||
{
|
||||
if (temp_part.order < system_store[i].mul_poly || system_store[i].mul_poly < 0)
|
||||
{
|
||||
system_store[i].type = 1;
|
||||
system_store[i].reductor = temp_part.order;
|
||||
system_store[i].mul_poly = temp_part.order;
|
||||
system_store[i].mul_mono = (temp_leader / temp_part).order;
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
auto setR_end = high_resolution_clock::now();
|
||||
setReducer_time += duration_cast<milliseconds>(setR_end - setR_start).count() / 1000.0;
|
||||
}
|
||||
|
||||
//高斯消元全过程
|
||||
void Gauss_elimination(int degree)
|
||||
{
|
||||
//强清空ps_temp中的稀疏多项式
|
||||
for (int i = 0; i < reductors; i++)
|
||||
{
|
||||
ps_temp.sparse_data[i].clear(1);
|
||||
}
|
||||
//清空矩阵，为高斯消元准备好空间
|
||||
matrix.clear();
|
||||
//导入被消元多项式
|
||||
auto append_start = high_resolution_clock::now();
|
||||
int ps_rec = ps.count_polys; //当前存储的多项式数量
|
||||
int count_reducible = 0;
|
||||
|
||||
//导入诱导多项式
|
||||
int begin_index = ps.count_polys;
|
||||
parallel_for(0, induce_polys.polys, [&](int i) {
|
||||
ps.insert(induce_polys[i].first, induce_polys[i].second, begin_index + i);
|
||||
});
|
||||
ps.count_polys += induce_polys.polys;
|
||||
begin_index += induce_polys.polys;
|
||||
|
||||
//导入交叉多项式
|
||||
cross_polys.remove_null(cross_polys.polys); //首先清除已被加入消元子的交叉多项式
|
||||
parallel_for(0, cross_polys.polys, [&](int i) {
|
||||
ps.insert(cross_polys[i].first, cross_polys[i].second, begin_index + i);
|
||||
});
|
||||
ps.count_polys += cross_polys.polys;
|
||||
|
||||
//计算被消元多项式总数
|
||||
count_reducible = induce_polys.polys + cross_polys.polys;
|
||||
auto append_end = high_resolution_clock::now();
|
||||
append_time += duration_cast<milliseconds>(append_end - append_start).count() / 1000.0;
|
||||
cout << "当前系统中存储多项式数量：" << ps.count_polys << endl;
|
||||
cout << "诱导多项式数量：" << induce_polys.polys << endl;
|
||||
cout << "交叉多项式数量：" << cross_polys.polys << endl;
|
||||
matrix.add_elim_poly(ps, ps_rec, ps_rec + count_reducible);
|
||||
cout << "实际导入多项式数量：" << count_reducible << endl;
|
||||
induce_polys.clear();
|
||||
cross_polys.clear();
|
||||
|
||||
//循环导入消元子进行消元
|
||||
auto gauss_start = high_resolution_clock::now();
|
||||
int max_leader = matrix.get_max_leader();
|
||||
for (int i = max_leader + 1; i > 0;)
|
||||
{
|
||||
//导入消元子
|
||||
int reductor_polys = 0;
|
||||
if (i % reductors)
|
||||
{
|
||||
reductor_polys = matrix.add_reductor(system_store, ps_temp, i - i % reductors, i);
|
||||
i -= i % reductors;
|
||||
}
|
||||
else
|
||||
{
|
||||
reductor_polys = matrix.add_reductor(system_store, ps_temp, i - reductors, i);
|
||||
i -= reductors;
|
||||
}
|
||||
matrix.Gauss_elimination();
|
||||
}
|
||||
auto gauss_end = high_resolution_clock::now();
|
||||
gauss_time += duration_cast<milliseconds>(gauss_end - gauss_start).count() / 1000.0;
|
||||
cout << "本段高斯消元时间：" << duration_cast<milliseconds>(gauss_end - gauss_start).count() / 1000.0 << "秒" << endl << endl;
|
||||
ps.remove_zero(); //去掉消元过程产生的零多项式
|
||||
//将高斯消元结果导入存储
|
||||
for (int i = ps_rec; i < ps.count_polys; i++)
|
||||
system_store.Store_Polynomial(&ps[i], i);
|
||||
}
|
||||
|
||||
void Reduce_Groebner()
|
||||
{
|
||||
cout << "下面得到既约Groebner基：" << endl;
|
||||
for (int d0 = 0; d0 <= D; d0++)
|
||||
{
|
||||
//强清空ps_temp中的稀疏多项式
|
||||
for (int i = 0; i < reductors; i++)
|
||||
{
|
||||
ps_temp.sparse_data[i].clear(1);
|
||||
}
|
||||
cout << "正在化简" << d0 << "次的多项式" << endl;
|
||||
//使用d0次消元子进行消元
|
||||
set_reducer(d0, 1);
|
||||
//清空矩阵，为高斯消元准备好空间
|
||||
matrix.clear();
|
||||
//导入系统中存储的多项式
|
||||
int matrix_rows = matrix.add_elim_poly(ps, 0, ps.count_polys, d0);
|
||||
//循环导入消元子进行消元
|
||||
auto gauss_start = high_resolution_clock::now();
|
||||
int max_leader = matrix.get_max_leader();
|
||||
for (int i = max_leader + 1; i > 0;)
|
||||
{
|
||||
//导入消元子
|
||||
int reductor_polys = 0;
|
||||
if (i % reductors)
|
||||
{
|
||||
reductor_polys = matrix.add_reductor(system_store, ps_temp, i - i % reductors, i, 2);
|
||||
i -= i % reductors;
|
||||
}
|
||||
else
|
||||
{
|
||||
reductor_polys = matrix.add_reductor(system_store, ps_temp, i - reductors, i, 2);
|
||||
i -= reductors;
|
||||
}
|
||||
//int reductor_polys = matrix.add_reductor(system_store, ps_temp, max(i - reductors, 0), i, 2);
|
||||
matrix.Gauss_elimination();
|
||||
}
|
||||
auto gauss_end = high_resolution_clock::now();
|
||||
gauss_time += duration_cast<milliseconds>(gauss_end - gauss_start).count() / 1000.0;
|
||||
cout << "用时：" << duration_cast<milliseconds>(gauss_end - gauss_start).count() / 1000.0 << "秒" << endl;
|
||||
ps.remove_zero(); //去掉消元过程产生的零多项式
|
||||
//将高斯消元结果导入存储
|
||||
system_store.clear();
|
||||
for (int i = 0; i < ps.count_polys; i++)
|
||||
system_store.Store_Polynomial(&ps[i], i);
|
||||
//继续简化本轮高斯消元结果
|
||||
set_reducer(d0, 1);
|
||||
int count_reduced = 0;
|
||||
//记录矩阵中实际存在的最大非零多项式
|
||||
max_leader = -1;
|
||||
for (int i = 0; i < matrix_rows; i++)
|
||||
{
|
||||
if (matrix.E[i] != NULL)
|
||||
{
|
||||
if (matrix.E[i]->leader_order > max_leader)
|
||||
{
|
||||
max_leader = matrix.E[i]->leader_order;
|
||||
}
|
||||
}
|
||||
}
|
||||
//重新记录系统中存储的多项式数量，为了内存复用，后续计算中间结果均存在ps对象中
|
||||
int polys_rec = ps.count_polys;
|
||||
for (int i = 1; i <= max_leader; i++)
|
||||
{
|
||||
if (system_store[i].type == 1)
|
||||
{
|
||||
ps.append(system_store[i].mul_poly, system_store[i].mul_mono); //存入ps对象中
|
||||
system_store[i].poly = &ps[ps.count_polys - 1];
|
||||
}
|
||||
if (system_store[i].type == 0)
|
||||
{
|
||||
for (int j = i - 1; j > 0; j--)
|
||||
{
|
||||
if (system_store[j].type >= 0 && system_store[i].poly->get(j))
|
||||
{
|
||||
*system_store[i].poly ^= (*system_store[j].poly);
|
||||
}
|
||||
}
|
||||
//count_reduced++;
|
||||
//if (count_reduced >= polys_rec)
|
||||
// break;
|
||||
}
|
||||
}
|
||||
ps.count_polys = polys_rec; //将后续计算中间结果强制清空
|
||||
}
|
||||
|
||||
}
|
||||
//计算Groebner基
|
||||
void Compute_Groebner()
|
||||
{
|
||||
int floor_4 = 0;
|
||||
for (d = 0; d <= D; d++)
|
||||
{
|
||||
cout << "正在处理次数为" << d << "的多项式\n";
|
||||
int polys_rec = 0;
|
||||
int polys_input = 0;
|
||||
//输入多项式
|
||||
for (int i = 0; i < ps_input.count_polys; i++)
|
||||
{
|
||||
if (ps_input.data[i].degree == d)
|
||||
{
|
||||
ps.append(ps_input.data[i]); //考虑到多项式编号问题，本处需要进行复制
|
||||
polys_input++;
|
||||
}
|
||||
}
|
||||
|
||||
//将输入多项式写入矩阵并进行高斯消元
|
||||
matrix.clear();
|
||||
int input_rows = matrix.add_elim_poly(ps, ps.count_polys - polys_input, ps.count_polys);
|
||||
if (input_rows > 0)
|
||||
{
|
||||
auto gauss_start = high_resolution_clock::now();
|
||||
//循环导入消元子进行消元
|
||||
for (int i = grevlex_degree_interval[d]; i > 0;)
|
||||
{
|
||||
//导入消元子
|
||||
int reductor_polys = 0;
|
||||
if (i % reductors)
|
||||
{
|
||||
reductor_polys = matrix.add_reductor(system_store, ps_temp, i - i % reductors, i);
|
||||
i -= i % reductors;
|
||||
}
|
||||
else
|
||||
{
|
||||
reductor_polys = matrix.add_reductor(system_store, ps_temp, i - reductors, i);
|
||||
i -= reductors;
|
||||
}
|
||||
matrix.Gauss_elimination();
|
||||
}
|
||||
auto gauss_end = high_resolution_clock::now();
|
||||
gauss_time += duration_cast<milliseconds>(gauss_end - gauss_start).count() / 1000.0;
|
||||
ps.remove_zero();
|
||||
}
|
||||
//将高斯消元结果导入存储
|
||||
for (int i = 0; i < input_rows; i++)
|
||||
system_store.Store_Polynomial(matrix.get_elim_poly(i), i);
|
||||
|
||||
while (polys_rec < ps.count_polys)
|
||||
{
|
||||
//生成诱导多项式
|
||||
auto make_start = high_resolution_clock::now();
|
||||
parallel_for(0, 240, [&](int index) {
|
||||
int insert_index = index * PART_ROWS;
|
||||
int diff = ps.count_polys - polys_rec;
|
||||
for (int i = polys_rec + index * diff / 240; i < polys_rec + (index + 1) * diff / 240; i++)
|
||||
{
|
||||
if (ps[i].degree + 1 <= d) //需增加筛选条件
|
||||
{
|
||||
Monomial m = ps[i].get_leader();
|
||||
for (int j = 0; j < m.degree; j++)
|
||||
{
|
||||
induce_polys.insert(ps[i].leader_order, N - m.factor[j], insert_index); //一次单项式项序顺序为 N-变量下标
|
||||
insert_index++;
|
||||
}
|
||||
}
|
||||
}
|
||||
});
|
||||
induce_polys.remove_null();
|
||||
//生成交叉多项式
|
||||
parallel_for(0, __PARALLEL, [&](int index) {
|
||||
int insert_index = index;
|
||||
for (int i = index; i < ps.count_polys; i += __PARALLEL)
|
||||
{
|
||||
memset(cross_mono_disable[index], 0, grevlex_degree_interval[d - ps[i].degree] * sizeof(bool));
|
||||
int& leader = ps[i].leader_order;
|
||||
if (ps[i].degree == d)
|
||||
continue;
|
||||
Monomial a;
|
||||
a.from_order(leader);
|
||||
//生成(f,m)形式的交叉多项式，其中m次数为1
|
||||
for (int j = 1; j < grevlex_degree_interval[d - ps[i].degree]; j++)
|
||||
{
|
||||
if (cross_mono_disable[index][j])
|
||||
continue;
|
||||
Monomial m_j;
|
||||
m_j.from_order(j);
|
||||
if (!coprime(a, m_j)) //a与m_j不互质，存在类诱导多项式的情形
|
||||
continue;
|
||||
bool ischoose = false;
|
||||
for (int k = 0; k < (1 << ps[i].degree); k++)
|
||||
{
|
||||
int& part = part_varis[leader][k];
|
||||
int& product = mult_res[part][j];
|
||||
if (system_store[product].type == 0 && max(i, system_store[product].id) >= polys_rec)
|
||||
{
|
||||
ischoose = true;
|
||||
Monomial b;
|
||||
b.from_order(product);
|
||||
Monomial ab = a * b;
|
||||
Monomial c;
|
||||
for (int r = 1; r < (1 << ab.degree) - 1; r++)
|
||||
{
|
||||
c.from_order(part_varis[ab.order][r]);
|
||||
if (c.order == a.order || c.order == b.order)
|
||||
continue;
|
||||
if (system_store[c.order].type == 0)
|
||||
{
|
||||
Monomial ac = a * c;
|
||||
Monomial bc = b * c;
|
||||
bool ac1 = !coprime(a, c);
|
||||
bool bc1 = !coprime(b, c);
|
||||
bool ac2 = !isfactor(c, b);
|
||||
bool bc2 = !isfactor(c, a);
|
||||
int id_a = i;
|
||||
int id_b = system_store[b.order].id;
|
||||
int id_c = system_store[c.order].id;
|
||||
auto tab = make_tuple(true, ab.degree, true, max(id_a, id_b), min(id_a, id_b));
|
||||
auto tac = make_tuple(ac1, ac.degree, ac2, max(id_a, id_c), min(id_a, id_c));
|
||||
auto tbc = make_tuple(bc1, bc.degree, bc2, max(id_b, id_c), min(id_b, id_c));
|
||||
if (tab > tac && tab > tbc)
|
||||
{
|
||||
ischoose = false;
|
||||
break;
|
||||
}
|
||||
}
|
||||
}
|
||||
if (ischoose)
|
||||
break;
|
||||
}
|
||||
}
|
||||
if (ischoose)
|
||||
{
|
||||
cross_polys.insert(leader, j, insert_index);
|
||||
for (int k = 1; k < grevlex_degree_interval[d - ps[i].degree - m_j.degree]; k++)
|
||||
{
|
||||
cross_mono_disable[index][mult_res[j][k]] = true;
|
||||
}
|
||||
insert_index += __PARALLEL;
|
||||
}
|
||||
}
|
||||
}
|
||||
});
|
||||
|
||||
/*
|
||||
parallel_for(0, 240, [&](int index) {
|
||||
int insert_index = index * PART_ROWS;
|
||||
//for (int i = index * ps.count_polys / 240; i < (index + 1) * ps.count_polys / 240; i++)
|
||||
for (int i = index; i < ps.count_polys; i+=240)
|
||||
{
|
||||
Monomial a = ps[i].get_leader();
|
||||
for (int j = max(polys_rec, i + 1); j < ps.count_polys; j++)
|
||||
{
|
||||
Monomial b = ps[j].get_leader();
|
||||
|
||||
//if (coprime(a, b))
|
||||
// continue;
|
||||
|
||||
Monomial ab = a * b;
|
||||
if (ab.degree > d)
|
||||
continue;
|
||||
|
||||
Monomial f_temp = ab / a; //本质为lcm(f0,g0)/f0
|
||||
Monomial g_temp = ab / b;
|
||||
|
||||
//筛选该交叉多项式是不是更简单的交叉多项式的倍数多项式
|
||||
bool ischoose = true;
|
||||
char part_varis[D + 1] = { 0 };
|
||||
Monomial c;
|
||||
for (int k = 0; k < (1 << ab.degree) - 1; k++)
|
||||
{
|
||||
int part_degree = get_part_varis(ab.factor, part_varis, k);
|
||||
c.from_factor(part_varis, part_degree);
|
||||
if (c.order == a.order || c.order == b.order)
|
||||
continue;
|
||||
if (system_store[c.order].type == 0)
|
||||
{
|
||||
Monomial ac = a * c;
|
||||
Monomial bc = b * c;
|
||||
bool ac1 = !coprime(a, c);
|
||||
bool bc1 = !coprime(b, c);
|
||||
bool ac2 = !isfactor(c, b);
|
||||
bool bc2 = !isfactor(c, a);
|
||||
int id_a = i;
|
||||
int id_b = j;
|
||||
int id_c = system_store[c.order].id;
|
||||
auto tab = make_tuple(true, ab.degree, true, max(id_a, id_b), min(id_a, id_b));
|
||||
auto tac = make_tuple(ac1, ac.degree, ac2, max(id_a, id_c), min(id_a, id_c));
|
||||
auto tbc = make_tuple(bc1, bc.degree, bc2, max(id_b, id_c), min(id_b, id_c));
|
||||
if (tab > tac && tab > tbc)
|
||||
{
|
||||
ischoose = false;
|
||||
break;
|
||||
}
|
||||
}
|
||||
}
|
||||
if (ischoose)
|
||||
{
|
||||
cross_polys.insert(ps[i].leader_order, f_temp.order, insert_index);
|
||||
cross_polys.insert(ps[j].leader_order, g_temp.order, insert_index + 1);
|
||||
insert_index += 2;
|
||||
}
|
||||
}
|
||||
}
|
||||
});
|
||||
*/
|
||||
//为cross_polys去重并记录剩余的多项式数量
|
||||
cross_polys.remove_null();
|
||||
sort(cross_polys.data, cross_polys.data + cross_polys.polys);
|
||||
cross_polys.polys = unique(cross_polys.data, cross_polys.data + cross_polys.polys) - cross_polys.data;
|
||||
|
||||
auto make_end = high_resolution_clock::now();
|
||||
make_time += duration_cast<milliseconds>(make_end - make_start).count() / 1000.0;
|
||||
//设置消元子
|
||||
set_reducer(d);
|
||||
|
||||
//记录当前系统中多项式数量
|
||||
polys_rec = ps.count_polys;
|
||||
//高斯消元
|
||||
Gauss_elimination(d);
|
||||
}
|
||||
}
|
||||
//得到既约Groebner基
|
||||
Reduce_Groebner();
|
||||
}
|
||||
|
||||
|
||||
int main()
|
||||
{
|
||||
cout<<sizeof(Polynomial)<<endl;
|
||||
cout<<ps.data[0].mono.data<<endl;
|
||||
cout<<ps.data[1].mono.data<<endl;
|
||||
ofstream logFile("log.txt");
|
||||
tbb::task_scheduler_init mt(__PARALLEL);
|
||||
//cout.rdbuf(logFile.rdbuf());
|
||||
//初始化
|
||||
mult_count = 0;
|
||||
Initial();
|
||||
|
||||
//打开文件读取多项式
|
||||
ifstream input_file(FILENAME);
|
||||
ofstream fout("Groebner.txt");
|
||||
ps_input.input(input_file);
|
||||
//ps_input.print();
|
||||
cout << ps_input.count_polys << endl;
|
||||
ps_input.set_leader_degree(M - 1);
|
||||
ps_input.sort_data();
|
||||
|
||||
auto start = high_resolution_clock::now();
|
||||
Compute_Groebner();
|
||||
auto end = high_resolution_clock::now();
|
||||
|
||||
double used_time = duration_cast<milliseconds>(end - start).count() / 1000.0;
|
||||
double other_time = used_time - gauss_time - append_time - make_time - setReducer_time;
|
||||
ps.print(fout);
|
||||
cout << setw(24) << left << "程序运行用时：" << used_time << "秒" << endl;
|
||||
cout << setw(24) << left << "高斯消元用时：" << gauss_time << "秒" << "\t占比：" << gauss_time * 100 / used_time << '%' << endl;
|
||||
cout << setw(24) << left << "导入多项式用时：" << append_time << "秒" << "\t占比：" << append_time * 100 / used_time << '%' << endl;
|
||||
cout << setw(24) << left << "生成乘积多项式用时：" << make_time << "秒" << "\t占比：" << make_time * 100 / used_time << '%' << endl;
|
||||
cout << setw(24) << left << "设置消元子用时：" << setReducer_time << "秒" << "\t占比：" << setReducer_time * 100 / used_time << '%' << endl;
|
||||
cout << setw(24) << left << "其他工作用时：" << other_time << "秒" << "\t占比：" << other_time * 100 / used_time << '%' << endl;
|
||||
logFile.close();
|
||||
}
|
||||
|
||||
Loading…
Reference in New Issue