[mathicgb] 115/393: Added an option to output all F4 matrices with more than the specified number of non-zero entries. The matrix format is as for Lachartre and Fayssal's code, except each file contains the 4 sub-matrices instead of one big matrix. There is now also a matrix action in addition to gb and siggb. It currently just reads a matrix and writes it out - I'm stopping today at the point where I got the written matrix to actually match the read matrix byte-for-byte. :) This also involved a reorganization of someo of the command line interface code which is a bit nicer now.
Doug Torrance
dtorrance-guest at moszumanska.debian.org
Fri Apr 3 15:58:44 UTC 2015
This is an automated email from the git hooks/post-receive script.
dtorrance-guest pushed a commit to branch upstream
in repository mathicgb.
commit 4753243529af1e995d97710afeaf9bf10c694108
Author: Bjarke Hammersholt Roune <bjarkehr.code at gmail.com>
Date: Fri Nov 16 20:15:34 2012 +0100
Added an option to output all F4 matrices with more than the specified number of non-zero entries. The matrix format is as for Lachartre and Fayssal's code, except each file contains the 4 sub-matrices instead of one big matrix. There is now also a matrix action in addition to gb and siggb. It currently just reads a matrix and writes it out - I'm stopping today at the point where I got the written matrix to actually match the read matrix byte-for-byte. :) This also involved a reorgani [...]
---
.gitignore | 2 +
Makefile.am | 6 +-
build/vs12/mathicgb-exe/mathicgb-exe.vcxproj | 4 +
.../vs12/mathicgb-exe/mathicgb-exe.vcxproj.filters | 16 ++-
src/cli/CommonParams.cpp | 77 ++++++++++++
src/cli/CommonParams.hpp | 41 ++++++
src/cli/GBAction.cpp | 137 +++++++++++++++++++++
src/cli/GBAction.hpp | 38 ++++++
src/cli/GBCommonParams.cpp | 83 +++++++++++++
src/cli/GBCommonParams.hpp | 23 ++++
src/cli/GBMain.cpp | 2 +
src/cli/MatrixAction.cpp | 64 ++++++++++
src/cli/MatrixAction.hpp | 31 +++++
src/cli/SigGBAction.cpp | 125 +++++++++++++++++++
src/cli/SigGBAction.hpp | 36 ++++++
src/mathicgb/BuchbergerAlg.cpp | 28 ++---
src/mathicgb/BuchbergerAlg.hpp | 6 +-
src/mathicgb/F4MatrixReducer.cpp | 4 +-
src/mathicgb/F4Reducer.cpp | 126 ++++++++++---------
src/mathicgb/F4Reducer.hpp | 44 ++++---
src/mathicgb/QuadMatrix.cpp | 60 ++++++++-
src/mathicgb/QuadMatrix.hpp | 18 +++
src/mathicgb/Reducer.cpp | 11 +-
src/mathicgb/Reducer.hpp | 6 +-
src/mathicgb/SparseMatrix.cpp | 120 +++++++++++++++++-
src/mathicgb/SparseMatrix.hpp | 11 +-
src/test/F4MatrixReducer.cpp | 8 +-
src/test/SparseMatrix.cpp | 2 +-
src/test/gb-test.cpp | 7 +-
29 files changed, 1007 insertions(+), 129 deletions(-)
diff --git a/.gitignore b/.gitignore
index 424e2c1..3e9d79e 100755
--- a/.gitignore
+++ b/.gitignore
@@ -1,6 +1,8 @@
# MathicGB generated files
*.gb
*.stats
+*.mat
+*.qmat
# Compiled Object files
*.slo
diff --git a/Makefile.am b/Makefile.am
index c3ac285..86cd0f9 100755
--- a/Makefile.am
+++ b/Makefile.am
@@ -61,7 +61,11 @@ libmathicgb_la_SOURCES = src/mathicgb/BjarkeGeobucket2.cpp \
src/mathicgb/QuadMatrix.cpp src/mathicgb/F4MatrixReducer.cpp \
src/mathicgb/F4MatrixReducer.hpp src/mathicgb/MonomialMap.hpp \
src/mathicgb/RawVector.hpp src/mathicgb/Atomic.hpp \
- src/mathicgb/FixedSizeMonomialMap.hpp
+ src/mathicgb/FixedSizeMonomialMap.hpp src/cli/CommonParams.hpp \
+ src/cli/CommonParams.cpp src/cli/GBAction.hpp src/cli/GBAction.cpp \
+ src/cli/GBCommonParams.hpp src/cli/GBCommonParams.cpp \
+ src/cli/MatrixAction.cpp src/cli/MatrixAction.hpp \
+ src/cli/SigGBAction.hpp src/cli/SigGBAction.cpp
# When making a distribution file, Automake knows to include all files
# that are necessary to build the project. EXTRA_DIST specifies files
diff --git a/build/vs12/mathicgb-exe/mathicgb-exe.vcxproj b/build/vs12/mathicgb-exe/mathicgb-exe.vcxproj
index dc8095d..426958b 100755
--- a/build/vs12/mathicgb-exe/mathicgb-exe.vcxproj
+++ b/build/vs12/mathicgb-exe/mathicgb-exe.vcxproj
@@ -43,9 +43,11 @@
</ProjectConfiguration>
</ItemGroup>
<ItemGroup>
+ <ClCompile Include="..\..\..\src\cli\CommonParams.cpp" />
<ClCompile Include="..\..\..\src\cli\GBAction.cpp" />
<ClCompile Include="..\..\..\src\cli\GBCommonParams.cpp" />
<ClCompile Include="..\..\..\src\cli\GBMain.cpp" />
+ <ClCompile Include="..\..\..\src\cli\MatrixAction.cpp" />
<ClCompile Include="..\..\..\src\cli\SigGBAction.cpp" />
</ItemGroup>
<ItemGroup>
@@ -60,8 +62,10 @@
</ProjectReference>
</ItemGroup>
<ItemGroup>
+ <ClInclude Include="..\..\..\src\cli\CommonParams.hpp" />
<ClInclude Include="..\..\..\src\cli\GBAction.hpp" />
<ClInclude Include="..\..\..\src\cli\GBCommonParams.hpp" />
+ <ClInclude Include="..\..\..\src\cli\MatrixAction.hpp" />
<ClInclude Include="..\..\..\src\cli\SigGBAction.hpp" />
</ItemGroup>
<PropertyGroup Label="Globals">
diff --git a/build/vs12/mathicgb-exe/mathicgb-exe.vcxproj.filters b/build/vs12/mathicgb-exe/mathicgb-exe.vcxproj.filters
index c8ac080..0a02513 100755
--- a/build/vs12/mathicgb-exe/mathicgb-exe.vcxproj.filters
+++ b/build/vs12/mathicgb-exe/mathicgb-exe.vcxproj.filters
@@ -27,16 +27,28 @@
<ClCompile Include="..\..\..\src\cli\SigGBAction.cpp">
<Filter>Source Files</Filter>
</ClCompile>
+ <ClCompile Include="..\..\..\src\cli\MatrixAction.cpp">
+ <Filter>Source Files</Filter>
+ </ClCompile>
+ <ClCompile Include="..\..\..\src\cli\CommonParams.cpp">
+ <Filter>Source Files</Filter>
+ </ClCompile>
</ItemGroup>
<ItemGroup>
<ClInclude Include="..\..\..\src\cli\GBCommonParams.hpp">
<Filter>Header Files</Filter>
</ClInclude>
<ClInclude Include="..\..\..\src\cli\GBAction.hpp">
- <Filter>Source Files</Filter>
+ <Filter>Header Files</Filter>
</ClInclude>
<ClInclude Include="..\..\..\src\cli\SigGBAction.hpp">
- <Filter>Source Files</Filter>
+ <Filter>Header Files</Filter>
+ </ClInclude>
+ <ClInclude Include="..\..\..\src\cli\MatrixAction.hpp">
+ <Filter>Header Files</Filter>
+ </ClInclude>
+ <ClInclude Include="..\..\..\src\cli\CommonParams.hpp">
+ <Filter>Header Files</Filter>
</ClInclude>
</ItemGroup>
</Project>
\ No newline at end of file
diff --git a/src/cli/CommonParams.cpp b/src/cli/CommonParams.cpp
new file mode 100644
index 0000000..fecbaec
--- /dev/null
+++ b/src/cli/CommonParams.cpp
@@ -0,0 +1,77 @@
+#include "mathicgb/stdinc.h"
+#include "CommonParams.hpp"
+
+CommonParams::CommonParams():
+ mInputFile("inputFile",
+ "The file to read input from.",
+ ""),
+
+ mTracingLevel("tracingLevel",
+ "How much information to print out about what the program does. No "
+ "information is shown if the value is zero. Higher values "
+ "result in more information.",
+ 0),
+
+ mThreadCount("threadCount",
+ "Specifies how many threads to use at a time.",
+ 1)
+{
+}
+
+void CommonParams::directOptions(
+ std::vector<std::string> tokens,
+ mathic::CliParser& parser
+) {
+ if (tokens.size() == 1)
+ mInputFile.processArgument(tokens.back());
+ if (tokens.size() > 1)
+ mathic::reportError("Too many direct options.");
+}
+
+void CommonParams::pushBackParameters(
+ std::vector<mathic::CliParameter*>& parameters
+) {
+ parameters.push_back(&mInputFile);
+ parameters.push_back(&mTracingLevel);
+ parameters.push_back(&mThreadCount);
+}
+
+void CommonParams::perform() {
+ tracingLevel = mTracingLevel.value();
+
+ // delete the old init object first to make the new one take control.
+ mTbbInit.reset();
+ std::unique_ptr<tbb::task_scheduler_init> mTbbInit;
+
+ mTbbInit = make_unique<tbb::task_scheduler_init>(
+ mThreadCount.value() == 0 ?
+ tbb::task_scheduler_init::automatic :
+ mThreadCount.value()
+ );
+}
+
+void CommonParams::registerFileNameExtension(std::string extension) {
+ MATHICGB_ASSERT(!extension.empty());
+ mExtensions.push_back(std::move(extension));
+}
+
+std::string CommonParams::inputFileNameStem() {
+ const auto& str = mInputFile.value();
+ const auto toStrip = inputFileNameExtension();
+ MATHICGB_ASSERT
+ (toStrip.size() < str.size() || (toStrip.empty() && str.empty()));
+ return str.substr(0, str.size() - toStrip.size());
+}
+
+std::string CommonParams::inputFileNameExtension() {
+ const auto& str = mInputFile.value();
+ const auto end = mExtensions.end();
+ for (auto it = mExtensions.begin(); it != end; ++it) {
+ if (
+ str.size() >= it->size() &&
+ str.substr(str.size() - it->size(), it->size()) == *it
+ )
+ return *it;
+ }
+ return std::string();
+}
diff --git a/src/cli/CommonParams.hpp b/src/cli/CommonParams.hpp
new file mode 100644
index 0000000..b54ab63
--- /dev/null
+++ b/src/cli/CommonParams.hpp
@@ -0,0 +1,41 @@
+#ifndef MATHICGB_COMMON_PARAMS_GUARD
+#define MATHICGB_COMMON_PARAMS_GUARD
+
+#include <mathic.h>
+#include <tbb/tbb.h>
+#include <vector>
+
+class CommonParams {
+public:
+ CommonParams();
+
+ void directOptions
+ (std::vector<std::string> tokens, mathic::CliParser& parser);
+
+ void pushBackParameters(std::vector<mathic::CliParameter*>& parameters);
+
+ /// Takes appropriate action depending on the parameters. For example this
+ /// will set the number of threads in tbb.
+ void perform();
+
+ /// If called with string X, then X will be considered an extension
+ /// for a file name instead of part of the file name.
+ void registerFileNameExtension(std::string extensions);
+
+ /// Returns the stem of the input file name, with any registered extensions
+ /// stripped off.
+ std::string inputFileNameStem();
+
+ /// Returns the registered extension of the input file name, if any.
+ std::string inputFileNameExtension();
+
+ mathic::StringParameter mInputFile;
+ mathic::IntegerParameter mTracingLevel;
+ mathic::IntegerParameter mThreadCount;
+
+private:
+ std::vector<std::string> mExtensions; // to recognize file type
+ std::unique_ptr<tbb::task_scheduler_init> mTbbInit; // to set thread count
+};
+
+#endif
diff --git a/src/cli/GBAction.cpp b/src/cli/GBAction.cpp
new file mode 100644
index 0000000..4b7f4e8
--- /dev/null
+++ b/src/cli/GBAction.cpp
@@ -0,0 +1,137 @@
+#include "mathicgb/stdinc.h"
+#include "GBAction.hpp"
+
+#include "mathicgb/BuchbergerAlg.hpp"
+#include "mathicgb/Ideal.hpp"
+#include "mathicgb/io-util.hpp"
+#include "mathicgb/F4Reducer.hpp"
+#include <fstream>
+#include <iostream>
+
+GBAction::GBAction():
+ mAutoTailReduce("autoTailReduce",
+ "Reduce the non-leading terms of all polynomials whenever an element "
+ "is inserted into the basis. Only relevant to the "
+ "classic Buchberger algorithm.",
+ false),
+
+ mAutoTopReduce("autoTopReduce",
+ "Reduce any basis element whose lead term becomes reducible "
+ "by a different basis element. Only relevant to the "
+ "classic Buchberger algorithm.",
+ true),
+
+ //mTermOrder("termOrder",
+ // "The term order to compute a Grobner basis with respect to. This is "
+ // "currently actually a free module term order, but that should be changed.",
+ // 4),
+
+ mSPairGroupSize("sPairGroupSize",
+ "Specifies how many S-pair to reduce at one time. A value of 0 "
+ "indicates not to group S-pairs together. Only currently relevant "
+ "for the classic Buchberger algorithm.",
+ 0),
+
+ mMinMatrixToStore("storeMatrices",
+ "If using a matrix-based reducer, store the matrices that are generated in "
+ "files named X-1.mat, X-2.mat and so on where X is the project name. Only "
+ "matrices with at least as many entries as the parameter are stored. "
+ "A value of 0 indicates not to store any matrices.",
+ 0)
+
+{
+ std::ostringstream orderOut;
+ FreeModuleOrder::displayOrderTypes(orderOut);
+ //mTermOrder.appendToDescription(orderOut.str());
+}
+
+void GBAction::directOptions(
+ std::vector<std::string> tokens,
+ mic::CliParser& parser
+) {
+ mParams.directOptions(tokens, parser);
+}
+
+void GBAction::performAction() {
+ mParams.perform();
+ mGBParams.perform();
+ const std::string projectName = mParams.inputFileNameStem();
+
+ // read input
+ std::unique_ptr<Ideal> ideal;
+ {
+ const std::string inputIdealFile = projectName + ".ideal";
+ std::ifstream inputFile(inputIdealFile.c_str());
+ if (inputFile.fail())
+ mic::reportError("Could not read input file \"" + inputIdealFile + '\n');
+ ideal = Ideal::parse(inputFile);
+ }
+ std::unique_ptr<PolyRing const> ring(&(ideal->ring()));
+
+ // run algorithm
+ const auto reducerType = Reducer::reducerType(mGBParams.mReducer.value());
+ std::unique_ptr<Reducer> reducer;
+ if (reducerType != Reducer::Reducer_F4)
+ reducer = Reducer::makeReducer(reducerType, *ring);
+ else {
+ auto f4Reducer = make_unique<F4Reducer>(*ring);
+ if (mMinMatrixToStore.value() > 0)
+ f4Reducer->writeMatricesTo(projectName, mMinMatrixToStore);
+ reducer = std::move(f4Reducer);
+ }
+
+ BuchbergerAlg alg(
+ *ideal,
+ 4 /*mModuleOrder.value()*/ , // todo: alg should not take a *module* order
+ *reducer,
+ mGBParams.mDivisorLookup.value(),
+ mGBParams.mPreferSparseReducers.value(),
+ mGBParams.mSPairQueue.value());
+ alg.setBreakAfter(mGBParams.mBreakAfter.value());
+ alg.setPrintInterval(mGBParams.mPrintInterval.value());
+ alg.setSPairGroupSize(mSPairGroupSize.value());
+ alg.setReducerMemoryQuantum(mGBParams.mMemoryQuantum.value());
+ alg.setUseAutoTopReduction(mAutoTopReduce.value());
+ alg.setUseAutoTailReduction(mAutoTailReduce.value());
+
+ alg.computeGrobnerBasis();
+ alg.printStats(std::cerr);
+
+ /*
+ std::ofstream statsOut((mProjectName.value() + ".stats").c_str());
+ alg.printStats(statsOut);
+ {
+ std::string basisFileName = mProjectName.value() + ".gb";
+ FILE* basisOut = std::fopen(basisFileName.c_str(), "w");
+ output(basisOut, alg.basis());
+ }
+ */
+}
+
+const char* GBAction::staticName() {
+ return "gb";
+}
+
+const char* GBAction::name() const {
+ return staticName();
+}
+
+const char* GBAction::description() const {
+ return "Compute a Grobner basis. "
+ "The project name is an optional direct parameter.";
+}
+
+const char* GBAction::shortDescription() const {
+ return "Compute a Grobner basis.";
+}
+
+void GBAction::pushBackParameters(
+ std::vector<mic::CliParameter*>& parameters
+) {
+ mParams.pushBackParameters(parameters);
+ mGBParams.pushBackParameters(parameters);
+ parameters.push_back(&mAutoTailReduce);
+ parameters.push_back(&mAutoTopReduce);
+ parameters.push_back(&mSPairGroupSize);
+ parameters.push_back(&mMinMatrixToStore);
+}
diff --git a/src/cli/GBAction.hpp b/src/cli/GBAction.hpp
new file mode 100644
index 0000000..83849ff
--- /dev/null
+++ b/src/cli/GBAction.hpp
@@ -0,0 +1,38 @@
+#ifndef MATHICGB_G_B_ACTION_GUARD
+#define MATHICGB_G_B_ACTION_GUARD
+
+#include "GBCommonParams.hpp"
+#include "CommonParams.hpp"
+#include <mathic.h>
+
+/// Calculates a classic Grobner basis using Buchberger's algorithm
+class GBAction : public mathic::Action {
+public:
+ GBAction();
+
+ virtual void directOptions(
+ std::vector<std::string> tokens,
+ mic::CliParser& parser
+ );
+
+ virtual void performAction();
+
+ static const char* staticName();
+
+ virtual const char* name() const;
+ virtual const char* description() const;
+ virtual const char* shortDescription() const;
+
+ virtual void pushBackParameters(std::vector<mic::CliParameter*>& parameters);
+
+private:
+ CommonParams mParams;
+ GBCommonParams mGBParams;
+ mathic::BoolParameter mAutoTailReduce;
+ mathic::BoolParameter mAutoTopReduce;
+ //mic::IntegerParameter mTermOrder;
+ mathic::IntegerParameter mSPairGroupSize;
+ mathic::IntegerParameter mMinMatrixToStore;
+};
+
+#endif
diff --git a/src/cli/GBCommonParams.cpp b/src/cli/GBCommonParams.cpp
new file mode 100755
index 0000000..54117f2
--- /dev/null
+++ b/src/cli/GBCommonParams.cpp
@@ -0,0 +1,83 @@
+#include "mathicgb/stdinc.h"
+#include "GBCommonParams.hpp"
+
+#include "mathicgb/MTArray.hpp"
+#include "mathicgb/PolyReducer.hpp"
+#include "mathicgb/DivisorLookup.hpp"
+
+GBCommonParams::GBCommonParams():
+ mPreferSparseReducers("preferSparseReducers",
+ "If true, always use the sparsest reducer in polynomial reduction. "
+ "This option impacts both classic and signature constrained "
+ "polynomial reduction. Ties are broken by taking the oldest reducer. "
+ "If this option is false, the oldest reducer is always used.",
+ true),
+
+ mSPairQueue("spairQueue",
+ "The priority queue used to order S-pairs.\n"
+ " 0 tournament tree in front of triangle\n"
+ " 1 heap in front of triangle\n"
+ " 2 tournament tree\n"
+ " 3 heap\n",
+ 0),
+
+ mBreakAfter("breakAfter",
+ "Stop the computation after this many elements have been added to "
+ "the basis. The computation runs uninterrupted if the value is zero.",
+ 0),
+
+ mPrintInterval("printInterval",
+ "Print information about the computation every time this many S-pair "
+ "reductions have been performed. Do not print information like this "
+ "during the computation if the value is zero.",
+ 0),
+
+ mMonomialTable("monomialTable",
+ "The kind of monomial table data structure to use.\n",
+ 2),
+
+ mDivisorLookup("divisorLookup",
+ "The divisor lookup data structure to use.\n",
+ 2),
+
+ mReducer("reducer",
+ "The data structure to use for polynomial reduction.\n",
+ 4),
+
+ mMemoryQuantum("memoryQuantumForReducer",
+ "Specifies how many items to allocate memory for at a time for the reducer.",
+ 1024 * 1024)
+{
+ {
+ std::ostringstream reducerOut;
+ Reducer::displayReducerTypes(reducerOut);
+ mReducer.appendToDescription(reducerOut.str());
+ }
+ {
+ std::ostringstream divisorLookupOut;
+ DivisorLookup::displayDivisorLookupTypes(divisorLookupOut);
+ mDivisorLookup.appendToDescription(divisorLookupOut.str());
+ }
+ {
+ std::ostringstream monomialTableOut;
+ MonomialTableArray::displayMTTypes(monomialTableOut);
+ mMonomialTable.appendToDescription(monomialTableOut.str());
+ }
+}
+
+void GBCommonParams::pushBackParameters(
+ std::vector<mathic::CliParameter*>& parameters
+) {
+ parameters.push_back(&mPreferSparseReducers);
+ parameters.push_back(&mSPairQueue);
+ parameters.push_back(&mBreakAfter);
+ parameters.push_back(&mPrintInterval);
+ parameters.push_back(&mMonomialTable);
+ parameters.push_back(&mDivisorLookup);
+ parameters.push_back(&mReducer);
+ parameters.push_back(&mMemoryQuantum);
+}
+
+void GBCommonParams::perform() {
+ // currently there is nothing to do
+}
diff --git a/src/cli/GBCommonParams.hpp b/src/cli/GBCommonParams.hpp
new file mode 100755
index 0000000..c7f353c
--- /dev/null
+++ b/src/cli/GBCommonParams.hpp
@@ -0,0 +1,23 @@
+#ifndef MATHICGB_GB_COMMON_PARAMS_GUARD
+#define MATHICGB_GB_COMMON_PARAMS_GUARD
+
+#include <mathic.h>
+
+class GBCommonParams {
+public:
+ GBCommonParams();
+
+ void pushBackParameters(std::vector<mathic::CliParameter*>& parameters);
+ void perform();
+
+ mathic::BoolParameter mPreferSparseReducers;
+ mathic::IntegerParameter mSPairQueue;
+ mathic::IntegerParameter mBreakAfter;
+ mathic::IntegerParameter mPrintInterval;
+ mathic::IntegerParameter mMonomialTable;
+ mathic::IntegerParameter mDivisorLookup;
+ mathic::IntegerParameter mReducer;
+ mathic::IntegerParameter mMemoryQuantum;
+};
+
+#endif
diff --git a/src/cli/GBMain.cpp b/src/cli/GBMain.cpp
index 6725956..3c10e58 100755
--- a/src/cli/GBMain.cpp
+++ b/src/cli/GBMain.cpp
@@ -2,6 +2,7 @@
#include "GBAction.hpp"
#include "SigGBAction.hpp"
+#include "MatrixAction.hpp"
#include <mathic.h>
#include <cctype>
#include <iostream>
@@ -12,6 +13,7 @@ int main(int argc, char **argv) {
mathic::CliParser parser;
parser.registerAction<SigGBAction>();
parser.registerAction<GBAction>();
+ parser.registerAction<MatrixAction>();
parser.registerAction<mathic::HelpAction>();
std::vector<std::string> commandLine(argv, argv + argc);
diff --git a/src/cli/MatrixAction.cpp b/src/cli/MatrixAction.cpp
new file mode 100644
index 0000000..c5545e2
--- /dev/null
+++ b/src/cli/MatrixAction.cpp
@@ -0,0 +1,64 @@
+#include "mathicgb/stdinc.h"
+#include "MatrixAction.hpp"
+
+#include "mathicgb/QuadMatrix.hpp"
+#include "mathicgb/SparseMatrix.hpp"
+#include <fstream>
+#include <iostream>
+
+MatrixAction::MatrixAction() {
+ mParams.registerFileNameExtension(".mat");
+}
+
+void MatrixAction::directOptions(
+ std::vector<std::string> tokens,
+ mic::CliParser& parser
+) {
+ mParams.directOptions(tokens, parser);
+}
+
+void MatrixAction::performAction() {
+ mParams.perform();
+
+ QuadMatrix matrix;
+
+ // @todo: fix leak of file on exception
+
+ // read matrix
+ SparseMatrix::Scalar modulus;
+ {
+ FILE* file = fopen((mParams.inputFileNameStem() + ".mat").c_str(), "rb");
+ modulus = matrix.read(file);
+ fclose(file);
+ }
+
+ // write matrix
+ {
+ FILE* file = fopen((mParams.inputFileNameStem() + ".out").c_str(), "wb");
+ matrix.write(modulus, file);
+ fclose(file);
+ }
+}
+
+const char* MatrixAction::staticName() {
+ return "matrix";
+}
+
+const char* MatrixAction::name() const {
+ return staticName();
+}
+
+const char* MatrixAction::description() const {
+ return "Perform matrix computations. "
+ "The name of the matrix file is an optional direct parameter.";
+}
+
+const char* MatrixAction::shortDescription() const {
+ return "Perform matrix computations.";
+}
+
+void MatrixAction::pushBackParameters(
+ std::vector<mic::CliParameter*>& parameters
+) {
+ mParams.pushBackParameters(parameters);
+}
diff --git a/src/cli/MatrixAction.hpp b/src/cli/MatrixAction.hpp
new file mode 100644
index 0000000..7156454
--- /dev/null
+++ b/src/cli/MatrixAction.hpp
@@ -0,0 +1,31 @@
+#ifndef MATHICGB_MATRIX_ACTION_GUARD
+#define MATHICGB_MATRIX_ACTION_GUARD
+
+#include "CommonParams.hpp"
+#include <mathic.h>
+
+/// Performs computations on matrices.
+class MatrixAction : public mathic::Action {
+public:
+ MatrixAction();
+
+ virtual void directOptions(
+ std::vector<std::string> tokens,
+ mic::CliParser& parser
+ );
+
+ virtual void performAction();
+
+ static const char* staticName();
+
+ virtual const char* name() const;
+ virtual const char* description() const;
+ virtual const char* shortDescription() const;
+
+ virtual void pushBackParameters(std::vector<mic::CliParameter*>& parameters);
+
+private:
+ CommonParams mParams;
+};
+
+#endif
diff --git a/src/cli/SigGBAction.cpp b/src/cli/SigGBAction.cpp
new file mode 100644
index 0000000..6bb2bb0
--- /dev/null
+++ b/src/cli/SigGBAction.cpp
@@ -0,0 +1,125 @@
+#include "mathicgb/stdinc.h"
+#include "SigGBAction.hpp"
+
+#include "mathicgb/Ideal.hpp"
+#include "mathicgb/SignatureGB.hpp"
+#include "mathicgb/io-util.hpp"
+#include <fstream>
+#include <iostream>
+
+SigGBAction::SigGBAction():
+ mUseSingularCriterionEarly("earlySingularCriterion",
+ "Apply the singular S-pair elimination criterion before queueing "
+ "that S-pair. Otherwise, the criterion is only checked just before "
+ "the S-pair would otherwise cause a polynomial reduction to occur. "
+ "This criterion is only relevant to the signature Buchberger "
+ "algorithm.",
+ false),
+
+ mPostponeKoszul("postponeKoszul",
+ "Postpone the construction of Koszul syzygy signatures.",
+ true),
+
+ mUseBaseDivisors("useBaseDivisors",
+ "Use high ratio and low ratio base divisors to eliminate "
+ "S-spairs quickly based on signature.",
+ true),
+
+ mModuleOrder("moduleOrder",
+ "The free module term order.\n",
+ 4)
+{
+ std::ostringstream orderOut;
+ FreeModuleOrder::displayOrderTypes(orderOut);
+ mModuleOrder.appendToDescription(orderOut.str());
+}
+
+void SigGBAction::directOptions(
+ std::vector<std::string> tokens,
+ mic::CliParser& parser
+) {
+ mParams.directOptions(tokens, parser);
+}
+
+void SigGBAction::performAction() {
+ mParams.perform();
+ mGBParams.perform();
+
+ // read input file
+ std::unique_ptr<Ideal> ideal;
+ {
+ std::string const inputIdealFile = mParams.inputFileNameStem() + ".ideal";
+ std::ifstream inputFile(inputIdealFile.c_str());
+ if (inputFile.fail())
+ mic::reportError("Could not read input file \"" + inputIdealFile + '\n');
+ ideal = Ideal::parse(inputFile);
+ }
+ std::unique_ptr<PolyRing const> ring(&(ideal->ring()));
+
+ SignatureGB alg(
+ *ideal,
+ mModuleOrder.value(),
+ Reducer::reducerType(mGBParams.mReducer.value()),
+ mGBParams.mDivisorLookup.value(),
+ mGBParams.mMonomialTable.value(),
+ mPostponeKoszul.value(),
+ mUseBaseDivisors.value(),
+ mGBParams.mPreferSparseReducers.value(),
+ mUseSingularCriterionEarly.value(),
+ mGBParams.mSPairQueue.value());
+ alg.setBreakAfter(mGBParams.mBreakAfter.value());
+ alg.setPrintInterval(mGBParams.mPrintInterval.value());
+ alg.computeGrobnerBasis();
+
+ // print statistics
+ alg.displayStats(std::cout);
+ alg.displayPaperStats(std::cout);
+ {
+ std::ofstream statsOut((mParams.inputFileNameStem() + ".stats").c_str());
+ alg.displayStats(statsOut);
+ alg.displayPaperStats(statsOut);
+ }
+
+ // print basis
+ {
+ std::ofstream ogb((mParams.inputFileNameStem() + ".gb").c_str());
+ ogb << "-- gb: ----\n";
+ alg.getGB()->display(ogb);
+ }
+
+ // print syzygy basis
+ {
+ std::ofstream syzygyOut((mParams.inputFileNameStem() + ".syz").c_str());
+ syzygyOut << "-- syz: ----\n";
+ alg.getSyzTable()->display(syzygyOut, 1);
+ syzygyOut << std::endl;
+ }
+}
+
+const char* SigGBAction::staticName() {
+ return "siggb";
+}
+
+const char* SigGBAction::name() const {
+ return staticName();
+}
+
+const char* SigGBAction::description() const {
+ return "Compute a signature Grobner basis. "
+ "The project name is an optional direct parameter.";
+}
+
+const char* SigGBAction::shortDescription() const {
+ return "Compute a signature Grobner basis";
+}
+
+void SigGBAction::pushBackParameters(
+ std::vector<mic::CliParameter*>& parameters
+) {
+ mParams.pushBackParameters(parameters);
+ mGBParams.pushBackParameters(parameters);
+ parameters.push_back(&mUseSingularCriterionEarly);
+ parameters.push_back(&mPostponeKoszul);
+ parameters.push_back(&mUseBaseDivisors);
+ parameters.push_back(&mModuleOrder);
+}
diff --git a/src/cli/SigGBAction.hpp b/src/cli/SigGBAction.hpp
new file mode 100644
index 0000000..b63b409
--- /dev/null
+++ b/src/cli/SigGBAction.hpp
@@ -0,0 +1,36 @@
+#ifndef MATHICGB_SIG_G_B_ACTION_GUARD
+#define MATHICGB_SIG_G_B_ACTION_GUARD
+
+#include "GBCommonParams.hpp"
+#include "CommonParams.hpp"
+#include <mathic.h>
+
+class SigGBAction : public mic::Action {
+public:
+ SigGBAction();
+
+ virtual void directOptions(
+ std::vector<std::string> tokens,
+ mic::CliParser& parser
+ );
+
+ virtual void performAction();
+
+ static const char* staticName();
+
+ virtual const char* name() const;
+ virtual const char* description() const;
+ virtual const char* shortDescription() const;
+
+ virtual void pushBackParameters(std::vector<mic::CliParameter*>& parameters);
+
+private:
+ CommonParams mParams;
+ GBCommonParams mGBParams;
+ mic::BoolParameter mUseSingularCriterionEarly;
+ mic::BoolParameter mPostponeKoszul;
+ mic::BoolParameter mUseBaseDivisors;
+ mic::IntegerParameter mModuleOrder;
+};
+
+#endif
diff --git a/src/mathicgb/BuchbergerAlg.cpp b/src/mathicgb/BuchbergerAlg.cpp
index 32e45c6..9386bd1 100755
--- a/src/mathicgb/BuchbergerAlg.cpp
+++ b/src/mathicgb/BuchbergerAlg.cpp
@@ -7,7 +7,7 @@
BuchbergerAlg::BuchbergerAlg(
const Ideal& ideal,
FreeModuleOrderType orderType,
- Reducer::ReducerType reducerType,
+ Reducer& reducer,
int divisorLookupType,
bool preferSparseReducers,
size_t queueType
@@ -19,7 +19,7 @@ BuchbergerAlg::BuchbergerAlg(
mUseAutoTailReduction(false),
mRing(*ideal.getPolyRing()),
mOrder(FreeModuleOrder::makeOrder(orderType, &ideal)),
- mReducer(Reducer::makeReducer(reducerType, *ideal.getPolyRing())),
+ mReducer(reducer),
mBasis(mRing, *mOrder, DivisorLookup::makeFactory(
*ideal.getPolyRing(),
divisorLookupType)->create(preferSparseReducers, true)
@@ -44,7 +44,7 @@ void BuchbergerAlg::insertPolys
if ((*it)->isZero())
continue;
if (mBasis.divisor((*it)->getLeadMonomial()) != static_cast<size_t>(-1)) {
- *it = mReducer->classicReduce(**it, mBasis);
+ *it = mReducer.classicReduce(**it, mBasis);
if ((*it)->isZero())
continue;
}
@@ -90,7 +90,7 @@ void BuchbergerAlg::insertPolys
MATHICGB_ASSERT(toRetire.empty());
// reduce everything in toReduce
- mReducer->classicReducePolySet(toReduce, mBasis, toInsert);
+ mReducer.classicReducePolySet(toReduce, mBasis, toInsert);
toReduce.clear();
}
@@ -140,7 +140,7 @@ void BuchbergerAlg::insertReducedPoly(
if (toReduce.empty()) // if first iteration
reduced = std::move(polyToInsert);
else {
- reduced = mReducer->classicReduce(*toReduce.back(), mBasis);
+ reduced = mReducer.classicReduce(*toReduce.back(), mBasis);
if (tracingLevel > 20) {
if (reduced->isZero()) {
std::cerr << "auto-top-reduce cascade: "
@@ -236,7 +236,7 @@ void BuchbergerAlg::step() {
<< p.second << ")" << std::endl;
}
std::unique_ptr<Poly> reduced
- (mReducer->classicReduceSPoly
+ (mReducer.classicReduceSPoly
(mBasis.poly(p.first), mBasis.poly(p.second), mBasis));
if (!reduced->isZero()) {
insertReducedPoly(std::move(reduced));
@@ -263,7 +263,7 @@ void BuchbergerAlg::step() {
if (spairGroup.empty())
return; // no more s-pairs
std::vector<std::unique_ptr<Poly> > reduced;
- mReducer->classicReduceSPolySet(spairGroup, mBasis, reduced);
+ mReducer.classicReduceSPolySet(spairGroup, mBasis, reduced);
insertPolys(reduced);
/*
@@ -271,7 +271,7 @@ void BuchbergerAlg::step() {
auto p = std::move(*it);
MATHICGB_ASSERT(!p->isZero());
if (it != reduced.begin())
- p = mReducer->classicReduce(*p, mBasis);
+ p = mReducer.classicReduce(*p, mBasis);
if (!p->isZero())
insertReducedPoly(std::move(p));
}
@@ -290,7 +290,7 @@ void BuchbergerAlg::autoTailReduce() {
if (mBasis.usedAsReducerCount(i) < 1000)
continue;
mBasis.replaceSameLeadTerm
- (i, mReducer->classicTailReduce(mBasis.poly(i), mBasis));
+ (i, mReducer.classicTailReduce(mBasis.poly(i), mBasis));
}
}
@@ -298,13 +298,13 @@ size_t BuchbergerAlg::getMemoryUse() const {
return
mBasis.getMemoryUse() +
mRing.getMonomialPool().getMemoryUse() +
- mReducer->getMemoryUse() +
+ mReducer.getMemoryUse() +
mSPairs.getMemoryUse();
}
void BuchbergerAlg::printStats(std::ostream& out) const {
out << " term order: " << mBasis.order().description() << '\n';
- out << " reduction type: " << mReducer->description() << '\n';
+ out << " reduction type: " << mReducer.description() << '\n';
out << " divisor tab type: " << mBasis.divisorLookup().getName() << '\n';
out << " S-pair queue type: " << mSPairs.name() << '\n';
out << " total compute time: " << mTimer.getMilliseconds()/1000.0 << " seconds " << '\n';
@@ -376,7 +376,7 @@ void BuchbergerAlg::printStats(std::ostream& out) const {
value << mic::ColumnPrinter::commafy(reductions) << '\n';
extra << '\n';
- Reducer::Stats reducerStats = mReducer->sigStats();
+ Reducer::Stats reducerStats = mReducer.sigStats();
SPairs::Stats sPairStats = mSPairs.stats();
unsigned long long const primeElim = sPairStats.relativelyPrimeHits;
@@ -418,7 +418,7 @@ void BuchbergerAlg::printStats(std::ostream& out) const {
value << mic::ColumnPrinter::commafy(longestReduction) << '\n';
extra << '\n';
- Reducer::Stats classicRedStats = mReducer->classicStats();
+ Reducer::Stats classicRedStats = mReducer.classicStats();
const unsigned long long clReductions = classicRedStats.reductions;
name << "Classic reductions:\n";
value << mic::ColumnPrinter::commafy(clReductions) << '\n';
@@ -509,7 +509,7 @@ void BuchbergerAlg::printMemoryUse(std::ostream& out) const
extra << mic::ColumnPrinter::percent(sPairMem, total) << '\n';
}
{ // Reducer
- const size_t reducerMem = mReducer->getMemoryUse();
+ const size_t reducerMem = mReducer.getMemoryUse();
name << "Reducer:\n";
value << mic::ColumnPrinter::bytesInUnit(reducerMem) << '\n';
extra << mic::ColumnPrinter::percent(reducerMem, total) << '\n';
diff --git a/src/mathicgb/BuchbergerAlg.hpp b/src/mathicgb/BuchbergerAlg.hpp
index 91cce45..712d90a 100755
--- a/src/mathicgb/BuchbergerAlg.hpp
+++ b/src/mathicgb/BuchbergerAlg.hpp
@@ -16,7 +16,7 @@ public:
BuchbergerAlg(
const Ideal& ideal,
FreeModuleOrderType orderType,
- Reducer::ReducerType reducerType,
+ Reducer& reducer,
int divisorLookupType,
bool preferSparseReducers,
size_t queueType);
@@ -51,7 +51,7 @@ public:
}
void setReducerMemoryQuantum(size_t memoryQuantum) {
- mReducer->setMemoryQuantum(memoryQuantum);
+ mReducer.setMemoryQuantum(memoryQuantum);
}
void setUseAutoTopReduction(bool value) {
@@ -81,7 +81,7 @@ private:
const PolyRing& mRing;
std::unique_ptr<FreeModuleOrder> mOrder;
- std::unique_ptr<Reducer> mReducer;
+ Reducer& mReducer;
PolyBasis mBasis;
SPairs mSPairs;
mic::Timer mTimer;
diff --git a/src/mathicgb/F4MatrixReducer.cpp b/src/mathicgb/F4MatrixReducer.cpp
index a8a37a4..28376df 100755
--- a/src/mathicgb/F4MatrixReducer.cpp
+++ b/src/mathicgb/F4MatrixReducer.cpp
@@ -335,7 +335,7 @@ namespace {
}});
//std::cout << "done reducing that batch" << std::endl;
- reduced.clear(colCount);
+ reduced.clear();
std::sort(nextReducers.begin(), nextReducers.end());
for (size_t i = 0; i < nextReducers.size(); ++i) {
size_t const row = nextReducers[i].second;
@@ -361,7 +361,7 @@ namespace {
dense[row].takeModulus(modulus);
}});
- toReduce.clear(colCount);
+ toReduce.clear();
for (size_t row = 0; row < rowCount; ++row)
if (!dense[row].empty())
dense[row].appendTo(toReduce);
diff --git a/src/mathicgb/F4Reducer.cpp b/src/mathicgb/F4Reducer.cpp
index 6cef465..5941f48 100755
--- a/src/mathicgb/F4Reducer.cpp
+++ b/src/mathicgb/F4Reducer.cpp
@@ -3,17 +3,24 @@
#include "F4MatrixBuilder.hpp"
#include "F4MatrixReducer.hpp"
+#include "QuadMatrix.hpp"
#include <iostream>
+#include <limits>
-F4Reducer::F4Reducer(
- const PolyRing& ring,
- std::unique_ptr<Reducer> fallback
-):
- mFallback(std::move(fallback)),
+F4Reducer::F4Reducer(const PolyRing& ring):
+ mFallback(Reducer::makeReducer(Reducer::Reducer_BjarkeGeo, ring)),
mRing(ring),
- mThreadCount(1),
- mMemoryQuantum(0) {
+ mMemoryQuantum(0),
+ mStoreToFile(""),
+ mMinEntryCountForStore(0),
+ mMatrixSaveCount(0) {
}
+
+void F4Reducer::writeMatricesTo(std::string file, size_t minEntries) {
+ mStoreToFile = std::move(file);
+ mMinEntryCountForStore = minEntries;
+ mMatrixSaveCount = 0;
+}
std::unique_ptr<Poly> F4Reducer::classicReduce
(const Poly& poly, const PolyBasis& basis) {
@@ -41,60 +48,32 @@ std::unique_ptr<Poly> F4Reducer::classicTailReduce
return p;
}
-std::unique_ptr<Poly> F4Reducer::classicReduceSPoly
-(const Poly& a, const Poly& b, const PolyBasis& basis) {
+std::unique_ptr<Poly> F4Reducer::classicReduceSPoly(
+ const Poly& a,
+ const Poly& b,
+ const PolyBasis& basis
+) {
if (tracingLevel >= 2)
- std::cerr << "F4Reducer: Reducing single S-pair.\n";
-
-
- QuadMatrix qm;
- {
- F4MatrixBuilder builder(basis, mThreadCount);
- builder.addSPolynomialToMatrix(a, b);
- builder.buildMatrixAndClear(qm);
-
- // there has to be something to reduce
- MATHICGB_ASSERT(qm.bottomLeft.rowCount() > 0);
- }
-
- SparseMatrix reduced;
- {
- F4MatrixReducer red(basis.ring());
- reduced = red.reduce(qm);
- }
-
- auto p = make_unique<Poly>(basis.ring());
- if (reduced.rowCount() > 0) {
- MATHICGB_ASSERT(reduced.rowCount() == 1);
- reduced.rowToPolynomial(0, qm.rightColumnMonomials, *p);
- }
-
- for (auto it = qm.leftColumnMonomials.begin();
- it != qm.leftColumnMonomials.end(); ++it)
- mRing.freeMonomial(*it);
- for (auto it = qm.rightColumnMonomials.begin();
- it != qm.rightColumnMonomials.end(); ++it)
- mRing.freeMonomial(*it);
- return p;
+ std::cerr << "F4Reducer: "
+ "Using fall-back reducer for single classic S-pair reduction\n";
+ return mFallback->classicReduceSPoly(a, b, basis);
}
-void F4Reducer::classicReduceSPolySet
-(std::vector<std::pair<size_t, size_t> >& spairs,
- const PolyBasis& basis,
- std::vector<std::unique_ptr<Poly> >& reducedOut)
-{
+void F4Reducer::classicReduceSPolySet(
+ std::vector<std::pair<size_t, size_t> >& spairs,
+ const PolyBasis& basis,
+ std::vector<std::unique_ptr<Poly> >& reducedOut
+) {
if (spairs.size() <= 1) {
if (tracingLevel >= 2)
std::cerr << "F4Reducer: Using fall-back reducer for "
- << spairs.size() << " S-pairs.\n";
+ << spairs.size() << " S-pairs.\n";
mFallback->classicReduceSPolySet(spairs, basis, reducedOut);
return;
}
-
reducedOut.clear();
- if (spairs.empty())
- return;
+ MATHICGB_ASSERT(!spairs.empty());
if (tracingLevel >= 2)
std::cerr << "F4Reducer: Reducing " << spairs.size() << " S-polynomials.\n";
@@ -113,10 +92,11 @@ void F4Reducer::classicReduceSPolySet
// there has to be something to reduce
MATHICGB_ASSERT(qm.bottomLeft.rowCount() > 0);
}
+ saveMatrix(qm);
reduced = F4MatrixReducer(basis.ring()).reduce(qm);
monomials = std::move(qm.rightColumnMonomials);
- for (auto it = qm.leftColumnMonomials.begin();
- it != qm.leftColumnMonomials.end(); ++it)
+ const auto end = qm.leftColumnMonomials.end();
+ for (auto it = qm.leftColumnMonomials.begin(); it != end; ++it)
mRing.freeMonomial(*it);
}
@@ -129,7 +109,8 @@ void F4Reducer::classicReduceSPolySet
reduced.rowToPolynomial(row, monomials, *p);
reducedOut.push_back(std::move(p));
}
- for (auto it = monomials.begin(); it != monomials.end(); ++it)
+ const auto end = monomials.end();
+ for (auto it = monomials.begin(); it != end; ++it)
mRing.freeMonomial(*it);
}
@@ -147,9 +128,8 @@ void F4Reducer::classicReducePolySet
}
reducedOut.clear();
- if (polys.empty())
- return;
+ MATHICGB_ASSERT(!polys.empty());
if (tracingLevel >= 2)
std::cerr << "F4Reducer: Reducing " << polys.size() << " polynomials.\n";
@@ -166,6 +146,7 @@ void F4Reducer::classicReducePolySet
// there has to be something to reduce
MATHICGB_ASSERT(qm.bottomLeft.rowCount() > 0);
}
+ saveMatrix(qm);
reduced = F4MatrixReducer(basis.ring()).reduce(qm);
monomials = std::move(qm.rightColumnMonomials);
for (auto it = qm.leftColumnMonomials.begin();
@@ -182,15 +163,20 @@ void F4Reducer::classicReducePolySet
reduced.rowToPolynomial(row, monomials, *p);
reducedOut.push_back(std::move(p));
}
- for (auto it = monomials.begin(); it != monomials.end(); ++it)
+ const auto end = monomials.end();
+ for (auto it = monomials.begin(); it != end; ++it)
mRing.freeMonomial(*it);
}
-Poly* F4Reducer::regularReduce
-(const_monomial sig,
- const_monomial multiple,
- size_t basisElement,
- const GroebnerBasis& basis) {
+Poly* F4Reducer::regularReduce(
+ const_monomial sig,
+ const_monomial multiple,
+ size_t basisElement,
+ const GroebnerBasis& basis
+) {
+ if (tracingLevel >= 2)
+ std::cerr <<
+ "F4Reducer: Using fall-back reducer for single regular reduction\n";
Poly* p = mFallback->regularReduce(sig, multiple, basisElement, basis);
mSigStats = mFallback->sigStats();
mClassicStats = mFallback->classicStats();
@@ -208,3 +194,21 @@ std::string F4Reducer::description() const {
size_t F4Reducer::getMemoryUse() const {
return 0; // @todo: implement
}
+
+void F4Reducer::saveMatrix(const QuadMatrix& matrix) {
+ const auto entryCount = matrix.entryCount();
+ if (mMinEntryCountForStore > entryCount)
+ return;
+ ++mMatrixSaveCount;
+ std::ostringstream fileName;
+ fileName << mStoreToFile << '-' << mMatrixSaveCount << ".mat";
+ if (tracingLevel > 2) {
+ std::cerr << "F4Reducer: Saving "
+ << mathic::ColumnPrinter::commafy(entryCount)
+ << " entry matrix to " << fileName.str() << '\n';
+ }
+ FILE* file = fopen(fileName.str().c_str(), "wb");
+ // @todo: fix leak of file on exception
+ matrix.write(static_cast<SparseMatrix::Scalar>(mRing.charac()), file);
+ fclose(file);
+}
diff --git a/src/mathicgb/F4Reducer.hpp b/src/mathicgb/F4Reducer.hpp
index a4eba3f..c86676b 100755
--- a/src/mathicgb/F4Reducer.hpp
+++ b/src/mathicgb/F4Reducer.hpp
@@ -3,36 +3,46 @@
#include "Reducer.hpp"
#include "PolyRing.hpp"
+#include <string>
+class QuadMatrix;
class F4Reducer : public Reducer {
public:
- F4Reducer(const PolyRing& ring,
- std::unique_ptr<Reducer> fallback);
+ F4Reducer(const PolyRing& ring);
+
+ /// Store all future matrices to file-1.mat, file-2.mat and so on.
+ /// Matrices with less than minEntries non-zero entries are not stored.
+ /// If file is an empty string then no matrices are stored. If this method
+ /// is never called then no matrices are stored.
+ void writeMatricesTo(std::string file, size_t minEntries);
virtual std::unique_ptr<Poly> classicReduce
- (const Poly& poly, const PolyBasis& basis);
+ (const Poly& poly, const PolyBasis& basis);
virtual std::unique_ptr<Poly> classicTailReduce
- (const Poly& poly, const PolyBasis& basis);
+ (const Poly& poly, const PolyBasis& basis);
virtual std::unique_ptr<Poly> classicReduceSPoly
- (const Poly& a, const Poly& b, const PolyBasis& basis);
+ (const Poly& a, const Poly& b, const PolyBasis& basis);
- virtual void classicReduceSPolySet
- (std::vector<std::pair<size_t, size_t> >& spairs,
- const PolyBasis& basis,
- std::vector<std::unique_ptr<Poly> >& reducedOut);
+ virtual void classicReduceSPolySet(
+ std::vector<std::pair<size_t, size_t> >& spairs,
+ const PolyBasis& basis,
+ std::vector<std::unique_ptr<Poly> >& reducedOut
+ );
- virtual void classicReducePolySet
- (const std::vector<std::unique_ptr<Poly> >& polys,
- const PolyBasis& basis,
- std::vector<std::unique_ptr<Poly> >& reducedOut);
+ virtual void classicReducePolySet(
+ const std::vector<std::unique_ptr<Poly> >& polys,
+ const PolyBasis& basis,
+ std::vector<std::unique_ptr<Poly> >& reducedOut
+ );
virtual Poly* regularReduce(
const_monomial sig,
const_monomial multiple,
size_t basisElement,
- const GroebnerBasis& basis);
+ const GroebnerBasis& basis
+ );
virtual void setMemoryQuantum(size_t quantum);
@@ -40,10 +50,14 @@ public:
virtual size_t getMemoryUse() const;
private:
+ void saveMatrix(const QuadMatrix& matrix);
+
std::unique_ptr<Reducer> mFallback;
const PolyRing& mRing;
- int mThreadCount;
size_t mMemoryQuantum;
+ std::string mStoreToFile; /// stem of file names to save matrices to
+ size_t mMinEntryCountForStore; /// don't save matrices with fewer entries
+ size_t mMatrixSaveCount; // how many matrices have been saved
};
#endif
diff --git a/src/mathicgb/QuadMatrix.cpp b/src/mathicgb/QuadMatrix.cpp
index edfe89e..b2315cd 100755
--- a/src/mathicgb/QuadMatrix.cpp
+++ b/src/mathicgb/QuadMatrix.cpp
@@ -28,9 +28,6 @@ bool QuadMatrix::debugAssertValid() const {
void QuadMatrix::print(std::ostream& out) const {
MATHICGB_ASSERT(debugAssertValid());
- // @todo: fix the code duplication here from QuadMatrixBuilder's
- // string code.
-
typedef SparseMatrix::ColIndex ColIndex;
mathic::ColumnPrinter printer;
printer.addColumn(true, "", "");
@@ -65,6 +62,24 @@ void QuadMatrix::print(std::ostream& out) const {
out << printer;
}
+size_t QuadMatrix::rowCount() const {
+ return topLeft.rowCount() + bottomLeft.rowCount();
+}
+
+size_t QuadMatrix::computeLeftColCount() const {
+ return std::max(topLeft.computeColCount(), bottomLeft.computeColCount());
+}
+
+size_t QuadMatrix::computeRightColCount() const {
+ return std::max(topRight.computeColCount(), bottomRight.computeColCount());
+}
+
+size_t QuadMatrix::entryCount() const {
+ return
+ topLeft.entryCount() + topRight.entryCount() +
+ bottomLeft.entryCount() + bottomRight.entryCount();
+}
+
std::string QuadMatrix::toString() const {
std::ostringstream out;
print(out);
@@ -131,9 +146,9 @@ void QuadMatrix::printSizes(std::ostream& out) const {
pr[2] << line << "/\n";
out << '\n' << pr
- << "Total memory: " << memoryUse() << " ("
- << ColPr::percent(memoryUseTrimmed(), memoryUse())
- << " written to)\n";
+ << "Total memory: " << ColPr::bytesInUnit(memoryUse()) << " ("
+ << ColPr::percent(memoryUseTrimmed(), memoryUse())
+ << " in use)\n";
}
QuadMatrix QuadMatrix::toCanonical() const {
@@ -300,3 +315,36 @@ void QuadMatrix::sortColumnsLeftRightParallel() {
}
});
}
+
+void QuadMatrix::write(
+ const SparseMatrix::Scalar modulus,
+ FILE* file
+) const {
+ MATHICGB_ASSERT(file != 0);
+ topLeft.write(modulus, file);
+ topRight.write(modulus, file);
+ bottomLeft.write(modulus, file);
+ bottomRight.write(modulus, file);
+}
+
+SparseMatrix::Scalar QuadMatrix::read(FILE* file) {
+ MATHICGB_ASSERT(file != 0);
+
+ leftColumnMonomials.clear();
+ rightColumnMonomials.clear();
+ ring = 0;
+
+ const auto topLeftModulus = topLeft.read(file);
+ const auto topRightModulus = topRight.read(file);
+ const auto bottomLeftModulus = bottomLeft.read(file);
+ const auto bottomRightModulus = bottomRight.read(file);
+
+ // todo: this should throw some kind of invalid format exception instead of
+ // these asserts.
+ MATHICGB_ASSERT(topLeftModulus == topRightModulus);
+ MATHICGB_ASSERT(bottomLeftModulus == topRightModulus);
+ MATHICGB_ASSERT(bottomRightModulus == topRightModulus);
+ MATHICGB_ASSERT(debugAssertValid());
+
+ return topLeftModulus;
+}
diff --git a/src/mathicgb/QuadMatrix.hpp b/src/mathicgb/QuadMatrix.hpp
index e9d858a..4e3d84f 100755
--- a/src/mathicgb/QuadMatrix.hpp
+++ b/src/mathicgb/QuadMatrix.hpp
@@ -55,6 +55,24 @@ public:
/// Shows whole matrix in a string. Useful for debugging.
std::string toString() const;
+ /// Return the combined number of non-zero entries.
+ size_t entryCount() const;
+
+ /// Return the combined number of left and right columns.
+ size_t rowCount() const;
+
+ /// Return the number of left columns.
+ size_t computeLeftColCount() const;
+
+ /// Return the number of right columns.
+ size_t computeRightColCount() const;
+
+ void write(SparseMatrix::Scalar modulus, FILE* file) const;
+
+ /// Read a matrix from file into *this. Return the modulus from file.
+ /// This method clears the column monomials and the ring pointer.
+ SparseMatrix::Scalar read(FILE* file);
+
/// Sort the left columns to be in decreasing order according to the monomial
/// order from the ring. The operation is done in parallel.
void sortColumnsLeftRightParallel();
diff --git a/src/mathicgb/Reducer.cpp b/src/mathicgb/Reducer.cpp
index cc6db0d..b465f79 100755
--- a/src/mathicgb/Reducer.cpp
+++ b/src/mathicgb/Reducer.cpp
@@ -60,19 +60,15 @@ std::unique_ptr<Reducer> Reducer::makeReducerNullOnUnknown(
return std::unique_ptr<Reducer>(new BjarkeGeobucket2(&ring));
case Reducer_TournamentTree:
return std::unique_ptr<Reducer>(new TournamentReducer(ring));
- //return std::unique_ptr<Reducer>
- //(new ReducerPack<mic::TourTree>(ring));
case Reducer_HashTourTree:
return std::unique_ptr<Reducer>(new HashTourReducer(ring));
-
case Reducer_TourTree_NoDedup:
return std::unique_ptr<Reducer>(new ReducerNoDedup<mic::TourTree>(ring));
case Reducer_TourTree_Dedup:
return std::unique_ptr<Reducer>(new ReducerDedup<mic::TourTree>(ring));
case Reducer_TourTree_Hashed:
return std::unique_ptr<Reducer>(new ReducerHash<mic::TourTree>(ring));
- //break;
case Reducer_TourTree_NoDedup_Packed:
return std::unique_ptr<Reducer>(new ReducerPack<mic::TourTree>(ring));
case Reducer_TourTree_Dedup_Packed:
@@ -86,7 +82,6 @@ std::unique_ptr<Reducer> Reducer::makeReducerNullOnUnknown(
return std::unique_ptr<Reducer>(new ReducerDedup<mic::Heap>(ring));
case Reducer_Heap_Hashed:
return std::unique_ptr<Reducer>(new ReducerHash<mic::Heap>(ring));
- //break;
case Reducer_Heap_NoDedup_Packed:
return std::unique_ptr<Reducer>(new ReducerPack<mic::Heap>(ring));
case Reducer_Heap_Dedup_Packed:
@@ -108,11 +103,7 @@ std::unique_ptr<Reducer> Reducer::makeReducerNullOnUnknown(
return std::unique_ptr<Reducer>(new ReducerHashPack<mic::Geobucket>(ring));
case Reducer_F4:
- {
- std::unique_ptr<Reducer> fallback = makeReducer(Reducer_BjarkeGeo, ring);
- return std::unique_ptr<Reducer>
- (new F4Reducer(ring, std::move(fallback)));
- }
+ return make_unique<F4Reducer>(ring);
default:
break;
diff --git a/src/mathicgb/Reducer.hpp b/src/mathicgb/Reducer.hpp
index 675e0f8..2de1983 100755
--- a/src/mathicgb/Reducer.hpp
+++ b/src/mathicgb/Reducer.hpp
@@ -102,13 +102,13 @@ public:
};
static std::unique_ptr<Reducer> makeReducer
- (ReducerType t, PolyRing const& ring);
+ (ReducerType t, const PolyRing& ring);
static std::unique_ptr<Reducer> makeReducerNullOnUnknown
- (ReducerType t, PolyRing const& ring);
+ (ReducerType t, const PolyRing& ring);
static ReducerType reducerType(int typ);
- static void displayReducerTypes(std::ostream &o);
+ static void displayReducerTypes(std::ostream& o);
// ***** Obtaining statistics about the reduction process
diff --git a/src/mathicgb/SparseMatrix.cpp b/src/mathicgb/SparseMatrix.cpp
index 3cef908..9a24622 100755
--- a/src/mathicgb/SparseMatrix.cpp
+++ b/src/mathicgb/SparseMatrix.cpp
@@ -61,7 +61,7 @@ void SparseMatrix::sortRowsByIncreasingPivots() {
std::sort(order.begin(), order.end());
// construct ordered with pivot columns in increaing order
- ordered.clear(lColCount);
+ ordered.clear();
for (size_t i = 0; i < lRowCount; ++i) {
const auto row = order[i].second;
const auto end = rowEnd(row);
@@ -165,7 +165,7 @@ SparseMatrix::ColIndex SparseMatrix::computeColCount() const {
return colCount;
}
-void SparseMatrix::clear(const ColIndex newColCount) {
+void SparseMatrix::clear() {
Block* block = &mBlock;
while (block != 0) {
delete[] block->mColIndices.releaseMemory();
@@ -271,12 +271,13 @@ void SparseMatrix::reserveFreeEntries(const size_t freeCount) {
std::distance(mRows.back().mIndicesEnd, mBlock.mColIndices.end())
);
+ // @todo: fix memory leaks on exception
+
auto oldBlock = new Block(std::move(mBlock));
MATHICGB_ASSERT(mBlock.mColIndices.begin() == 0);
MATHICGB_ASSERT(mBlock.mScalars.begin() == 0);
MATHICGB_ASSERT(mBlock.mHasNoRows);
MATHICGB_ASSERT(mBlock.mPreviousBlock == 0);
- mBlock.mPreviousBlock = oldBlock;
{
const auto begin = new ColIndex[count];
@@ -296,11 +297,20 @@ void SparseMatrix::reserveFreeEntries(const size_t freeCount) {
(oldBlock->mColIndices.begin(), oldBlock->mColIndices.end());
mBlock.mScalars.rawAssign
(oldBlock->mScalars.begin(), oldBlock->mScalars.end());
+ delete oldBlock; // no reason to keep it around
} else {
mBlock.mColIndices.rawAssign
(mRows.back().mIndicesEnd, oldBlock->mColIndices.end());
mBlock.mScalars.rawAssign
(mRows.back().mScalarsEnd, oldBlock->mScalars.end());
+
+ // remove the pending entries from old block so that counting the number
+ // of entries will give the correct result in future.
+ oldBlock->mColIndices.resize
+ (std::distance(oldBlock->mColIndices.begin(), mRows.back().mIndicesEnd));
+ oldBlock->mScalars.resize
+ (std::distance(oldBlock->mScalars.begin(), mRows.back().mScalarsEnd));
+ mBlock.mPreviousBlock = oldBlock;
}
}
@@ -350,3 +360,107 @@ std::ostream& operator<<(std::ostream& out, const SparseMatrix& matrix) {
matrix.print(out);
return out;
}
+
+namespace {
+ template<class T>
+ T readOne(FILE* file) {
+ T t;
+ fread(&t, sizeof(T), 1, file);
+ return t;
+ }
+
+ template<class T>
+ void writeOne(const T& t, FILE* file) {
+ fwrite(&t, sizeof(T), 1, file);
+ }
+
+ template<class T>
+ void writeMany(const std::vector<T>& v, FILE* file) {
+ if (v.empty())
+ return;
+ fwrite(&v[0], sizeof(T), v.size(), file);
+ }
+
+ template<class T>
+ void readMany(FILE* file, size_t count, std::vector<T>& v) {
+ size_t const originalSize = v.size();
+ v.resize(originalSize + count);
+ fread(&v[originalSize], sizeof(T), count, file);
+ }
+
+}
+
+void SparseMatrix::write(const Scalar modulus, FILE* file) const {
+ const auto storedRowCount = rowCount();
+
+ writeOne(static_cast<uint32>(storedRowCount), file);
+ writeOne(static_cast<uint32>(computeColCount()), file);
+ writeOne(static_cast<uint32>(modulus), file);
+ writeOne(static_cast<uint64>(entryCount()), file);
+
+ // write scalars
+ for (SparseMatrix::RowIndex row = 0; row < storedRowCount; ++row)
+ fwrite(&rowBegin(row).index(), sizeof(uint16), entryCountInRow(row), file);
+
+ // write indices
+ for (SparseMatrix::RowIndex row = 0; row < storedRowCount; ++row)
+ fwrite(&rowBegin(row).index(), sizeof(uint32), entryCountInRow(row), file);
+
+ std::vector<uint32> entryCounts;
+ for (SparseMatrix::RowIndex row = 0; row < storedRowCount; ++row)
+ entryCounts.push_back(entryCountInRow(row));
+ writeMany<uint32>(entryCounts, file);
+}
+
+SparseMatrix::Scalar SparseMatrix::read(FILE* file) {
+ MATHICGB_ASSERT(file != 0);
+
+ const auto rowCount = readOne<uint32>(file);
+ const auto colCount = readOne<uint32>(file);
+ const auto modulus = readOne<uint32>(file);
+ const auto entryCount = readOne<uint64>(file);
+
+ // Allocate memory to hold the matrix in one block.
+ clear();
+ reserveFreeEntries(entryCount);
+ mRows.reserve(rowCount);
+ MATHICGB_ASSERT(mBlock.mPreviousBlock == 0); // only one block
+
+ // @todo: we can read directly into the block. Do that.
+
+ // Read scalars.
+ {
+ mBlock.mScalars.resize(entryCount);
+ std::vector<uint16> scalars;
+ readMany(file, entryCount, scalars);
+ std::copy(scalars.begin(), scalars.end(), mBlock.mScalars.begin());
+ }
+
+ // Read column indices.
+ {
+ mBlock.mColIndices.resize(entryCount);
+ std::vector<uint32> indices;
+ readMany(file, entryCount, indices);
+ std::copy(indices.begin(), indices.end(), mBlock.mColIndices.begin());
+ }
+
+ // Read where rows begin and end.
+ {
+ std::vector<uint32> sizes;
+ readMany(file, rowCount, sizes);
+ uint32 runningOffset = 0;
+ for (auto it = sizes.begin(); it != sizes.end(); ++it) {
+ Row row;
+ row.mIndicesBegin = mBlock.mColIndices.begin() + runningOffset;
+ row.mScalarsBegin = mBlock.mScalars.begin() + runningOffset;
+ runningOffset += *it;
+ row.mIndicesEnd = mBlock.mColIndices.begin() + runningOffset;
+ row.mScalarsEnd = mBlock.mScalars.begin() + runningOffset;
+ mRows.push_back(row);
+ }
+ MATHICGB_ASSERT(runningOffset == entryCount);
+ }
+
+ MATHICGB_ASSERT(mBlock.mPreviousBlock == 0); // still only one block
+ return modulus;
+}
diff --git a/src/mathicgb/SparseMatrix.hpp b/src/mathicgb/SparseMatrix.hpp
index 22557fd..9980ba8 100755
--- a/src/mathicgb/SparseMatrix.hpp
+++ b/src/mathicgb/SparseMatrix.hpp
@@ -73,7 +73,9 @@ public:
SparseMatrix& operator=(const SparseMatrix&);
void swap(SparseMatrix& matrix);
- void clear(ColIndex newColCount = 0);
+
+ // Removes all rows from *this.
+ void clear();
/// Appends the rows from matrix to this object. Avoids most of the copies
/// that would otherwise be required for a big matrix insert by taking
@@ -220,6 +222,12 @@ public:
/// slow and it makes a copy internally.
void sortRowsByIncreasingPivots();
+ // Write *this and modulus to file.
+ void write(Scalar modulus, FILE* file) const;
+
+ // Set *this to a matrix read from file and return the modulus from the file.
+ Scalar read(FILE* file);
+
/// Iterates through the entries in a row.
class ConstRowIterator {
public:
@@ -269,7 +277,6 @@ public:
const Scalar* mScalarIt;
};
-
private:
MATHICGB_NO_INLINE void growEntryCapacity();
diff --git a/src/test/F4MatrixReducer.cpp b/src/test/F4MatrixReducer.cpp
index ec55c1c..527d911 100755
--- a/src/test/F4MatrixReducer.cpp
+++ b/src/test/F4MatrixReducer.cpp
@@ -29,7 +29,7 @@ TEST(F4MatrixReducer, Reduce) {
}
// top left
- m.topLeft.clear(4);
+ m.topLeft.clear();
m.topLeft.appendEntry(0, 1);
m.topLeft.appendEntry(1, 2);
m.topLeft.appendEntry(3, 3);
@@ -44,7 +44,7 @@ TEST(F4MatrixReducer, Reduce) {
m.topLeft.rowDone();
// top right
- m.topRight.clear(5);
+ m.topRight.clear();
m.topRight.appendEntry(2,8);
m.topRight.rowDone();
m.topRight.appendEntry(3,9);
@@ -54,7 +54,7 @@ TEST(F4MatrixReducer, Reduce) {
m.topRight.rowDone();
// bottom left
- m.bottomLeft.clear(4);
+ m.bottomLeft.clear();
m.bottomLeft.rowDone();
m.bottomLeft.appendEntry(1, 9);
m.bottomLeft.rowDone();
@@ -72,7 +72,7 @@ TEST(F4MatrixReducer, Reduce) {
m.bottomLeft.rowDone();
// bottom right
- m.bottomRight.clear(5);
+ m.bottomRight.clear();
m.bottomRight.rowDone();
m.bottomRight.appendEntry(1, 2);
m.bottomRight.appendEntry(3, 11);
diff --git a/src/test/SparseMatrix.cpp b/src/test/SparseMatrix.cpp
index ab3f299..d140553 100755
--- a/src/test/SparseMatrix.cpp
+++ b/src/test/SparseMatrix.cpp
@@ -66,7 +66,7 @@ TEST(SparseMatrix, toRow) {
monomials.push_back(it.getMonomial());
SparseMatrix mat(5);
- mat.clear(6);
+ mat.clear();
mat.rowDone();
mat.appendEntry(0,10);
mat.rowDone();
diff --git a/src/test/gb-test.cpp b/src/test/gb-test.cpp
index 474803e..3579d14 100755
--- a/src/test/gb-test.cpp
+++ b/src/test/gb-test.cpp
@@ -281,8 +281,10 @@ spairQueue reducerType divLookup monTable buchberger postponeKoszul useBaseDivis
tbb::task_scheduler_init scheduler(threadCount);
if (buchberger) {
+ const auto reducer = Reducer::makeReducer
+ (Reducer::reducerType(reducerType), I->ring());
BuchbergerAlg alg
- (*I, freeModuleOrder, Reducer::reducerType(reducerType),
+ (*I, freeModuleOrder, *reducer,
divLookup, preferSparseReducers, spairQueue);
alg.setUseAutoTopReduction(autoTopReduce);
alg.setUseAutoTailReduction(autoTailReduce);
@@ -296,7 +298,8 @@ spairQueue reducerType divLookup monTable buchberger postponeKoszul useBaseDivis
} else {
SignatureGB basis
(*I, freeModuleOrder, Reducer::reducerType(reducerType),
- divLookup, monTable, postponeKoszul, useBaseDivisors, preferSparseReducers, useSingularCriterionEarly, spairQueue);
+ divLookup, monTable, postponeKoszul, useBaseDivisors,
+ preferSparseReducers, useSingularCriterionEarly, spairQueue);
basis.computeGrobnerBasis();
EXPECT_EQ(sigBasisStr, toString(basis.getGB(), 1))
<< reducerType << ' ' << divLookup << ' '
--
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/debian-science/packages/mathicgb.git
More information about the debian-science-commits
mailing list