diff --git a/CMakeLists.txt b/CMakeLists.txt index 1461cef9..817fda4d 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -15,6 +15,15 @@ set(CMAKE_CXX_STANDARD 17) include(FortranCInterface) FortranCInterface_VERIFY(CXX) + +FortranCInterface_HEADER(${CMAKE_BINARY_DIR}/include/fortran_interface.h + MACRO_NAMESPACE "FC_" + SYMBOL_NAMESPACE "FC_") + +# directories +set(DIRECTORIES uno ${CMAKE_BINARY_DIR}/include) + + set(CMAKE_CXX_FLAGS "-Wall -Wextra -Wnon-virtual-dtor -pedantic -Wunused-value -Wconversion") set(CMAKE_CXX_FLAGS_DEBUG "-pg") set(CMAKE_CXX_FLAGS_RELEASE "-O3 -DNDEBUG") # disable asserts @@ -28,8 +37,6 @@ set(CMAKE_MODULE_PATH ${CMAKE_MODULE_PATH} ${CMAKE_CURRENT_SOURCE_DIR}/cmake ${C option(WITH_GTEST "Enable GoogleTest" OFF) message(STATUS "GoogleTest: WITH_GTEST=${WITH_GTEST}") -# directories -set(DIRECTORIES uno) # source files file(GLOB UNO_SOURCE_FILES diff --git a/uno/solvers/BQPD/BQPDSolver.cpp b/uno/solvers/BQPD/BQPDSolver.cpp index 3382b61c..a888e48a 100644 --- a/uno/solvers/BQPD/BQPDSolver.cpp +++ b/uno/solvers/BQPD/BQPDSolver.cpp @@ -12,6 +12,11 @@ #include "tools/Infinity.hpp" #include "tools/Logger.hpp" #include "options/Options.hpp" +#include "fortran_interface.h" + +#define WSC FC_GLOBAL(wsc,WSC) +#define KKTALPHAC FC_GLOBAL(kktalphac,KKTALPHAC) +#define BQPD FC_GLOBAL(bqpd,BQPD) namespace uno { #define BIG 1e30 @@ -20,15 +25,15 @@ namespace uno { // fortran common block used in bqpd/bqpd.f extern struct { int kk, ll, kkk, lll, mxws, mxlws; - } wsc_; + } WSC; // fortran common for inertia correction in wdotd extern struct { double alpha; - } kktalphac_; + } KKTALPHAC; extern void - bqpd_(const int* n, const int* m, int* k, int* kmax, double* a, int* la, double* x, double* bl, double* bu, double* f, double* fmin, double* g, + BQPD(const int* n, const int* m, int* k, int* kmax, double* a, int* la, double* x, double* bl, double* bu, double* f, double* fmin, double* g, double* r, double* w, double* e, int* ls, double* alp, int* lp, int* mlp, int* peq, double* ws, int* lws, const int* mode, int* ifail, int* info, int* iprint, int* nout); } @@ -95,11 +100,11 @@ namespace uno { const WarmstartInformation& warmstart_information) { // initialize wsc_ common block (Hessian & workspace for BQPD) // setting the common block here ensures that several instances of BQPD can run simultaneously - wsc_.kk = static_cast(this->number_hessian_nonzeros); - wsc_.ll = static_cast(this->size_hessian_sparsity); - wsc_.mxws = static_cast(this->size_hessian_workspace); - wsc_.mxlws = static_cast(this->size_hessian_sparsity_workspace); - kktalphac_.alpha = 0; // inertia control + WSC.kk = static_cast(this->number_hessian_nonzeros); + WSC.ll = static_cast(this->size_hessian_sparsity); + WSC.mxws = static_cast(this->size_hessian_workspace); + WSC.mxlws = static_cast(this->size_hessian_sparsity_workspace); + KKTALPHAC.alpha = 0; // inertia control if (this->print_subproblem) { DEBUG << "objective gradient: " << linear_objective; @@ -142,7 +147,7 @@ namespace uno { const int mode_integer = static_cast(mode); // solve the LP/QP - bqpd_(&n, &m, &this->k, &this->kmax, this->jacobian.data(), this->jacobian_sparsity.data(), direction.primals.data(), this->lb.data(), + BQPD(&n, &m, &this->k, &this->kmax, this->jacobian.data(), this->jacobian_sparsity.data(), direction.primals.data(), this->lb.data(), this->ub.data(), &direction.subproblem_objective, &this->fmin, this->gradient_solution.data(), this->residuals.data(), this->w.data(), this->e.data(), this->active_set.data(), this->alp.data(), this->lp.data(), &this->mlp, &this->peq_solution, this->hessian_values.data(), this->hessian_sparsity.data(), &mode_integer, &this->ifail, this->info.data(), &this->iprint, &this->nout); diff --git a/uno/solvers/BQPD/BQPDSolver.hpp b/uno/solvers/BQPD/BQPDSolver.hpp index 5d8b9382..0304f517 100644 --- a/uno/solvers/BQPD/BQPDSolver.hpp +++ b/uno/solvers/BQPD/BQPDSolver.hpp @@ -4,6 +4,7 @@ #ifndef UNO_BQPDSOLVER_H #define UNO_BQPDSOLVER_H +#include #include #include "ingredients/subproblems/SubproblemStatus.hpp" #include "linear_algebra/Vector.hpp" diff --git a/uno/solvers/MA57/MA57Solver.cpp b/uno/solvers/MA57/MA57Solver.cpp index 5f5b0d53..6e482d48 100644 --- a/uno/solvers/MA57/MA57Solver.cpp +++ b/uno/solvers/MA57/MA57Solver.cpp @@ -6,23 +6,30 @@ #include "linear_algebra/SymmetricMatrix.hpp" #include "linear_algebra/Vector.hpp" #include "tools/Logger.hpp" +#include "fortran_interface.h" + +#define MA57ID FC_GLOBAL(ma57id, MA57ID) +#define MA57AD FC_GLOBAL(ma57ad, MA57AD) +#define MA57BD FC_GLOBAL(ma57bd, MA57BD) +#define MA57CD FC_GLOBAL(ma57cd, MA57CD) +#define MA57DD FC_GLOBAL(ma57dd, MA57DD) namespace uno { extern "C" { // MA57 // default values of controlling parameters - void ma57id_(double cntl[], int icntl[]); + void MA57ID(double cntl[], int icntl[]); // symbolic factorization - void ma57ad_(const int* n, const int* ne, const int irn[], const int jcn[], const int* lkeep, int keep[], int iwork[], int icntl[], int info[], double + void MA57AD(const int* n, const int* ne, const int irn[], const int jcn[], const int* lkeep, int keep[], int iwork[], int icntl[], int info[], double rinfo[]); // numerical factorization - void ma57bd_(const int* n, int* ne, const double a[], /* out */ double fact[], const int* lfact, /* out */ int ifact[], const int* lifact, + void MA57BD(const int* n, int* ne, const double a[], /* out */ double fact[], const int* lfact, /* out */ int ifact[], const int* lifact, const int* lkeep, const int keep[], int iwork[], int icntl[], double cntl[], /* out */ int info[], /* out */ double rinfo[]); // linear system solve without iterative refinement - void ma57cd_(const int* job, const int* n, double fact[], int* lfact, int ifact[], int* lifact, const int* nrhs, double rhs[], const int* lrhs, double + void MA57CD(const int* job, const int* n, double fact[], int* lfact, int ifact[], int* lifact, const int* nrhs, double rhs[], const int* lrhs, double work[], int* lwork, int iwork[], int icntl[], int info[]); // linear system solve with iterative refinement - void ma57dd_(const int* job, const int* n, int* ne, const double a[], const int irn[], const int jcn[], double fact[], int* lfact, int ifact[], int* + void MA57DD(const int* job, const int* n, int* ne, const double a[], const int irn[], const int jcn[], double fact[], int* lfact, int ifact[], int* lifact, const double rhs[], double x[], double resid[], double work[], int iwork[], int icntl[], double cntl[], int info[], double rinfo[]); } @@ -36,7 +43,7 @@ namespace uno { this->row_indices.reserve(number_nonzeros); this->column_indices.reserve(number_nonzeros); // set the default values of the controlling parameters - ma57id_(this->cntl.data(), this->icntl.data()); + MA57ID(this->cntl.data(), this->icntl.data()); // suppress warning messages this->icntl[4] = 0; // iterative refinement enabled @@ -61,7 +68,7 @@ namespace uno { const int nnz = static_cast(matrix.number_nonzeros()); // symbolic factorization - ma57ad_(/* const */ &n, + MA57AD(/* const */ &n, /* const */ &nnz, /* const */ this->row_indices.data(), /* const */ this->column_indices.data(), @@ -95,7 +102,7 @@ namespace uno { int nnz = static_cast(matrix.number_nonzeros()); // numerical factorization - ma57bd_(&n, + MA57BD(&n, &nnz, /* const */ matrix.data_pointer(), /* out */ this->fact.data(), @@ -116,7 +123,7 @@ namespace uno { // solve the linear system if (this->use_iterative_refinement) { - ma57dd_(&this->job, &n, &nnz, matrix.data_pointer(), this->row_indices.data(), this->column_indices.data(), + MA57DD(&this->job, &n, &nnz, matrix.data_pointer(), this->row_indices.data(), this->column_indices.data(), this->fact.data(), &this->factorization.lfact, this->ifact.data(), &this->factorization.lifact, rhs.data(), result.data(), this->residuals.data(), this->work.data(), this->iwork.data(), this->icntl.data(), this->cntl.data(), this->info.data(), this->rinfo.data()); @@ -125,7 +132,7 @@ namespace uno { // copy rhs into result (overwritten by MA57) result = rhs; - ma57cd_(&this->job, &n, this->fact.data(), &this->factorization.lfact, this->ifact.data(), + MA57CD(&this->job, &n, this->fact.data(), &this->factorization.lfact, this->ifact.data(), &this->factorization.lifact, &this->nrhs, result.data(), &lrhs, this->work.data(), &this->lwork, this->iwork.data(), this->icntl.data(), this->info.data()); }