Merge remote-tracking branch 'upstream/main' into trilinos

Author: yib0liu

CMakeLists.txt

@@ -171,7 +171,7 @@ target_include_directories(polysolve PUBLIC ${PROJECT_SOURCE_DIR}/src)
 ################################################################################
 
 if(POLYSOLVE_LARGE_INDEX)
-    target_compile_definitions(polysolve_linear PUBLIC -DPOLYSOLVE_LARGE_INDEX)
+    target_compile_definitions(polysolve_linear PUBLIC POLYSOLVE_LARGE_INDEX)
 endif()
 
 target_compile_definitions(polysolve_linear PRIVATE POLYSOLVE_LINEAR_SPEC="${PROJECT_SOURCE_DIR}/linear-solver-spec.json")
@@ -200,7 +200,7 @@ if(POLYSOLVE_WITH_ACCELERATE)
     find_package(BLAS REQUIRED)
     find_package(LAPACK REQUIRED)
     target_link_libraries(polysolve_linear PRIVATE BLAS::BLAS LAPACK::LAPACK)
-    target_compile_definitions(polysolve_linear PRIVATE -DPOLYSOLVE_WITH_ACCELERATE)
+    target_compile_definitions(polysolve_linear PUBLIC POLYSOLVE_WITH_ACCELERATE)
 endif()
 
 include(eigen)
@@ -222,7 +222,7 @@ target_link_libraries(polysolve_linear PUBLIC jse::jse)
 if(POLYSOLVE_WITH_HYPRE AND NOT POLYSOLVE_WITH_TRILINOS)
     include(hypre)
     target_link_libraries(polysolve_linear PUBLIC HYPRE::HYPRE)
-    target_compile_definitions(polysolve_linear PRIVATE -DPOLYSOLVE_WITH_HYPRE)
+    target_compile_definitions(polysolve_linear PUBLIC POLYSOLVE_WITH_HYPRE)
     if(HYPRE_WITH_MPI)
         target_compile_definitions(polysolve_linear PUBLIC HYPRE_WITH_MPI)
     endif()
@@ -232,22 +232,22 @@ endif()
 if(POLYSOLVE_WITH_CHOLMOD)
     include(suitesparse)
     target_link_libraries(polysolve_linear PRIVATE SuiteSparse::CHOLMOD)
-    target_compile_definitions(polysolve_linear PRIVATE -DPOLYSOLVE_WITH_CHOLMOD)
+    target_compile_definitions(polysolve_linear PUBLIC POLYSOLVE_WITH_CHOLMOD)
 endif()
 
 # MKL library
 if(POLYSOLVE_WITH_MKL)
     include(mkl)
     target_link_libraries(polysolve_linear PRIVATE mkl::mkl)
-    target_compile_definitions(polysolve_linear PRIVATE -DPOLYSOLVE_WITH_MKL)
+    target_compile_definitions(polysolve_linear PUBLIC POLYSOLVE_WITH_MKL)
 endif()
 
 # Pardiso solver
 if(POLYSOLVE_WITH_PARDISO)
     include(pardiso)
     if(TARGET Pardiso::Pardiso)
         target_link_libraries(polysolve_linear PRIVATE Pardiso::Pardiso)
-        target_compile_definitions(polysolve_linear PRIVATE -DPOLYSOLVE_WITH_PARDISO)
+        target_compile_definitions(polysolve_linear PUBLIC POLYSOLVE_WITH_PARDISO)
     else()
         message(WARNING "Pardiso not found, solver will not be available.")
     endif()
@@ -257,15 +257,15 @@ endif()
 if(POLYSOLVE_WITH_UMFPACK)
     include(suitesparse)
     target_link_libraries(polysolve_linear PRIVATE SuiteSparse::UMFPACK)
-    target_compile_definitions(polysolve_linear PRIVATE -DPOLYSOLVE_WITH_UMFPACK)
+    target_compile_definitions(polysolve_linear PUBLIC POLYSOLVE_WITH_UMFPACK)
 endif()
 
 # SuperLU solver
 if(POLYSOLVE_WITH_SUPERLU)
     include(superlu)
     if(TARGET SuperLU::SuperLU)
         target_link_libraries(polysolve_linear PRIVATE SuperLU::SuperLU)
-        target_compile_definitions(polysolve_linear PRIVATE -DPOLYSOLVE_WITH_SUPERLU)
+        target_compile_definitions(polysolve_linear PUBLIC POLYSOLVE_WITH_SUPERLU)
     else()
         message(WARNING "SuperLU not found, solver will not be available.")
     endif()
@@ -275,14 +275,14 @@ endif()
 if(POLYSOLVE_WITH_AMGCL)
     include(amgcl)
     target_link_libraries(polysolve_linear PUBLIC amgcl::amgcl)
