Skip to content

Commit

Permalink
Merge branch 'main' into drop-cppoptlib
Browse files Browse the repository at this point in the history
  • Loading branch information
@zfergus
zfergus committed Feb 8, 2024
2 parents 4454147 + 811c9a7 commit ba80533
Show file tree
Hide file tree
Showing 6 changed files with 38 additions and 23 deletions.
22 changes: 11 additions & 11 deletions CMakeLists.txt
View file
Original file line number Diff line number Diff line change
Expand Up @@ -169,7 +169,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")
Expand Down Expand Up @@ -198,7 +198,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)
Expand All @@ -220,7 +220,7 @@ target_link_libraries(polysolve_linear PUBLIC jse::jse)
if(POLYSOLVE_WITH_HYPRE)
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()
Expand All @@ -230,22 +230,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()
Expand All @@ -255,15 +255,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()
Expand All @@ -273,22 +273,22 @@ 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()

# cuSolver solvers
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()
Expand Down
8 changes: 8 additions & 0 deletions nonlinear-solver-spec.json
View file
Original file line number Diff line number Diff line change
Expand Up @@ -694,6 +694,7 @@
"optional": [
"f_delta",
"f_delta_step_tol",
"derivative_along_delta_x_tol",
"apply_gradient_fd",
"gradient_fd_eps"
],
Expand All @@ -712,6 +713,13 @@
"type": "int",
"doc": "Dangerous Option: Quit the optimization if the solver reduces the energy by less than f_delta for consecutive f_delta_step_tol steps."
},
{
"pointer": "/advanced/derivative_along_delta_x_tol",
"default": 0,
"min": 0,
"type": "float",
"doc": "Quit the optimization if the directional derivative along the descent direction is smaller than this tolerance."
},
{
"pointer": "/advanced/apply_gradient_fd",
"default": "None",
Expand Down
2 changes: 1 addition & 1 deletion src/polysolve/linear/Solver.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -83,7 +83,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;
Expand Down
7 changes: 7 additions & 0 deletions src/polysolve/nonlinear/Criteria.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -22,6 +22,7 @@ namespace polysolve::nonlinear
gradNorm = 0;
firstGradNorm = 0;
fDeltaCount = 0;
xDeltaDotGrad = 0;
}

void Criteria::print(std::ostream &os) const
Expand All @@ -30,6 +31,7 @@ namespace polysolve::nonlinear
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;
}
Expand All @@ -53,6 +55,11 @@ namespace polysolve::nonlinear
{
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;
}

Expand Down
1 change: 1 addition & 0 deletions src/polysolve/nonlinear/Criteria.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -34,6 +34,7 @@ namespace polysolve::nonlinear
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();
Expand Down
21 changes: 10 additions & 11 deletions src/polysolve/nonlinear/Solver.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -202,12 +202,14 @@ namespace polysolve::nonlinear
m_stop.fDelta = solver_params["advanced"]["f_delta"];
m_stop.gradNorm = solver_params["grad_norm"];
m_stop.firstGradNorm = solver_params["first_grad_norm_tol"];
m_stop.xDeltaDotGrad = -solver_params["advanced"]["derivative_along_delta_x_tol"].get<double>();

// Make these relative to the characteristic length
m_stop.xDelta *= characteristic_length;
m_stop.fDelta *= characteristic_length;
m_stop.gradNorm *= characteristic_length;
m_stop.firstGradNorm *= characteristic_length;
// m_stop.xDeltaDotGrad *= characteristic_length;

m_stop.iterations = solver_params["max_iterations"];
allow_out_of_iterations = solver_params["allow_out_of_iterations"];
Expand Down Expand Up @@ -337,8 +339,9 @@ namespace polysolve::nonlinear
// --- Update direction --------------------------------------------

const bool ok = compute_update_direction(objFunc, x, grad, delta_x);
m_current.xDeltaDotGrad = delta_x.dot(grad);

if (!ok || (m_strategies[m_descent_strategy]->is_direction_descent() && m_current.gradNorm != 0 && delta_x.dot(grad) >= 0))
if (!ok || (m_strategies[m_descent_strategy]->is_direction_descent() && m_current.gradNorm != 0 && m_current.xDeltaDotGrad >= 0))
{
const std::string current_name = descent_strategy_name();

Expand All @@ -350,18 +353,14 @@ namespace polysolve::nonlinear
m_status = Status::NotDescentDirection;
log_and_throw_error(m_logger, "[{}][{}] direction is not a descent direction on last strategy (‖Δx‖={:g}; ‖g‖={:g}; Δx⋅g={:g}≥0); stopping",
current_name, m_line_search->name(),
delta_x.norm(), compute_grad_norm(x, grad), delta_x.dot(grad));
delta_x.norm(), compute_grad_norm(x, grad), m_current.xDeltaDotGrad);
}

if (ok) // if ok, then the direction is not a descent direction
{
m_logger.debug(
"[{}][{}] direction is not a descent direction (‖Δx‖={:g}; ‖g‖={:g}; Δx⋅g={:g}≥0); reverting to {}",
current_name, m_line_search->name(),
delta_x.norm(), compute_grad_norm(x, grad), delta_x.dot(grad), descent_strategy_name());
}

m_status = Status::Continue;
m_logger.debug(
"[{}][{}] direction is not a descent direction (‖Δx‖={:g}; ‖g‖={:g}; Δx⋅g={:g}≥0); reverting to {}",
current_name, m_line_search->name(),
delta_x.norm(), compute_grad_norm(x, grad), m_current.xDeltaDotGrad, descent_strategy_name());
m_status = Status::NotDescentDirection;
continue;
}

Expand Down

0 comments on commit ba80533

Please sign in to comment.