diff --git a/BaseLib/TemplateLogogFormatterSuppressedGCC.h b/BaseLib/TemplateLogogFormatterSuppressedGCC.h
new file mode 100644
index 0000000000000000000000000000000000000000..6fa00e00405d731a581f6e9e9e5041285624368e
--- /dev/null
+++ b/BaseLib/TemplateLogogFormatterSuppressedGCC.h
@@ -0,0 +1,40 @@
+/**
+ * Copyright (c) 2012, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/project/license
+ *
+ *
+ * \file TemplateLogogFormatterSuppressedGCC.h
+ *
+ * Created on 2012-12-06 by Norihiro Watanabe
+ */
+
+#ifndef TEMPLATELOGOGFORMATTERSUPPRESSEDGCC_H_
+#define TEMPLATELOGOGFORMATTERSUPPRESSEDGCC_H_
+
+// ** INCLUDES **
+#include "logog/include/logog.hpp"
+
+namespace BaseLib {
+
+/**
+ * \brief TemplateLogogFormatterSuppressedGCC strips topics given as a template
+ * parameter from logog::FormatterGCC.
+ * See http://johnwbyrd.github.com/logog/customformatting.html for details.
+ **/
+template <int T_SUPPPRESS_TOPIC_FLAG>
+class TemplateLogogFormatterSuppressedGCC : public logog::FormatterGCC
+{
+
+ virtual TOPIC_FLAGS GetTopicFlags( const logog::Topic &topic )
+ {
+ return ( logog::Formatter::GetTopicFlags( topic ) &
+ ~( T_SUPPPRESS_TOPIC_FLAG ));
+ }
+
+};
+
+#endif // TEMPLATELOGOGFORMATTERSUPPRESSEDGCC_H_
+
+} // namespace BaseLib
diff --git a/CMakeLists.txt b/CMakeLists.txt
index 384eb58798739aa1409a429cf7fd2b6f9c3b4f15..08dae53b6abe21e9aaa5341860a308b5c1292b83 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -55,6 +55,9 @@ OPTION(OGS_BUILD_UTILS "Should the utilities programms be built?" OFF)
OPTION(OGS_NO_EXTERNAL_LIBS "Builds OGS without any external dependencies." OFF)
+# Linear solvers
+OPTION(OGS_USE_LIS "Use Lis" OFF)
+
# Logging
OPTION(OGS_DISABLE_LOGGING "Disables all logog messages." OFF)
@@ -98,6 +101,10 @@ INCLUDE_DIRECTORIES( ${CMAKE_SOURCE_DIR}/ThirdParty )
INCLUDE_DIRECTORIES( ${CMAKE_SOURCE_DIR}/ThirdParty/quickcheck )
INCLUDE_DIRECTORIES( ${CMAKE_BINARY_DIR}/ThirdParty/zlib )
+IF(OGS_USE_LIS)
+ ADD_DEFINITIONS(-DUSE_LIS)
+ENDIF()
+
ADD_SUBDIRECTORY( BaseLib )
# TODO This is a hack but we have to make sure that Boost is built first
ADD_DEPENDENCIES(BaseLib Boost)
diff --git a/MathLib/CMakeLists.txt b/MathLib/CMakeLists.txt
index f6e6b904d90f458a98bd38fd575b55ff24e3103b..1e42e687a7934791f5246e507954ae6cf83539a2 100644
--- a/MathLib/CMakeLists.txt
+++ b/MathLib/CMakeLists.txt
@@ -20,6 +20,16 @@ SET ( SOURCES ${SOURCES} ${SOURCES_LINALG_SOLVERS})
GET_SOURCE_FILES(SOURCES_LINALG_PRECOND LinAlg/Preconditioner)
SET ( SOURCES ${SOURCES} ${SOURCES_LINALG_PRECOND})
+GET_SOURCE_FILES(SOURCES_LINALG_LEQS LinAlg/SystemOfLinearEquations)
+IF (NOT OGS_USE_LIS)
+ LIST(REMOVE_ITEM SOURCES_LINALG_LEQS
+ LinAlg/SystemOfLinearEquations/LisOption.h
+ LinAlg/SystemOfLinearEquations/LisOption.cpp
+ LinAlg/SystemOfLinearEquations/LisLinearSystem.h
+ LinAlg/SystemOfLinearEquations/LisLinearSystem.cpp
+ )
+ENDIF ()
+SET ( SOURCES ${SOURCES} ${SOURCES_LINALG_LEQS})
IF (METIS_FOUND)
GET_SOURCE_FILES(SOURCES_LINALG_SPARSE_NESTEDDISSECTION LinAlg/Sparse/NestedDissectionPermutation)
@@ -38,8 +48,17 @@ IF(METIS_FOUND)
)
ENDIF()
+IF (LIS_FOUND)
+ INCLUDE_DIRECTORIES(${LIS_INCLUDE_DIR})
+ENDIF()
+
+
# Create the library
ADD_LIBRARY( MathLib STATIC ${SOURCES} )
SET_TARGET_PROPERTIES(MathLib PROPERTIES LINKER_LANGUAGE CXX)
+IF (LIS_FOUND)
+ TARGET_LINK_LIBRARIES( MathLib ${LIS_LIBRARIES} )
+ENDIF()
+
diff --git a/MathLib/LinAlg/Sparse/Sparsity.h b/MathLib/LinAlg/Sparse/Sparsity.h
new file mode 100644
index 0000000000000000000000000000000000000000..63d5ca099938b4c4735034e7a9535063b3315a74
--- /dev/null
+++ b/MathLib/LinAlg/Sparse/Sparsity.h
@@ -0,0 +1,29 @@
+/**
+ * Copyright (c) 2012, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/project/license
+ *
+ *
+ * \file Sparsity.h
+ *
+ * Created on 2012-06-25 by Norihiro Watanabe
+ */
+
+#ifndef SPARSITY_H_
+#define SPARSITY_H_
+
+#include <vector>
+#include <set>
+
+namespace MathLib
+{
+
+/**
+ * \brief Row-major sparse pattern
+ */
+typedef std::vector<std::set<std::size_t> > RowMajorSparsity;
+
+} //MathLib
+
+#endif // SPARSITY_H_
diff --git a/MathLib/LinAlg/SystemOfLinearEquations/ISystemOfLinearEquations.h b/MathLib/LinAlg/SystemOfLinearEquations/ISystemOfLinearEquations.h
new file mode 100644
index 0000000000000000000000000000000000000000..750002822ad378d8db3446b6ab918013cf951b13
--- /dev/null
+++ b/MathLib/LinAlg/SystemOfLinearEquations/ISystemOfLinearEquations.h
@@ -0,0 +1,202 @@
+/**
+ * Copyright (c) 2012, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/project/license
+ *
+ *
+ * \file ISystemOfLinearEquations.h
+ *
+ * Created on 2012-06-25 by Norihiro Watanabe
+ */
+
+#ifndef ISYSTEMOFLINEAREQUATIONS_H_
+#define ISYSTEMOFLINEAREQUATIONS_H_
+
+#include <iostream>
+#include <vector>
+#include <boost/property_tree/ptree.hpp>
+#include "MathLib/LinAlg/Sparse/Sparsity.h"
+
+namespace MathLib
+{
+
+/**
+ * \brief Interface to any system of linear equations Ax=b
+ *
+ * This class defines unified interface to a system of linear equations
+ * including linear solvers.
+ *
+ */
+class ISystemOfLinearEquations
+{
+public:
+ /// Index representing invalid
+ const static std::size_t index_npos = -1;
+
+ ///
+ virtual ~ISystemOfLinearEquations()
+ {
+ }
+
+ /**
+ * set properties
+ *
+ * @param option
+ */
+ virtual void setOption(const boost::property_tree::ptree &option) = 0;
+
+ /**
+ * initialize this system with zero
+ *
+ * This function initialize the system by setting zero to all entries in
+ * a matrix, a RHS vector and solution vector. Dimension of the system is
+ * not changed.
+ */
+ virtual void setZero() = 0;
+
+ /// return dimension of this equation
+ virtual std::size_t getDimension() const = 0;
+
+ /**
+ * set an entry in a matrix
+ *
+ * @param rowId
+ * @param colId
+ * @param v
+ */
+ virtual void setMatEntry(std::size_t rowId, std::size_t colId, double v) = 0;
+
+ /**
+ * add a value into a matrix
+ *
+ * @param rowId
+ * @param colId
+ * @param v
+ */
+ virtual void addMatEntry(std::size_t rowId, std::size_t colId, double v) = 0;
+
+ /**
+ * add a sub-matrix into a matrix at given row and column positions
+ *
+ * @param vec_row_pos
+ * A vector of global row index. The number of entries of vec_row_pos have
+ * to be identical to the number of rows of the local matrix.
+ * @param vec_col_pos
+ * A vector of global column index. The number of vec_col_pos have to be
+ * identical to the number of columns of the local matrix.
+ * @param sub_matrix A sub-matrix
+ * @param fac scaling parameter
+ */
+ template <class T_DENSE_MATRIX>
+ void addSubMat(const std::vector<std::size_t> &vec_row_pos,
+ const std::vector<std::size_t> &vec_col_pos,
+ const T_DENSE_MATRIX &sub_matrix, double fac = 1.0);
+
+ /**
+ * add a sub-matrix into a matrix
+ *
+ * @param vec_row_col_pos
+ * A vector of global row and column index. The number of entries of
+ * vec_row_col_pos have to be identical to the number of rows and columns
+ * of the local matrix.
+ * @param sub_matrix A sub-matrix
+ * @param fac scaling parameter
+ */
+ template <class T_DENSE_MATRIX>
+ void addSubMat(const std::vector<std::size_t> &vec_row_col_pos,
+ const T_DENSE_MATRIX &sub_matrix, double fac = 1.0);
+
+ /**
+ * get RHS entry
+ *
+ * @param rowId
+ * @return
+ */
+ virtual double getRHSVec(std::size_t rowId) const = 0;
+
+ /**
+ * set RHS entry
+ *
+ * @param rowId
+ * @param v
+ */
+ virtual void setRHSVec(std::size_t rowId, double v) = 0;
+
+ /**
+ * add RHS entry
+ *
+ * @param rowId
+ * @param v
+ */
+ virtual void addRHSVec(std::size_t rowId, double v) = 0;
+
+ /**
+ * add a sub vector to RHS
+ *
+ * @param vec_row_pos
+ * A vector of global row index. The number of entries of vec_row_pos
+ * have to be identical to the number of rows of the sub vector.
+ * @param sub_vector
+ * Pointer to a sub-vector. Its length must be the same as the length of vec_row_pos.
+ * @param fac
+ * Scaling parameter. Default value is 1.
+ */
+ inline virtual void addSubRHS(const std::vector<std::size_t> &vec_row_pos,
+ const double* sub_vector, double fac = 1.0);
+
+ /**
+ * add a sub vector to RHS
+ *
+ * @param vec_row_pos
+ * A vector of global row index. The number of entries of vec_row_pos
+ * have to be identical to the number of rows of the sub vector.
+ * @param sub_vector A sub-vector
+ * @param fac Scaling parameter
+ */
+ template <class T_DENSE_VECTOR>
+ void addSubRHS(const std::vector<std::size_t> &vec_row_pos,
+ const T_DENSE_VECTOR &sub_vector, double fac = 1.0);
+
+ /**
+ * get an entry in a solution vector
+ * @param rowId
+ * @return
+ */
+ virtual double getSolVec(std::size_t rowId) = 0;
+
+ /**
+ * set an entry in a solution vector
+ * @param rowId
+ * @param v
+ */
+ virtual void setSolVec(std::size_t rowId, double v) = 0;
+
+ /**
+ * set prescribed value
+ * @param row_id
+ * @param x
+ */
+ virtual void setKnownSolution(std::size_t row_id, double x) = 0;
+
+ /**
+ * set prescribed value
+ * @param vec_id A vector of global row index
+ * @param vec_x A vector of prescribed values
+ */
+ virtual void setKnownSolution(const std::vector<std::size_t> &vec_id,
+ const std::vector<double> &vec_x) = 0;
+
+ /// solve the linear system
+ virtual void solve() = 0;
+
+ /// print out the equation for debugging
+ virtual void printout(std::ostream &os = std::cout) const = 0;
+
+};
+
+}
+
+#include "ISystemOfLinearEquations.tpp"
+
+#endif //ISYSTEMOFLINEAREQUATIONS_H_
diff --git a/MathLib/LinAlg/SystemOfLinearEquations/ISystemOfLinearEquations.tpp b/MathLib/LinAlg/SystemOfLinearEquations/ISystemOfLinearEquations.tpp
new file mode 100644
index 0000000000000000000000000000000000000000..627ed4b4e45b13810897f88eb6bd69889aeb6a6e
--- /dev/null
+++ b/MathLib/LinAlg/SystemOfLinearEquations/ISystemOfLinearEquations.tpp
@@ -0,0 +1,63 @@
+/**
+ * Copyright (c) 2012, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/project/license
+ *
+ *
+ * \file ISystemOfLinearEquations.tpp
+ *
+ * Created on 2012-06-25 by Norihiro Watanabe
+ */
+
+#ifndef ISYSTEMOFLINEAREQUATIONS_TPP_
+#define ISYSTEMOFLINEAREQUATIONS_TPP_
+
+namespace MathLib
+{
+
+template <class T_DENSE_MATRIX>
+void ISystemOfLinearEquations::addSubMat(const std::vector<std::size_t> &vec_row_pos, const std::vector<std::size_t> &vec_col_pos, const T_DENSE_MATRIX &sub_matrix, double fkt)
+{
+ const std::size_t n_rows = vec_row_pos.size();
+ const std::size_t n_cols = vec_col_pos.size();
+ for (std::size_t i=0; i<n_rows; i++) {
+ const std::size_t rowId = vec_row_pos[i];
+ if (rowId==index_npos) continue;
+ for (std::size_t j=0; j<n_cols; j++) {
+ const std::size_t colId = vec_col_pos[j];
+ if (colId==index_npos) continue;
+ addMatEntry(rowId, colId, fkt*sub_matrix(i,j));
+ }
+ }
+}
+
+template <class T_DENSE_MATRIX>
+void ISystemOfLinearEquations::addSubMat(const std::vector<std::size_t> &vec_pos, const T_DENSE_MATRIX &sub_matrix, double fkt)
+{
+ addSubMat(vec_pos, vec_pos, sub_matrix, fkt);
+}
+
+inline void ISystemOfLinearEquations::addSubRHS(const std::vector<std::size_t> &vec_row_pos, const double *sub_vector, double fkt)
+{
+ for (std::size_t i=0; i<vec_row_pos.size(); i++) {
+ const std::size_t rowId = vec_row_pos[i];
+ if (rowId==index_npos) continue;
+ addRHSVec(rowId, sub_vector[i]*fkt);
+ }
+}
+
+template <class T_DENSE_VECTOR>
+void ISystemOfLinearEquations::addSubRHS(const std::vector<std::size_t> &vec_row_pos, const T_DENSE_VECTOR &sub_vector, double fkt)
+{
+ for (std::size_t i=0; i<vec_row_pos.size(); i++) {
+ const std::size_t rowId = vec_row_pos[i];
+ if (rowId==index_npos) continue;
+ addRHSVec(rowId, sub_vector(i)*fkt);
+ }
+}
+
+}
+
+#endif //ISYSTEMOFLINEAREQUATIONS_TPP_
+
diff --git a/MathLib/LinAlg/SystemOfLinearEquations/LisLinearSystem.cpp b/MathLib/LinAlg/SystemOfLinearEquations/LisLinearSystem.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..3c24746b71f07a2e618c7decd74539f3dfc0e583
--- /dev/null
+++ b/MathLib/LinAlg/SystemOfLinearEquations/LisLinearSystem.cpp
@@ -0,0 +1,189 @@
+/**
+ * Copyright (c) 2012, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/project/license
+ *
+ *
+ * \file LisLinearSystem.cpp
+ *
+ * Created on 2012-06-25 by Norihiro Watanabe
+ */
+
+#include "LisLinearSystem.h"
+
+#include <string>
+#include <sstream>
+#ifdef _OPENMP
+#include <omp.h>
+#endif
+#include "logog/include/logog.hpp"
+
+namespace MathLib
+{
+
+bool LisLinearSystem::checkError(int err)
+{
+ bool ok = (err == LIS_SUCCESS);
+ if (!ok) {
+ ERR("***ERROR: Lis error code = %d", err);
+ }
+ return ok;
+}
+
+LisLinearSystem::LisLinearSystem(std::size_t dimension)
+: _dim(dimension), _max_diag_coeff(.0)
+{
+ int ierr = 0;
+ ierr = lis_matrix_create(0, &_AA); checkError(ierr);
+ ierr = lis_matrix_set_size(_AA, 0, dimension); checkError(ierr);
+ ierr = lis_vector_duplicate(_AA, &_bb); checkError(ierr);
+ ierr = lis_vector_duplicate(_AA, &_xx); checkError(ierr);
+}
+
+LisLinearSystem::~LisLinearSystem()
+{
+ int ierr = 0;
+ ierr = lis_matrix_destroy(_AA); checkError(ierr);
+ ierr = lis_vector_destroy(_bb); checkError(ierr);
+ ierr = lis_vector_destroy(_xx); checkError(ierr);
+}
+
+void LisLinearSystem::setOption(const boost::property_tree::ptree &option)
+{
+ using boost::property_tree::ptree;
+
+ boost::optional<ptree> ptSolver = option.get_child("LinearSolver");
+ if (!ptSolver)
+ return;
+
+ boost::optional<std::string> solver_type = ptSolver->get_optional<std::string>("solver_type");
+ if (solver_type) {
+ _option.solver_type = _option.getSolverType(*solver_type);
+ }
+ boost::optional<std::string> precon_type = ptSolver->get_optional<std::string>("precon_type");
+ if (precon_type) {
+ _option.precon_type = _option.getPreconType(*precon_type);
+ }
+ boost::optional<std::string> matrix_type = ptSolver->get_optional<std::string>("matrix_type");
+ if (matrix_type) {
+ _option.matrix_type = _option.getMatrixType(*matrix_type);
+ }
+ boost::optional<double> error_tolerance = ptSolver->get_optional<double>("error_tolerance");
+ if (error_tolerance) {
+ _option.error_tolerance = *error_tolerance;
+ }
+ boost::optional<int> max_iteration_step = ptSolver->get_optional<int>("max_iteration_step");
+ if (max_iteration_step) {
+ _option.max_iterations = *max_iteration_step;
+ }
+}
+
+void LisLinearSystem::setZero()
+{
+ int ierr = 0;
+ // A matrix has to be removed and created every time because Lis doesn't provide a
+ // function to set matrix entries to zero
+ ierr = lis_matrix_destroy(_AA); checkError(ierr);
+ ierr = lis_matrix_create(0, &_AA); checkError(ierr);
+ ierr = lis_matrix_set_size(_AA, 0, getDimension()); checkError(ierr);
+ // set zero RHS, x
+ ierr = lis_vector_set_all(0.0, _bb); checkError(ierr);
+ ierr = lis_vector_set_all(0.0, _xx); checkError(ierr);
+
+ _max_diag_coeff = .0;
+}
+
+void LisLinearSystem::applyKnownSolution()
+{
+ //Use penalty parameter
+ const double penalty_scaling = 1e+10;
+ const double penalty = _max_diag_coeff * penalty_scaling;
+ INFO("-> max. absolute value of diagonal entries = %e", _max_diag_coeff);
+ INFO("-> penalty scaling = %e", penalty_scaling);
+ const std::size_t n_bc = _vec_knownX_id.size();
+ for (std::size_t i_bc=0; i_bc<n_bc; i_bc++) {
+ const std::size_t rowId = _vec_knownX_id[i_bc];
+ const double x = _vec_knownX_x[i_bc];
+ //A(k, k) = penalty
+ setMatEntry(rowId, rowId, penalty);
+ //b(k) = x*penalty
+ setRHSVec(rowId, x*penalty);
+ }
+}
+
+void LisLinearSystem::solve()
+{
+ INFO("------------------------------------------------------------------");
+ INFO("*** LIS solver computation");
+#ifdef _OPENMP
+ INFO("-> max number of threads = %d", omp_get_num_procs());
+ INFO("-> number of threads = %d", omp_get_max_threads());
+#endif
+
+ applyKnownSolution();
+ // assemble a matrix
+ int ierr = 0;
+ ierr = lis_matrix_set_type(_AA, static_cast<int>(_option.matrix_type)); checkError(ierr);
+ ierr = lis_matrix_assemble(_AA); checkError(ierr);
+
+ // configure option
+ std::string solver_options;
+ if (_option.solver_precon_arg.empty()) {
+ std::stringstream ss;
+ ss << "-i " << static_cast<int>(_option.solver_type);
+ ss << " -p " << static_cast<int>(_option.precon_type);
+ if (!_option.extra_arg.empty())
+ ss << " " << _option.extra_arg;
+ solver_options = ss.str();
+ } else {
+ solver_options = _option.solver_precon_arg;
+ }
+ std::string tol_option;
+ {
+ std::stringstream ss;
+ ss << "-tol " << _option.error_tolerance;
+ ss << " -maxiter " << _option.max_iterations;
+ ss << " -initx_zeros 0"; //0: use given x as initial guess, 1: x0=0
+#ifdef _OPENMP
+ const int nthreads = omp_get_max_threads();
+ ss << " -omp_num_threads " << nthreads;
+#endif
+ tol_option = ss.str();
+ }
+
+ // Create solver
+ LIS_SOLVER solver;
+ ierr = lis_solver_create(&solver); checkError(ierr);
+ ierr = lis_solver_set_option(const_cast<char*>(solver_options.c_str()), solver); checkError(ierr);
+ ierr = lis_solver_set_option(const_cast<char*>(tol_option.c_str()), solver); checkError(ierr);
+ ierr = lis_solver_set_option(const_cast<char*>("-print mem"), solver); checkError(ierr);
+
+ // solve
+ INFO("-> solve");
+ ierr = lis_solve(_AA, _bb, _xx, solver); checkError(ierr);
+ checkError(ierr);
+ //lis_output(_AA, _bb, _xx, LIS_FMT_MM, "/home/localadmin/tasks/20120814_ogs6test/matrix1.txt");
+
+ int iter = 0;
+ double resid = 0.0;
+ ierr = lis_solver_get_iters(solver, &iter); checkError(ierr);
+ ierr = lis_solver_get_residualnorm(solver, &resid); checkError(ierr);
+ INFO("\t iteration: %d/%ld\n", iter, _option.max_iterations);
+ INFO("\t residuals: %e\n", resid);
+
+ // Clear solver
+ ierr = lis_solver_destroy(solver); checkError(ierr);
+ INFO("------------------------------------------------------------------");
+}
+
+void LisLinearSystem::printout(std::ostream &os) const
+{
+ os << "#A=" << std::endl;
+ os << "#x=" << std::endl;
+ lis_vector_print(_xx);
+ os << "#b=" << std::endl;
+ lis_vector_print(_bb);
+}
+
+} //MathLib
diff --git a/MathLib/LinAlg/SystemOfLinearEquations/LisLinearSystem.h b/MathLib/LinAlg/SystemOfLinearEquations/LisLinearSystem.h
new file mode 100644
index 0000000000000000000000000000000000000000..8e038b3e43ba82010e397b3a413383676e6d95c3
--- /dev/null
+++ b/MathLib/LinAlg/SystemOfLinearEquations/LisLinearSystem.h
@@ -0,0 +1,183 @@
+/**
+ * Copyright (c) 2012, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/project/license
+ *
+ *
+ * \file LisLinearSystem.h
+ *
+ * Created on 2012-06-25 by Norihiro Watanabe
+ */
+
+#ifndef LISLINEARSYSTEM_H_
+#define LISLINEARSYSTEM_H_
+
+#include <iostream>
+#include <string>
+#include <vector>
+#include <boost/property_tree/ptree.hpp>
+
+#include "lis.h"
+
+#include "LisOption.h"
+
+namespace MathLib
+{
+
+/**
+ * \brief Linear system using Lis solver (http://www.ssisc.org/lis/)
+ *
+ * This class utilizes Lis solver.
+ */
+class LisLinearSystem
+{
+public:
+ //---------------------------------------------------------------
+ // realization of ISystemOfLinearEquations
+ //---------------------------------------------------------------
+
+ /**
+ * Create a linear system using Lis
+ *
+ * @param length System dimension
+ */
+ LisLinearSystem(std::size_t length);
+
+ /**
+ *
+ */
+ virtual ~LisLinearSystem();
+
+ /**
+ * configure linear solvers
+ * @param option
+ */
+ virtual void setOption(const boost::property_tree::ptree &option);
+
+ /// return the system dimension
+ virtual std::size_t getDimension() const { return _dim; };
+
+ /// reset this equation
+ virtual void setZero();
+
+ /// set entry in A
+ virtual void setMatEntry(std::size_t rowId, std::size_t colId, double v)
+ {
+ if (rowId==colId)
+ _max_diag_coeff = std::max(_max_diag_coeff, std::abs(v));
+ lis_matrix_set_value(LIS_INS_VALUE, rowId, colId, v, _AA);
+ }
+
+ /// add value into A
+ virtual void addMatEntry(std::size_t rowId, std::size_t colId, double v)
+ {
+ if (rowId==colId)
+ _max_diag_coeff = std::max(_max_diag_coeff, std::abs(v));
+ lis_matrix_set_value(LIS_ADD_VALUE, rowId, colId, v, _AA);
+ }
+
+ /// get RHS entry
+ virtual double getRHSVec(std::size_t rowId) const
+ {
+ double v;
+ lis_vector_get_value(_bb, rowId, &v);
+ return v;
+ }
+
+ /// set RHS entry
+ virtual void setRHSVec(std::size_t rowId, double v)
+ {
+ lis_vector_set_value(LIS_INS_VALUE, rowId, v, _bb);
+ }
+
+ /// add RHS entry
+ virtual void addRHSVec(std::size_t rowId, double v)
+ {
+ lis_vector_set_value(LIS_ADD_VALUE, rowId, v, _bb);
+ }
+
+ /// get an entry in a solution vector
+ virtual double getSolVec(std::size_t rowId)
+ {
+ double v;
+ lis_vector_get_value(_xx, rowId, &v);
+ return v;
+ }
+
+ /// set a solution vector
+ virtual void setSolVec(std::size_t rowId, double v)
+ {
+ lis_vector_set_value(LIS_INS_VALUE, rowId, v, _xx);
+ }
+
+ /// set prescribed value
+ virtual void setKnownSolution(std::size_t rowId, double x)
+ {
+ _vec_knownX_id.push_back(rowId);
+ _vec_knownX_x.push_back(x);
+ }
+
+ /// set prescribed values
+ virtual void setKnownSolution(const std::vector<std::size_t> &vec_id, const std::vector<double> &vec_x)
+ {
+ _vec_knownX_id.insert(_vec_knownX_id.end(), vec_id.begin(), vec_id.end());
+ _vec_knownX_x.insert(_vec_knownX_x.end(), vec_x.begin(), vec_x.end());
+ }
+
+ /// solve this equation
+ virtual void solve();
+
+ /// printout this equation for debugging
+ virtual void printout(std::ostream &os=std::cout) const;
+
+ //---------------------------------------------------------------
+ // specific to this class
+ //---------------------------------------------------------------
+ /**
+ * configure linear solvers
+ * @param option
+ */
+ void setOption(const LisOption &option)
+ {
+ _option = option;
+ }
+
+ /**
+ * get linear solver options
+ * @return
+ */
+ LisOption &getOption()
+ {
+ return _option;
+ }
+
+ /**
+ * get a solution vector
+ * @param x
+ */
+ void getSolVec(double* x)
+ {
+ lis_vector_gather(_xx, x);
+ }
+
+private:
+ bool checkError(int err);
+ void applyKnownSolution();
+
+private:
+ const std::size_t _dim;
+ double _max_diag_coeff;
+ LisOption _option;
+ LIS_MATRIX _AA;
+ LIS_VECTOR _bb;
+ LIS_VECTOR _xx;
+ std::vector<std::size_t> _vec_knownX_id;
+ std::vector<double> _vec_knownX_x;
+};
+
+
+} // MathLib
+
+#endif //LISLINEARSYSTEM_H_
+
diff --git a/MathLib/LinAlg/SystemOfLinearEquations/LisOption.cpp b/MathLib/LinAlg/SystemOfLinearEquations/LisOption.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..8ee5d504b0932d7ef41b2be903afc0a6f692dd8f
--- /dev/null
+++ b/MathLib/LinAlg/SystemOfLinearEquations/LisOption.cpp
@@ -0,0 +1,100 @@
+/**
+ * Copyright (c) 2012, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/project/license
+ *
+ *
+ * \file LisOption.cpp
+ *
+ * Created on 2012-12-05 by Norihiro Watanabe
+ */
+
+#include "LisOption.h"
+
+namespace MathLib
+{
+
+LisOption::LisOption()
+{
+ solver_type = SolverType::CG;
+ precon_type = PreconType::NONE;
+ matrix_type = MatrixType::CRS;
+ max_iterations = 500;
+ error_tolerance = 1.e-10;
+}
+
+LisOption::SolverType LisOption::getSolverType(const std::string &solver_name)
+{
+#define RETURN_SOLVER_ENUM_IF_SAME_STRING(str, TypeName) \
+ if (#TypeName==str) return SolverType::TypeName;
+
+ RETURN_SOLVER_ENUM_IF_SAME_STRING(solver_name, CG);
+ RETURN_SOLVER_ENUM_IF_SAME_STRING(solver_name, BiCG);
+ RETURN_SOLVER_ENUM_IF_SAME_STRING(solver_name, CGS);
+ RETURN_SOLVER_ENUM_IF_SAME_STRING(solver_name, BiCGSTAB);
+ RETURN_SOLVER_ENUM_IF_SAME_STRING(solver_name, BiCGSTABl);
+ RETURN_SOLVER_ENUM_IF_SAME_STRING(solver_name, GPBiCG);
+ RETURN_SOLVER_ENUM_IF_SAME_STRING(solver_name, TFQMR);
+ RETURN_SOLVER_ENUM_IF_SAME_STRING(solver_name, Orthomin);
+ RETURN_SOLVER_ENUM_IF_SAME_STRING(solver_name, GMRES);
+ RETURN_SOLVER_ENUM_IF_SAME_STRING(solver_name, Jacobi);
+ RETURN_SOLVER_ENUM_IF_SAME_STRING(solver_name, GaussSeidel);
+ RETURN_SOLVER_ENUM_IF_SAME_STRING(solver_name, SOR);
+ RETURN_SOLVER_ENUM_IF_SAME_STRING(solver_name, BiCGSafe);
+ RETURN_SOLVER_ENUM_IF_SAME_STRING(solver_name, CR);
+ RETURN_SOLVER_ENUM_IF_SAME_STRING(solver_name, BiCR);
+ RETURN_SOLVER_ENUM_IF_SAME_STRING(solver_name, CRS);
+ RETURN_SOLVER_ENUM_IF_SAME_STRING(solver_name, BiCRSTAB);
+ RETURN_SOLVER_ENUM_IF_SAME_STRING(solver_name, GPBiCR);
+ RETURN_SOLVER_ENUM_IF_SAME_STRING(solver_name, BiCRSafe);
+ RETURN_SOLVER_ENUM_IF_SAME_STRING(solver_name, FGMRESm);
+ RETURN_SOLVER_ENUM_IF_SAME_STRING(solver_name, IDRs);
+ RETURN_SOLVER_ENUM_IF_SAME_STRING(solver_name, MINRES);
+
+ return SolverType::INVALID;
+#undef RETURN_SOLVER_ENUM_IF_SAME_STRING
+}
+
+LisOption::PreconType LisOption::getPreconType(const std::string &precon_name)
+{
+#define RETURN_PRECOM_ENUM_IF_SAME_STRING(str, TypeName) \
+ if (#TypeName==str) return PreconType::TypeName;
+
+ RETURN_PRECOM_ENUM_IF_SAME_STRING(precon_name, NONE);
+ RETURN_PRECOM_ENUM_IF_SAME_STRING(precon_name, JACOBI);
+ RETURN_PRECOM_ENUM_IF_SAME_STRING(precon_name, ILU);
+ RETURN_PRECOM_ENUM_IF_SAME_STRING(precon_name, SSOR);
+ RETURN_PRECOM_ENUM_IF_SAME_STRING(precon_name, Hybrid);
+ RETURN_PRECOM_ENUM_IF_SAME_STRING(precon_name, IplusS);
+ RETURN_PRECOM_ENUM_IF_SAME_STRING(precon_name, SAINV);
+ RETURN_PRECOM_ENUM_IF_SAME_STRING(precon_name, SAAMG);
+ RETURN_PRECOM_ENUM_IF_SAME_STRING(precon_name, CroutILU);
+ RETURN_PRECOM_ENUM_IF_SAME_STRING(precon_name, ILUT);
+
+ return PreconType::NONE;
+#undef RETURN_PRECOM_ENUM_IF_SAME_STRING
+}
+
+LisOption::MatrixType LisOption::getMatrixType(const std::string &matrix_name)
+{
+#define RETURN_MATRIX_ENUM_IF_SAME_STRING(str, TypeName) \
+ if (#TypeName==str) return MatrixType::TypeName;
+
+ RETURN_MATRIX_ENUM_IF_SAME_STRING(matrix_name, CRS);
+ RETURN_MATRIX_ENUM_IF_SAME_STRING(matrix_name, CCS);
+ RETURN_MATRIX_ENUM_IF_SAME_STRING(matrix_name, MSR);
+ RETURN_MATRIX_ENUM_IF_SAME_STRING(matrix_name, DIA);
+ RETURN_MATRIX_ENUM_IF_SAME_STRING(matrix_name, ELL);
+ RETURN_MATRIX_ENUM_IF_SAME_STRING(matrix_name, JDS);
+ RETURN_MATRIX_ENUM_IF_SAME_STRING(matrix_name, BSR);
+ RETURN_MATRIX_ENUM_IF_SAME_STRING(matrix_name, BSC);
+ RETURN_MATRIX_ENUM_IF_SAME_STRING(matrix_name, VBR);
+ RETURN_MATRIX_ENUM_IF_SAME_STRING(matrix_name, COO);
+ RETURN_MATRIX_ENUM_IF_SAME_STRING(matrix_name, DNS);
+
+ return MatrixType::CRS;
+#undef RETURN_MATRIX_ENUM_IF_SAME_STRING
+}
+
+} //MathLib
diff --git a/MathLib/LinAlg/SystemOfLinearEquations/LisOption.h b/MathLib/LinAlg/SystemOfLinearEquations/LisOption.h
new file mode 100644
index 0000000000000000000000000000000000000000..704b38f80ee8b6d4087ed341cd93e86b493c9793
--- /dev/null
+++ b/MathLib/LinAlg/SystemOfLinearEquations/LisOption.h
@@ -0,0 +1,143 @@
+/**
+ * Copyright (c) 2012, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/project/license
+ *
+ *
+ * \file LisOption.h
+ *
+ * Created on 2012-06-25 by Norihiro Watanabe
+ */
+
+#ifndef LIS_OPTION_H_
+#define LIS_OPTION_H_
+
+#include <string>
+
+namespace MathLib
+{
+
+/**
+ * \brief Option for Lis solver
+ */
+struct LisOption
+{
+ /// Solver type
+ enum class SolverType : int
+ {
+ INVALID = 0,
+ CG = 1,
+ BiCG = 2,
+ CGS = 3,
+ BiCGSTAB = 4,
+ BiCGSTABl = 5,
+ GPBiCG = 6,
+ TFQMR = 7,
+ Orthomin = 8,
+ GMRES = 9,
+ Jacobi = 10,
+ GaussSeidel = 11,
+ SOR = 12,
+ BiCGSafe = 13,
+ CR = 14,
+ BiCR = 15,
+ CRS = 16,
+ BiCRSTAB = 17,
+ GPBiCR = 18,
+ BiCRSafe = 19,
+ FGMRESm = 20,
+ IDRs = 21,
+ MINRES = 22
+ };
+
+ /// Preconditioner type
+ enum class PreconType : int
+ {
+ NONE = 0,
+ JACOBI = 1,
+ ILU = 2,
+ SSOR = 3,
+ Hybrid = 4,
+ IplusS = 5,
+ SAINV = 6,
+ SAAMG = 7,
+ CroutILU = 8,
+ ILUT = 9
+ };
+
+ /// Matrix type
+ enum class MatrixType : int
+ {
+ CRS = 1,
+ CCS = 2,
+ MSR = 3,
+ DIA = 4,
+ ELL = 5,
+ JDS = 6,
+ BSR = 7,
+ BSC = 8,
+ VBR = 9,
+ COO = 10,
+ DNS = 11
+ };
+
+ /// Linear solver type
+ SolverType solver_type;
+ /// Preconditioner type
+ PreconType precon_type;
+ /// Matrix type
+ MatrixType matrix_type;
+ /// Maximum iteration count
+ long max_iterations;
+ /// Error tolerance
+ double error_tolerance;
+ /// Extra option
+ std::string extra_arg;
+ /// Arguments for solver and preconditioner. This variable is always preferred
+ /// to other variables.
+ std::string solver_precon_arg;
+
+ /**
+ * Constructor
+ *
+ * Default options are CG, no preconditioner, iteration count 500 and
+ * tolerance 1e-10. Default matrix storage type is CRS.
+ */
+ LisOption();
+
+ /// Destructor
+ ~LisOption() {};
+
+ /**
+ * return a linear solver type from the solver name
+ *
+ * @param solver_name
+ * @return a linear solver type
+ * If there is no solver type matched with the given name, INVALID
+ * is returned.
+ */
+ static SolverType getSolverType(const std::string &solver_name);
+
+ /**
+ * return a preconditioner type from the name
+ *
+ * @param precon_name
+ * @return a preconditioner type
+ * If there is no preconditioner type matched with the given name, NONE
+ * is returned.
+ */
+ static PreconType getPreconType(const std::string &precon_name);
+
+ /**
+ * return a matrix type
+ * @param matrix_name
+ * @return
+ * If there is no matrix type matched with the given name, CRS
+ * is returned.
+ */
+ static MatrixType getMatrixType(const std::string &matrix_name);
+};
+
+}
+#endif //LIS_OPTION_H_
diff --git a/Tests/CMakeLists.txt b/Tests/CMakeLists.txt
index be2627dea1dfa4bc84baf0cd6f7b19a23cbff050..faee29a1ed283ab9c0b2ba96b1b92e4af193bd76 100644
--- a/Tests/CMakeLists.txt
+++ b/Tests/CMakeLists.txt
@@ -7,13 +7,19 @@ ENDIF()
APPEND_SOURCE_FILES(TEST_SOURCES)
APPEND_SOURCE_FILES(TEST_SOURCES BaseLib)
APPEND_SOURCE_FILES(TEST_SOURCES GeoLib)
+APPEND_SOURCE_FILES(TEST_SOURCES MathLib)
INCLUDE_DIRECTORIES(
${CMAKE_SOURCE_DIR}/Tests/gtest
${CMAKE_SOURCE_DIR}/BaseLib
${CMAKE_SOURCE_DIR}/GeoLib
+ ${CMAKE_SOURCE_DIR}/MathLib
)
+IF (LIS_FOUND)
+ INCLUDE_DIRECTORIES(${LIS_INCLUDE_DIR})
+ENDIF()
+
ADD_EXECUTABLE (testrunner testrunner.cpp ${TEST_SOURCES})
SET_TARGET_PROPERTIES(testrunner PROPERTIES FOLDER Testing)
@@ -21,6 +27,8 @@ TARGET_LINK_LIBRARIES(testrunner
GTest
BaseLib
GeoLib
+ MathLib
+ logog
${Boost_LIBRARIES}
${CMAKE_THREAD_LIBS_INIT}
)
diff --git a/Tests/MathLib/TestLinearSystem.cpp b/Tests/MathLib/TestLinearSystem.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..dea1c61aa005a29bc10346f5933b80a0ae17ecc0
--- /dev/null
+++ b/Tests/MathLib/TestLinearSystem.cpp
@@ -0,0 +1,127 @@
+/**
+ * Copyright (c) 2012, OpenGeoSys Community (http://www.opengeosys.com)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.com/LICENSE.txt
+ *
+ *
+ * \file TestLinearSystem.cpp
+ *
+ * Created on 2012-08-03 by Norihiro Watanabe
+ */
+
+#include <gtest/gtest.h>
+#include <boost/property_tree/ptree.hpp>
+
+#include "MathLib/LinAlg/Sparse/Sparsity.h"
+#ifdef USE_LIS
+#include "MathLib/LinAlg/SystemOfLinearEquations/LisLinearSystem.h"
+#endif
+
+namespace
+{
+
+inline void ASSERT_DOUBLE_ARRAY_EQ(const double* Expected, const double* Actual, size_t N, double epsilon=1.0e-8) {
+ for (size_t i=0; i<N; i++) \
+ ASSERT_NEAR(Expected[i], Actual[i], epsilon);
+}
+
+struct Example1
+{
+ std::vector<double> mat;
+ std::vector<size_t> list_dirichlet_bc_id;
+ std::vector<double> list_dirichlet_bc_value;
+ static const size_t dim_eqs = 9;
+ std::vector<double> exH;
+ MathLib::RowMajorSparsity sparse;
+
+ Example1()
+ {
+ double d_mat[] = {
+ 6.66667e-012, -1.66667e-012, 0, -1.66667e-012, -3.33333e-012, 0, 0, 0, 0,
+ -1.66667e-012, 1.33333e-011, -1.66667e-012, -3.33333e-012, -3.33333e-012, -3.33333e-012, 0, 0, 0,
+ 0, -1.66667e-012, 6.66667e-012, 0, -3.33333e-012, -1.66667e-012, 0, 0, 0,
+ -1.66667e-012, -3.33333e-012, 0, 1.33333e-011, -3.33333e-012, 0, -1.66667e-012, -3.33333e-012, 0,
+ -3.33333e-012, -3.33333e-012, -3.33333e-012, -3.33333e-012, 2.66667e-011, -3.33333e-012, -3.33333e-012, -3.33333e-012, -3.33333e-012,
+ 0, -3.33333e-012, -1.66667e-012, 0, -3.33333e-012, 1.33333e-011, 0, -3.33333e-012, -1.66667e-012,
+ 0, 0, 0, -1.66667e-012, -3.33333e-012, 0, 6.66667e-012, -1.66667e-012, 0,
+ 0, 0, 0, -3.33333e-012, -3.33333e-012, -3.33333e-012, -1.66667e-012, 1.33333e-011, -1.66667e-012,
+ 0, 0, 0, 0, -3.33333e-012, -1.66667e-012, 0, -1.66667e-012, 6.66667e-012
+ };
+ mat.assign(d_mat, d_mat+dim_eqs*dim_eqs);
+ size_t int_dirichlet_bc_id[] = {2,5,8,0,3,6};
+ list_dirichlet_bc_id.assign(int_dirichlet_bc_id, int_dirichlet_bc_id+6);
+ list_dirichlet_bc_value.resize(6);
+ fill(list_dirichlet_bc_value.begin(), list_dirichlet_bc_value.begin()+3, .0);
+ fill(list_dirichlet_bc_value.begin()+3, list_dirichlet_bc_value.end(), 1.0);
+ exH.resize(9);
+ for (size_t i=0; i<9; i++) {
+ if (i%3==0) exH[i] = 1.0;
+ if (i%3==1) exH[i] = 0.5;
+ if (i%3==2) exH[i] = 0.;
+ }
+ sparse.resize(dim_eqs);
+ for (size_t i=0; i<dim_eqs; i++) {
+ for (size_t j=0; j<dim_eqs; j++) {
+ if (mat[i*dim_eqs+j]!=.0)
+ sparse[i].insert(j);
+ }
+ }
+
+ }
+};
+}
+
+#ifdef USE_LIS
+TEST(Math, LinearSystemLis)
+{
+ // set a problem
+ Example1 ex1;
+
+ // create a linear system
+ MathLib::LisLinearSystem eqs(ex1.dim_eqs);
+
+ // construct
+ for (size_t i=0; i<ex1.dim_eqs; i++) {
+ for (size_t j=0; j<ex1.dim_eqs; j++) {
+ double v = ex1.mat[i*ex1.dim_eqs+j];
+ if (v!=.0)
+ eqs.addMatEntry(i, j, v);
+ }
+ }
+
+ // apply BC
+ eqs.setKnownSolution(ex1.list_dirichlet_bc_id, ex1.list_dirichlet_bc_value);
+
+ // set solver options using Boost property tree
+ boost::property_tree::ptree t_root;
+ boost::property_tree::ptree t_solver;
+ t_solver.put("solver_type", "CG");
+ t_solver.put("precon_type", "NONE");
+ t_solver.put("matrix_type", "CCS");
+ t_solver.put("error_tolerance", 1e-15);
+ t_solver.put("max_iteration_step", 1000);
+ t_root.put_child("LinearSolver", t_solver);
+ eqs.setOption(t_root);
+
+ // check if the option was correctly parsed
+ MathLib::LisOption &lisOption = eqs.getOption();
+ ASSERT_EQ(MathLib::LisOption::SolverType::CG, lisOption.solver_type);
+ ASSERT_EQ(MathLib::LisOption::PreconType::NONE, lisOption.precon_type);
+ ASSERT_EQ(MathLib::LisOption::MatrixType::CCS, lisOption.matrix_type);
+ ASSERT_EQ(1e-15, lisOption.error_tolerance);
+ ASSERT_EQ(1000, lisOption.max_iterations);
+
+ // solve
+ eqs.solve();
+
+// eqs.printout();
+
+ // check solution
+ std::vector<double> vec_x(eqs.getDimension());
+ eqs.getSolVec(&vec_x[0]);
+ ASSERT_DOUBLE_ARRAY_EQ(&ex1.exH[0], &vec_x[0], ex1.dim_eqs, 1.e-5);
+}
+
+#endif
+
diff --git a/Tests/testrunner.cpp b/Tests/testrunner.cpp
index 3d3cf4f5fb858fff5ce1c51c12e98516eeb6126b..1aa85ef15b26b6e592b18f4c59fe4b82e8f1a75f 100644
--- a/Tests/testrunner.cpp
+++ b/Tests/testrunner.cpp
@@ -12,10 +12,43 @@
// ** INCLUDES **
#include "gtest/gtest.h"
+#include "logog/include/logog.hpp"
+#ifdef USE_LIS
+#include "lis.h"
+#endif
+#include "BaseLib/TemplateLogogFormatterSuppressedGCC.h"
/// Implementation of the googletest testrunner
int main(int argc, char* argv[])
{
- testing::InitGoogleTest ( &argc, argv );
- return RUN_ALL_TESTS();
+ int ret = 0;
+ LOGOG_INITIALIZE();
+ {
+ logog::Cout out;
+ BaseLib::TemplateLogogFormatterSuppressedGCC<TOPIC_LEVEL_FLAG | TOPIC_FILE_NAME_FLAG | TOPIC_LINE_NUMBER_FLAG> custom_format;
+ out.SetFormatter(custom_format);
+
+ try {
+ // initialize libraries which will be used while testing
+#ifdef USE_LIS
+ lis_initialize(&argc, &argv);
+#endif
+ // start google test
+ testing::InitGoogleTest ( &argc, argv );
+ ret = RUN_ALL_TESTS();
+ } catch (char* e) {
+ ERR(e);
+ } catch (std::exception& e) {
+ ERR(e.what());
+ } catch (...) {
+ ERR("Unknown exception occurred!");
+ }
+ // finalize libraries
+#ifdef USE_LIS
+ lis_finalize();
+#endif
+ } // make sure no logog objects exist when LOGOG_SHUTDOWN() is called.
+ LOGOG_SHUTDOWN();
+
+ return ret;
}
diff --git a/scripts/cmake/Find.cmake b/scripts/cmake/Find.cmake
index 4feea6704f4d5ca2c681a5bf9a4e3b3ca65ec3a5..fdff334950dbfd024ed2971aaedd3ca9c096e14f 100644
--- a/scripts/cmake/Find.cmake
+++ b/scripts/cmake/Find.cmake
@@ -138,6 +138,11 @@ IF(Shapelib_FOUND)
ADD_DEFINITIONS(-DShapelib_FOUND)
ENDIF() # Shapelib_FOUND
+## lis ##
+IF(OGS_USE_LIS)
+ FIND_PACKAGE( LIS REQUIRED )
+ENDIF()
+
########################
### Find other stuff ###
########################
diff --git a/scripts/cmake/Functions.cmake b/scripts/cmake/Functions.cmake
index a759e11b3f369da64464e7f4d764d50b0b295fd5..8c1c0b833a77284541b6ec7c08767dc4f6d46ce9 100644
--- a/scripts/cmake/Functions.cmake
+++ b/scripts/cmake/Functions.cmake
@@ -15,8 +15,9 @@ MACRO(GET_SOURCE_FILES SOURCE_FILES)
ENDIF()
# Get all files in the directory
- FILE(GLOB GET_SOURCE_FILES_HEADERS ${DIR}/*.h)
- FILE(GLOB GET_SOURCE_FILES_SOURCES ${DIR}/*.cpp)
+ FILE(GLOB GET_SOURCE_FILES_HEADERS RELATIVE ${CMAKE_CURRENT_SOURCE_DIR} ${DIR}/*.h)
+ FILE(GLOB GET_SOURCE_FILES_SOURCES RELATIVE ${CMAKE_CURRENT_SOURCE_DIR} ${DIR}/*.tpp)
+ FILE(GLOB GET_SOURCE_FILES_SOURCES RELATIVE ${CMAKE_CURRENT_SOURCE_DIR} ${DIR}/*.cpp)
SET(${SOURCE_FILES} ${GET_SOURCE_FILES_HEADERS} ${GET_SOURCE_FILES_SOURCES})
@@ -54,4 +55,4 @@ MACRO(ADD_GOOGLE_TESTS executable)
# message ("Adding test: ${test_name}")
ENDFOREACH(hit)
ENDFOREACH()
-ENDMACRO()
\ No newline at end of file
+ENDMACRO()
diff --git a/scripts/cmake/cmake/FindLIS.cmake b/scripts/cmake/cmake/FindLIS.cmake
new file mode 100644
index 0000000000000000000000000000000000000000..4deb12a127f3ce1160d59360c321749ce6eb6c75
--- /dev/null
+++ b/scripts/cmake/cmake/FindLIS.cmake
@@ -0,0 +1,34 @@
+# - Try to find LIS
+# Once done, this will define
+#
+# LIS_FOUND
+# LIS_INCLUDE_DIRS
+# LIS_LIBRARIES
+
+set(LIS_ROOT_DIR
+ "${LIS_ROOT_DIR}"
+ CACHE
+ PATH
+ "Directory to search for Lis library")
+
+find_path( LIS_INCLUDE_DIR
+ NAMES lis.h
+ HINTS
+ ${LIS_ROOT_DIR}/include
+ /usr/include/lis
+ $ENV{HOME}/include/
+ )
+
+find_library(LIS_LIBRARY
+ NAMES lis
+ HINTS
+ ${LIS_ROOT_DIR}/lib
+ /usr/lib
+ $ENV{HOME}/lib/
+ )
+
+SET(LIS_LIBRARIES ${LIS_LIBRARY})
+
+INCLUDE(FindPackageHandleStandardArgs)
+FIND_PACKAGE_HANDLE_STANDARD_ARGS(LIS DEFAULT_MSG LIS_LIBRARY LIS_INCLUDE_DIR)
+