-    target_compile_definitions(polysolve_linear PRIVATE -DPOLYSOLVE_WITH_AMGCL)
+    target_compile_definitions(polysolve_linear PUBLIC POLYSOLVE_WITH_AMGCL)
 endif()
 
 # Spectra (MPL 2.0)
 if(POLYSOLVE_WITH_SPECTRA)
     include(spectra)
     target_link_libraries(polysolve_linear PUBLIC Spectra::Spectra)
-    target_compile_definitions(polysolve_linear PRIVATE -DPOLYSOLVE_WITH_SPECTRA)
+    target_compile_definitions(polysolve_linear PUBLIC POLYSOLVE_WITH_SPECTRA)
 endif()
 
 # Trilinos
@@ -297,7 +297,7 @@ if(POLYSOLVE_WITH_CUSOLVER)
     include(cusolverdn)
     if(TARGET CUDA::cusolver)
         target_link_libraries(polysolve_linear PUBLIC CUDA::cusolver)
-        target_compile_definitions(polysolve_linear PRIVATE -DPOLYSOLVE_WITH_CUSOLVER)
+        target_compile_definitions(polysolve_linear PUBLIC POLYSOLVE_WITH_CUSOLVER)
     else()
         message(WARNING "cuSOLVER not found, solver will not be available.")
     endif()
@@ -320,17 +320,13 @@ target_link_libraries(polysolve_linear PRIVATE polysolve::warnings)
 # polysolve::linear
 target_link_libraries(polysolve PUBLIC polysolve::linear)
 
-# CppNumericalSolvers
-include(cppoptlib)
-target_link_libraries(polysolve PUBLIC cppoptlib)
-
 # LBFGSpp
 include(LBFGSpp)
-target_link_libraries(polysolve PUBLIC LBFGSpp::LBFGSpp)
+target_link_libraries(polysolve PRIVATE LBFGSpp::LBFGSpp)
 
 # finite-diff (include this after eigen)
 include(finite-diff)
-target_link_libraries(polysolve PUBLIC finitediff::finitediff)
+target_link_libraries(polysolve PRIVATE finitediff::finitediff)
 
 # Sanitizers
 if(POLYSOLVE_WITH_SANITIZERS)

cmake/recipes/cppoptlib.cmake

@@ -0,0 +1,11 @@
+coverage:
+  status:
+    project:
+      default:
+        target: auto
+        threshold: 1%
+    patch:
+      default:
+        target: 75%
+        threshold: 5%
+        only_pulls: true

codecov.yml

@@ -1,5 +1,4 @@
-[
-    {
+[{
         "pointer": "/",
         "default": null,
         "type": "object",
@@ -562,10 +561,8 @@
         "type": "string",
         "options": [
             "Armijo",
-            "ArmijoAlt",
             "RobustArmijo",
             "Backtracking",
-            "MoreThuente",
             "None"
         ],
         "doc": "Line-search type"

nonlinear-solver-spec.json

@@ -113,7 +113,7 @@ namespace polysolve::linear
                 }
             }
             if (!accepted_solver.empty())
-                logger.info("Solver {} is the highest priority availble solver; using it.", accepted_solver);
+                logger.info("Solver {} is the highest priority available solver; using it.", accepted_solver);
             else
                 logger.warn("No valid solver found in the list of specified solvers!");
             params[lin_solver_ptr] = accepted_solver;

src/polysolve/linear/Solver.cpp

@@ -1,12 +1,14 @@
 set(SOURCES
 	BoxConstraintSolver.cpp
 	BoxConstraintSolver.hpp
-	Solver.hpp
-	Solver.cpp
-	Problem.hpp
-	Problem.cpp
-	PostStepData.hpp
+	Criteria.cpp
+	Criteria.hpp
 	PostStepData.cpp
+	PostStepData.hpp
+	Problem.cpp
+	Problem.hpp
+	Solver.cpp
+	Solver.hpp
 )
 
 source_group(TREE "${CMAKE_CURRENT_SOURCE_DIR}" PREFIX "Source Files" FILES ${SOURCES})

src/polysolve/nonlinear/CMakeLists.txt

@@ -0,0 +1,112 @@
+// Source: https://github.com/PatWie/CppNumericalSolvers/blob/7eddf28fa5a8872a956d3c8666055cac2f5a535d/include/cppoptlib/meta.h
+// License: MIT
+#include "Criteria.hpp"
+
+namespace polysolve::nonlinear
+{
+    bool is_converged_status(const Status s)
+    {
+        return s == Status::XDeltaTolerance || s == Status::FDeltaTolerance || s == Status::GradNormTolerance;
+    }
+
+    Criteria::Criteria()
+    {
+        reset();
+    }
+
+    void Criteria::reset()
+    {
+        iterations = 0;
+        xDelta = 0;
+        fDelta = 0;
+        gradNorm = 0;
+        firstGradNorm = 0;
+        fDeltaCount = 0;
+        xDeltaDotGrad = 0;
+    }
+
+    void Criteria::print(std::ostream &os) const
+    {
+        os << "Iterations: " << iterations << std::endl;
+        os << "xDelta:     " << xDelta << std::endl;
+        os << "fDelta:     " << fDelta << std::endl;
+        os << "GradNorm:   " << gradNorm << std::endl;
+        os << "xDeltaDotGrad: " << xDeltaDotGrad << std::endl;
+        os << "FirstGradNorm: " << firstGradNorm << std::endl;
+        os << "fDeltaCount: " << fDeltaCount << std::endl;
+    }
+
+    Status checkConvergence(const Criteria &stop, const Criteria &current)
+    {
+        if (stop.iterations > 0 && current.iterations > stop.iterations)
+        {
+            return Status::IterationLimit;
+        }
+        if (stop.xDelta > 0 && current.xDelta < stop.xDelta)
+        {
+            return Status::XDeltaTolerance;
+        }
+        if (stop.fDelta > 0 && current.fDelta < stop.fDelta && current.fDeltaCount >= stop.fDeltaCount)
+        {
+            return Status::FDeltaTolerance;
+        }
+        const double stopGradNorm = current.iterations == 0 ? stop.firstGradNorm : stop.gradNorm;
+        if (stopGradNorm > 0 && current.gradNorm < stopGradNorm)
+        {
+            return Status::GradNormTolerance;
+        }
+        // Δx⋅∇f ≥ 0 means that the search direction is not a descent direction
+        if (stop.xDeltaDotGrad < 0 && current.xDeltaDotGrad > stop.xDeltaDotGrad)
+        {
+            return Status::NotDescentDirection;
+        }
+        return Status::Continue;
+    }
+
+    std::ostream &operator<<(std::ostream &os, const Status &s)
+    {
+        switch (s)
+        {
+        case Status::NotStarted:
+            os << "Solver not started.";
+            break;
+        case Status::Continue:
+            os << "Convergence criteria not reached.";
+            break;
+        case Status::IterationLimit:
+            os << "Iteration limit reached.";
+            break;
+        case Status::XDeltaTolerance:
+            os << "Change in parameter vector too small.";
+            break;
+        case Status::FDeltaTolerance:
+            os << "Change in cost function value too small.";
+            break;
+        case Status::GradNormTolerance:
+            os << "Gradient vector norm too small.";
+            break;
+        case Status::ObjectiveCustomStop:
+            os << "Objective function specified to stop.";
+            break;
+        case Status::NanEncountered:
+            os << "Objective or gradient function returned NaN.";
+            break;
+        case Status::NotDescentDirection:
+            os << "Search direction not a descent direction.";
+            break;
+        case Status::LineSearchFailed:
+            os << "Line search failed.";
+            break;
+        default:
+            os << "Unknown status.";
+            break;
+        }
+        return os;
+    }
+
+    std::ostream &operator<<(std::ostream &os, const Criteria &c)
+    {
+        c.print(os);
+        return os;
+    }
+} // namespace polysolve::nonlinear

src/polysolve/nonlinear/Criteria.cpp

@@ -0,0 +1,52 @@
+#pragma once
+
+#include <cstddef>
+#include <iostream>
+
+namespace polysolve::nonlinear
+{
+    // Source: https://github.com/PatWie/CppNumericalSolvers/blob/7eddf28fa5a8872a956d3c8666055cac2f5a535d/include/cppoptlib/meta.h
+    // License: MIT
+
+    enum class Status
+    {
+        NotStarted = -1, ///< The solver has not been started
+        Continue = 0,    ///< The solver should continue
+        // Success cases
+        IterationLimit,      ///< The maximum number of iterations has been reached
+        XDeltaTolerance,     ///< The change in the parameter vector is below the tolerance
+        FDeltaTolerance,     ///< The change in the cost function is below the tolerance
+        GradNormTolerance,   ///< The norm of the gradient vector is below the tolerance
+        ObjectiveCustomStop, ///< The objective function specified to stop
+        // Failure cases
+        NanEncountered,      ///< The objective function returned NaN
+        NotDescentDirection, ///< The search direction is not a descent direction
+        LineSearchFailed,    ///< The line search failed
+    };
+
+    bool is_converged_status(const Status s);
+
+    class Criteria
+    {
+    public:
+        size_t iterations;    ///< Maximum number of iterations
+        double xDelta;        ///< Minimum change in parameter vector
+        double fDelta;        ///< Minimum change in cost function
+        double gradNorm;      ///< Minimum norm of gradient vector
+        double firstGradNorm; ///< Initial norm of gradient vector
+        double xDeltaDotGrad; ///< Dot product of parameter vector and gradient vector
+        unsigned fDeltaCount; ///< Number of steps where fDelta is satisfied
+
+        Criteria();
+
+        void reset();
+
+        void print(std::ostream &os) const;
+    };
+
+    Status checkConvergence(const Criteria &stop, const Criteria &current);
+
+    std::ostream &operator<<(std::ostream &os, const Status &s);
+
+    std::ostream &operator<<(std::ostream &os, const Criteria &c);
+} // namespace polysolve::nonlinear