Merge pull request #213 from zoq/report-callback

Optimizer Summary Report

.appveyor.yml

@@ -0,0 +1,57 @@
+environment:
+  ARMADILLO_DOWNLOAD: "http://ftp.fau.de/macports/distfiles/armadillo/armadillo-8.400.0.tar.xz"
+  BLAS_LIBRARY: "%APPVEYOR_BUILD_FOLDER%/OpenBLAS.0.2.14.1/lib/native/lib/x64/libopenblas.dll.a"
+  BLAS_LIBRARY_DLL: "%APPVEYOR_BUILD_FOLDER%/OpenBLAS.0.2.14.1/lib/native/lib/x64/libopenblas.dll"
+
+  matrix:
+    - APPVEYOR_BUILD_WORKER_IMAGE: Visual Studio 2015
+      VSVER: Visual Studio 14 2015 Win64
+      MSBUILD: C:\Program Files (x86)\MSBuild\14.0\bin\MSBuild.exe
+
+    - APPVEYOR_BUILD_WORKER_IMAGE: Visual Studio 2017
+      VSVER: Visual Studio 15 2017 Win64
+      MSBUILD: C:\Program Files (x86)\Microsoft Visual Studio\2017\Community\MSBuild\15.0\Bin\MSBuild.exe
+
+    - APPVEYOR_BUILD_WORKER_IMAGE: Visual Studio 2019
+      VSVER: Visual Studio 16 2019
+      MSBUILD: C:\Program Files (x86)\Microsoft Visual Studio\2019\Community\MSBuild\Current\Bin\MSBuild.exe
+
+configuration: Release
+
+install:
+  - ps: nuget install OpenBLAS -o "${env:APPVEYOR_BUILD_FOLDER}"
+
+build_script:
+  # First, download and build Armadillo.
+  - cd ..
+  - appveyor DownloadFile %ARMADILLO_DOWNLOAD% -FileName armadillo.tar.xz
+  - 7z x armadillo.tar.xz -so -txz | 7z x -si -ttar > nul
+  - cd armadillo-8.400.0 && mkdir build && cd build
+  - >
+    cmake -G "%VSVER%"
+    -DBLAS_LIBRARY:FILEPATH=%BLAS_LIBRARY%
+    -DLAPACK_LIBRARY:FILEPATH=%BLAS_LIBRARY%
+    -DCMAKE_PREFIX:FILEPATH="%APPVEYOR_BUILD_FOLDER%/armadillo"
+    -DBUILD_SHARED_LIBS=OFF
+    -DCMAKE_BUILD_TYPE=Release ..
+  - >
+    "%MSBUILD%" "armadillo.sln"
+    /m /verbosity:quiet /p:Configuration=Release;Platform=x64
+  - cd ../..
+
+  # Now build ensmallen.
+  - cd ensmallen && mkdir build && cd build
+  - >
+    cmake -G "%VSVER%"
+    -DARMADILLO_INCLUDE_DIR=%APPVEYOR_BUILD_FOLDER%/../armadillo-8.400.0/include/
+    -DARMADILLO_LIBRARIES=%BLAS_LIBRARY%
+    -DLAPACK_LIBRARY=%BLAS_LIBRARY%
+    -DBLAS_LIBRARY=%BLAS_LIBRARY%
+    -DCMAKE_BUILD_TYPE=Release ..
+  - >
+    "%MSBUILD%" "ensmallen.sln"
+    /m /verbosity:minimal /nologo /p:BuildInParallel=true
+
+  # Run tests after copying libraries.
+  - ps: cp C:\projects\ensmallen\OpenBLAS.0.2.14.1\lib\native\bin\x64\*.* C:\projects\ensmallen\build\
+  - ctest -C Release -V --output-on-failure .

.travis.yml

@@ -1,26 +1,32 @@
-sudo: required
+os: linux
 dist: trusty
 language: cpp
 
 env:
-  - ARMADILLO=latest
+  - ARMADILLO=latest SANITY_HISTORY=perform
   - ARMADILLO=minimum
 
-before_install:
+stages:
+  - test
+  - name: sanity
+    if: type = pull_request AND env(SANITY_HISTORY) = "perform"
+
+jobs:
+  include:
+    - stage: sanity
+      name: "HISTORY.md Check"
+      script: sh ./scripts/history-update-check.sh
+
+script:
   - sudo apt-get update
   - sudo apt-get install -y --allow-unauthenticated libopenblas-dev liblapack-dev g++ xz-utils
   - if [ $ARMADILLO == "latest" ]; then
       curl https://ftp.fau.de/macports/distfiles/armadillo/`curl https://ftp.fau.de/macports/distfiles/armadillo/ -- | grep '.tar.xz' | sed 's/^.*<a href="\(armadillo-[0-9]*.[0-9]*.[0-9]*.tar.xz\)".*$/\1/' | tail -1` | tar xvJ && cd armadillo*;
-    fi
-  - if [ $ARMADILLO == "minimum" ]; then
+    else
       curl https://ftp.fau.de/macports/distfiles/armadillo/armadillo-8.400.0.tar.xz | tar -xvJ && cd armadillo*;
     fi
   - cmake . && make && sudo make install && cd ..
-
-install:
   - mkdir build && cd build && cmake -DCMAKE_CXX_FLAGS="-Werror" -DCMAKE_C_FLAGS="-Werror" .. && make -j2
-
-script:
   - CTEST_OUTPUT_ON_FAILURE=1 travis_wait 30 ctest -j2
 
 notifications:

CMake/FindArmadillo.cmake

@@ -366,4 +366,10 @@ mark_as_advanced(
   ARMADILLO_INCLUDE_DIR
   ARMADILLO_LIBRARIES)
 
+if (ARMADILLO_FOUND AND NOT TARGET Armadillo:Armadillo)
+  add_library(Armadillo::Armadillo INTERFACE IMPORTED)
+  set_target_properties(Armadillo::Armadillo PROPERTIES INTERFACE_INCLUDE_DIRECTORIES "${ARMADILLO_INCLUDE_DIR}"
+      INTERFACE_LINK_LIBRARIES "${ARMADILLO_LIBRARIES}")
+endif()
+
 #======================

CMake/ensmallen-config.cmake.in

@@ -0,0 +1,4 @@
+@PACKAGE_INIT@
+
+include(${CMAKE_CURRENT_LIST_DIR}/@TARGETS_EXPORT_NAME@.cmake)
+check_required_components(ensmallen)

CMakeLists.txt

@@ -1,63 +1,92 @@
 # ensmallen CMake configuration.  This project has no configurable options---it
 # just installs the headers to the install location, and optionally builds the
 # test program.
-cmake_minimum_required(VERSION 2.8.10)
-project(ensmallen C CXX)
+cmake_minimum_required(VERSION 3.3.2)
+project(ensmallen 
+        LANGUAGES C CXX)
 
+# Configurable options for CMake.
 option(USE_OPENMP "If available, use OpenMP for parallelization." ON)
 option(BUILD_TESTS "Build tests." ON)
 
 set(CMAKE_MODULE_PATH ${CMAKE_MODULE_PATH} "${CMAKE_SOURCE_DIR}/CMake")
 
-# Ensure we have C++11 features.  Since we support CMake < 3.1, this needs a
-# little bit of special handling.
-if ((${CMAKE_MAJOR_VERSION} LESS 3 OR
-    (${CMAKE_MAJOR_VERSION} EQUAL 3 AND ${CMAKE_MINOR_VERSION} LESS 1))
-    AND NOT FORCE_CXX11)
-  # Older versions of CMake do not support target_compile_features(), so we have
-  # to use something kind of hacky.
-  include(CMake/CXX11.cmake)
-  check_for_cxx11_compiler(HAS_CXX11)
-  if(NOT HAS_CXX11)
-    message(FATAL_ERROR "No C++11 compiler available!")
-  endif()
-  enable_cxx11()
-else()
-  # set required standard to c++11
-  set(CMAKE_CXX_STANDARD 11)
-  set(CMAKE_CXX_STANDARD_REQUIRED ON)
-endif ()
+# Set required C++ standard to C++11.
+set(CMAKE_CXX_STANDARD 11)
+set(CMAKE_CXX_STANDARD_REQUIRED ON)
 
-# Detect OpenMP support in a compiler. If the compiler supports OpenMP, flags
-# to compile with OpenMP are returned and added for compilation.
-if (USE_OPENMP)
-  find_package(OpenMP)
-endif ()
+# Extract version from sources.
+set(ENSMALLEN_VERSION_FILE_NAME "${PROJECT_SOURCE_DIR}/include/ensmallen_bits/ens_version.hpp")
 
-if (OPENMP_FOUND)
-  set(CMAKE_C_FLAGS "${CMAKE_C_FLAGS} ${OpenMP_C_FLAGS}")
-  set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} ${OpenMP_CXX_FLAGS}")
-else ()
-  # Disable warnings for all the unknown OpenMP pragmas.
-  set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -Wno-unknown-pragmas")
-endif ()
-
-# Set the CFLAGS and CXXFLAGS depending on the options the user specified.
-if(CMAKE_COMPILER_IS_GNUCC OR "${CMAKE_CXX_COMPILER_ID}" STREQUAL "Clang")
-
-  set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -Wall -Wpedantic -Wunused-parameter")
-  set(CMAKE_C_FLAGS "${CMAKE_C_FLAGS} -Wall -Wpedantic -Wunused-parameter")
+if(NOT EXISTS ${ENSMALLEN_VERSION_FILE_NAME})
+  message(FATAL_ERROR "Can't read ${ENSMALLEN_VERSION_FILE_NAME}")
 endif()
 
-# The only dependency we need is Armadillo.
-#
-# We keep the minimum version in sync with mlpack, otherwise we could have
-# irritating compatibility issues.
-find_package(Armadillo 8.400.0 REQUIRED)
-include_directories(BEFORE "${ARMADILLO_INCLUDE_DIR}")
-include_directories(BEFORE "${CMAKE_SOURCE_DIR}/include/")
+file(READ ${ENSMALLEN_VERSION_FILE_NAME} ENSMALLEN_VERSION_FILE_CONTENTS)
+string(REGEX REPLACE ".*#define ENS_VERSION_MAJOR ([0-9]+).*" "\\1" ENSMALLEN_VERSION_MAJOR "${ENSMALLEN_VERSION_FILE_CONTENTS}")
+string(REGEX REPLACE ".*#define ENS_VERSION_MINOR ([0-9]+).*" "\\1" ENSMALLEN_VERSION_MINOR "${ENSMALLEN_VERSION_FILE_CONTENTS}")
+string(REGEX REPLACE ".*#define ENS_VERSION_PATCH ([0-9]+).*" "\\1" ENSMALLEN_VERSION_PATCH "${ENSMALLEN_VERSION_FILE_CONTENTS}")
 
-# Install the headers to the correct location.
+message(STATUS "Configuring ensmallen ${ENSMALLEN_VERSION_MAJOR}.${ENSMALLEN_VERSION_MINOR}.${ENSMALLEN_VERSION_PATCH}")
+set(VERSION "${ENSMALLEN_VERSION_MAJOR}.${ENSMALLEN_VERSION_MINOR}.${ENSMALLEN_VERSION_PATCH}")
+
+# Create library target.
+add_library(ensmallen INTERFACE)
+target_include_directories(ensmallen INTERFACE
+  $<BUILD_INTERFACE:${PROJECT_SOURCE_DIR}/include>
+  $<INSTALL_INTERFACE:include>)
+
+# Set warning flags for target.
+if(MSVC)
+  target_compile_options(ensmallen INTERFACE $<BUILD_INTERFACE:/Wall>)
+else()
+  target_compile_options(ensmallen INTERFACE $<BUILD_INTERFACE:-Wall -Wpedantic -Wunused-parameter>)
+endif()
+
+# Find OpenMP and link it.
+if(USE_OPENMP)
+  if(NOT TARGET OpenMP::OpenMP_CXX)
+    find_package(Threads REQUIRED)
+    add_library(OpenMP::OpenMP_CXX IMPORTED INTERFACE)
+    set_property(TARGET OpenMP::OpenMP_CXX
+        PROPERTY INTERFACE_COMPILE_OPTIONS ${OpenMP_CXX_FLAGS})
+    # Only works if the same flag is passed to the linker; use CMake 3.9+ otherwise (Intel, AppleClang).
+    set_property(TARGET OpenMP::OpenMP_CXX
+        PROPERTY INTERFACE_LINK_LIBRARIES ${OpenMP_CXX_FLAGS} Threads::Threads)
+  endif()
+  target_link_libraries(ensmallen INTERFACE OpenMP::OpenMP_CXX)
+endif()
+
+# Find Armadillo and link it.
+find_package(Armadillo 8.400.0 REQUIRED)
+target_link_libraries(ensmallen INTERFACE Armadillo::Armadillo)
+
+# Set helper variables for creating the version, config and target files.
+include(CMakePackageConfigHelpers)
+set(ENSMALLEN_CMAKE_DIR "lib/cmake/ensmallen" CACHE STRING
+      "Installation directory for cmake files, relative to ${CMAKE_INSTALL_PREFIX}.")
+set(VERSION_CONFIG "${PROJECT_BINARY_DIR}/ensmallen-config-version.cmake")
+set(PROJECT_CONFIG "${PROJECT_BINARY_DIR}/ensmallen-config.cmake")
+set(TARGETS_EXPORT_NAME ensmallen-targets)
+
+# Generate the version, config and target files into the build directory.
+write_basic_package_version_file(${VERSION_CONFIG}
+    VERSION ${VERSION}
+    COMPATIBILITY AnyNewerVersion)
+configure_package_config_file(${PROJECT_SOURCE_DIR}/CMake/ensmallen-config.cmake.in
+    ${PROJECT_CONFIG}
+    INSTALL_DESTINATION ${ENSMALLEN_CMAKE_DIR})
+export(TARGETS ensmallen NAMESPACE ensmallen::
+    FILE ${PROJECT_BINARY_DIR}/${TARGETS_EXPORT_NAME}.cmake)
+
+# Install version, config and target files.
+install(FILES ${PROJECT_CONFIG} ${VERSION_CONFIG}
+    DESTINATION ${ENSMALLEN_CMAKE_DIR})
+install(EXPORT ${TARGETS_EXPORT_NAME} DESTINATION ${ENSMALLEN_CMAKE_DIR}
+    NAMESPACE ensmallen::)
+
+# Export the targets and install the header files.
+install(TARGETS ensmallen EXPORT ${TARGETS_EXPORT_NAME} DESTINATION lib)
 install(DIRECTORY "${CMAKE_SOURCE_DIR}/include/ensmallen_bits"
         DESTINATION "${CMAKE_INSTALL_PREFIX}/include"
         PATTERN "*~" EXCLUDE
@@ -65,8 +94,8 @@ install(DIRECTORY "${CMAKE_SOURCE_DIR}/include/ensmallen_bits"
 install(FILES ${CMAKE_SOURCE_DIR}/include/ensmallen.hpp
         DESTINATION "${CMAKE_INSTALL_PREFIX}/include")
 
+# Enable testing and build tests.
 enable_testing()
-
 if (BUILD_TESTS)
   add_subdirectory(tests)
 endif()

CONTRIBUTING.md

@@ -108,8 +108,8 @@ $ cd ensmallen
 
 # - or -
 
-$ wget http://ensmallen.org/files/ensmallen-2.11.1.tar.gz
-$ tar -xvzpf ensmallen-2.11.1.tar.gz
+$ wget http://ensmallen.org/files/ensmallen-2.14.2.tar.gz
+$ tar -xvzpf ensmallen-2.14.2.tar.gz
 $ cd ensmallen-latest
 ```
 

COPYRIGHT.txt

@@ -38,6 +38,7 @@ Copyright:
   Copyright 2019, Rahul Ganesh Prabhu
   Copyright 2019, Roberto Hueso <robertohueso96@gmail.com>
   Copyright 2019, Sayan Goswami <sayan.goswami.106@gmail.com>
+  Copyright 2020, Joe Dinius <josephwdinius@gmail.com>
 
 License: BSD-3-clause
   All rights reserved.

HISTORY.md

@@ -1,3 +1,111 @@
+### ensmallen ?.??.?: "???"
+###### ????-??-??
+ * Make a few tests more robust
+   ([#228](https://github.com/mlpack/ensmallen/pull/228)).
+
+ * Add release date to version information. ([#226](https://github.com/mlpack/ensmallen/pull/226))
+ 
+ * Fix typo in release script
+   ([#236](https://github.com/mlpack/ensmallen/pull/236)).
+
+### ensmallen 2.14.2: "No Direction Home"
+###### 2020-08-31
+ * Fix implementation of fonesca fleming problem function f1 and f2 
+   type usage and negative signs. ([#223](https://github.com/mlpack/ensmallen/pull/223))
+
+### ensmallen 2.14.1: "No Direction Home"
+###### 2020-08-19
+ * Fix release script (remove hardcoded information, trim leading whitespaces
+   introduced by `wc -l` in MacOS)
+    ([#216](https://github.com/mlpack/ensmallen/pull/216),
+     [#220](https://github.com/mlpack/ensmallen/pull/220)).
+
+ * Adjust tolerance for AugLagrangian convergence based on element type
+   ([#217](https://github.com/mlpack/ensmallen/pull/217)).
+
+### ensmallen 2.14.0: "No Direction Home"
+###### 2020-08-10
+ * Add NSGA2 optimizer for multi-objective functions
+    ([#149](https://github.com/mlpack/ensmallen/pull/149)).
+
+ * Update automatic website update release script
+   ([#207](https://github.com/mlpack/ensmallen/pull/207)).
+
+ * Clarify and fix documentation for constrained optimizers
+   ([#201](https://github.com/mlpack/ensmallen/pull/201)).
+
+ * Fix L-BFGS convergence when starting from a minimum
+   ([#201](https://github.com/mlpack/ensmallen/pull/201)).
+
+* Add optimizer summary report callback 
+   ([#213](https://github.com/mlpack/ensmallen/pull/213)).
+
+### ensmallen 2.13.0: "Automatically Automated Automation"
+###### 2020-07-15
+ * Fix CMake package export
+    ([#198](https://github.com/mlpack/ensmallen/pull/198)).
+
+ * Allow early stop callback to accept a lambda function
+   ([#165](https://github.com/mlpack/ensmallen/pull/165)).
+
+### ensmallen 2.12.1: "Stir Crazy"
+###### 2020-04-20
+ * Fix total number of epochs and time estimation for ProgressBar callback
+    ([#181](https://github.com/mlpack/ensmallen/pull/181)).
+
+ * Handle SpSubview_col and SpSubview_row in Armadillo 9.870
+    ([#194](https://github.com/mlpack/ensmallen/pull/194)).
+
+ * Minor documentation fixes
+    ([#197](https://github.com/mlpack/ensmallen/pull/197)).
+
+### ensmallen 2.12.0: "Stir Crazy"
+###### 2020-03-28
+ * Correction in the formulation of sigma in CMA-ES
+    ([#183](https://github.com/mlpack/ensmallen/pull/183)).
+
+ * Remove deprecated methods from PrimalDualSolver implementation
+    ([#185](https://github.com/mlpack/ensmallen/pull/185).
+
+ * Update logo ([#186](https://github.com/mlpack/ensmallen/pull/186)).
+
+### ensmallen 2.11.5: "The Poster Session Is Full"
+###### 2020-03-11
+  * Change "mathematical optimization" term to "numerical optimization" in the
+    documentation ([#177](https://github.com/mlpack/ensmallen/pull/177)).
+
+### ensmallen 2.11.4: "The Poster Session Is Full"
+###### 2020-03-03
+  * Require new HISTORY.md entry for each PR.
+    ([#171](https://github.com/mlpack/ensmallen/pull/171),
+     [#172](https://github.com/mlpack/ensmallen/pull/172),
+     [#175](https://github.com/mlpack/ensmallen/pull/175)).
+
+  * Update/fix example documentation
+    ([#174](https://github.com/mlpack/ensmallen/pull/174)).
+
+### ensmallen 2.11.3: "The Poster Session Is Full"
+###### 2020-02-19
+  * Prevent spurious compiler warnings
+    ([#161](https://github.com/mlpack/ensmallen/pull/161)).
+
+  * Fix minor memory leaks
+    ([#167](https://github.com/mlpack/ensmallen/pull/167)).
+
+  * Revamp CMake configuration
+    ([#152](https://github.com/mlpack/ensmallen/pull/152)).
+
+### ensmallen 2.11.2: "The Poster Session Is Full"
+###### 2020-01-16
+  * Allow callback instantiation for SGD based optimizer
+    ([#138](https://github.com/mlpack/ensmallen/pull/155)).
+
+  * Minor test stability fixes on i386
+    ([#156](https://github.com/mlpack/ensmallen/pull/156)).
+
+  * Fix Lookahead MaxIterations() check.
+    ([#159](https://github.com/mlpack/ensmallen/pull/159)).
+
 ### ensmallen 2.11.1: "The Poster Session Is Full"
 ###### 2019-12-28
   * Fix Lookahead Synchronization period type

README.md

@@ -1,10 +1,14 @@
-**ensmallen** is a C++ header-only library for mathematical optimization.
+<h2 align="center">
+  <a href="http://ensmallen.org/"><img src="http://ensmallen.org/img/ensmallen_text.svg" style="background-color:rgba(0,0,0,0);" height=230 alt="ensmallen: a C++ header-only library for numerical optimization"></a>
+</h2>
+
+**ensmallen** is a C++ header-only library for numerical optimization.
 
 Documentation and downloads: http://ensmallen.org
 
 ensmallen provides a simple set of abstractions for writing an objective
 function to optimize. It also provides a large set of standard and cutting-edge
-optimizers that can be used for virtually any mathematical optimization task.
+optimizers that can be used for virtually any numerical optimization task.
 These include full-batch gradient descent techniques, small-batch techniques,
 gradient-free optimizers, and constrained optimization.
 
@@ -16,14 +20,21 @@ gradient-free optimizers, and constrained optimization.
 * OpenBLAS or Intel MKL or LAPACK (see Armadillo site for details)
 
 
+### Installation
+
+ensmallen can be installed with CMake 3.3 or later.
+If CMake is not already available on your system, it can be obtained from https://cmake.org
+
+If you are using an older system such as RHEL 7 or CentOS 7,
+an updated version of CMake is also available via the EPEL repository via the `cmake3` package.
+
+
 ### License
 
-Unless stated otherwise, the source code for **ensmallen**
-is licensed under the 3-clause BSD license (the "License").
-A copy of the License is included in the "LICENSE.txt" file.
-You may also obtain a copy of the License at
-http://opensource.org/licenses/BSD-3-Clause
-
+Unless stated otherwise, the source code for **ensmallen** is licensed under the
+3-clause BSD license (the "License").  A copy of the License is included in the
+"LICENSE.txt" file.  You may also obtain a copy of the License at
+http://opensource.org/licenses/BSD-3-Clause .
 
 ### Citation
 
@@ -35,6 +46,26 @@ the library.
   [ensmallen: a flexible C++ library for efficient function optimization](http://www.ensmallen.org/files/ensmallen_2018.pdf).
   Workshop on Systems for ML and Open Source Software at NIPS 2018.
 
+```
+@article{DBLP:journals/corr/abs-1810-09361,
+  author    = {Shikhar Bhardwaj and
+               Ryan R. Curtin and
+               Marcus Edel and
+               Yannis Mentekidis and
+               Conrad Sanderson},
+  title     = {ensmallen: a flexible {C++} library for efficient function optimization},
+  journal   = {CoRR},
+  volume    = {abs/1810.09361},
+  doi       = {10.5281/zenodo.2008650},
+  year      = {2018},
+  url       = {http://arxiv.org/abs/1810.09361},
+  archivePrefix = {arXiv},
+  eprint    = {1810.09361},
+  timestamp = {Wed, 31 Oct 2018 14:24:29 +0100},
+  biburl    = {https://dblp.org/rec/bib/journals/corr/abs-1810-09361},
+  bibsource = {dblp computer science bibliography, https://dblp.org}
+}
+```
 
 ### Developers and Contributors
 

doc/callbacks.md

@@ -23,7 +23,7 @@ MomentumSGD optimizer(0.01, 32, 100000, 1e-5, true, MomentumUpdate(0.5));
 optimizer.Optimize(f, coordinates, PrintLoss());
 ```
 
-</details>  
+</details>
 
 Passing multiple callbacks is just the same as passing a single callback:
 
@@ -42,7 +42,7 @@ MomentumSGD optimizer(0.01, 32, 100000, 1e-5, true, MomentumUpdate(0.5));
 optimizer.Optimize(f, coordinates, PrintLoss(), EarlyStopAtMinLoss());
 ```
 
-</details>  
+</details>
 
 It is also possible to pass a callback instantiation that allows accessing of
 internal callback parameters at a later state:
@@ -84,12 +84,20 @@ has been made.
 
  * `EarlyStopAtMinLoss()`
  * `EarlyStopAtMinLoss(`_`patience`_`)`
+ * `EarlyStopAtMinLoss(`_`func`_`)`
+ * `EarlyStopAtMinLoss(`_`func`_`,`_`patience`_`)`
 
 #### Attributes
 
 | **type** | **name** | **description** | **default** |
 |----------|----------|-----------------|-------------|
 | `size_t` | **`patience`** | The number of epochs to wait after the minimum loss has been reached. | `10` |
+| `std::function<double(const arma::mat&)>` | **`func`** | A callback to return immediate loss evaluated by the function. | |
+
+Note that for the `func` argument above, if a
+[different matrix type](#alternate-matrix-types) is desired, instead of using
+the class `EarlyStopAtMinLoss`, the class `EarlyStopAtMinLossType<MatType>`
+should be used.
 
 #### Examples:
 
@@ -104,6 +112,35 @@ RosenbrockFunction f;
 arma::mat coordinates = f.GetInitialPoint();
 optimizer.Optimize(f, coordinates, EarlyStopAtMinLoss());
 ```
+Another example of using lambda in the constructor.
+
+```c++
+// Generate random training data and labels.
+arma::mat trainingData(5, 100, arma::fill::randu);
+arma::Row<size_t> trainingLabels =
+    arma::randi<arma::Row<size_t>>(100, arma::distr_param(0, 1));
+// Generate a validation set.
+arma::mat validationData(5, 100, arma::fill::randu);
+arma::Row<size_t> validationLabels =
+    arma::randi<arma::Row<size_t>>(100, arma::distr_param(0, 1));
+
+// Create a LogisticRegressionFunction for both the training and validation data.
+LogisticRegressionFunction lrfTrain(trainingData, trainingLabels);
+LogisticRegressionFunction lrfValidation(validationData, validationLabels);
+
+// Create a callback that will terminate when the validation loss starts to
+// increase.
+EarlyStopAtMinLoss cb(
+    [&](const arma::mat& coordinates)
+    {
+      // You could also, e.g., print the validation loss here to watch it converge.
+      return lrfValidation.Evaluate(coordinates);
+    });
+    
+arma::mat coordinates = lrfTrain.GetInitialPoint();
+SMORMS3 smorms3;
+smorms3.Optimize(lrfTrain, coordinates, cb);
+```
 
 </details>
 
@@ -171,6 +208,95 @@ optimizer.Optimize(f, coordinates, ProgressBar());
 
 </details>
 
+### Report 
+
+Callback that prints a optimizer report to stdout or a specified output stream.
+
+#### Constructors
+
+ * `Report()`
+ * `Report(`_`iterationsPercentage`_`)`
+ * `Report(`_`iterationsPercentage, output`_`)`
+ * `Report(`_`iterationsPercentage, output, outputMatrixSize`_`)`
+
+#### Attributes
+
+| **type** | **name** | **description** | **default** |
+|----------|----------|-----------------|-------------|
+| `double` | **`iterationsPercentage`** | The number of iterations to report in percent, between [0, 1]. | `0.1` |
+| `std::ostream` | **`output`** | Ostream which receives output from this object. | `stdout` |
+| `size_t` | **`outputMatrixSize`** | The number of values to output for the function coordinates. | `4` |
+
+#### Examples:
+
+<details open>
+<summary>Click to collapse/expand example code.
+</summary>
+
+```c++
+AdaDelta optimizer(1.0, 1, 0.99, 1e-8, 1000, 1e-9, true);
+
+RosenbrockFunction f;
+arma::mat coordinates = f.GetInitialPoint();
+optimizer.Optimize(f, coordinates, Report(0.1));
+```
+
+<details open>
+<summary>Click to collapse/expand example output.
+</summary>
+
+```
+Optimization Report
+--------------------------------------------------------------------------------
+
+Initial Coordinates:
+  -1.2000   1.0000
+
+Final coordinates:
+  -1.0490   1.1070
+
+iter          loss          loss change   |gradient|    step size     total time
+0             24.2          0             233           1             4.27e-05
+100           8.6           15.6          104           1             0.000215
+200           5.26          3.35          48.7          1             0.000373
+300           4.49          0.767         23.4          1             0.000533
+400           4.31          0.181         11.3          1             0.000689
+500           4.27          0.0431        5.4           1             0.000846
+600           4.26          0.012         2.86          1             0.00101
+700           4.25          0.00734       2.09          1             0.00117
+800           4.24          0.00971       1.95          1             0.00132
+900           4.22          0.0146        1.91          1             0.00148
+
+--------------------------------------------------------------------------------
+
+Version:
+ensmallen:                    2.13.0 (Automatically Automated Automation)
+armadillo:                    9.900.1 (Nocturnal Misbehaviour)
+
+Function:
+Number of functions:          1
+Coordinates rows:             2
+Coordinates columns:          1
+
+Loss:
+Initial                       24.2
+Final                         4.2
+Change                        20
+
+Optimizer:
+Maximum iterations:           1000
+Reached maximum iterations:   true
+Batchsize:                    1
+Iterations:                   1000
+Number of epochs:             1001
+Initial step size:            1
+Final step size:              1
+Coordinates max. norm:        233
+Evaluate calls:               1000
+Gradient calls:               1000
+Time (in seconds):            0.00163
+```
+
 ### StoreBestCoordinates
 
 Callback that stores the model parameter after every epoch if the objective
@@ -205,14 +331,14 @@ StoreBestCoordinates<arma::mat> cb;
 optimizer.Optimize(f, coordinates, cb);
 
 std::cout << "The optimized model found by AdaDelta has the "
-      << "parameters " << cb.BestCoordinatest();
+      << "parameters " << cb.BestCoordinates();
 ```
 
 </details>
 
 ## Callback States
 
-Callbacks are called at different states during the optimization process:
+Callbacks are called at several states during the optimization process:
 
 * At the beginning and end of the optimization process.
 * After any call to `Evaluate()` and `EvaluateConstraint`.

doc/function_types.md

@@ -841,23 +841,72 @@ int main()
 
 </details>
 
+
+## Multi-objective functions
+
+A multi-objective optimizer does not return just one set of coordinates at the
+minimum of all objective functions, but instead finds a *front* or *frontier* of
+possible coordinates that are Pareto-optimal (that is, no individual objective
+function's value can be reduced without increasing at least one other
+objective function).
+
+In order to optimize a multi-objective function with ensmallen, a `std::tuple<>`
+containing multiple `ArbitraryFunctionType`s ([see here](#arbitrary-functions))
+should be passed to a multi-objective optimizer's `Optimize()` function.
+
+An example below simultaneously optimizes the generalized Rosenbrock function
+in 6 dimensions and the Wood function using [NSGA2](#nsga2).
+
+<details open>
+<summary>Click to collapse/expand example code.
+</summary>
+
+```c++
+GeneralizedRosenbrockFunction rf(6);
+WoodFunction wf;
+std::tuple<GeneralizedRosenbrockFunction, WoodFunction> objectives(rf, wf);
+
+// Create an initial point (a random point in 6 dimensions).
+arma::mat coordinates(6, 1, arma::fill::randu);
+
+// `coordinates` will be set to the coordinates on the best front that minimize the
+// sum of objective functions, and `bestFrontSum` will be the sum of all objectives
+// at that coordinate set.
+NSGA2 nsga;
+double bestFrontSum = nsga.Optimize(objectives, coordinates);
+
+// Set `bestFront` to contain all of the coordinates on the best front.
+std::vector<arma::mat> bestFront = optimizer.Front();
+}
+```
+
+</details>
+
+*Note*: all multi-objective function optimizers have both the function `Optimize()` to find the
+best front, and also the function `Front()` to return all sets of coordinates that are on the
+front.
+
+The following optimizers can be used with multi-objective functions:
+- [NSGA2](#nsga2)
+
 ## Constrained functions
 
 A constrained function is an objective function `f(x)` that is also subject to
 some constraints on `x`.  (For instance, perhaps a constraint could be that `x`
 is a positive semidefinite matrix.)  ensmallen is able to handle differentiable
 objective functions of this type---so, `f'(x)` must also be computable.  Given
-some set of constraints c_0(x), ... c_M(x), we can re-express our constrained
+some set of constraints `c_0(x)`, ..., `c_M(x)`, we can re-express our constrained
 objective function as
 
 ```
 f_C(x) = f(x) + c_0(x) + ... + c_M(x)
 ```
 
-where the constraint `c_i(x)` is `DBL_MAX` if it is not satisfied, and
-otherwise takes some real value.  For a "hard constraint", we can simply take
-`c_i(x) = 0` when it is satisfied.  But allowing `c_i(x)` to return anything
-allows us to handle "soft" constraints also.
+where the (soft) constraint `c_i(x)` is a positive value if it is not satisfied, and
+`0` if it is satisfied.  The soft constraint `c_i(x)` should take some value
+representing how far from a feasible solution `x` is.  It should be
+differentiable, since ensmallen's constrained optimizers will use the gradient
+of the constraint to find a feasible solution.
 
 In order to optimize a constrained function with ensmallen, a class
 implementing the API below is required.
@@ -880,16 +929,14 @@ class ConstrainedFunctionType
   size_t NumConstraints();
 
   // Evaluate constraint i at the parameters x.  If the constraint is
-  // unsatisfied, DBL_MAX should be returned.  If the constraint is satisfied,
-  // any real value can be returned.  The optimizer will add this value to its
-  // overall objective that it is trying to minimize.  (So, a hard constraint
-  // can just return 0 if it's satisfied.)
+  // unsatisfied, a value greater than 0 should be returned.  If the constraint
+  // is satisfied, 0 should be returned.  The optimizer will add this value to
+  // its overall objective that it is trying to minimize.
   double EvaluateConstraint(const size_t i, const arma::mat& x);
 
   // Evaluate the gradient of constraint i at the parameters x, storing the
-  // result in the given matrix g.  If this is a hard constraint you can set
-  // the gradient to 0.  If the constraint is not satisfied, it could be
-  // helpful to set the gradient in such a way that the gradient points in the
+  // result in the given matrix g.  If the constraint is not satisfied, the
+  // gradient should be set in such a way that the gradient points in the
   // direction where the constraint would be satisfied.
   void GradientConstraint(const size_t i, const arma::mat& x, arma::mat& g);
 };
@@ -1036,12 +1083,12 @@ int main()
   // use the PrimalDualSolver to solve it.
   // ens::PrimalDualSolver could be replaced with ens::LRSDP or other ensmallen
   // SDP solvers.
-  PrimalDualSolver<SDP<arma::sp_mat>> solver(sdp);
+  PrimalDualSolver solver;
   arma::mat X, Z;
   arma::vec ysparse, ydense;
   // ysparse, ydense, and Z hold the primal and dual variables found during the
   // optimization.
-  const double obj = solver.Optimize(X, ysparse, ydense, Z);
+  const double obj = solver.Optimize(sdp, X, ysparse, ydense, Z);
 
   std::cout << "SDP optimized with objective " << obj << "." << std::endl;
 }
@@ -1098,14 +1145,14 @@ class SquaredFunction
 
 int main()
 {
-  // The minimum is at x = [0 0 0].  Our initial point is chosen to be 
+  // The minimum is at x = [0 0 0].  Our initial point is chosen to be
   // [1.0, -1.0, 1.0].
   arma::fmat x("1.0 -1.0 1.0");
 
   // Create simulated annealing optimizer with default options.
   // The ens::SA<> type can be replaced with any suitable ensmallen optimizer
   // that is able to handle arbitrary functions.
-  ens::L_BFGS<> optimizer;
+  ens::L_BFGS optimizer;
   SquaredFunction f; // Create function to be optimized.
   optimizer.Optimize(f, x); // The optimizer will infer arma::fmat!
 

doc/optimizers.md

@@ -85,7 +85,7 @@ gradients.
 | `double` | **`stepSize`** | Step size for each iteration. | `1.0` |
 | `size_t` | **`batchSize`**| Number of points to process in one step. | `32` |
 | `double` | **`rho`** | Smoothing constant. Corresponding to fraction of gradient to keep at each time step. | `0.95` |
-| `double` | **`epsilon`** | Value used to initialise the mean squared gradient parameter. | `1e-6` |
+| `double` | **`epsilon`** | Value used to initialize the mean squared gradient parameter. | `1e-6` |
 | `size_t` | **`maxIterations`** | Maximum number of iterations allowed (0 means no limit). | `100000` |
 | `double` | **`tolerance`** | Maximum absolute tolerance to terminate algorithm. | `1e-5` |
 | `bool` | **`shuffle`** | If true, the function order is shuffled; otherwise, each function is visited in linear order. | `true` |
@@ -141,7 +141,7 @@ parameters.
 |----------|----------|-----------------|-------------|
 | `double` | **`stepSize`** | Step size for each iteration. | `0.01` |
 | `size_t` | **`batchSize`** | Number of points to process in one step. | `32` |
-| `double` | **`epsilon`** | Value used to initialise the mean squared gradient parameter. | `1e-8` |
+| `double` | **`epsilon`** | Value used to initialize the mean squared gradient parameter. | `1e-8` |
 | `size_t` | **`maxIterations`** | Maximum number of iterations allowed (0 means no limit). | `100000` |
 | `double` | **`tolerance`** | Maximum absolute tolerance to terminate algorithm. | `tolerance` |
 | `bool` | **`shuffle`** | If true, the function order is shuffled; otherwise, each function is visited in linear order. | `true` |
@@ -179,7 +179,7 @@ optimizer.Optimize(f, coordinates);
 
 *An optimizer for [differentiable separable functions](#differentiable-separable-functions).*
 
-Adam is an an algorithm for first-order gradient-based optimization of
+Adam is an algorithm for first-order gradient-based optimization of
 stochastic objective functions, based on adaptive estimates of lower-order
 moments.
 
@@ -424,7 +424,7 @@ optimizer uses [L-BFGS](#l-bfgs).
 
 #### Constructors
 
- * `AugLagrangian(`_`maxIterations, penaltyThresholdFactor sigmaUpdateFactor`_`)`
+ * `AugLagrangian(`_`maxIterations, penaltyThresholdFactor, sigmaUpdateFactor`_`)`
 
 #### Attributes
 
@@ -554,11 +554,11 @@ RosenbrockFunction f;
 arma::mat coordinates = f.GetInitialPoint();
 
 // Big-Batch SGD with the adaptive stepsize policy.
-BBS_BB optimizer(batchSize, 0.01, 0.1, 8000, 1e-4);
+BBS_BB optimizer(10, 0.01, 0.1, 8000, 1e-4);
 optimizer.Optimize(f, coordinates);
 
 // Big-Batch SGD with backtracking line search.
-BBS_Armijo optimizer2(batchSize, 0.01, 0.1, 8000, 1e-4);
+BBS_Armijo optimizer2(10, 0.01, 0.1, 8000, 1e-4);
 optimizer2.Optimize(f, coordinates);
 ```
 
@@ -629,11 +629,11 @@ RosenbrockFunction f;
 arma::mat coordinates = f.GetInitialPoint();
 
 // CMAES with the FullSelection policy.
-CMAES<> optimizer(0, -1, 1, 32, 200, 0.1e-4);
+CMAES<> optimizer(0, -1, 1, 32, 200, 1e-4);
 optimizer.Optimize(f, coordinates);
 
 // CMAES with the RandomSelection policy.
-ApproxCMAES<> approxOptimizer(batchSize, 0.01, 0.1, 8000, 1e-4);
+ApproxCMAES<> approxOptimizer(0, -1, 1. 32, 200, 1e-4);
 approxOptimizer.Optimize(f, coordinates);
 ```
 
@@ -766,7 +766,7 @@ Eve is a stochastic gradient based optimization method with locally and globally
 | `double` | **`beta2`** | Exponential decay rate for the weighted infinity norm estimates. | `0.999` |
 | `double` | **`beta3`** | Exponential decay rate for relative change. | `0.999` |
 | `double` | **`epsilon`** | Value used to initialize the mean squared gradient parameter. | `1e-8` |
-| `double` | **`clip`** | Clipping range to avoid extreme valus. | `10` |
+| `double` | **`clip`** | Clipping range to avoid extreme values. | `10` |
 | `size_t` | **`max_iterations`** | Maximum number of iterations allowed (0 means no limit). | `100000` |
 | `double` | **`tolerance`** | Maximum absolute tolerance to terminate algorithm. | `1e-5` |
 | `bool` | **`shuffle`** | If true, the function order is shuffled; otherwise, each function is visited in linear order. | `true` |
@@ -1167,7 +1167,7 @@ proximalOptimizer.Optimize(f, coordinates);
 
 *An optimizer for [differentiable functions](#differentiable-functions)*
 
-L-BFGS is an optimization algorithm in the family of quasi-Newton methods that approximates the Broyden-Fletcher-Goldfarb-Shanno (BFGS) algorithm using a limited amount of computer memory.  
+L-BFGS is an optimization algorithm in the family of quasi-Newton methods that approximates the Broyden-Fletcher-Goldfarb-Shanno (BFGS) algorithm using a limited amount of computer memory.
 
 #### Constructors
 
@@ -1563,6 +1563,73 @@ optimizer.Optimize(f, coordinates);
  * [SGD in Wikipedia](https://en.wikipedia.org/wiki/Stochastic_gradient_descent)
  * [Differentiable separable functions](#differentiable-separable-functions)
 
+## NSGA2
+
+*An optimizer for arbitrary multi-objective functions.*
+
+NSGA2 (Non-dominated Sorting Genetic Algorithm - II) is a multi-objective
+optimization algorithm. The algorithm works by generating a candidate population
+from a fixed starting point. At each stage of optimization, a new population of
+children is generated. This new population along with its predecessor is sorted
+using non-domination as the metric. Following this, the population is further
+segregated into fronts. A new population is generated from these fronts having
+size equal to that of the starting population.
+
+#### Constructors
+
+ * `NSGA2()`
+ * `NSGA2(`_`populationSize, maxGenerations, crossoverProb, mutationProb, mutationStrength, epsilon, lowerBound, upperBound`_`)`
+
+#### Attributes
+
+| **type** | **name** | **description** | **default** |
+|----------|----------|-----------------|-------------|
+| `size_t` | **`populationSize`** | The number of candidates in the population. This should be at least 4 in size and a multiple of 4. | `100` |
+| `size_t` | **`maxGenerations`** | The maximum number of generations allowed for NSGA2. | `2000` |
+| `double` | **`crossoverProb`** | Probability that a crossover will occur. | `0.6` |
+| `double` | **`mutationProb`** | Probability that a weight will get mutated. | `0.3` |
+| `double` | **`mutationStrength`** | The range of mutation noise to be added. This range is between 0 and mutationStrength. | `0.001` |
+| `double` | **`epsilon`** | The value used internally to evaluate approximate equality in crowding distance based sorting. | `1e-6` |
+| `double`, `arma::vec` | **`lowerBound`** | Lower bound of the coordinates on the coordinates of the whole population during the search process. | `0` |
+| `double`, `arma::vec` | **`upperBound`** | Lower bound of the coordinates on the coordinates of the whole population during the search process. | `1` |
+
+Note that the parameters `lowerBound` and `upperBound` are overloaded. Data types of `double` or `arma::mat` may be used. If they are initialized as single values of `double`, then the same value of the bound applies to all the axes, resulting in an initialization following a uniform distribution in a hypercube. If they are initialized as matrices of `arma::mat`, then the value of `lowerBound[i]` applies to axis `[i]`; similarly, for values in `upperBound`. This results in an initialization following a uniform distribution in a hyperrectangle within the specified bounds.
+
+Attributes of the optimizer may also be changed via the member methods
+`PopulationSize()`, `MaxGenerations()`, `CrossoverRate()`, `MutationProbability()`, `MutationStrength()`, `Epsilon()`, `LowerBound()` and `UpperBound()`.
+
+#### Examples:
+
+<details open>
+<summary>Click to collapse/expand example code.
+</summary>
+
+```c++
+SchafferFunctionN1<arma::mat> SCH;
+arma::vec lowerBound("-1000 -1000");
+arma::vec upperBound("1000 1000");
+NSGA2 opt(20, 5000, 0.5, 0.5, 1e-3, 1e-6, lowerBound, upperBound);
+
+typedef decltype(SCH.objectiveA) ObjectiveTypeA;
+typedef decltype(SCH.objectiveB) ObjectiveTypeB;
+
+arma::mat coords = SCH.GetInitialPoint();
+std::tuple<ObjectiveTypeA, ObjectiveTypeB> objectives = SCH.GetObjectives();
+
+// obj will contain the minimum sum of objectiveA and objectiveB found on the best front.
+double obj = opt.Optimize(objectives, coords);
+// Now obtain the best front.
+std::vector<arma::mat> bestFront = opt.Front();
+```
+
+</details>
+
+#### See also:
+
+ * [NSGA-II Algorithm](https://www.iitk.ac.in/kangal/Deb_NSGA-II.pdf)
+ * [Multi-objective Functions in Wikipedia](https://en.wikipedia.org/wiki/Test_functions_for_optimization#Test_functions_for_multi-objective_optimization)
+  * [Multi-objective functions](#multi-objective-functions)
+
 ## OptimisticAdam
 
 *An optimizer for [differentiable separable functions](#differentiable-separable-functions).*
@@ -1571,7 +1638,7 @@ OptimisticAdam is an optimizer which implements the Optimistic Adam algorithm
 which uses Optmistic Mirror Descent with the Adam Optimizer.  It addresses the
 problem of limit cycling while training GANs (generative adversarial networks).
 It uses OMD to achieve faster regret rates in solving the zero sum game of
-training a GAN. It consistently achieves a smaller KL divergnce with~ respect to
+training a GAN. It consistently achieves a smaller KL divergence with~ respect to
 the true underlying data distribution.  The implementation here can be used with
 any differentiable separable function, not just GAN training.
 
@@ -1914,7 +1981,7 @@ optimizer.Optimize(f, coordinates);
 
 QHAdam is an optimizer that uses quasi-hyperbolic descent with the Adam
 optimizer.  This replaces the moment estimators of Adam with quasi-hyperbolic
-terms, and different values of the `v1` and `v2` parameters are equivalent to
+terms, and various values of the `v1` and `v2` parameters are equivalent to
 the following other optimizers:
 
  * When `v1 = v2 = 1`, `QHAdam` is equivalent to `Adam`.
@@ -1992,7 +2059,7 @@ RMSProp utilizes the magnitude of recent gradients to normalize the gradients.
 | `double` | **`stepSize`** | Step size for each iteration. | `0.01` |
 | `size_t` | **`batchSize`** | Number of points to process in each step. | `32` |
 | `double` | **`alpha`** | Smoothing constant, similar to that used in AdaDelta and momentum methods. | `0.99` |
-| `double` | **`epsilon`** | Value used to initialise the mean squared gradient parameter. | `1e-8` |
+| `double` | **`epsilon`** | Value used to initialize the mean squared gradient parameter. | `1e-8` |
 | `size_t` | **`maxIterations`** | Maximum number of iterations allowed (0 means no limit). | `100000` |
 | `double` | **`tolerance`** | Maximum absolute tolerance to terminate algorithm. |
 | `bool` | **`shuffle`** | If true, the function order is shuffled; otherwise, each function is visited in linear order. | `true` |
@@ -2068,7 +2135,7 @@ shorter type `SA<>` may be used instead of the equivalent
 
 | **type** | **name** | **description** | **default** |
 |----------|----------|-----------------|-------------|
-| `CoolingScheduleType` | **`coolingSchedule`** | Instantiated cooling schedule (default ExponentialSchedule). | **n/a** |
+| `CoolingScheduleType` | **`coolingSchedule`** | Instantiated cooling schedule (default ExponentialSchedule). | **CoolingScheduleType()** |
 | `size_t` | **`maxIterations`** | Maximum number of iterations allowed (0 indicates no limit). | `1000000` |
 | `double` | **`initT`** | Initial temperature. | `10000.0` |
 | `size_t` | **`initMoves`** | Number of initial iterations without changing temperature. | `1000` |
@@ -2509,7 +2576,7 @@ optimizer.Optimize(f, coordinates);
 *An optimizer for [differentiable separable functions](#differentiable-separable-functions).*
 
 SMORMS3 is a hybrid of RMSprop, which is trying to estimate a safe and optimal
-distance based on curvature or perhaps just normalizing the stepsize in the
+distance based on curvature or perhaps just normalizing the step-size in the
 parameter space.
 
 #### Constructors
@@ -2525,7 +2592,7 @@ parameter space.
 |----------|----------|-----------------|-------------|
 | `double` | **`stepSize`** | Step size for each iteration. | `0.001` |
 | `size_t` | **`batchSize`** | Number of points to process at each step. | `32` |
-| `double` | **`epsilon`** | Value used to initialise the mean squared gradient parameter. | `1e-16` |
+| `double` | **`epsilon`** | Value used to initialize the mean squared gradient parameter. | `1e-16` |
 | `size_t` | **`maxIterations`** | Maximum number of iterations allowed (0 means no limit). | `100000` |
 | `double` | **`tolerance`** | Maximum absolute tolerance to terminate algorithm. | `1e-5` |
 | `bool` | **`shuffle`** | If true, the mini-batch order is shuffled; otherwise, each mini-batch is visited in linear order. | `true` |
@@ -2732,7 +2799,7 @@ the projection of Adam steps on the gradient subspace.
 | `size_t` | **`batchSize`** | Number of points to process at each step. | `32` |
 | `double` | **`beta1`** | Exponential decay rate for the first moment estimates. | `0.9` |
 | `double` | **`beta2`** | Exponential decay rate for the weighted infinity norm estimates. | `0.999` |
-| `double` | **`epsilon`** | Value used to initialise the mean squared gradient parameter. | `1e-16` |
+| `double` | **`epsilon`** | Value used to initialize the mean squared gradient parameter. | `1e-16` |
 | `size_t` | **`maxIterations`** | Maximum number of iterations allowed (0 means no limit). | `100000` |
 | `double` | **`tolerance`** | Maximum absolute tolerance to terminate algorithm. | `1e-5` |
 | `bool` | **`shuffle`** | If true, the mini-batch order is shuffled; otherwise, each mini-batch is visited in linear order. | `true` |

include/ensmallen.hpp

@@ -33,20 +33,24 @@
   #error "need Armadillo version 8.400 or later"
 #endif
 
-#include <cmath>
-#include <cstdlib>
-#include <cstdio>
-#include <cstring>
 #include <cctype>
-#include <climits>
 #include <cfloat>
+#include <climits>
+#include <cmath>
 #include <cstdint>
+#include <cstdio>
+#include <cstdlib>
+#include <cstring>
+#include <iostream>
+#include <map>
+#include <set>
+#include <limits>
+#include <sstream>
 #include <stdexcept>
+#include <string>
 #include <tuple>
 #include <utility>
-#include <iostream>
-#include <string>
-#include <sstream>
+#include <vector>
 
 // On Visual Studio, disable C4519 (default arguments for function templates)
 // since it's by default an error, which doesn't even make any sense because
@@ -67,6 +71,7 @@
 #include "ensmallen_bits/callbacks/early_stop_at_min_loss.hpp"
 #include "ensmallen_bits/callbacks/print_loss.hpp"
 #include "ensmallen_bits/callbacks/progress_bar.hpp"
+#include "ensmallen_bits/callbacks/report.hpp"
 #include "ensmallen_bits/callbacks/store_best_coordinates.hpp"
 #include "ensmallen_bits/callbacks/timer_stop.hpp"
 
@@ -94,6 +99,7 @@
 #include "ensmallen_bits/katyusha/katyusha.hpp"
 #include "ensmallen_bits/lbfgs/lbfgs.hpp"
 #include "ensmallen_bits/lookahead/lookahead.hpp"
+#include "ensmallen_bits/nsga2/nsga2.hpp"
 #include "ensmallen_bits/padam/padam.hpp"
 #include "ensmallen_bits/parallel_sgd/parallel_sgd.hpp"
 #include "ensmallen_bits/pso/pso.hpp"

include/ensmallen_bits/ada_bound/ada_bound.hpp

@@ -113,8 +113,9 @@ class AdaBoundType
            MatType& iterate,
            CallbackTypes&&... callbacks)
   {
-    return optimizer.template Optimize<DecomposableFunctionType, MatType, GradType,
-        CallbackTypes...>(function, iterate, callbacks...);
+    return optimizer.template Optimize<DecomposableFunctionType, MatType,
+        GradType, CallbackTypes...>(function, iterate,
+        std::forward<CallbackTypes>(callbacks)...);
   }
 
   //! Forward the MatType as GradType.

include/ensmallen_bits/ada_delta/ada_delta.hpp

@@ -105,7 +105,8 @@ class AdaDelta
            CallbackTypes&&... callbacks)
   {
     return optimizer.Optimize<SeparableFunctionType, MatType, GradType,
-        CallbackTypes...>(function, iterate, callbacks...);
+        CallbackTypes...>(function, iterate,
+        std::forward<CallbackTypes>(callbacks)...);
   }
 
   //! Forward the MatType as GradType.

include/ensmallen_bits/ada_grad/ada_grad.hpp

@@ -101,7 +101,8 @@ class AdaGrad
            CallbackTypes&&... callbacks)
   {
     return optimizer.Optimize<SeparableFunctionType, MatType, GradType,
-        CallbackTypes...>(function, iterate, callbacks...);
+        CallbackTypes...>(function, iterate,
+        std::forward<CallbackTypes>(callbacks)...);
   }
 
   //! Forward the MatType as GradType.

include/ensmallen_bits/adam/adam.hpp

@@ -128,7 +128,7 @@ class AdamType
   {
     return optimizer.template Optimize<
         SeparableFunctionType, MatType, GradType, CallbackTypes...>(
-        function, iterate, callbacks...);
+        function, iterate, std::forward<CallbackTypes>(callbacks)...);
   }
 
   //! Forward the MatType as GradType.

include/ensmallen_bits/aug_lagrangian/aug_lagrangian_impl.hpp

@@ -27,7 +27,8 @@ inline AugLagrangian::AugLagrangian(const size_t maxIterations,
     penaltyThresholdFactor(penaltyThresholdFactor),
     sigmaUpdateFactor(sigmaUpdateFactor),
     lbfgs(lbfgs),
-    terminate(false)
+    terminate(false),
+    sigma(0.0)
 {
 }
 
@@ -107,6 +108,10 @@ AugLagrangian::Optimize(
   // Track the last objective to compare for convergence.
   ElemType lastObjective = function.Evaluate(coordinates);
 
+  // Convergence tolerance---depends on the epsilon of the type we are using for
+  // optimization.
+  ElemType tolerance = 1e3 * std::numeric_limits<ElemType>::epsilon();
+
   // Then, calculate the current penalty.
   ElemType penalty = 0;
   for (size_t i = 0; i < function.NumConstraints(); i++)
@@ -134,6 +139,7 @@ AugLagrangian::Optimize(
     if (!lbfgs.Optimize(augfunc, coordinates, callbacks...))
       Info << "L-BFGS reported an error during optimization."
           << std::endl;
+    Info << "Done with L-BFGS: " << coordinates << "\n";
 
     const ElemType objective = function.Evaluate(coordinates);
 
@@ -142,7 +148,7 @@ AugLagrangian::Optimize(
 
     // Check if we are done with the entire optimization (the threshold we are
     // comparing with is arbitrary).
-    if (std::abs(lastObjective - objective) < 1e-10 &&
+    if (std::abs(lastObjective - objective) < tolerance &&
         augfunc.Sigma() > 500000)
     {
       lambda = std::move(augfunc.Lambda());
@@ -196,6 +202,13 @@ AugLagrangian::Optimize(
       // We multiply sigma by a constant value.
       augfunc.Sigma() *= sigmaUpdateFactor;
       Info << "Updated sigma to " << augfunc.Sigma() << "." << std::endl;
+      if (augfunc.Sigma() >= std::numeric_limits<ElemType>::max() / 2.0)
+      {
+        Warn << "AugLagrangian::Optimize(): sigma too large for element type; "
+            << "terminating." << std::endl;
+        Callback::EndOptimization(*this, function, coordinates, callbacks...);
+        return false;
+      }
     }
 
     terminate |= Callback::StepTaken(*this, function, coordinates,

include/ensmallen_bits/bigbatch_sgd/bigbatch_sgd.hpp

@@ -101,6 +101,12 @@ class BigBatchSGD
               const double tolerance = 1e-5,
               const bool shuffle = true,
               const bool exactObjective = false);
+
+  /**
+   * Clean any memory associated with the BigBatchSGD object.
+   */
+  ~BigBatchSGD();
+
   /**
    * Optimize the given function using big-batch SGD.  The given starting point
    * will be modified to store the finishing point of the algorithm, and the

include/ensmallen_bits/bigbatch_sgd/bigbatch_sgd_impl.hpp

@@ -38,6 +38,12 @@ BigBatchSGD<UpdatePolicyType>::BigBatchSGD(
     updatePolicy(UpdatePolicyType())
 { /* Nothing to do. */ }
 
+template<typename UpdatePolicyType>
+BigBatchSGD<UpdatePolicyType>::~BigBatchSGD()
+{
+  instUpdatePolicy.Clean();
+}
+
 //! Optimize the function (minimize).
 template<typename UpdatePolicyType>
 template<typename SeparableFunctionType,

include/ensmallen_bits/callbacks/callbacks.hpp

@@ -257,6 +257,7 @@ class Callback
   {
     // This will return immediately once a callback returns true.
     bool result = false;
+    (void)(objective);  // prevent spurious compiler warnings
     (void)std::initializer_list<bool>{ result =
         result || Callback::EvaluateFunction(callbacks, optimizer, function,
         coordinates, objective)... };
@@ -332,6 +333,8 @@ class Callback
   {
     // This will return immediately once a callback returns true.
     bool result = false;
+    (void)(constraint);  // prevent spurious compiler warnings
+    (void)(constraintValue);
     (void)std::initializer_list<bool>{ result =
         result || Callback::EvaluateConstraintFunction(callbacks, optimizer,
             function, coordinates, constraint, constraintValue)... };
@@ -478,6 +481,7 @@ class Callback
   {
     // This will return immediately once a callback returns true.
     bool result = false;
+    (void)(constraint);  // prevent spurious compiler warnings
     (void)std::initializer_list<bool>{ result =
         result || Callback::GradientConstraintFunction(callbacks, optimizer,
         function, coordinates, constraint, gradient)... };
@@ -509,6 +513,7 @@ class Callback
   {
     // This will return immediately once a callback returns true.
     bool result = false;
+    (void)(objective);  // prevent spurious compiler warnings
     (void)std::initializer_list<bool>{ result =
         result || Callback::EvaluateFunction(callbacks, optimizer, function,
         coordinates, objective)... };
@@ -584,6 +589,8 @@ class Callback
   {
     // This will return immediately once a callback returns true.
     bool result = false;
+    (void)(epoch);  // prevent spurious compiler warnings
+    (void)(objective);
     (void)std::initializer_list<bool>{ result =
         result || Callback::BeginEpochFunction(callbacks, optimizer, function,
         coordinates, epoch, objective)... };
@@ -672,6 +679,8 @@ class Callback
   {
     // This will return immediately once a callback returns true.
     bool result = false;
+    (void)(epoch);  // prevent spurious compiler warnings
+    (void)(objective);
     (void)std::initializer_list<bool>{ result =
         result || Callback::EndEpochFunction(callbacks, optimizer, function,
         coordinates, epoch, objective)... };

include/ensmallen_bits/callbacks/early_stop_at_min_loss.hpp

@@ -1,6 +1,7 @@
 /**
  * @file early_stop_at_min_loss.hpp
  * @author Marcus Edel
+ * @author Omar Shrit
  *
  * Implementation of the early stop at minimum loss callback function.
  *
@@ -12,13 +13,16 @@
 #ifndef ENSMALLEN_CALLBACKS_EARLY_STOP_AT_MIN_LOSS_HPP
 #define ENSMALLEN_CALLBACKS_EARLY_STOP_AT_MIN_LOSS_HPP
 
+#include <functional>
+
 namespace ens {
 
 /**
  * Early stopping to terminate the optimization process early if the loss stops
  * decreasing.
  */
-class EarlyStopAtMinLoss
+template<typename MatType = arma::mat>
+class EarlyStopAtMinLossType
 {
  public:
   /**
@@ -28,12 +32,33 @@ class EarlyStopAtMinLoss
    * @param patienceIn The number of epochs to wait after the minimum loss has
    *    been reached or no improvement has been made (Default: 10).
    */
-  EarlyStopAtMinLoss(const size_t patienceIn = 10) :
-      patience(patienceIn),
+  EarlyStopAtMinLossType<MatType>(const size_t patienceIn = 10) :
+      callbackUsed(false), 
+      patience(patienceIn), 
       bestObjective(std::numeric_limits<double>::max()),
       steps(0)
   { /* Nothing to do here */ }
 
+  /**
+   * Set up the early stop at min loss class, which keeps track of the minimum
+   * loss and stops the optimization process if the loss stops decreasing.
+   *
+   * @param func, callback to return immediate loss evaluated by the function
+   * @param patienceIn The number of epochs to wait after the minimum loss has
+   *    been reached or no improvement has been made (Default: 10).
+   */
+  EarlyStopAtMinLossType<MatType>(
+      std::function<double(const MatType&)> func,
+      const size_t patienceIn = 10)
+    : callbackUsed(true), 
+      patience(patienceIn), 
+      bestObjective(std::numeric_limits<double>::max()),
+      steps(0), 
+      localFunc(func)
+  {
+    // Nothing to do here
+  }
+
   /**
    * Callback function called at the end of a pass over the data.
    *
@@ -43,13 +68,18 @@ class EarlyStopAtMinLoss
    * @param epoch The index of the current epoch.
    * @param objective Objective value of the current point.
    */
-  template<typename OptimizerType, typename FunctionType, typename MatType>
+  template<typename OptimizerType, typename FunctionType>
   bool EndEpoch(OptimizerType& /* optimizer */,
                 FunctionType& /* function */,
-                const MatType& /* coordinates */,
+                const MatType& coordinates,
                 const size_t /* epoch */,
-                const double objective)
+                double objective)
   {
+    if (callbackUsed)
+    {
+      objective = localFunc(coordinates);
+    } 
+
     if (objective < bestObjective)
     {
       steps = 0;
@@ -68,6 +98,9 @@ class EarlyStopAtMinLoss
   }
 
  private:
+  //! False if the first constructor is called, true if the user passed a lambda. 
+  bool callbackUsed;
+
   //! The number of epochs to wait before terminating the optimization process.
   size_t patience;
 
@@ -76,8 +109,19 @@ class EarlyStopAtMinLoss
 
   //! Locally-stored number of steps since the loss improved.
   size_t steps;
+
+  //! Function to call at the end of the epoch.
+  std::function<double(const MatType&)> localFunc;
 };
 
+/*
+ * Note that the using definition is temporary, this definition should
+ * be removed when releasing ensmallen 3.0
+ * The renaming of the class is only to avoid a major version bump
+ * because if the template type added to this class
+ */
+using EarlyStopAtMinLoss = EarlyStopAtMinLossType<arma::mat>;
+
 } // namespace ens
 
 #endif

include/ensmallen_bits/callbacks/progress_bar.hpp

@@ -81,17 +81,9 @@ class ProgressBar
     if (function.NumFunctions() % optimizer.BatchSize() > 0)
       epochSize++;
 
-    if (!optimizer.MaxIterations())
-    {
-      Warn << "Maximum number of iterations not defined (no limit),"
-           << " no progress bar shown." << std::endl;
-    }
-    else
-    {
-      epochs = optimizer.MaxIterations() / epochSize;
-      if (optimizer.MaxIterations() % epochSize > 0)
-        epochs++;
-    }
+    epochs = optimizer.MaxIterations() / function.NumFunctions();
+    if (optimizer.MaxIterations() % function.NumFunctions() > 0)
+      epochs++;
 
     stepTimer.tic();
   }
@@ -138,8 +130,12 @@ class ProgressBar
   {
     if (newEpoch)
     {
-      output << "Epoch " << epoch << "/" << epochs << "\n";
-      output.flush();
+      output << "Epoch " << epoch;
+      if (epochs > 0)
+      {
+        output << "/" << epochs;
+      }
+      output << '\n';
       newEpoch = false;
     }
 
@@ -161,8 +157,8 @@ class ProgressBar
       }
     }
 
-    output << "] " << progress << "% - ETA: " << (size_t) stepTimer.toc() *
-        (epochSize - step + 1) % 60 << "s - loss: " <<
+    output << "] " << progress << "% - ETA: " << (size_t) (stepTimer.toc() *
+        (epochSize - step + 1)) % 60 << "s - loss: " <<
         objective / (double) step <<  "\r";
     output.flush();
 

include/ensmallen_bits/callbacks/report.hpp

@@ -0,0 +1,590 @@
+/**
+ * @file report.hpp
+ * @author Marcus Edel
+ *
+ * Implementation of a simple report callback function.
+ *
+ * ensmallen is free software; you may redistribute it and/or modify it under
+ * the terms of the 3-clause BSD license.  You should have received a copy of
+ * the 3-clause BSD license along with ensmallen.  If not, see
+ * http://www.opensource.org/licenses/BSD-3-Clause for more information.
+ */
+#ifndef ENSMALLEN_CALLBACKS_REPORT_HPP
+#define ENSMALLEN_CALLBACKS_REPORT_HPP
+
+#include <ensmallen_bits/function.hpp>
+#include <iomanip>
+
+namespace ens {
+
+/**
+ * A simple optimization report.
+ */
+class Report
+{
+ public:
+  /**
+   * Set up the report callback class with the given output stream.
+   *
+   * @param iterationsPercentageIn The number of iterations to report in
+   *     percent, between [0, 1]).
+   * @param outputIn Ostream which receives output from this object.
+   * @param outputMatrixSizeIn The number of values to output for the function
+   *     coordinates.
+   */
+  Report(const double iterationsPercentageIn = 0.1,
+         std::ostream& outputIn = arma::get_cout_stream(),
+         const size_t outputMatrixSizeIn = 4) :
+      iterationsPercentage(iterationsPercentageIn),
+      output(outputIn),
+      outputMatrixSize(outputMatrixSizeIn),
+      objective(0),
+      gradientNorm(0),
+      hasGradient(false),
+      hasEndEpoch(false),
+      gradientCalls(0),
+      evaluateCalls(0),
+      epochCalls(0)
+  { /* Nothing to do here. */ }
+
+  /**
+   * Callback function called at the begin of the optimization process.
+   *
+   * @param optimizer The optimizer used to update the function.
+   * @param function Function to optimize.
+   * @param coordinates Starting point.
+   */
+  template<typename OptimizerType, typename FunctionType, typename MatType>
+  void BeginOptimization(OptimizerType& /* optimizer */,
+                         FunctionType& /* function */,
+                         MatType& coordinates)
+  {
+    initialCoordinates = coordinates;
+    optimizationTimer.tic();
+  }
+
+  /**
+   * Callback function called at the begin of the optimization process.
+   *
+   * @param optimizer The optimizer used to update the function.
+   * @param function Function to optimize.
+   * @param coordinates Starting point.
+   */
+  template<typename OptimizerType, typename FunctionType, typename MatType>
+  void EndOptimization(OptimizerType& optimizer,
+                       FunctionType& function,
+                       MatType& coordinates)
+  {
+    output << "Optimization Report" << std::endl;
+    output << std::string(80, '-') << std::endl << std::endl;
+
+    std::streamsize streamPrecision = output.precision(4);
+
+    if (coordinates.n_rows > outputMatrixSize ||
+        coordinates.n_cols > outputMatrixSize)
+    {
+      output << "Initial coordinates: " << std::endl;
+      TruncatePrint(initialCoordinates, outputMatrixSize);
+      output << std::endl << "Final coordinates: " << std::endl;
+      TruncatePrint(coordinates, outputMatrixSize);
+    }
+    else
+    {
+      output << "Initial Coordinates:" << std::endl << initialCoordinates.t();
+      output << std::endl << "Final coordinates:" << std::endl
+          << coordinates.t() << std::endl;
+    }
+
+    PrettyPrintElement("iter");
+    PrettyPrintElement("loss");
+    PrettyPrintElement("loss change");
+
+    if (hasGradient)
+      PrettyPrintElement("|gradient|");
+
+    if (!stepsizes.empty())
+      PrettyPrintElement("step size");
+
+    PrettyPrintElement("total time");
+    output << std::endl;
+
+    size_t iterationStep = objectives.size() / (iterationsPercentage * 100);
+    if (iterationStep <= 0)
+      iterationStep = 1;
+
+    for (size_t i = 0; i < objectives.size(); i += iterationStep)
+    {
+      PrettyPrintElement(i);
+      PrettyPrintElement(objectives[i]);
+      PrettyPrintElement(
+          i > 0 ? objectives[i - iterationStep] - objectives[i] : 0);
+
+      if (hasGradient)
+        PrettyPrintElement(gradientsNorm[i]);
+
+      if (!stepsizes.empty())
+        PrettyPrintElement(stepsizes[i]);
+
+      PrettyPrintElement(timings[i]);
+      output << std::endl;
+    }
+
+    output << std::endl << std::string(80, '-') << std::endl << std::endl;
+    output << "Version:" << std::endl;
+    PrettyPrintElement("ensmallen:", 30);
+    output << ens::version::as_string() << std::endl;
+    PrettyPrintElement("armadillo:", 30);
+    output << arma::arma_version::as_string() << std::endl << std::endl;
+
+    output << "Function:" << std::endl;
+    std::stringstream functionStream;
+
+    PrintNumFunctions(function, functionStream);
+    if (functionStream.rdbuf()->in_avail() > 0)
+      output << functionStream.str();
+
+    PrettyPrintElement("Coordinates rows:", 30);
+    output << coordinates.n_rows << std::endl;
+    PrettyPrintElement("Coordinates columns:", 30);
+    output << coordinates.n_cols << std::endl;
+    output << std::endl;
+
+    output << "Loss:" << std::endl;
+    PrettyPrintElement("Initial", 30);
+    output << objectives[0] << std::endl;
+    PrettyPrintElement("Final", 30);
+    output << objectives[objectives.size() - 1] << std::endl;
+    PrettyPrintElement("Change", 30);
+    output << objectives[0] - objectives[objectives.size() - 1] << std::endl;
+
+    output << std::endl << "Optimizer:" << std::endl;
+    std::stringstream optimizerStream;
+
+    PrintMaxIterations(optimizer, optimizerStream);
+    PrintBatchSize(optimizer, optimizerStream);
+    if (functionStream.rdbuf()->in_avail() > 0)
+      output << optimizerStream.str();
+
+    PrettyPrintElement("Iterations:", 30);
+    output << objectives.size() << std::endl;
+
+    if (epochCalls > 0)
+    {
+      PrettyPrintElement("Number of epochs:", 30);
+      output << epochCalls << std::endl;
+    }
+
+    if (!stepsizes.empty())
+    {
+      PrettyPrintElement("Initial step size:", 30);
+      output << stepsizes.front() << std::endl;
+
+      PrettyPrintElement("Final step size:", 30);
+      output << stepsizes.back() << std::endl;
+    }
+
+    if (hasGradient)
+    {
+      PrettyPrintElement("Coordinates max. norm:", 30);
+      output << *std::max_element(std::begin(gradientsNorm),
+          std::end(gradientsNorm)) << std::endl;
+    }
+
+    PrettyPrintElement("Evaluate calls:", 30);
+    output << evaluateCalls << std::endl;
+
+    if (hasGradient)
+    {
+      PrettyPrintElement("Gradient calls:", 30);
+      output << gradientCalls << std::endl;
+    }
+
+    PrettyPrintElement("Time (in seconds):", 30);
+    output << timings[timings.size() - 1] << std::endl;
+
+    // Restore precision.
+    output.precision(streamPrecision);
+  }
+
+  /**
+   * Callback function called at the beginning of a pass over the data.
+   *
+   * @param optimizer The optimizer used to update the function.
+   * @param function Function to optimize.
+   * @param coordinates Starting point.
+   * @param epoch The index of the current epoch.
+   * @param objective Objective value of the current point.
+   */
+  template<typename OptimizerType, typename FunctionType, typename MatType>
+  void BeginEpoch(OptimizerType& /* optimizer */,
+                  FunctionType& /* function */,
+                  const MatType& /* coordinates */,
+                  const size_t /* epoch */,
+                  const double /* objective */)
+  {
+    epochCalls++;
+  }
+
+  /**
+   * Callback function called at the end of a pass over the data.
+   *
+   * @param optimizer The optimizer used to update the function.
+   * @param function Function to optimize.
+   * @param coordinates Starting point.
+   * @param epoch The index of the current epoch.
+   * @param objective Objective value of the current point.
+   */
+  template<typename OptimizerType, typename FunctionType, typename MatType>
+  void EndEpoch(OptimizerType& optimizer,
+                FunctionType& /* function */,
+                const MatType& /* coordinates */,
+                const size_t /* epoch */,
+                const double objective)
+  {
+    // In case StepTaken() has been called first we clear the existing data.
+    if (!hasEndEpoch)
+    {
+      hasEndEpoch = true;
+
+      objectives.clear();
+      timings.clear();
+      gradientsNorm.clear();
+      stepsizes.clear();
+    }
+
+    objectives.push_back(objective);
+    timings.push_back(optimizationTimer.toc());
+
+    if (hasGradient)
+      gradientsNorm.push_back(gradientNorm);
+
+    SaveStepSize(optimizer);
+  }
+
+  /**
+   * Callback function called once a step is taken.
+   *
+   * @param optimizer The optimizer used to update the function.
+   * @param function Function to optimize.
+   * @param coordinates Starting point.
+   * @param objective Objective value of the current point.
+   */
+  template<typename OptimizerType, typename FunctionType, typename MatType>
+  void StepTaken(OptimizerType& optimizer,
+                 FunctionType& /* function */,
+                 const MatType& /* coordinates */)
+  {
+    if (!hasEndEpoch)
+    {
+      objectives.push_back(objective);
+      timings.push_back(optimizationTimer.toc());
+
+      if (hasGradient)
+        gradientsNorm.push_back(gradientNorm);
+
+      SaveStepSize(optimizer);
+    }
+  }
+
+  /**
+   * Callback function called at any call to Evaluate().
+   *
+   * @param optimizer The optimizer used to update the function.
+   * @param function Function to optimize.
+   * @param coordinates Starting point.
+   * @param objectiveIn Objective value of the current point.
+   */
+  template<typename OptimizerType, typename FunctionType, typename MatType>
+  void Evaluate(OptimizerType& /* optimizer */,
+                FunctionType& /* function */,
+                const MatType& /* coordinates */,
+                const double objectiveIn)
+  {
+    objective = objectiveIn;
+    evaluateCalls++;
+  }
+
+  /**
+    * Callback function called at any call to EvaluateConstraint().
+    *
+    * @param optimizer The optimizer used to update the function.
+    * @param function Function to optimize.
+    * @param coordinates Starting point.
+    * @param constraint The index of the constraint;
+    * @param objectiveIn Objective value of the current point.
+    */
+  template<typename OptimizerType, typename FunctionType, typename MatType>
+  void EvaluateConstraint(OptimizerType& /* optimizer */,
+                          FunctionType& /* function */,
+                          const MatType& /* coordinates */,
+                          const size_t /* constraint */,
+                          const double objectiveIn)
+  {
+    objective += objectiveIn;
+    evaluateCalls++;
+  }
+
+  /**
+   * Callback function called at any call to Gradient().
+   *
+   * @param optimizer The optimizer used to update the function.
+   * @param function Function to optimize.
+   * @param coordinates Starting point.
+   * @param gradientIn Matrix that holds the gradient.
+   */
+  template<typename OptimizerType, typename FunctionType, typename MatType>
+  void Gradient(OptimizerType& /* optimizer */,
+                FunctionType& /* function */,
+                const MatType& /* coordinates */,
+                const MatType& gradientIn)
+  {
+    hasGradient = true;
+    gradientNorm = arma::norm(gradientIn);
+    gradientCalls++;
+  }
+
+  /**
+   * Callback function called at any call to GradientConstraint().
+   *
+   * @param optimizer The optimizer used to update the function.
+   * @param function Function to optimize.
+   * @param coordinates Starting point.
+   * @param constraint The index of the constraint;
+   * @param gradient Matrix that holds the gradient;
+   */
+  template<typename OptimizerType, typename FunctionType, typename MatType>
+  void GradientConstraint(OptimizerType& optimizer,
+                          FunctionType& function,
+                          const MatType& coordinates,
+                          const size_t /* constraint */,
+                          const MatType& gradient)
+  {
+    Gradient(optimizer, function, coordinates, gradient);
+  }
+
+ private:
+  /**
+   * Helper function to print the number of function to the specified output
+   * stream.
+   *
+   * @param function The instantiated function that implements NumFunctions().
+   * @param stream The output stream.
+   */
+  template<typename FunctionType>
+  typename std::enable_if<
+      traits::HasNumFunctionsSignature<FunctionType>::value, void>::type
+  PrintNumFunctions(const FunctionType& function, std::stringstream& stream)
+  {
+    PrettyPrintElement(stream, "Number of functions:", 30);
+    stream << function.NumFunctions() << std::endl;
+  }
+
+  template<typename FunctionType>
+  typename std::enable_if<
+      !traits::HasNumFunctionsSignature<FunctionType>::value, void>::type
+  PrintNumFunctions(const FunctionType& /* function */,
+                    std::stringstream& /* stream */) { }
+
+  /**
+   * Helper function to output the max-iterations to the specified output
+   * stream.
+   *
+   * @param optimizer The instantiated optimizer that implements
+   *     MaxIterations().
+   * @param stream The output stream.
+   */
+  template<typename OptimizerType>
+  typename std::enable_if<
+      traits::HasMaxIterationsSignature<OptimizerType>::value, void>::type
+  PrintMaxIterations(const OptimizerType& optimizer, std::stringstream& stream)
+  {
+    PrettyPrintElement(stream, "Maximum iterations:", 30);
+    stream << optimizer.MaxIterations() << std::endl;
+
+    PrettyPrintElement(stream, "Reached maximum iterations:", 30);
+    stream << std::string(optimizer.MaxIterations() == objectives.size() ?
+        "true" : "false") << std::endl;
+  }
+
+  template<typename OptimizerType>
+  typename std::enable_if<
+      !traits::HasMaxIterationsSignature<OptimizerType>::value, void>::type
+  PrintMaxIterations(const OptimizerType& /* optimizer */,
+                     std::stringstream& /* stream */) { }
+
+  /**
+   * Helper function to output the batch-size to the specified output stream.
+   *
+   * @param optimizer The instantiated optimizer that implements BatchSize().
+   * @param stream The output stream.
+   */
+  template<typename OptimizerType>
+  typename std::enable_if<traits::HasBatchSizeSignature<OptimizerType>::value,
+      void>::type
+  PrintBatchSize(const OptimizerType& optimizer, std::stringstream& stream)
+  {
+    PrettyPrintElement(stream, "Batch size:", 30);
+    stream << optimizer.BatchSize() << std::endl;
+  }
+
+  template<typename OptimizerType>
+  typename std::enable_if<!traits::HasBatchSizeSignature<OptimizerType>::value,
+      void>::type
+  PrintBatchSize(const OptimizerType& /* optimizer */,
+                 std::stringstream& /* stream */) { }
+
+  /**
+   * Output formatted data.
+   *
+   * @param out Output stream.
+   * @param data The data to print on the given stream.
+   * @param width The width of the the formatted output data.
+   */
+  template<typename T>
+  void PrettyPrintElement(std::ostream& out,
+                          const T& data,
+                          const size_t width = 14)
+  {
+    out << std::left << std::setw(width) << std::setfill(' ')
+        << std::setprecision(3) << data;
+  }
+
+  /**
+   * Output formatted data.
+   *
+   * @param data The data to print on the given stream.
+   * @param width The width of the the formatted output data.
+   */
+  template<typename T>
+  void PrettyPrintElement(const T& data, const size_t width = 14)
+  {
+    PrettyPrintElement(output, data, width);
+  }
+
+  /**
+   * Outputs the given matrix in a truncated format. For example, the matrix:
+   *
+   * 1 2 3 4 5
+   * 6 7 8 9 10
+   * 11 12 13 14
+   * 15 16 17 18
+   *
+   * will be truncated to:
+   *
+   * 1 2 ... 5
+   * 6 7 ... 10
+   * ...
+   * 15 16 ... 18
+   *
+   * @param data The data to print on the given stream in a truncated format.
+   * @param size The number of elements per column/row.
+   */
+  template<typename T>
+  void TruncatePrint(const T& data, const size_t size)
+  {
+    // We can't directly output the result of submat or use .print, because
+    // both introduce a new line at the end, so we iterate over the elements.
+    for (size_t c = 0, n = 0; c < data.n_cols; ++c)
+    {
+      // Skip to the last column.
+      if (c >= (size - 1))
+      {
+        output << "..." << std::endl;
+        n = (data.n_cols - 2) * data.n_rows - 1;
+      }
+
+      for (size_t r = 0; r < data.n_rows; ++r)
+      {
+        // Check if need to skip to the last row.
+        if (r < (size - 1))
+        {
+          output << std::fixed;
+
+          // Add space for positive value, to align with negative values.
+          if (data(n) >= 0)
+            output << " ";
+
+          output << data(n++) << " ";
+        }
+        else
+        {
+          n = (c + 1) * data.n_rows - 1;
+          output << " ... " << data(n) << std::endl;
+          break;
+        }
+      }
+
+      if (c >= (size - 1))
+        break;
+    }
+  }
+
+  /**
+   * Helper function to store the step-size.
+   *
+   * @param optimizer The instantiated optimzer that implements StepSize().
+   */
+  template<typename OptimizerType>
+  typename std::enable_if<traits::HasStepSizeSignature<OptimizerType>::value,
+      void>::type
+  SaveStepSize(const OptimizerType& optimizer)
+  {
+    stepsizes.push_back(optimizer.StepSize());
+  }
+
+  template<typename OptimizerType>
+  typename std::enable_if<!traits::HasStepSizeSignature<OptimizerType>::value,
+      void>::type
+  SaveStepSize(const OptimizerType& /* optimizer */) { }
+
+  //! The number of iterations to print in percent.
+  double iterationsPercentage;
+
+  //! The output stream that all data is to be sent to; example: std::cout.
+  std::ostream& output;
+
+  //! The number of values to print for the function coordinates.
+  size_t outputMatrixSize;
+  //! The initial coordinates.
+  arma::mat initialCoordinates;
+
+  //! Gradient norm storage.
+  std::vector<double> gradientsNorm;
+
+  //! Objective storage.
+  std::vector<double> objectives;
+
+  //! Timing storage.
+  std::vector<double> timings;
+
+  //! Step-size storage.
+  std::vector<double> stepsizes;
+
+  //! Objective over the current epoch.
+  double objective;
+
+  //! Locally-stored gradient norm for a single step.
+  double gradientNorm;
+
+  //! Whether Gradient() was called.
+  bool hasGradient;
+
+  //! Whether EndEpoch() was called.
+  bool hasEndEpoch;
+
+  //! The number of Gradient() calls.
+  size_t gradientCalls;
+
+  //! The number of Evaluate() calls.
+  size_t evaluateCalls;
+
+  //! The number of BeginEpoch() calls.
+  size_t epochCalls;
+
+  //! Locally-stored optimization step timer object.
+  arma::wall_clock optimizationTimer;
+};
+
+} // namespace ens
+
+#endif

include/ensmallen_bits/callbacks/store_best_coordinates.hpp

@@ -55,7 +55,7 @@ class StoreBestCoordinates
   //! Get the best coordinates.
   ModelMatType const& BestCoordinates() const { return bestCoordinates; }
   //! Modify the best coordinates.
-  ModelMatType& BestCoordinatesl() { return bestCoordinates; }
+  ModelMatType& BestCoordinates() { return bestCoordinates; }
 
   //! Get the best objective.
   double const& BestObjective() const { return bestObjective; }

include/ensmallen_bits/cmaes/cmaes_impl.hpp

@@ -77,7 +77,7 @@ typename MatType::elem_type CMAES<SelectionPolicyType>::Optimize(
   const double muEffective = 1 / arma::accu(arma::pow(w, 2));
 
   // Step size control parameters.
-  BaseMatType sigma(3, 1); // sigma is vector-shaped.
+  BaseMatType sigma(2, 1); // sigma is vector-shaped.
   sigma(0) = 0.3 * (upperBound - lowerBound);
   const double cs = (muEffective + 2) / (iterate.n_elem + muEffective + 5);
   const double ds = 1 + cs + 2 * std::max(std::sqrt((muEffective - 1) /
@@ -97,7 +97,7 @@ typename MatType::elem_type CMAES<SelectionPolicyType>::Optimize(
       muEffective) / (std::pow(iterate.n_elem + 2, 2) +
       alphaMu * muEffective / 2));
 
-  std::vector<BaseMatType> mPosition(3, BaseMatType(iterate.n_rows,
+  std::vector<BaseMatType> mPosition(2, BaseMatType(iterate.n_rows,
       iterate.n_cols));
   mPosition[0] = lowerBound + arma::randu<BaseMatType>(
       iterate.n_rows, iterate.n_cols) * (upperBound - lowerBound);
@@ -216,8 +216,7 @@ typename MatType::elem_type CMAES<SelectionPolicyType>::Optimize(
     }
 
     const ElemType psNorm = arma::norm(ps[idx1]);
-    sigma(idx1) = sigma(idx0) * std::pow(
-        std::exp(cs / ds * psNorm / enn - 1), 0.3);
+    sigma(idx1) = sigma(idx0) * std::exp(cs / ds * ( psNorm / enn - 1));
 
     // Update covariance matrix.
     if ((psNorm / sqrt(1 - std::pow(1 - cs, 2 * i))) < h)

include/ensmallen_bits/ens_version.hpp

@@ -15,13 +15,17 @@
 #define ENS_VERSION_MAJOR 2
 // The minor version is two digits so regular numerical comparisons of versions
 // work right.  The first minor version of a release is always 10.
-#define ENS_VERSION_MINOR 11
-#define ENS_VERSION_PATCH 1
+#define ENS_VERSION_MINOR 14
+#define ENS_VERSION_PATCH 2
 // If this is a release candidate, it will be reflected in the version name
 // (i.e. the version name will be "RC1", "RC2", etc.).  Otherwise the version
 // name will typically be a seemingly arbitrary set of words that does not
 // contain the capitalized string "RC".
-#define ENS_VERSION_NAME "The Poster Session Is Full"
+#define ENS_VERSION_NAME  "No Direction Home"
+// Incorporate the date the version was released.
+#define ENS_VERSION_YEAR  "2020"
+#define ENS_VERSION_MONTH "09"
+#define ENS_VERSION_DAY   "05"
 
 namespace ens {
 
@@ -41,6 +45,14 @@ struct version
 
     return ss.str();
   }
+    
+  static inline std::string date()
+  {
+    std::stringstream ss;
+    ss << ENS_VERSION_YEAR << '-' << ENS_VERSION_MONTH << '-' << ENS_VERSION_DAY;
+
+    return ss.str();
+  }
 };
 
 } // namespace ens

include/ensmallen_bits/ftml/ftml.hpp

@@ -105,7 +105,8 @@ class FTML
            CallbackTypes&&... callbacks)
   {
     return optimizer.Optimize<SeparableFunctionType, MatType, GradType,
-        CallbackTypes...>(function, iterate, callbacks...);
+        CallbackTypes...>(function, iterate,
+        std::forward<CallbackTypes>(callbacks)...);
   }
 
   //! Forward the MatType as GradType.

include/ensmallen_bits/function/arma_traits.hpp

@@ -76,6 +76,31 @@ struct MatTypeTraits<arma::SpSubview<eT>>
       "or a matrix alias instead!");
 };
 
+
+#if ((ARMA_VERSION_MAJOR >= 10) || \
+    ((ARMA_VERSION_MAJOR == 9) && (ARMA_VERSION_MINOR >= 869)))
+
+// Armadillo 9.869+ has SpSubview_col and SpSubview_row
+
+template<typename eT>
+struct MatTypeTraits<arma::SpSubview_col<eT>>
+{
+  static_assert(sizeof(arma::SpSubview_col<eT>) == 0,
+      "Armadillo subviews cannot be passed to Optimize()!  Create a matrix "
+      "or a matrix alias instead!");
+};
+
+template<typename eT>
+struct MatTypeTraits<arma::SpSubview_row<eT>>
+{
+  static_assert(sizeof(arma::SpSubview_row<eT>) == 0,
+      "Armadillo subviews cannot be passed to Optimize()!  Create a matrix "
+      "or a matrix alias instead!");
+};
+
+#endif
+
+
 template<typename eT>
 struct MatTypeTraits<arma::Cube<eT>>
 {

include/ensmallen_bits/function/traits.hpp

@@ -45,6 +45,8 @@ ENS_HAS_EXACT_METHOD_FORM(MaxIterations, HasMaxIterations)
 ENS_HAS_EXACT_METHOD_FORM(ResetPolicy, HasResetPolicy)
 //! Detect an BatchSize() method.
 ENS_HAS_EXACT_METHOD_FORM(BatchSize, HasBatchSize)
+//! Detect an StepSize() method.
+ENS_HAS_EXACT_METHOD_FORM(StepSize, HasStepSize)
 
 template<typename MatType, typename GradType>
 struct TypedForms
@@ -391,6 +393,22 @@ struct HasBatchSizeSignature
       HasBatchSize<OptimizerType, BatchSizeConstForm>::value;
 };
 
+//! Utility struct, check if size_t StepSize() const or size_t StepSize()
+//! exists.
+template<typename OptimizerType>
+struct HasStepSizeSignature
+{
+  template<typename C>
+  using StepSizeConstForm = double(C::*)(void) const;
+
+  template<typename C>
+  using StepSizeForm = double(C::*)(void);
+
+  const static bool value =
+      HasStepSize<OptimizerType, StepSizeForm>::value ||
+      HasStepSize<OptimizerType, StepSizeConstForm>::value;
+};
+
 //! Utility struct, check if size_t MaxIterations() const exists.
 template<typename OptimizerType>
 struct HasMaxIterationsSignature

include/ensmallen_bits/lbfgs/lbfgs_impl.hpp

@@ -391,7 +391,7 @@ L_BFGS::Optimize(FunctionType& function,
     //
     // But don't do this on the first iteration to ensure we always take at
     // least one descent step.
-    if (itNum > 0 && (arma::norm(gradient, 2) < minGradientNorm))
+    if (arma::norm(gradient, 2) < minGradientNorm)
     {
       Info << "L-BFGS gradient norm too small (terminating successfully)."
           << std::endl;

include/ensmallen_bits/lookahead/lookahead.hpp

@@ -108,6 +108,11 @@ class Lookahead
             const bool resetPolicy = false,
             const bool exactObjective = false);
 
+  /**
+   * Clean any memory associated with the Lookahead object.
+   */
+  ~Lookahead();
+
   /**
    * Optimize the given function using Lookahead. The given starting point will
    * be modified to store the finishing point of the algorithm, and the final
@@ -188,6 +193,31 @@ class Lookahead
   bool& ExactObjective() { return exactObjective; }
 
  private:
+  /**
+   * Set the maximum number of iterations if the given optimizer implements
+   * MaxIterations().
+   *
+   * @param optimizer Optimizer to check for MaxIterations().
+   * @param k The number of iterations.
+   */
+  template<typename OptimizerType>
+  static typename std::enable_if<traits::HasMaxIterationsSignature<
+      OptimizerType>::value, void>::type
+  SetMaxIterations(OptimizerType& optimizer, const size_t k)
+  {
+    optimizer.MaxIterations() = k;
+  }
+
+  template<typename OptimizerType>
+  static typename std::enable_if<!traits::HasMaxIterationsSignature<
+      OptimizerType>::value, void>::type
+  SetMaxIterations(const OptimizerType& /* optimizer */, const size_t /* k */)
+  {
+    Warn << "The base optimizer does not have a definition of "
+        << "MaxIterations(), the base optimizer will have its configuration "
+        << "unchanged.";
+  }
+
   //! The base optimizer for the forward step.
   BaseOptimizerType baseOptimizer;
 

include/ensmallen_bits/lookahead/lookahead_impl.hpp

@@ -60,6 +60,12 @@ inline Lookahead<BaseOptimizerType, DecayPolicyType>::Lookahead(
     isInitialized(false)
 { /* Nothing to do. */ }
 
+template<typename BaseOptimizerType, typename DecayPolicyType>
+inline Lookahead<BaseOptimizerType, DecayPolicyType>::~Lookahead()
+{
+  instDecayPolicy.Clean();
+}
+
 //! Optimize the function (minimize).
 template<typename BaseOptimizerType, typename DecayPolicyType>
 template<typename SeparableFunctionType,
@@ -98,16 +104,7 @@ Lookahead<BaseOptimizerType, DecayPolicyType>::Optimize(
 
   // Check if the optimizer implements HasMaxIterations() and override the
   // parameter with k.
-  if (traits::HasMaxIterationsSignature<BaseOptimizerType>::value)
-  {
-    baseOptimizer.MaxIterations() = k;
-  }
-  else
-  {
-    Warn << "The base optimizer does not have a definition of "
-        << "MaxIterations(), the base optimizer will have its configuration "
-        << "unchanged.";
-  }
+  SetMaxIterations(baseOptimizer, k);
 
   // Check if the optimizer implements ResetPolicy() and override the reset
   // policy.

include/ensmallen_bits/nsga2/nsga2.hpp

@@ -0,0 +1,349 @@
+/**
+ * @file nsga2.hpp
+ * @author Sayan Goswami
+ *
+ * NSGA-II is a multi-objective optimization algorithm, widely used in
+ * many real-world applications. NSGA-II generates offsprings using
+ * crossover and mutation and then selects the next generation according
+ * to non-dominated-sorting and crowding distance comparison.
+ *
+ * ensmallen is free software; you may redistribute it and/or modify it under
+ * the terms of the 3-clause BSD license.  You should have received a copy of
+ * the 3-clause BSD license along with ensmallen.  If not, see
+ * http://www.opensource.org/licenses/BSD-3-Clause for more information.
+ */
+
+#ifndef ENSMALLEN_NSGA2_NSGA2_HPP
+#define ENSMALLEN_NSGA2_NSGA2_HPP
+
+namespace ens {
+
+/**
+ * NSGA-II (Non-dominated Sorting Genetic Algorithm - II) is a multi-objective
+ * optimization algorithm. This class implements the NSGA-II algorithm.
+ *
+ * The algorithm works by generating a candidate population from a fixed
+ * starting point. At each stage of optimization, a new population of children
+ * is generated. This new population along with its predecessor is sorted using
+ * non-domination as the metric. Following this, the population is further
+ * segregated in fronts. A new population is generated from these fronts having
+ * size equal to that of the starting population.
+ *
+ * During evolution, two parents are randomly chosen using binary tournament
+ * selection. A pair of children are generated by crossing over these two
+ * candidates followed by mutation.
+ *
+ * The best front (Pareto optimal) is returned by the Optimize() method.
+ *
+ * For more information, see the following:
+ *
+ * @code
+ * @article{10.1109/4235.996017,
+ *   author = {Deb, K. and Pratap, A. and Agarwal, S. and Meyarivan, T.},
+ *   title = {A Fast and Elitist Multiobjective Genetic Algorithm: NSGA-II},
+ *   year = {2002},
+ *   url = {https://doi.org/10.1109/4235.996017},
+ *   journal = {Trans. Evol. Comp}}
+ * @endcode
+ *
+ * NSGA-II can optimize arbitrary multi-objective functions. For more details,
+ * see the documentation on function types included with this distribution or
+ * on the ensmallen website.
+ */
+class NSGA2 {
+ public:
+  /**
+   * Constructor for the NSGA2 optimizer.
+   *
+   * The default values provided over here are not necessarily suitable for a
+   * given function. Therefore it is highly recommended to adjust the
+   * parameters according to the problem.
+   *
+   * @param populationSize The number of candidates in the population.
+   *     This should be atleast 4 in size and a multiple of 4.
+   * @param maxGenerations The maximum number of generations allowed for NSGA-II.
+   * @param crossoverProb The probability that a crossover will occur.
+   * @param mutationProb The probability that a mutation will occur.
+   * @param mutationStrength The strength of the mutation.
+   * @param epsilon The minimum difference required to distinguish between
+   *     candidate solutions.
+   * @param lowerBound Lower bound of the coordinates of the initial population.
+   * @param upperBound Upper bound of the coordinates of the initial population.
+   */
+  NSGA2(const size_t populationSize = 100,
+        const size_t maxGenerations = 2000,
+        const double crossoverProb = 0.6,
+        const double mutationProb = 0.3,
+        const double mutationStrength = 1e-3,
+        const double epsilon = 1e-6,
+        const arma::vec& lowerBound = arma::zeros(1, 1),
+        const arma::vec& upperBound = arma::ones(1, 1));
+
+  /**
+   * Constructor for the NSGA2 optimizer. This constructor provides an overload
+   * to use `lowerBound` and `upperBound` of type double.
+   *
+   * The default values provided over here are not necessarily suitable for a
+   * given function. Therefore it is highly recommended to adjust the
+   * parameters according to the problem.
+   *
+   * @param populationSize The number of candidates in the population.
+   *     This should be atleast 4 in size and a multiple of 4.
+   * @param maxGenerations The maximum number of generations allowed for NSGA-II.
+   * @param crossoverProb The probability that a crossover will occur.
+   * @param mutationProb The probability that a mutation will occur.
+   * @param mutationStrength The strength of the mutation.
+   * @param epsilon The minimum difference required to distinguish between
+   *     candidate solutions.
+   * @param lowerBound Lower bound of the coordinates of the initial population.
+   * @param upperBound Upper bound of the coordinates of the initial population.
+   */
+  NSGA2(const size_t populationSize = 100,
+        const size_t maxGenerations = 2000,
+        const double crossoverProb = 0.6,
+        const double mutationProb = 0.3,
+        const double mutationStrength = 1e-3,
+        const double epsilon = 1e-6,
+        const double lowerBound = 0,
+        const double upperBound = 1);
+
+  /**
+   * Optimize a set of objectives. The initial population is generated using the
+   * starting point. The output is the best generated front.
+   *
+   * @tparam ArbitraryFunctionType std::tuple of multiple objectives.
+   * @tparam MatType Type of matrix to optimize.
+   * @tparam CallbackTypes Types of callback functions.
+   * @param objectives Vector of objective functions to optimize for.
+   * @param iterate Starting point.
+   * @param callbacks Callback functions.
+   * @return MatType::elem_type The minimum of the accumulated sum over the
+   *     objective values in the best front.
+   */
+  template<typename MatType,
+           typename... ArbitraryFunctionType,
+           typename... CallbackTypes>
+ typename MatType::elem_type Optimize(
+     std::tuple<ArbitraryFunctionType...>& objectives,
+     MatType& iterate,
+     CallbackTypes&&... callbacks);
+
+  //! Get the population size.
+  size_t PopulationSize() const { return populationSize; }
+  //! Modify the population size.
+  size_t& PopulationSize() { return populationSize; }
+
+  //! Get the maximum number of generations.
+  size_t MaxGenerations() const { return maxGenerations; }
+  //! Modify the maximum number of generations.
+  size_t& MaxGenerations() { return maxGenerations; }
+
+  //! Get the crossover rate.
+  double CrossoverRate() const { return crossoverProb; }
+  //! Modify the crossover rate.
+  double& CrossoverRate() { return crossoverProb; }
+
+  //! Get the mutation probability.
+  double MutationProbability() const { return mutationProb; }
+  //! Modify the mutation probability.
+  double& MutationProbability() { return mutationProb; }
+
+  //! Get the mutation strength.
+  double MutationStrength() const { return mutationStrength; }
+  //! Modify the mutation strength.
+  double& MutationStrength() { return mutationStrength; }
+
+  //! Get the tolerance.
+  double Epsilon() const { return epsilon; }
+  //! Modify the tolerance.
+  double& Epsilon() { return epsilon; }
+
+  //! Retrieve value of lowerBound.
+  const arma::vec& LowerBound() const { return lowerBound; }
+  //! Modify value of lowerBound.
+  arma::vec& LowerBound() { return lowerBound; }
+
+  //! Retrieve value of upperBound.
+  const arma::vec& UpperBound() const { return upperBound; }
+  //! Modify value of upperBound.
+  arma::vec& UpperBound() { return upperBound; }
+
+  //! Retrieve the best front (the Pareto frontier).  This returns an empty vector until `Optimize()`
+  //! has been called.
+  const std::vector<arma::mat>& Front() const { return bestFront; }
+
+ private:
+  /**
+   * Evaluate objectives for the elite population.
+   *
+   * @tparam ArbitraryFunctionType std::tuple of multiple function types.
+   * @tparam MatType Type of matrix to optimize.
+   * @param population The elite population.
+   * @param objectives The set of objectives.
+   * @param calculatedObjectives Vector to store calculated objectives.
+   */
+  template<std::size_t I = 0,
+           typename MatType,
+           typename ...ArbitraryFunctionType>
+  typename std::enable_if<I == sizeof...(ArbitraryFunctionType), void>::type
+  EvaluateObjectives(std::vector<MatType>&,
+                     std::tuple<ArbitraryFunctionType...>&,
+                     std::vector<arma::Col<double> >&);
+
+  template<std::size_t I = 0,
+           typename MatType,
+           typename ...ArbitraryFunctionType>
+  typename std::enable_if<I < sizeof...(ArbitraryFunctionType), void>::type
+  EvaluateObjectives(std::vector<MatType>& population,
+                     std::tuple<ArbitraryFunctionType...>& objectives,
+                     std::vector<arma::Col<double> >& calculatedObjectives);
+
+  /**
+   * Reproduce candidates from the elite population to generate a new
+   * population.
+   *
+   * @tparam MatType Type of matrix to optimize.
+   * @param population The elite population.
+   * @param objectives The set of objectives.
+   * @param lowerBound Lower bound of the coordinates of the initial population.
+   * @param upperBound Upper bound of the coordinates of the initial population.
+   */
+  template<typename MatType>
+  void BinaryTournamentSelection(std::vector<MatType>& population,
+                                 const arma::vec& lowerBound,
+                                 const arma::vec& upperBound);
+
+  /**
+   * Crossover two parents to create a pair of new children.
+   *
+   * @tparam MatType Type of matrix to optimize.
+   * @param childA A newly generated candidate.
+   * @param childB Another newly generated candidate.
+   * @param parentA First parent from elite population.
+   * @param parentB Second parent from elite population.
+   */
+  template<typename MatType>
+  void Crossover(MatType& childA,
+                 MatType& childB,
+                 const MatType& parentA,
+                 const MatType& parentB);
+
+  /**
+   * Mutate the coordinates for a candidate.
+   *
+   * @tparam MatType Type of matrix to optimize.
+   * @param child The candidate whose coordinates are being modified.
+   * @param objectives The set of objectives.
+   * @param lowerBound Lower bound of the coordinates of the initial population.
+   * @param upperBound Upper bound of the coordinates of the initial population.
+   */
+  template<typename MatType>
+  void Mutate(MatType& child,
+              const arma::vec& lowerBound,
+              const arma::vec& upperBound);
+
+  /**
+   * Sort the candidate population using their domination count and the set of
+   * dominated nodes.
+   *
+   * @tparam MatType Type of matrix to optimize.
+   * @param fronts The population is sorted into these Pareto fronts. The first
+   *     front is the best, the second worse and so on.
+   * @param ranks The assigned ranks, used for crowding distance based sorting.
+   * @param calculatedObjectives The previously calculated objectives.
+   */
+  template<typename MatType>
+  void FastNonDominatedSort(
+      std::vector<std::vector<size_t> >& fronts,
+      std::vector<size_t>& ranks,
+      std::vector<arma::Col<typename MatType::elem_type> >& calculatedObjectives);
+
+  /**
+   * Operator to check if one candidate Pareto-dominates the other.
+   *
+   * A candidate is said to dominate the other if it is at least as good as the
+   * other candidate for all the objectives and there exists at least one
+   * objective for which it is strictly better than the other candidate.
+   *
+   * @tparam MatType Type of matrix to optimize.
+   * @param calculatedObjectives The previously calculated objectives.
+   * @param candidateP The candidate being compared from the elite population.
+   * @param candidateQ The candidate being compared against.
+   * @return true if candidateP Pareto dominates candidateQ, otherwise, false.
+   */
+  template<typename MatType>
+  bool Dominates(
+      std::vector<arma::Col<typename MatType::elem_type> >& calculatedObjectives,
+      size_t candidateP,
+      size_t candidateQ);
+
+  /**
+   * Assigns crowding distance metric for sorting.
+   *
+   * @param front The previously generated Pareto fronts.
+   * @param objectives The set of objectives.
+   * @param crowdingDistance The previously calculated objectives.
+   */
+  void CrowdingDistanceAssignment(const std::vector<size_t>& front,
+                                  std::vector<double>& crowdingDistance);
+
+  /**
+   * The operator used in the crowding distance based sorting.
+   *
+   * If a candidates has a lower rank then it is preferred.
+   * Otherwise, if the ranks are equal then the candidate with the larger
+   * crowding distance is preferred.
+   *
+   * @param idxP The index of the first cadidate from the elite population being
+   *     sorted.
+   * @param idxQ The index of the second cadidate from the elite population
+   *     being sorted.
+   * @param ranks The previously calculated ranks.
+   * @param crowdingDistance The previously calculated objectives.
+   * @return true if the first candidate is preferred, otherwise, false.
+   */
+  bool CrowdingOperator(size_t idxP,
+                        size_t idxQ,
+                        const std::vector<size_t>& ranks,
+                        const std::vector<double>& crowdingDistance);
+
+  //! The number of objectives being optimised for.
+  size_t numObjectives;
+
+  //! The numbeer of variables used per objectives.
+  size_t numVariables;
+
+  //! The number of candidates in the population.
+  size_t populationSize;
+
+  //! Maximum number of generations before termination criteria is met.
+  size_t maxGenerations;
+
+  //! Probability that crossover will occur.
+  double crossoverProb;
+
+  //! Probability that mutation will occur.
+  double mutationProb;
+
+  //! Strength of the mutation.
+  double mutationStrength;
+
+  //! The tolerance for termination.
+  double epsilon;
+
+  //! Lower bound of the initial swarm.
+  arma::vec lowerBound;
+
+  //! Upper bound of the initial swarm.
+  arma::vec upperBound;
+
+  //! Best front, stored after Optimize() is called.
+  std::vector<arma::mat> bestFront;
+};
+
+} // namespace ens
+
+// Include implementation.
+#include "nsga2_impl.hpp"
+
+#endif

include/ensmallen_bits/nsga2/nsga2_impl.hpp

@@ -0,0 +1,438 @@
+/**
+ * @file nsga2_impl.hpp
+ * @author Sayan Goswami
+ *
+ * Implementation of the NSGA-II algorithm. Used for multi-objective
+ * optimization problems on arbitrary functions.
+ *
+ * ensmallen is free software; you may redistribute it and/or modify it under
+ * the terms of the 3-clause BSD license.  You should have received a copy of
+ * the 3-clause BSD license along with ensmallen.  If not, see
+ * http://www.opensource.org/licenses/BSD-3-Clause for more Information.
+ */
+
+#ifndef ENSMALLEN_NSGA2_NSGA2_IMPL_HPP
+#define ENSMALLEN_NSGA2_NSGA2_IMPL_HPP
+
+#include "nsga2.hpp"
+#include <assert.h>
+
+namespace ens {
+
+inline NSGA2::NSGA2(const size_t populationSize,
+                    const size_t maxGenerations,
+                    const double crossoverProb,
+                    const double mutationProb,
+                    const double mutationStrength,
+                    const double epsilon,
+                    const arma::vec& lowerBound,
+                    const arma::vec& upperBound) :
+    populationSize(populationSize),
+    maxGenerations(maxGenerations),
+    crossoverProb(crossoverProb),
+    mutationProb(mutationProb),
+    mutationStrength(mutationStrength),
+    epsilon(epsilon),
+    lowerBound(lowerBound),
+    upperBound(upperBound)
+{ /* Nothing to do here. */ }
+
+inline NSGA2::NSGA2(const size_t populationSize,
+                    const size_t maxGenerations,
+                    const double crossoverProb,
+                    const double mutationProb,
+                    const double mutationStrength,
+                    const double epsilon,
+                    const double lowerBound,
+                    const double upperBound) :
+    populationSize(populationSize),
+    maxGenerations(maxGenerations),
+    crossoverProb(crossoverProb),
+    mutationProb(mutationProb),
+    mutationStrength(mutationStrength),
+    epsilon(epsilon),
+    lowerBound(lowerBound * arma::ones(1, 1)),
+    upperBound(upperBound * arma::ones(1, 1))
+{ /* Nothing to do here. */ }
+
+//! Optimize the function.
+template<typename MatType,
+         typename... ArbitraryFunctionType,
+         typename... CallbackTypes>
+typename MatType::elem_type NSGA2::Optimize(
+    std::tuple<ArbitraryFunctionType...>& objectives,
+    MatType& iterate,
+    CallbackTypes&&... callbacks)
+{
+  // Make sure for evolution to work at least four candidates are present.
+  if (populationSize < 4 && populationSize % 4 != 0)
+  {
+    throw std::logic_error("NSGA2::Optimize(): population size should be at"
+        " least 4, and, a multiple of 4!");
+  }
+
+  // Check if lower bound is a vector of a single dimension.
+  if (lowerBound.n_rows == 1)
+    lowerBound = lowerBound(0, 0) * arma::ones(iterate.n_rows, iterate.n_cols);
+
+  // Check if lower bound is a vector of a single dimension.
+  if (upperBound.n_rows == 1)
+    upperBound = upperBound(0, 0) * arma::ones(iterate.n_rows, iterate.n_cols);
+
+  // Check the dimensions of lowerBound and upperBound.
+  assert(lowerBound.n_rows == iterate.n_rows && "The dimensions of "
+      "lowerBound are not the same as the dimensions of iterate.");
+  assert(upperBound.n_rows == iterate.n_rows && "The dimensions of "
+      "upperBound are not the same as the dimensions of iterate.");
+
+  // Convenience typedefs.
+  typedef typename MatType::elem_type ElemType;
+
+  numObjectives = sizeof...(ArbitraryFunctionType);
+  numVariables = iterate.n_rows;
+
+  // Cache calculated objectives.
+  std::vector<arma::Col<ElemType> > calculatedObjectives;
+  // Pre-allocate space for the calculated objectives.
+  calculatedObjectives.resize(populationSize);
+
+  // Population size reserved to 2 * populationSize + 1 to accommodate
+  // for the size of intermediate candidate population.
+  std::vector<MatType> population;
+  population.reserve(2 * populationSize + 1);
+
+  // Pareto fronts, initialized during non-dominated sorting.
+  std::vector<std::vector<size_t> > fronts;
+  // Initialised in CrowdingDistanceAssignment.
+  std::vector<double> crowdingDistance;
+  // Initialised during non-dominated sorting.
+  std::vector<size_t> ranks;
+
+  // Controls early termination of the optimization process.
+  bool terminate = false;
+
+  // Generate the population based on a uniform distribution around the given
+  // starting point.
+  for (size_t i = 0; i < populationSize; i++)
+  {
+    population.push_back(arma::randu<MatType>(iterate.n_rows,
+        iterate.n_cols) - 0.5 + iterate);
+  }
+
+  Info << "NSGA2 initialized successfully. Optimization started." << std::endl;
+
+  // Evaluate the fitness before optimization.
+  for (size_t i = 0; i < population.size(); i++)
+    calculatedObjectives[i] = arma::Col<ElemType>(numObjectives, arma::fill::zeros);
+  EvaluateObjectives(population, objectives, calculatedObjectives);
+
+  // Iterate until maximum number of generations is obtained.
+  terminate |= Callback::BeginOptimization(*this, objectives, iterate, callbacks...);
+
+  for (size_t generation = 1; generation <= maxGenerations && !terminate; generation++)
+  {
+    Info << "NSGA2: iteration " << generation << "." << std::endl;
+    terminate |= Callback::StepTaken(*this, objectives, iterate, callbacks...);
+
+    // Create new population of candidate from the present elite population.
+    // Have P_t, generate G_t using P_t.
+    BinaryTournamentSelection(population, lowerBound, upperBound);
+
+    // Evaluate the objectives for the new population.
+    calculatedObjectives.resize(population.size());
+    for (size_t i = 0; i < population.size(); i++)
+      calculatedObjectives[i] = arma::Col<ElemType>(numObjectives, arma::fill::zeros);
+    EvaluateObjectives(population, objectives, calculatedObjectives);
+
+    // Perform fast non dominated sort on P_t ∪ G_t.
+    ranks.resize(population.size());
+    FastNonDominatedSort<MatType>(fronts, ranks, calculatedObjectives);
+
+    // Perform crowding distance assignment.
+    crowdingDistance.resize(population.size());
+
+    for (size_t fNum = 0; fNum < fronts.size(); fNum++)
+    {
+      CrowdingDistanceAssignment(fronts[fNum], crowdingDistance);
+    }
+
+    // Sort based on crowding distance.
+    std::sort(population.begin(), population.end(),
+              [this, ranks, crowdingDistance, population](MatType candidateP,
+                                                          MatType candidateQ)
+              {
+                size_t idxP, idxQ;
+                for (size_t i = 0; i < population.size(); i++)
+                {
+                  if (arma::approx_equal(population[i], candidateP, "absdiff", epsilon))
+                    idxP = i;
+
+                  if (arma::approx_equal(population[i], candidateQ, "absdiff", epsilon))
+                    idxQ = i;
+                }
+
+                return CrowdingOperator(idxP, idxQ, ranks, crowdingDistance);
+              }
+    );
+
+    // Yield a new population P_{t+1} of size populationSize.
+    population.resize(populationSize);
+  }
+
+  // Set the candidates from the best front as the output.
+  std::vector<MatType> front;
+
+  for (size_t f: fronts[0])
+    front.push_back(population[f]);
+
+  // bestFront is stored, can be obtained by the Front() getter.
+  bestFront = front;
+
+  // Assign iterate to first element of the best front.
+  iterate = bestFront[0];
+
+  Callback::EndOptimization(*this, objectives, iterate, callbacks...);
+
+  ElemType performance = std::numeric_limits<ElemType>::max();
+
+  for(arma::Col<ElemType> objective: calculatedObjectives)
+    if (arma::accu(objective) < performance)
+      performance = arma::accu(objective);
+
+  return performance;
+}
+
+//! No objectives to evaluate.
+template<std::size_t I,
+         typename MatType,
+         typename ...ArbitraryFunctionType>
+typename std::enable_if<I == sizeof...(ArbitraryFunctionType), void>::type
+NSGA2::EvaluateObjectives(
+    std::vector<MatType>&,
+    std::tuple<ArbitraryFunctionType...>&,
+    std::vector<arma::Col<double> >&)
+{
+  // Nothing to do here.
+}
+
+//! Evaluate the objectives for the entire population.
+template<std::size_t I,
+         typename MatType,
+         typename ...ArbitraryFunctionType>
+typename std::enable_if<I < sizeof...(ArbitraryFunctionType), void>::type
+NSGA2::EvaluateObjectives(
+    std::vector<MatType>& population,
+    std::tuple<ArbitraryFunctionType...>& objectives,
+    std::vector<arma::Col<double> >& calculatedObjectives)
+{
+  for (size_t i = 0; i < populationSize; i++)
+  {
+    calculatedObjectives[i](I) = std::get<I>(objectives).Evaluate(population[i]);
+    EvaluateObjectives<I+1, MatType, ArbitraryFunctionType...>(population, objectives,
+                                                               calculatedObjectives);
+  }
+}
+
+//! Reproduce and generate new candidates.
+template<typename MatType>
+inline void NSGA2::BinaryTournamentSelection(std::vector<MatType>& population,
+                                             const arma::vec& lowerBound,
+                                             const arma::vec& upperBound)
+{
+  std::vector<MatType> children;
+
+  while (children.size() < population.size())
+  {
+    // Choose two random parents for reproduction from the elite population.
+    size_t indexA = arma::randi<size_t>(arma::distr_param(0, populationSize - 1));
+    size_t indexB = arma::randi<size_t>(arma::distr_param(0, populationSize - 1));
+
+    // Make sure that the parents differ.
+    if (indexA == indexB)
+    {
+      if (indexB < populationSize - 1)
+        indexB++;
+      else
+        indexB--;
+    }
+
+    // Initialize the children to the respective parents.
+    MatType childA = population[indexA], childB = population[indexB];
+
+    Crossover(childA, childB, population[indexA], population[indexB]);
+
+    Mutate(childA, lowerBound, upperBound);
+    Mutate(childB, lowerBound, upperBound);
+
+    // Add the children to the candidate population.
+    children.push_back(childA);
+    children.push_back(childB);
+  }
+
+  // Add the candidates to the elite population.
+  population.insert(std::end(population), std::begin(children), std::end(children));
+}
+
+//! Perform crossover of genes for the children.
+template<typename MatType>
+inline void NSGA2::Crossover(MatType& childA,
+                             MatType& childB,
+                             const MatType& parentA,
+                             const MatType& parentB)
+{
+  // Indices at which crossover is to occur.
+  const arma::umat idx = arma::randu<MatType>(childA.n_rows, childA.n_cols) < crossoverProb;
+
+  // Use traits from parentA for indices where idx is 1 and parentB otherwise.
+  childA = parentA % idx + parentB % (1 - idx);
+  // Use traits from parentB for indices where idx is 1 and parentA otherwise.
+  childB = parentA % (1 - idx) + parentA % idx;
+}
+
+//! Perform mutation of the candidates weights with some noise.
+template<typename MatType>
+inline void NSGA2::Mutate(MatType& child,
+                          const arma::vec& lowerBound,
+                          const arma::vec& upperBound)
+{
+  child += (arma::randu<MatType>(child.n_rows, child.n_cols) < mutationProb) %
+      (mutationStrength * arma::randn<MatType>(child.n_rows, child.n_cols));
+
+  // Constrain all genes to be between bounds.
+  for (size_t idx = 0; idx < numVariables; idx++)
+  {
+    if (child[idx] < lowerBound(idx))
+      child[idx] = lowerBound(idx);
+    else if (child[idx] > upperBound(idx))
+      child[idx] = upperBound(idx);
+  }
+}
+
+//! Sort population into Pareto fronts.
+template<typename MatType>
+inline void NSGA2::FastNonDominatedSort(
+    std::vector<std::vector<size_t> >& fronts,
+    std::vector<size_t>& ranks,
+    std::vector<arma::Col<typename MatType::elem_type> >& calculatedObjectives)
+{
+  std::map<size_t, size_t> dominationCount;
+  std::map<size_t, std::set<size_t> > dominated;
+
+  // Reset and initialize fronts.
+  fronts.clear();
+  fronts.push_back(std::vector<size_t>());
+
+  for (size_t p = 0; p < populationSize; p++)
+  {
+    dominated[p] = std::set<size_t>();
+    dominationCount[p] = 0;
+
+    for (size_t q = 0; q < populationSize; q++)
+    {
+      if (Dominates<MatType>(calculatedObjectives, p, q))
+        dominated[p].insert(q);
+      else if (Dominates<MatType>(calculatedObjectives, q, p))
+        dominationCount[p] += 1;
+    }
+
+    if (dominationCount[p] == 0)
+    {
+      ranks[p] = 0;
+      fronts[0].push_back(p);
+    }
+  }
+
+  size_t i = 0;
+
+  while (fronts[i].size() > 0)
+  {
+    std::vector<size_t> nextFront;
+
+    for (size_t p: fronts[i])
+    {
+      for (size_t q: dominated[p])
+      {
+        dominationCount[q]--;
+
+        if (dominationCount[q] == 0)
+        {
+          ranks[q] = i + 1;
+          nextFront.push_back(q);
+        }
+      }
+    }
+
+    i++;
+    fronts.push_back(nextFront);
+  }
+}
+
+//! Check if a candidate Pareto dominates another candidate.
+template<typename MatType>
+inline bool NSGA2::Dominates(
+    std::vector<arma::Col<typename MatType::elem_type> >& calculatedObjectives,
+    size_t candidateP,
+    size_t candidateQ)
+{
+  bool allBetterOrEqual = true;
+  bool atleastOneBetter = false;
+  size_t n_objectives = calculatedObjectives[0].n_elem;
+
+  for (size_t i = 0; i < n_objectives; i++)
+  {
+    // P is worse than Q for the i-th objective function.
+    if (calculatedObjectives[candidateP](i) > calculatedObjectives[candidateQ](i))
+      allBetterOrEqual = false;
+
+    // P is better than Q for the i-th objective function.
+    else if (calculatedObjectives[candidateP](i) < calculatedObjectives[candidateQ](i))
+      atleastOneBetter = true;
+  }
+
+  return allBetterOrEqual && atleastOneBetter;
+}
+
+//! Assign crowding distance to the population.
+inline void NSGA2::CrowdingDistanceAssignment(const std::vector<size_t>& front,
+                                              std::vector<double>& crowdingDistance)
+{
+  if (front.size() > 0)
+  {
+    for (size_t elem: front)
+      crowdingDistance[elem] = 0;
+
+    size_t fSize = front.size();
+
+    for (size_t m = 0; m < numObjectives; m++)
+    {
+      crowdingDistance[front[0]] = std::numeric_limits<double>::max();
+      crowdingDistance[front[fSize - 1]] = std::numeric_limits<double>::max();
+
+      for (size_t i = 1; i < fSize - 1 ; i++)
+      {
+        crowdingDistance[front[i]] += (crowdingDistance[front[i - 1]] -
+            crowdingDistance[front[i + 1]]) /
+            (std::numeric_limits<double>::max() -
+            std::numeric_limits<double>::min());
+      }
+    }
+  }
+}
+
+//! Comparator for crowding distance based sorting.
+inline bool NSGA2::CrowdingOperator(size_t idxP,
+                                    size_t idxQ,
+                                    const std::vector<size_t>& ranks,
+                                    const std::vector<double>& crowdingDistance)
+{
+  if (ranks[idxP] < ranks[idxQ])
+    return true;
+  else if (ranks[idxP] == ranks[idxQ] && crowdingDistance[idxP] > crowdingDistance[idxQ])
+    return true;
+
+  return false;
+}
+
+} // namespace ens
+
+#endif

include/ensmallen_bits/padam/padam.hpp

@@ -113,7 +113,7 @@ class Padam
   {
     return optimizer.template Optimize<
         SeparableFunctionType, MatType, GradType, CallbackTypes...>(
-        function, iterate, callbacks...);
+        function, iterate, std::forward<CallbackTypes>(callbacks)...);
   }
 
   //! Forward the MatType as GradType.

include/ensmallen_bits/problems/aug_lagrangian_test_functions_impl.hpp

@@ -129,14 +129,15 @@ inline typename MatType::elem_type GockenbachFunction::EvaluateConstraint(
   switch (index)
   {
     case 0: // g(x) = (x_3 - x_2 - x_1 - 1) = 0
-      constraint = (coordinates[2] - coordinates[1] - coordinates[0] - 1);
+      constraint = (coordinates[2] - coordinates[1] - coordinates[0] -
+          typename MatType::elem_type(1));
       break;
 
     case 1: // h(x) = (x_3 - x_1^2) >= 0
       // To deal with the inequality, the constraint will simply evaluate to 0
       // when h(x) >= 0.
-      constraint = std::min(0.0,
-          (coordinates[2] - std::pow(coordinates[0], 2)));
+      constraint = std::min(typename MatType::elem_type(0), (coordinates[2] -
+          std::pow(coordinates[0], typename MatType::elem_type(2))));
       break;
   }
 

include/ensmallen_bits/problems/fonseca_fleming_function.hpp

@@ -0,0 +1,103 @@
+/**
+ * @file fonseca_fleming_function_n1.hpp
+ * @author Sayan Goswami
+ *
+ * Implementation of Fonseca Fleming function.
+ *
+ * ensmallen is free software; you may redistribute it and/or modify it under
+ * the terms of the 3-clause BSD license.  You should have received a copy of
+ * the 3-clause BSD license along with ensmallen.  If not, see
+ * http://www.opensource.org/licenses/BSD-3-Clause for more information.
+ */
+
+#ifndef ENSMALLEN_PROBLEMS_FONSECA_FLEEMING_FUNCTION_HPP
+#define ENSMALLEN_PROBLEMS_FONSECA_FLEEMING_FUNCTION_HPP
+
+#include <tuple>
+
+namespace ens {
+namespace test {
+
+/**
+ * The Fonseca Fleming function N.1 is defined by
+ *
+ * \f[
+ * f_{1}\left(\boldsymbol{x}\right) = 1 - \exp \left[-\sum_{i=1}^{3} \left(x_{i} - \frac{1}{\sqrt{n}} \right)^{2} \right] \\
+ * f_{2}\left(\boldsymbol{x}\right) = 1 - \exp \left[-\sum_{i=1}^{3} \left(x_{i} + \frac{1}{\sqrt{n}} \right)^{2} \right] \\
+ * \f]
+ *
+ * The optimal solutions to this multi-objective function lie in the
+ * range [-1/sqrt(3), 1/sqrt(3)].
+ *
+ * @tparam arma::mat Type of matrix to optimize.
+ */
+template<typename MatType = arma::mat>
+class FonsecaFlemingFunction
+{
+ private:
+  size_t numObjectives;
+  size_t numVariables;
+
+ public:
+  FonsecaFlemingFunction() : numObjectives(2), numVariables(3)
+  {/* Nothing to do here. */}
+
+  /**
+   * Evaluate the objectives with the given coordinate.
+   *
+   * @param coords The function coordinates.
+   * @return arma::Col<typename MatType::elem_type>
+   */
+  arma::Col<typename MatType::elem_type> Evaluate(const MatType& coords)
+  {
+    // Convenience typedef.
+    typedef typename MatType::elem_type ElemType;
+
+    arma::Col<ElemType> objectives(numObjectives);
+
+    objectives(0) = objectiveA.Evaluate(coords);
+    objectives(1) = objectiveB.Evaluate(coords);
+
+    return objectives;
+  }
+
+  //! Get the starting point.
+  MatType GetInitialPoint()
+  {
+    return arma::vec(numVariables, 1, arma::fill::zeros);
+  }
+
+  struct ObjectiveA
+  {
+    typename MatType::elem_type Evaluate(const MatType& coords)
+    {
+        return 1.0 - exp(
+             -pow(static_cast<double>(coords[0]) - 1.0 / sqrt(3.0), 2.0)
+             -pow(static_cast<double>(coords[1]) - 1.0 / sqrt(3.0), 2.0)
+             -pow(static_cast<double>(coords[2]) - 1.0 / sqrt(3.0), 2.0)
+        );
+    }
+  } objectiveA;
+
+  struct ObjectiveB
+  {
+    typename MatType::elem_type Evaluate(const MatType& coords)
+    {
+        return 1.0 - exp(
+            -pow(static_cast<double>(coords[0]) + 1.0 / sqrt(3.0), 2.0)
+            -pow(static_cast<double>(coords[1]) + 1.0 / sqrt(3.0), 2.0)
+            -pow(static_cast<double>(coords[2]) + 1.0 / sqrt(3.0), 2.0)
+        );
+    }
+  } objectiveB;
+
+  //! Get objective functions.
+  std::tuple<ObjectiveA, ObjectiveB> GetObjectives()
+  {
+    return std::make_tuple(objectiveA, objectiveB);
+  }
+};
+} // namespace test
+} // namespace ens
+
+#endif

include/ensmallen_bits/problems/problems.hpp

@@ -17,6 +17,7 @@
 #include "drop_wave_function.hpp"
 #include "easom_function.hpp"
 #include "eggholder_function.hpp"
+#include "fonseca_fleming_function.hpp"
 #include "fw_test_function.hpp"
 #include "generalized_rosenbrock_function.hpp"
 #include "goldstein_price_function.hpp"
@@ -30,6 +31,7 @@
 #include "rastrigin_function.hpp"
 #include "rosenbrock_function.hpp"
 #include "rosenbrock_wood_function.hpp"
+#include "schaffer_function_n1.hpp"
 #include "schaffer_function_n2.hpp"
 #include "schaffer_function_n4.hpp"
 #include "schwefel_function.hpp"

include/ensmallen_bits/problems/schaffer_function_n1.hpp

@@ -0,0 +1,94 @@
+/**
+ * @file schaffer_function_n1.hpp
+ * @author Sayan Goswami
+ *
+ * Implementation of Schaffer function N.1.
+ *
+ * ensmallen is free software; you may redistribute it and/or modify it under
+ * the terms of the 3-clause BSD license.  You should have received a copy of
+ * the 3-clause BSD license along with ensmallen.  If not, see
+ * http://www.opensource.org/licenses/BSD-3-Clause for more information.
+ */
+
+#ifndef ENSMALLEN_PROBLEMS_SCHAFFER_FUNCTION_N1_HPP
+#define ENSMALLEN_PROBLEMS_SCHAFFER_FUNCTION_N1_HPP
+
+namespace ens {
+namespace test {
+
+/**
+ * The Schaffer function N.1 is defined by
+ *
+ * \f[
+ * f_1(x) = x^2
+ * f_2(x) = (x-2)^2
+ * \f]
+ *
+ * The optimal solutions to this multi-objective function lie in the
+ * range [0, 2].
+ *
+ * @tparam arma::mat Type of matrix to optimize.
+ */
+template<typename MatType = arma::mat>
+class SchafferFunctionN1
+{
+ private:
+  size_t numObjectives;
+  size_t numVariables;
+
+ public:
+ //! Initialize the SchafferFunctionN1
+  SchafferFunctionN1() : numObjectives(2), numVariables(1)
+  {/* Nothing to do here. */}
+
+  /**
+   * Evaluate the objectives with the given coordinate.
+   *
+   * @param coords The function coordinates.
+   * @return arma::Col<typename MatType::elem_type>
+   */
+  arma::Col<typename MatType::elem_type> Evaluate(const MatType& coords)
+  {
+    // Convenience typedef.
+    typedef typename MatType::elem_type ElemType;
+
+    arma::Col<ElemType> objectives(numObjectives);
+
+    objectives(0) = std::pow(coords[0], 2);
+    objectives(1) = std::pow(coords[0] - 2, 2);
+
+    return objectives;
+  }
+
+  //! Get the starting point.
+  MatType GetInitialPoint()
+  {
+    return arma::vec(numVariables, 1, arma::fill::zeros);
+  }
+
+  struct ObjectiveA
+  {
+    typename MatType::elem_type Evaluate(const MatType& coords)
+    {
+        return std::pow(coords[0], 2);
+    }
+  } objectiveA;
+
+  struct ObjectiveB
+  {
+    typename MatType::elem_type Evaluate(const MatType& coords)
+    {
+        return std::pow(coords[0] - 2, 2);
+    }
+  } objectiveB;
+
+  //! Get objective functions.
+  std::tuple<ObjectiveA, ObjectiveB> GetObjectives()
+  {
+    return std::make_tuple(objectiveA, objectiveB);
+  }
+};
+} // namespace test
+} // namespace ens
+
+#endif

include/ensmallen_bits/problems/schaffer_function_n2.hpp

@@ -41,7 +41,7 @@ namespace test {
 class SchafferFunctionN2
 {
  public:
-  //! Initialize the SchafferFunctionN4.
+  //! Initialize the SchafferFunctionN2.
   SchafferFunctionN2();
 
   /**

include/ensmallen_bits/pso/pso.hpp

@@ -110,6 +110,14 @@ class PSOType
           initPolicy(initPolicy)
   { /* Nothing to do. */ }
 
+  /**
+   * Clean memory associated with the PSO object.
+   */
+  ~PSOType()
+  {
+    instUpdatePolicy.Clean();
+  }
+
   /**
    * Construct the particle swarm optimizer with the given function and
    * parameters. The defaults here are not necessarily good for the given

include/ensmallen_bits/pso/pso_impl.hpp

@@ -18,7 +18,9 @@
 #include <queue>
 
 namespace ens {
-/* After the velocity of each particle is updated at the end of each iteration
+
+/**
+ * After the velocity of each particle is updated at the end of each iteration
  * in PSO, the position of particle i (in iteration j) is updated as:
  *
  * \f[

include/ensmallen_bits/qhadam/qhadam.hpp

@@ -107,7 +107,8 @@ class QHAdam
            CallbackTypes&&... callbacks)
   {
     return optimizer.Optimize<SeparableFunctionType, MatType, GradType,
-        CallbackTypes...>(function, iterate, callbacks...);
+        CallbackTypes...>(function, iterate,
+        std::forward<CallbackTypes>(callbacks)...);
   }
 
   //! Forward the MatType as GradType.

include/ensmallen_bits/rmsprop/rmsprop.hpp

@@ -116,7 +116,8 @@ class RMSProp
            CallbackTypes&&... callbacks)
   {
     return optimizer.Optimize<SeparableFunctionType, MatType, GradType,
-        CallbackTypes...>(function, iterate, callbacks...);
+        CallbackTypes...>(function, iterate,
+        std::forward<CallbackTypes>(callbacks)...);
   }
 
   //! Forward the MatType as GradType.

include/ensmallen_bits/sa/sa.hpp

@@ -73,7 +73,7 @@ class SA
    * @param initMoveCoef Initial move size.
    * @param gain Proportional control in feedback move control.
    */
-  SA(CoolingScheduleType& coolingSchedule,
+  SA(const CoolingScheduleType& coolingSchedule = CoolingScheduleType(),
      const size_t maxIterations = 1000000,
      const double initT = 10000.,
      const size_t initMoves = 1000,
@@ -102,6 +102,11 @@ class SA
                                        MatType& iterate,
                                        CallbackTypes&&... callbacks);
 
+  //! Get the cooling schedule.
+  CoolingScheduleType CoolingSchedule() const { return coolingSchedule; }
+  //! Modify the cooling schedule.
+  CoolingScheduleType& CoolingSchedule() { return coolingSchedule; }
+
   //! Get the temperature.
   double Temperature() const { return temperature; }
   //! Modify the temperature.
@@ -139,7 +144,7 @@ class SA
 
  private:
   //! The cooling schedule being used.
-  CoolingScheduleType& coolingSchedule;
+  CoolingScheduleType coolingSchedule;
   //! The maximum number of iterations.
   size_t maxIterations;
   //! The current temperature.

include/ensmallen_bits/sa/sa_impl.hpp

@@ -18,7 +18,7 @@ namespace ens {
 
 template<typename CoolingScheduleType>
 SA<CoolingScheduleType>::SA(
-    CoolingScheduleType& coolingSchedule,
+    const CoolingScheduleType& coolingSchedule,
     const size_t maxIterations,
     const double initT,
     const size_t initMoves,

include/ensmallen_bits/sdp/lrsdp_function.hpp

@@ -58,6 +58,11 @@ class LRSDPFunction
                 const size_t numDenseConstraints,
                 const arma::Mat<typename SDPType::ElemType>& initialPoint);
 
+  /**
+   * Clean any memory associated with the LRSDPFunction.
+   */
+  ~LRSDPFunction();
+
   /**
    * Evaluate the objective function of the LRSDP (no constraints) at the given
    * coordinates.

include/ensmallen_bits/sdp/lrsdp_function_impl.hpp

@@ -49,6 +49,12 @@ LRSDPFunction<SDPType>::LRSDPFunction(
   }
 }
 
+template<typename SDPType>
+LRSDPFunction<SDPType>::~LRSDPFunction()
+{
+  rrt.Clean();
+}
+
 template<typename SDPType>
 template<typename MatType>
 typename MatType::elem_type LRSDPFunction<SDPType>::Evaluate(

include/ensmallen_bits/sdp/primal_dual.hpp

@@ -22,11 +22,7 @@ namespace ens {
  * PrimalDualSolver can optimize semidefinite programs.  For more details, see the
  * documentation on function types included with this distribution or on the
  * ensmallen website.
- *
- * @tparam DeprecatedSDPType Type of SDP to solve.  This parameter is deprecated
- *      and will be removed in ensmallen 2.10.0.
  */
-template<typename DeprecatedSDPType = SDP<arma::mat>>
 class PrimalDualSolver
 {
  public:
@@ -46,72 +42,6 @@ class PrimalDualSolver
                    const double primalInfeasTol = 1e-7,
                    const double dualInfeasTol = 1e-7);
 
-  /**
-   * Construct a new solver instance from a given SDP instance.  Uses a random,
-   * positive initialization point.
-   *
-   * This constructor is deprecated.  Use the constructor that does not take an
-   * SDPType and then call Optimize() with the SDP to be solved.
-   *
-   * This constructor will be removed in ensmallen 2.10.0.
-   *
-   * @param sdp Initialized SDP to be solved.
-   */
-  ens_deprecated PrimalDualSolver(const DeprecatedSDPType& sdp);
-
-  /**
-   * Construct a new solver instance with the given SDP instance and initial
-   * points for optimization.
-   *
-   * This constructor is deprecated.  Use the constructor that does not take an
-   * SDPType and then call Optimize() with the SDP to be solved.
-   *
-   * This constructor will be removed in ensmallen 2.10.0.
-   *
-   * @param sdp Initialized SDP to be solved.
-   * @param initialX Initial primal point for optimization.
-   * @param initialYSparse Initial y values for sparse constraints.
-   * @param initialYDense Initial y values for dense constraints.
-   * @param initialZ Initial dual point for optimization.
-   */
-  ens_deprecated PrimalDualSolver(const DeprecatedSDPType& sdp,
-                                  const arma::mat& initialX,
-                                  const arma::vec& initialYSparse,
-                                  const arma::vec& initialYDense,
-                                  const arma::mat& initialZ);
-
-  /**
-   * Optimize the stored SDP instance, storing the primal coordinates in X and
-   * returning the primal objective value.  Any initial point in X will be
-   * ignored.
-   *
-   * This function is deprecated and will be removed in ensmallen 2.10.0.  Use
-   * the overload of Optimize() that takes an SDP.
-   *
-   * @param X Matrix to store final primal coordinates for optimization.
-   * @return Primal objective value.
-   */
-  ens_deprecated double Optimize(arma::mat& X);
-
-  /**
-   * Optimize the stored SDP instance, storing the primal coordinates, dual
-   * coordinates, and sparse and dense y values into the given matrices and
-   * vectors.  The primal objective is returned.  Any initial setting of the
-   * given matrices will be ignored.
-   *
-   * This function is deprecated and will be rmeoved in ensmallen 2.10.0.  Use
-   * the overload of Optimize() that takes an SDP.
-   *
-   * @param X Matrix to store final primal coordinates into.
-   * @param ySparse Vector to store final sparse y values into.
-   * @param yDense Vector to store final dense y values into.
-   * @param Z Matrix to store final dual coordinates into.
-   */
-  ens_deprecated double Optimize(arma::mat& X,
-                                 arma::vec& ySparse,
-                                 arma::vec& yDense,
-                                 arma::mat& Z);
-
   /**
    * Optimize the given SDP with the given initial coordinates.  To get a set of
    * initial coordinates from an SDP class, consider calling
@@ -181,15 +111,6 @@ class PrimalDualSolver
   double& DualInfeasTol() { return dualInfeasTol; }
 
  private:
-  /**
-   * These are deprecated and will be removed in ensmallen 2.10.0.
-   */
-  DeprecatedSDPType deprecatedSDP;
-  arma::mat initialX;
-  arma::vec initialYSparse;
-  arma::vec initialYDense;
-  arma::mat initialZ;
-
   //! Maximum number of iterations to run. Set to 0 for no limit.
   size_t maxIterations;
 

include/ensmallen_bits/sdp/primal_dual_impl.hpp

@@ -33,52 +33,11 @@
 #include "lin_alg.hpp"
 
 namespace ens {
-
-template<typename DeprecatedSDPType>
-ens_deprecated
-PrimalDualSolver<DeprecatedSDPType>::PrimalDualSolver(
-    const DeprecatedSDPType& sdp) :
-    deprecatedSDP(sdp),
-    initialX(arma::eye<arma::mat>(sdp.N(), sdp.N())),
-    initialYSparse(arma::ones<arma::vec>(sdp.NumSparseConstraints())),
-    initialYDense(arma::ones<arma::vec>(sdp.NumDenseConstraints())),
-    initialZ(arma::eye<arma::mat>(sdp.N(), sdp.N())),
-    maxIterations(1000),
-    tau(0.99),
-    normXzTol(1e-7),
-    primalInfeasTol(1e-7),
-    dualInfeasTol(1e-7)
-{ /* Nothing to do. */ }
-
-template<typename DeprecatedSDPType>
-ens_deprecated
-PrimalDualSolver<DeprecatedSDPType>::PrimalDualSolver(
-    const DeprecatedSDPType& sdp,
-    const arma::mat& initialX,
-    const arma::vec& initialYSparse,
-    const arma::vec& initialYDense,
-    const arma::mat& initialZ) :
-    deprecatedSDP(sdp),
-    initialX(initialX),
-    initialYSparse(initialYSparse),
-    initialYDense(initialYDense),
-    initialZ(initialZ),
-    maxIterations(1000),
-    tau(0.99),
-    normXzTol(1e-7),
-    primalInfeasTol(1e-7),
-    dualInfeasTol(1e-7)
-{
-  // Nothing to do.
-}
-
-template<typename DeprecatedSDPType>
-PrimalDualSolver<DeprecatedSDPType>::PrimalDualSolver(
-    const size_t maxIterations,
-    const double tau,
-    const double normXzTol,
-    const double primalInfeasTol,
-    const double dualInfeasTol) :
+inline PrimalDualSolver::PrimalDualSolver(const size_t maxIterations,
+                                          const double tau,
+                                          const double normXzTol,
+                                          const double primalInfeasTol,
+                                          const double dualInfeasTol) :
     maxIterations(maxIterations),
     tau(tau),
     normXzTol(normXzTol),
@@ -227,46 +186,8 @@ SolveKKTSystem(const SparseConstraintType& aSparse,
   dsZ = rd - subTerm;
 }
 
-template<typename DeprecatedSDPType>
-ens_deprecated
-double PrimalDualSolver<DeprecatedSDPType>::Optimize(
-    arma::mat& coordinates)
-{
-  // Use the internally-held initial parameters.
-  coordinates = initialX;
-  arma::mat ySparse = initialYSparse;
-  arma::mat yDense = initialYDense;
-  arma::mat Z = initialZ;
-  return Optimize<DeprecatedSDPType, arma::mat>(deprecatedSDP, coordinates,
-      ySparse, yDense, Z);
-}
-
-template<typename DeprecatedSDPType>
-ens_deprecated
-double PrimalDualSolver<DeprecatedSDPType>::Optimize(
-    arma::mat& coordinates,
-    arma::vec& ySparse,
-    arma::vec& yDense,
-    arma::mat& z)
-{
-  // Initialize internally then call the other overload.
-  coordinates = initialX;
-  arma::mat ySparseMat = initialYSparse;
-  arma::mat yDenseMat = initialYDense;
-  z = initialZ;
-
-  const double result = Optimize<DeprecatedSDPType, arma::mat>(deprecatedSDP,
-      coordinates, ySparseMat, yDenseMat, z);
-
-  ySparse = ySparseMat.col(0);
-  yDense = yDenseMat.col(0);
-
-  return result;
-}
-
-template<typename DeprecatedSDPType>
 template<typename SDPType, typename MatType, typename... CallbackTypes>
-typename MatType::elem_type PrimalDualSolver<DeprecatedSDPType>::Optimize(
+typename MatType::elem_type PrimalDualSolver::Optimize(
     const SDPType& sdp,
     MatType& coordinates,
     CallbackTypes&&... callbacks)
@@ -279,9 +200,8 @@ typename MatType::elem_type PrimalDualSolver<DeprecatedSDPType>::Optimize(
   return Optimize(sdp, coordinates, ySparse, yDense, z, callbacks...);
 }
 
-template<typename DeprecatedSDPType>
 template<typename SDPType, typename MatType, typename... CallbackTypes>
-typename MatType::elem_type PrimalDualSolver<DeprecatedSDPType>::Optimize(
+typename MatType::elem_type PrimalDualSolver::Optimize(
     const SDPType& sdp,
     MatType& coordinates,
     MatType& ySparse,

include/ensmallen_bits/sgd/sgd.hpp

@@ -101,6 +101,11 @@ class SGD
       const bool resetPolicy = true,
       const bool exactObjective = false);
 
+  /**
+   * Clean any memory associated with the SGD object.
+   */
+  ~SGD();
+
   /**
    * Optimize the given function using stochastic gradient descent.  The given
    * starting point will be modified to store the finishing point of the

include/ensmallen_bits/sgd/sgd_impl.hpp

@@ -44,6 +44,14 @@ SGD<UpdatePolicyType, DecayPolicyType>::SGD(
     isInitialized(false)
 { /* Nothing to do. */ }
 
+template<typename UpdatePolicyType, typename DecayPolicyType>
+SGD<UpdatePolicyType, DecayPolicyType>::~SGD()
+{
+  // Clean decay and update policies, if they were initialized.
+  instDecayPolicy.Clean();
+  instUpdatePolicy.Clean();
+}
+
 //! Optimize the function (minimize).
 template<typename UpdatePolicyType, typename DecayPolicyType>
 template<typename SeparableFunctionType,

include/ensmallen_bits/smorms3/smorms3.hpp

@@ -101,7 +101,8 @@ class SMORMS3
     // TODO: disallow sp_mat
 
     return optimizer.Optimize<SeparableFunctionType, MatType, GradType,
-        CallbackTypes...>(function, iterate, callbacks...);
+        CallbackTypes...>(function, iterate,
+        std::forward<CallbackTypes>(callbacks)...);
   }
 
   //! Forward the MatType as GradType.

include/ensmallen_bits/spalera_sgd/spalera_sgd.hpp

@@ -112,6 +112,11 @@ class SPALeRASGD
              const bool resetPolicy = true,
              const bool exactObjective = false);
 
+  /**
+   * Clean any memory associated with the SPALeRA SGD object.
+   */
+  ~SPALeRASGD();
+
   /**
    * Optimize the given function using SPALeRA SGD.  The given starting point
    * will be modified to store the finishing point of the algorithm, and the

include/ensmallen_bits/spalera_sgd/spalera_sgd_impl.hpp

@@ -45,6 +45,13 @@ SPALeRASGD<DecayPolicyType>::SPALeRASGD(const double stepSize,
     isInitialized(false)
 { /* Nothing to do. */ }
 
+template<typename DecayPolicyType>
+SPALeRASGD<DecayPolicyType>::~SPALeRASGD()
+{
+  instUpdatePolicy.Clean();
+  instDecayPolicy.Clean();
+}
+
 //! Optimize the function (minimize).
 template<typename DecayPolicyType>
 template<typename SeparableFunctionType,

include/ensmallen_bits/svrg/svrg.hpp

@@ -120,6 +120,11 @@ class SVRGType
            const bool resetPolicy = true,
            const bool exactObjective = false);
 
+  /**
+   * Clean any memory associated with the SVRGType object.
+   */
+  ~SVRGType();
+
   /**
    * Optimize the given function using SVRG. The given starting point will be
    * modified to store the finishing point of the algorithm, and the final

include/ensmallen_bits/svrg/svrg_impl.hpp

@@ -42,6 +42,13 @@ SVRGType<UpdatePolicyType, DecayPolicyType>::SVRGType(
     isInitialized(false)
 { /* Nothing to do. */ }
 
+template<typename UpdatePolicyType, typename DecayPolicyType>
+SVRGType<UpdatePolicyType, DecayPolicyType>::~SVRGType()
+{
+  instUpdatePolicy.Clean();
+  instDecayPolicy.Clean();
+}
+
 //! Optimize the function (minimize).
 template<typename UpdatePolicyType, typename DecayPolicyType>
 template<typename SeparableFunctionType,

include/ensmallen_bits/swats/swats.hpp

@@ -105,7 +105,8 @@ class SWATS
            CallbackTypes&&... callbacks)
   {
     return optimizer.Optimize<SeparableFunctionType, MatType, GradType,
-        CallbackTypes...>(function, iterate, callbacks...);
+        CallbackTypes...>(function, iterate,
+        std::forward<CallbackTypes>(callbacks)...);
   }
 
   //! Forward the MatType as GradType.

include/ensmallen_bits/utility/arma_traits.hpp

@@ -98,6 +98,27 @@ struct IsArmaType<arma::SpSubview<eT> >
   const static bool value = true;
 };
 
+
+#if ((ARMA_VERSION_MAJOR >= 10) || \
+    ((ARMA_VERSION_MAJOR == 9) && (ARMA_VERSION_MINOR >= 869)))
+
+  // Armadillo 9.869+ has SpSubview_col and SpSubview_row
+
+  template<typename eT>
+  struct IsArmaType<arma::SpSubview_col<eT> >
+  {
+    const static bool value = true;
+  };
+
+  template<typename eT>
+  struct IsArmaType<arma::SpSubview_row<eT> >
+  {
+    const static bool value = true;
+  };
+
+#endif
+
+
 // template<>
 template<typename eT>
 struct IsArmaType<arma::Mat<eT> >
@@ -141,4 +162,4 @@ struct tuple_element<N, T0, T...> {
 
 } // namespace ens
 
-#endif
+#endif

include/ensmallen_bits/wn_grad/wn_grad.hpp

@@ -94,7 +94,8 @@ class WNGrad
            CallbackTypes&&... callbacks)
   {
     return optimizer.Optimize<SeparableFunctionType, MatType, GradType,
-        CallbackTypes...>(function, iterate, callbacks...);
+        CallbackTypes...>(function, iterate,
+        std::forward<CallbackTypes>(callbacks)...);
   }
 
   //! Forward the MatType as GradType.

rel/ensmallen-release.sh

@@ -1,68 +0,0 @@
-#!/usr/bin/env bash
-#
-# Release a new version of ensmallen.
-#
-# Arguments:
-#   $ ensmallen-release.sh <major> <minor> <patch> [<name>]
-#
-# This should be run from the root of the repository.
-#
-# Make sure to update HISTORY.md manually first!
-set -e
-
-if [ "$#" -lt 3 ]; then
-  echo "At least three arguments required!";
-  echo "$ ensmallen-release.sh <major> <minor> <patch> [<name>]";
-  exit 1;
-fi
-
-if [ "$#" -gt 4 ]; then
-  echo "Too many arguments!"
-  echo "$ ensmallen-release.sh <major> <minor> <patch> [<name>]";
-  exit 1;
-fi
-
-lines=`git diff | wc -l`;
-if [ "$lines" != "0" ]; then
-  echo "git diff returned a nonzero result!"
-  git diff
-  exit 1;
-fi
-
-MAJOR=$1;
-MINOR=$2;
-PATCH=$3;
-
-sed -i 's/ENS_VERSION_MAJOR[ ]*[0-9]*$/ENS_VERSION_MAJOR '$MAJOR'/' include/ensmallen_bits/ens_version.hpp;
-sed -i 's/ENS_VERSION_MINOR[ ]*[0-9]*$/ENS_VERSION_MINOR '$MINOR'/' include/ensmallen_bits/ens_version.hpp;
-sed -i 's/ENS_VERSION_PATCH[ ]*[0-9]*$/ENS_VERSION_PATCH '$PATCH'/' include/ensmallen_bits/ens_version.hpp;
-
-if [ "$#" -eq "4" ]; then
-  sed -i 's/ENS_VERSION_NAME[ ]*\".*\"$/ENS_VERSION_NAME \"'"$4"'\"/' include/ensmallen_bits/ens_version.hpp;
-fi
-
-# update CONTRIBUTING.md
-sed -i "s/ensmallen-[0-9]*\.[0-9]*\.[0-9]*/ensmallen-$MAJOR.$MINOR.$PATCH/g" CONTRIBUTING.md;
-
-git pull
-git add include/ensmallen_bits/ens_version.hpp;
-git add CONTRIBUTING.md
-git commit -m "Update and release version $MAJOR.$MINOR.$PATCH.";
-git tag $MAJOR.$MINOR.$PATCH;
-git push origin $MAJOR.$MINOR.$PATCH;
-git push origin master;
-
-git clone https://github.com/mlpack/ensmallen.org /tmp/ensmallen.org/;
-git archive --prefix=ensmallen-$MAJOR.$MINOR.$PATCH/ $MAJOR.$MINOR.$PATCH | gzip > /tmp/ensmallen.org/files/ensmallen-$MAJOR.$MINOR.$PATCH.tar.gz;
-cd /tmp/ensmallen.org/;
-git add files/ensmallen-$MAJOR.$MINOR.$PATCH.tar.gz;
-cd files/;
-rm ensmallen-latest.tar.gz;
-ln -s ensmallen-$MAJOR.$MINOR.$PATCH.tar.gz ensmallen-latest.tar.gz;
-cd ../
-git add files/ensmallen-latest.tar.gz;
-git commit -m "Release version $MAJOR.$MINOR.$PATCH.";
-git push origin;
-cd -
-
-rm -rf /tmp/ensmallen.org;

scripts/ensmallen-release.sh

@@ -0,0 +1,182 @@
+#!/usr/bin/env bash
+#
+# Release a new version of ensmallen.
+#
+# Arguments:
+#   $ ensmallen-release.sh <github username> <major> <minor> <patch> [<name>]
+#
+# This should be run from the root of the repository.
+set -e
+
+if [ "$#" -lt 4 ]; then
+  echo "At least four arguments required!";
+  echo "$ ensmallen-release.sh <github username> <major> <minor> <patch>" \
+      "[<name>]";
+  exit 1;
+fi
+
+if [ "$#" -gt 5 ]; then
+  echo "Too many arguments!";
+  echo "$ ensmallen-release.sh <github username> <major> <minor> <patch>" \
+      "[<name>]";
+  exit 1;
+fi
+
+# Make sure that the branch is clean.
+# Truncate leading whitespaces since wc -l on MacOS adds an extra \t.
+lines=`git diff | wc -l | sed -e 's/^\s*//g'`;
+if [ "$lines" != "0" ]; then
+  echo "git diff returned a nonzero result!";
+  echo "";
+  git diff;
+  exit 1;
+fi
+
+# Make sure that the mlpack repository exists.
+dest_remote_name=`git remote -v |\
+                  grep "mlpack/ensmallen (fetch)" |\
+                  head -1 |\
+                  awk -F' ' '{ print $1 }'`;
+
+if [ "a$dest_remote_name" == "a" ]; then
+  echo "No git remote found for https://github.com/mlpack/ensmallen!";
+  echo "Make sure that you've got the ensmallen repository as a remote, and" \
+      "that the master branch from that remote is checked out.";
+  echo "You can do this with a fresh repository via \`git clone" \
+      "https://github.com/mlpack/ensmallen\`.";
+  exit 1;
+fi
+
+# Also check that we're on the master branch, from the correct origin.
+current_branch=`git branch --no-color | grep '^\* ' | awk -F' ' '{ print $2 }'`;
+current_origin=`git rev-parse --abbrev-ref --symbolic-full-name @{u} |\
+                awk -F'/' '{ print $1 }'`;
+if [ "a$current_branch" != "amaster" ]; then
+  echo "Current branch is $current_branch."
+  echo "This script has to be run from the master branch."
+  exit 1;
+elif [ "a$current_origin" != "a$dest_remote_name" ]; then
+  echo "Current branch does not track from remote mlpack repository!";
+  echo "Instead, it tracks from $current_origin/master.";
+  echo "Make sure to check out a branch that tracks $dest_remote_name/master.";
+  exit 1;
+fi
+
+# Make sure 'gh' is installed.
+hub_output="`which hub`" || true;
+if [ "a$hub_output" == "a" ]; then
+  echo "The Hub command-line tool must be installed for this script to run" \
+      "successfully.";
+  echo "See https://hub.github.com for more details and installation" \
+      "instructions.";
+  echo "";
+  echo "(apt-get install hub on Debian and Ubuntu)";
+  echo "(brew install hub via Homebrew)";
+  exit 1;
+fi
+
+# Check git remotes: we need to make sure we have a fork to push to.
+github_user=$1;
+remote_name=`git remote -v |\
+             grep "$github_user/ensmallen (push)" |\
+             head -1 |\
+             awk -F' ' '{ print $1 }'`;
+if [ "a$remote_name" == "a" ]; then
+  echo "No git remote found for $github_user/ensmallen!";
+  echo "Adding remote '$github_user'.";
+  git remote add $github_user https://github.com/$github_user/ensmallen;
+  remote_name="$github_user";
+fi
+git fetch $github_user;
+
+# Make sure everything is up to date.
+git pull;
+
+# Make updates to files that will be needed for the release.
+MAJOR=$2;
+MINOR=$3;
+PATCH=$4;
+
+sed -i 's/ENS_VERSION_MAJOR[ ]*[0-9]*$/ENS_VERSION_MAJOR '$MAJOR'/' \
+    include/ensmallen_bits/ens_version.hpp;
+sed -i 's/ENS_VERSION_MINOR[ ]*[0-9]*$/ENS_VERSION_MINOR '$MINOR'/' \
+    include/ensmallen_bits/ens_version.hpp;
+sed -i 's/ENS_VERSION_PATCH[ ]*[0-9]*$/ENS_VERSION_PATCH '$PATCH'/' \
+    include/ensmallen_bits/ens_version.hpp;
+
+if [ "$#" -eq "5" ]; then
+  sed -i 's/ENS_VERSION_NAME[ ]*\".*\"$/ENS_VERSION_NAME \"'"$5"'\"/' \
+      include/ensmallen_bits/ens_version.hpp;
+fi
+
+# Update CONTRIBUTING.md.
+sed -i "s/ensmallen-[0-9]*\.[0-9]*\.[0-9]*/ensmallen-$MAJOR.$MINOR.$PATCH/g" \
+    CONTRIBUTING.md;
+
+# Update HISTORY.md with the release date and possibly name.
+version_name=`grep ENS_VERSION_NAME include/ensmallen_bits/ens_version.hpp |\
+              head -1 |\
+              sed 's/.*\"\(.*\)\"/\1/'`;
+year=`date +%Y`;
+month=`date +%m`;
+day=`date +%d`;
+new_line="ensmallen $MAJOR.$MINOR.$PATCH: \"$version_name\"";
+sed -i "s/### ensmallen ?.??.?: \"???\"/### $new_line/" HISTORY.md;
+sed -i "s/###### ????-??-??/###### $year-$month-$day/" HISTORY.md;
+
+# Update date in ens_version.hpp
+sed -i 's/ENS_VERSION_YEAR[ ]*\".*\"$/ENS_VERSION_YEAR \"'"$year"'\"/' \
+    include/ensmallen_bits/ens_version.hpp;
+sed -i 's/ENS_VERSION_MONTH[ ]*\".*\"$/ENS_VERSION_MONTH \"'"$month"'\"/' \
+    include/ensmallen_bits/ens_version.hpp;
+sed -i 's/ENS_VERSION_DAY[ ]*\".*\"$/ENS_VERSION_DAY \"'"$day"'\"/' \
+    include/ensmallen_bits/ens_version.hpp;
+
+# Now, we'll do all this on a new release branch.
+git checkout -b release-$MAJOR.$MINOR.$PATCH;
+
+git add include/ensmallen_bits/ens_version.hpp;
+git add CONTRIBUTING.md;
+git add HISTORY.md;
+git commit -m "Update and release version $MAJOR.$MINOR.$PATCH.";
+
+changelog_str=`cat HISTORY.md |\
+    awk '/^### /{f=0} /^### ensmallen '"$MAJOR"'.'"$MINOR"'.'"$PATCH"': "'"$version_name"'"/{f=1} f{print}' |\
+    grep -v '^#' |\
+    tr '\n' '!' |\
+    sed -e 's/!  [ ]*/ /g' |\
+    tr '!' '\n'`;
+echo "Changelog string:"
+echo "$changelog_str"
+
+# Add one more commit to create the new HISTORY block.
+echo "### ensmallen ?.??.?: \"???\"" > HISTORY.md.new;
+echo "###### ????-??-??" >> HISTORY.md.new;
+echo "" >> HISTORY.md.new;
+cat HISTORY.md >> HISTORY.md.new;
+mv HISTORY.md.new HISTORY.md;
+git add HISTORY.md;
+git commit -m "Add new block for next release to HISTORY.md.";
+
+# Push to new branch.
+git push --set-upstream $github_user release-$MAJOR.$MINOR.$PATCH;
+
+# Next, we have to actually open the PR for the release.  These lines would be
+# hard to wrap so they are longer than the length limit. :)
+hub pull-request \
+    -b mlpack:master \
+    -h $github_user:release-$MAJOR.$MINOR.$PATCH \
+    -m "Release version $MAJOR.$MINOR.$PATCH: \"$version_name\"" \
+    -m "This automatically-generated pull request adds the commits necessary to make the $MAJOR.$MINOR.$PATCH release." \
+    -m "Once the PR is merged, mlpack-bot will tag the release as HEAD~1 (so that it doesn't include the new HISTORY block) and publish it." \
+    -m "Or, well, hopefully that will happen someday." \
+    -m "When you merge this PR, be sure to merge it using a *rebase*." \
+    -m "### Changelog" \
+    -m "$changelog_str" \
+    -l "t: release"
+
+echo "";
+echo "Switching back to 'master' branch.";
+echo "If you want to access the release branch again, use \`git checkout " \
+    "release-$MAJOR.$MINOR.$PATCH\`.";
+exit 0;

scripts/history-update-check.sh

@@ -0,0 +1,20 @@
+#!/usr/bin/env bash
+#
+# Check each PR has an entry in HISTORY.md during a CI routine.
+#
+# Arguments:
+#   $ history-update-check.sh
+#
+# This should be run from the root of the repository.
+res=$(git diff origin/master --name-only | grep ^HISTORY.md | wc -l)
+echo "Files Changed:"
+git diff --name-only
+if [ $res -gt 0 ]; then
+    echo "HISTORY.md was updated with a change for the PR ..."
+else
+    echo "Please describe your PR changes in HISTORY.md ..."
+    echo "Exiting CI process ... "
+    exit 1
+fi
+
+exit 0

scripts/update-website-after-release.sh

@@ -0,0 +1,60 @@
+#!/bin/bash
+#
+# This script is used to update the website after a release is made.  Push
+# access to the ensmallen.org website is needed.  Generally, this script will be
+# run by mlpack-bot, so it never needs to be run by hand.
+#
+# Usage: update-website-after-release.sh <major> <minor> <patch>
+
+MAJOR=$1;
+MINOR=$2;
+PATCH=$3;
+
+# Make sure that the mlpack repository exists.
+dest_remote_name=`git remote -v |\
+                  grep "mlpack/ensmallen (fetch)" |\
+                  head -1 |\
+                  awk -F' ' '{ print $1 }'`;
+
+if [ "a$dest_remote_name" == "a" ]; then
+  echo "No git remote found for mlpack/ensmallen!";
+  echo "Make sure that you've got the ensmallen repository as a remote, and" \
+      "that the master branch from that remote is checked out.";
+  echo "You can do this with a fresh repository via \`git clone" \
+      "https://github.com/mlpack/ensmallen\`.";
+  exit 1;
+fi
+
+# Update the checked out repository, so that we can get the tags.
+git fetch $dest_remote_name;
+
+# Check out a copy of the ensmallen.org repository.
+git clone git@github.com:mlpack/ensmallen.org /tmp/ensmallen.org/;
+
+# Create the release file.
+git archive --prefix=ensmallen-$MAJOR.$MINOR.$PATCH/ $MAJOR.$MINOR.$PATCH |\
+    gzip > /tmp/ensmallen.org/files/ensmallen-$MAJOR.$MINOR.$PATCH.tar.gz;
+
+# Now update the website.
+wd=`pwd`;
+cd /tmp/ensmallen.org/;
+git add files/ensmallen-$MAJOR.$MINOR.$PATCH.tar.gz;
+
+# Update the link to the latest version.
+cd files/;
+rm ensmallen-latest.tar.gz;
+ln -s ensmallen-$MAJOR.$MINOR.$PATCH.tar.gz ensmallen-latest.tar.gz;
+cd ../
+
+# Update the index page.
+sed -i 's/\[ensmallen-[0-9]*\.[0-9]*\.[0-9]\.tar\.gz\](files\/ensmallen-[0-9]*\.[0-9]*\.[0-9]*\.tar.gz)/[ensmallen-'$MAJOR'.'$MINOR'.'$PATCH'.tar.gz](files\/ensmallen-'$MAJOR'.'$MINOR'.'$PATCH'.tar.gz)/' index.md
+
+git add files/ensmallen-latest.tar.gz;
+git add index.md;
+git commit -m "Release version $MAJOR.$MINOR.$PATCH.";
+
+# Finally, push, and we're done.
+git push origin;
+cd $wd;
+
+rm -rf /tmp/ensmallen.org;

tests/CMakeLists.txt

@@ -1,5 +1,4 @@
-project(ensmallen_tests CXX)
-
+# The tests that need to be compiled.
 set(ENSMALLEN_TESTS_SOURCES
     main.cpp
     ada_bound_test.cpp
@@ -26,6 +25,7 @@ set(ENSMALLEN_TESTS_SOURCES
     lrsdp_test.cpp
     momentum_sgd_test.cpp
     nesterov_momentum_sgd_test.cpp
+    nsga2_test.cpp
     parallel_sgd_test.cpp
     proximal_test.cpp
     pso_test.cpp
@@ -47,17 +47,16 @@ set(ENSMALLEN_TESTS_SOURCES
 )
 
 set(CMAKE_RUNTIME_OUTPUT_DIRECTORY ${CMAKE_BINARY_DIR})
-add_executable(${PROJECT_NAME} ${ENSMALLEN_TESTS_SOURCES})
-
-target_link_libraries(${PROJECT_NAME} ${ARMADILLO_LIBRARIES})
+add_executable(ensmallen_tests ${ENSMALLEN_TESTS_SOURCES})
+target_link_libraries(ensmallen_tests PRIVATE ensmallen)
 
 # Copy test data into place.
-add_custom_command(TARGET ${PROJECT_NAME}
+add_custom_command(TARGET ensmallen_tests
   POST_BUILD
   COMMAND ${CMAKE_COMMAND} -E copy_directory ${CMAKE_CURRENT_SOURCE_DIR}/data/
       ${CMAKE_BINARY_DIR}/data/
 )
 
 enable_testing()
-add_test(NAME ${PROJECT_NAME} COMMAND ${PROJECT_NAME}
-WORKING_DIRECTORY ${CMAKE_BINARY_DIR})
+add_test(NAME ensmallen_tests COMMAND ensmallen_tests
+WORKING_DIRECTORY ${CMAKE_BINARY_DIR})

tests/ada_delta_test.cpp

@@ -65,11 +65,20 @@ TEST_CASE("AdaDeltaLogisticRegressionTest", "[AdaDeltaTest]")
  */
 TEST_CASE("SimpleAdaDeltaTestFunctionFMat", "[AdaDeltaTest]")
 {
+  size_t trials = 3;
   SGDTestFunction f;
-  AdaDelta optimizer(1.0, 1, 0.05, 1e-6, 5000000, 1e-15, true, true);
+  arma::fmat coordinates;
 
-  arma::fmat coordinates = f.GetInitialPoint<arma::fmat>();
-  optimizer.Optimize(f, coordinates);
+  for (size_t i = 0; i < trials; ++i)
+  {
+    coordinates = f.GetInitialPoint<arma::fmat>();
+
+    AdaDelta optimizer(2.0, 1, 0.05, 1e-6, 5000000, 1e-8, true, true);
+    optimizer.Optimize(f, coordinates);
+
+    if (arma::max(arma::vectorise(arma::abs(coordinates))) < 0.01f)
+      break;
+  }
 
   REQUIRE(coordinates(0) == Approx(0.0f).margin(0.01));
   REQUIRE(coordinates(1) == Approx(0.0f).margin(0.01));

tests/ada_grad_test.cpp

@@ -63,11 +63,20 @@ TEST_CASE("AdaGradLogisticRegressionTest", "[AdaGradTest]")
  */
 TEST_CASE("SimpleAdaGradTestFunctionFMat", "[AdaGradTest]")
 {
+  size_t trials = 3;
   SGDTestFunction f;
-  AdaGrad optimizer(0.99, 1, 1e-8, 5000000, 1e-9, true);
+  arma::fmat coordinates;
 
-  arma::fmat coordinates = f.GetInitialPoint<arma::fmat>();
-  optimizer.Optimize(f, coordinates);
+  for (size_t i = 0; i < trials; ++i)
+  {
+    coordinates = f.GetInitialPoint<arma::fmat>();
+
+    AdaGrad optimizer(0.99, 1, 1e-8, 5000000, 1e-9, true);
+    optimizer.Optimize(f, coordinates);
+
+    if (arma::max(arma::vectorise(arma::abs(coordinates))) < 0.01f)
+      break;
+  }
 
   REQUIRE(coordinates(0) == Approx(0.0f).margin(0.01));
   REQUIRE(coordinates(1) == Approx(0.0f).margin(0.01));

tests/callbacks_test.cpp

@@ -152,6 +152,105 @@ void CallbacksFullFunctionTest(OptimizerType& optimizer,
   REQUIRE(cb.calledStepTaken == calledStepTaken);
 }
 
+template<typename OptimizerType>
+void CallbacksFullMultiobjectiveFunctionTest(OptimizerType& optimizer,
+                                             bool calledEvaluate,
+                                             bool calledGradient,
+                                             bool calledBeginEpoch,
+                                             bool calledEndEpoch,
+                                             bool calledBeginOptimization,
+                                             bool calledEndOptimization,
+                                             bool calledEvaluateConstraint,
+                                             bool calledGradientConstraint,
+                                             bool calledStepTaken)
+{
+  SchafferFunctionN1<arma::mat> SCH;
+
+  typedef decltype(SCH.objectiveA) ObjectiveTypeA;
+  typedef decltype(SCH.objectiveB) ObjectiveTypeB;
+
+  CompleteCallbackTestFunction cb;
+
+  arma::mat coordinates = SCH.GetInitialPoint();
+  std::tuple<ObjectiveTypeA, ObjectiveTypeB> objectives = SCH.GetObjectives();
+
+  optimizer.Optimize(objectives, coordinates, cb);
+
+  REQUIRE(cb.calledEvaluate == calledEvaluate);
+  REQUIRE(cb.calledGradient == calledGradient);
+  REQUIRE(cb.calledBeginEpoch == calledBeginEpoch);
+  REQUIRE(cb.calledEndEpoch == calledEndEpoch);
+  REQUIRE(cb.calledBeginOptimization == calledBeginOptimization);
+  REQUIRE(cb.calledEndOptimization == calledEndOptimization);
+  REQUIRE(cb.calledEvaluateConstraint == calledEvaluateConstraint);
+  REQUIRE(cb.calledGradientConstraint == calledGradientConstraint);
+  REQUIRE(cb.calledStepTaken == calledStepTaken);
+}
+
+template<typename OptimizerType>
+void EarlyStopCallbacksLambdaFunctionTest(OptimizerType& optimizer)
+{
+  arma::mat data, testData, shuffledData;
+  arma::Row<size_t> responses, testResponses, shuffledResponses;
+
+  LogisticRegressionTestData(data, testData, shuffledData,
+      responses, testResponses, shuffledResponses);
+
+  LogisticRegression<> lr(shuffledData, shuffledResponses, 0.5);
+  arma::mat coordinates = lr.GetInitialPoint();
+
+  EarlyStopAtMinLoss cb(
+      [&](const arma::mat& /* coordinates */)
+      {
+        return lr.ComputeAccuracy(testData, testResponses,
+                                  coordinates);
+      });
+
+  optimizer.Optimize(lr, coordinates, cb);
+}
+
+TEST_CASE("EarlyStopAtMinLossLambdaCallbackTest", "[CallbacksTest]")
+{
+  SMORMS3 smorms3;
+  EarlyStopCallbacksLambdaFunctionTest(smorms3);
+}
+
+TEST_CASE("EarlyStopAtMinLossCustomLambdaTest", "[CallbacksTest]")
+{
+  // Use the 50-dimensional Rosenbrock function.
+  GeneralizedRosenbrockFunction f(50);
+  // Start at some really large point.
+  arma::mat coordinates = f.GetInitialPoint();
+  coordinates.fill(100.0);
+
+  EarlyStopAtMinLoss cb(
+      [&](const arma::mat& coordinates)
+      {
+        // Terminate if any coordinate has a value less than 10.
+        double minValue = arma::abs(coordinates).min();
+        return (minValue < 10.0) ?
+          std::numeric_limits<double>::max() : minValue;
+      });
+
+  SMORMS3 smorms3;
+  smorms3.Optimize(f, coordinates, cb);
+
+  // Make sure that we did not get to the optimum.
+  for (size_t i = 0; i < coordinates.n_elem; ++i)
+    REQUIRE(std::abs(coordinates[i]) >= 3.0);
+}
+
+/**
+ * Make sure we invoke all callbacks (AdaBound).
+ */
+TEST_CASE("AdaBoundCallbacksFullFunctionTest", "[CallbacksTest]")
+{
+  AdaBound optimizer(0.001, 2, 0.1, 1e-3, 0.9, 0.999, 1e-8, 1000,
+      1e-3, false);
+  CallbacksFullFunctionTest(optimizer, true, true, true, true, true, true,
+      false, false, true);
+}
+
 /**
  * Make sure we invoke all callbacks (AdaDelta).
  */
@@ -272,6 +371,59 @@ TEST_CASE("KatyushaCallbacksFullFunctionTest", "[CallbacksTest]")
       false, false, true);
 }
 
+/**
+ * Make sure we invoke all callbacks (NSGA2).
+ */
+TEST_CASE("NSGA2CallbacksFullFunctionTest", "[CallbackTest]")
+{
+  arma::vec lowerBound = {-1000};
+  arma::vec upperBound = {1000};
+  NSGA2 optimizer(20, 5000, 0.5, 0.5, 1e-3, 1e-6, lowerBound, upperBound);
+  CallbacksFullMultiobjectiveFunctionTest(optimizer, false, false, false, false,
+      true, true, false, false, true);
+}
+
+/**
+ * Make sure we invoke all callbacks (Lookahead).
+ */
+TEST_CASE("LookaheadCallbacksFullFunctionTest", "[CallbacksTest]")
+{
+  Adam adam(0.001, 1, 0.9, 0.999, 1e-8, 100, 1e-10, false, true);
+  Lookahead<Adam> optimizer(adam, 0.5, 1000, 10, -10, NoDecay(),
+      false, true);
+  CallbacksFullFunctionTest(optimizer, true, true, true, true, true, true,
+      false, false, true);
+}
+
+/**
+ * Make sure we invoke all callbacks (Padam).
+ */
+TEST_CASE("PadamCallbacksFullFunctionTest", "[CallbacksTest]")
+{
+  Padam optimizer(1e-2, 1, 0.9, 0.99, 0.25, 1e-5, 1000);
+  CallbacksFullFunctionTest(optimizer, true, true, true, true, true, true,
+      false, false, true);
+}
+
+/**
+ * Make sure we invoke all callbacks (QHAdam).
+ */
+TEST_CASE("QHAdamCallbacksFullFunctionTest", "[CallbacksTest]")
+{
+  QHAdam optimizer(0.02, 2, 0.6, 0.9, 0.9, 0.999, 1e-8, 1000, 1e-7, true);
+  CallbacksFullFunctionTest(optimizer, true, true, true, true, true, true,
+      false, false, true);
+}
+
+/**
+ * Make sure we invoke all callbacks (RMSProp).
+ */
+TEST_CASE("RMSPropCallbacksFullFunctionTest", "[CallbacksTest]")
+{
+  RMSProp optimizer(1e-3, 1, 0.99, 1e-8, 1000, 1e-9, true);
+  CallbacksFullFunctionTest(optimizer, true, true, true, true, true, true,
+      false, false, true);
+}
 /**
  * Make sure we invoke all callbacks (SARAH).
  */
@@ -312,6 +464,16 @@ TEST_CASE("SGDRCallbacksFullFunctionTest", "[CallbacksTest]")
       false, false, true);
 }
 
+/**
+ * Make sure we invoke all callbacks (SMORMS3).
+ */
+TEST_CASE("SMORMS3CallbacksFullFunctionTest", "[CallbacksTest]")
+{
+  SMORMS3 optimizer(0.001, 1, 1e-16, 1000, 1e-9, true);
+  CallbacksFullFunctionTest(optimizer, true, true, true, true, true, true,
+      false, false, true);
+}
+
 /**
  * Make sure we invoke all callbacks (SPALeRASGD).
  */
@@ -342,6 +504,26 @@ TEST_CASE("SVRGCallbacksFullFunctionTest", "[CallbacksTest]")
       false, false, true);
 }
 
+/**
+ * Make sure we invoke all callbacks (SWATS).
+ */
+TEST_CASE("SWATSCallbacksFullFunctionTest", "[CallbacksTest]")
+{
+  SWATS optimizer(0.01, 10, 0.9, 0.999, 1e-6, 1000, 1e-9, true);
+  CallbacksFullFunctionTest(optimizer, true, true, true, true, true, true,
+      false, false, true);
+}
+
+/**
+ * Make sure we invoke all callbacks (WNGrad).
+ */
+TEST_CASE("WNGradCallbacksFullFunctionTest", "[CallbacksTest]")
+{
+  WNGrad optimizer(0.56, 1, 1000, 1e-9, true);
+  CallbacksFullFunctionTest(optimizer, true, true, true, true, true, true,
+      false, false, true);
+}
+
 /**
  * Make sure we invoke all callbacks (ParallelSGD).
  */
@@ -470,14 +652,106 @@ TEST_CASE("TimerStopCallbackTest", "[CallbacksTest]")
 
   // Instantiate the optimizer with a number of iterations that will take a
   // long time to finish.
-  StandardSGD s(0.0003, 1, 2000000000, -100, true);
-
+  Adam opt(0.5, 2, 0.7, 0.999, 1e-8, 2000000000, -100, false);
   arma::wall_clock timer;
+
   timer.tic();
-
   // The optimization process should return in one second.
-  s.Optimize(f, coordinates, TimerStop(0.5));
-
+  opt.Optimize(f, coordinates, TimerStop(0.5));
   // Add some time to account for the function to return.
   REQUIRE(timer.toc() < 2);
 }
+
+/**
+ * Make sure the ProgressBar callback will show the progress on the specified
+ * output stream if the MaxIterations parameter of the optimizer is 0.
+ */
+TEST_CASE("ProgressBarCallbackNoMaxIterationsTest", "[CallbacksTest]")
+{
+  SGDTestFunction f;
+  arma::mat coordinates = f.GetInitialPoint();
+
+  StandardSGD s(0.0003, 1, 0, DBL_MAX, true);
+
+  std::stringstream stream;
+  s.Optimize(f, coordinates, ProgressBar(10, stream));
+
+  REQUIRE(stream.str().length() > 0);
+}
+
+/**
+ * Make sure the ProgressBar callback will show the progress on the specified
+ * output stream with the correct epoch number if the MaxIterations parameter
+ * of the optimizer is 0.
+ */
+TEST_CASE("ProgressBarCallbackNoMaxIterationsEpochTest", "[CallbacksTest]")
+{
+  SGDTestFunction f;
+  arma::mat coordinates = f.GetInitialPoint();
+
+  StandardSGD s(0.0003, 1, 0, DBL_MAX, true);
+
+  std::stringstream stream;
+  s.Optimize(f, coordinates, ProgressBar(10, stream));
+  REQUIRE(stream.str().find("Epoch 1") != std::string::npos);
+  REQUIRE(stream.str().find("Epoch 1/") == std::string::npos);
+}
+
+/**
+ * Make sure the ProgressBar callback will show the progress on the specified
+ * output stream with the correct epoch number if the MaxIterations parameter
+ * of the optimizer is not equal to 0.
+ */
+TEST_CASE("ProgressBarCallbackEpochTest", "[CallbacksTest]")
+{
+  SGDTestFunction f;
+  arma::mat coordinates = f.GetInitialPoint();
+
+  StandardSGD s(0.0003, 1, 1, 1e-9, true);
+
+  std::stringstream stream;
+  s.Optimize(f, coordinates, ProgressBar(10, stream));
+  REQUIRE(stream.str().find("Epoch 1/1") != std::string::npos);
+}
+
+/**
+ * Make sure the Report callback will show the report on the specified
+ * output stream.
+ */
+TEST_CASE("ReportCallbackTest", "[CallbacksTest]")
+{
+  std::stringstream stream;
+
+  SGDTestFunction f0;
+  StandardSGD s(0.0003, 1, 10000, 1e-9, true);
+
+  arma::mat coordinates = f0.GetInitialPoint();
+  s.Optimize(f0, coordinates, Report(0.1, stream));
+  REQUIRE(stream.str().length() > 0);
+
+  stream.str("");
+  RosenbrockWoodFunction f1;
+  L_BFGS lbfgs;
+  lbfgs.MaxIterations() = 100;
+
+  coordinates = f1.GetInitialPoint();
+  lbfgs.Optimize(f1, coordinates, Report(0.1, stream));
+  REQUIRE(stream.str().length() > 0);
+
+  stream.str("");
+  SchafferFunctionN2 f2;
+  CNE cne;
+  cne.MaxGenerations() = 100;
+
+  coordinates = f2.GetInitialPoint();
+  cne.Optimize(f2, coordinates, Report(0.1, stream));
+  REQUIRE(stream.str().length() > 0);
+
+  stream.str("");
+  AugLagrangianTestFunction f3;
+  AugLagrangian aug;
+
+  coordinates = f3.GetInitialPoint();
+  aug.Optimize(f3, coordinates, Report(0.1, stream));
+  REQUIRE(stream.str().length() > 0);
+}

tests/lbfgs_test.cpp

@@ -47,11 +47,11 @@ TEST_CASE("RosenbrockFunctionFloatTest", "[LBFGSTest]")
   arma::fmat coords = f.GetInitialPoint<arma::fvec>();
   lbfgs.Optimize(f, coords);
 
-  double finalValue = f.Evaluate(coords);
+  float finalValue = f.Evaluate(coords);
 
-  REQUIRE(finalValue == Approx(0.0).margin(1e-5));
-  REQUIRE(coords(0) == Approx(1.0).epsilon(1e-7));
-  REQUIRE(coords(1) == Approx(1.0).epsilon(1e-7));
+  REQUIRE(finalValue == Approx(0.0f).margin(1e-3));
+  REQUIRE(coords(0) == Approx(1.0f).epsilon(1e-4));
+  REQUIRE(coords(1) == Approx(1.0f).epsilon(1e-4));
 }
 
 /**

tests/nsga2_test.cpp

@@ -0,0 +1,189 @@
+/**
+ * @file nsga2_test.cpp
+ * @author Sayan Goswami
+ *
+ * ensmallen is free software; you may redistribute it and/or modify it under
+ * the terms of the 3-clause BSD license.  You should have received a copy of
+ * the 3-clause BSD license along with ensmallen.  If not, see
+ * http://www.opensource.org/licenses/BSD-3-Clause for more information.
+ */
+
+#include <ensmallen.hpp>
+#include "catch.hpp"
+#include "test_function_tools.hpp"
+
+using namespace ens;
+using namespace ens::test;
+using namespace std;
+
+/**
+ * Checks if low <= value <= high. Used by NSGA2FonsecaFlemingTest.
+ *
+ * @param value The value being checked.
+ * @param low The lower bound.
+ * @param high The upper bound.
+ * @return true if value lies in the range [low, high].
+ * @return false if value does not lie in the range [low, high].
+ */
+bool IsInBounds(const double& value, const double& low, const double& high)
+{
+  return !(value < low) && !(high < value);
+}
+
+/**
+ * Optimize for the Schaffer N.1 function using NSGA-II optimizer.
+ */
+TEST_CASE("NSGA2SchafferN1Test", "[NSGA2Test]")
+{
+  SchafferFunctionN1<arma::mat> SCH;
+  const double lowerBound = -1000;
+  const double upperBound = 1000;
+
+  NSGA2 opt(20, 5000, 0.5, 0.5, 1e-3, 1e-6, lowerBound, upperBound);
+
+  typedef decltype(SCH.objectiveA) ObjectiveTypeA;
+  typedef decltype(SCH.objectiveB) ObjectiveTypeB;
+
+  arma::mat coords = SCH.GetInitialPoint();
+  std::tuple<ObjectiveTypeA, ObjectiveTypeB> objectives = SCH.GetObjectives();
+
+  opt.Optimize(objectives, coords);
+  std::vector<arma::mat> bestFront = opt.Front();
+
+  bool allInRange = true;
+
+  for (arma::mat solution: bestFront)
+  {
+    double val = arma::as_scalar(solution);
+
+    if (val < 0.0 || val > 2.0)
+    {
+      allInRange = false;
+      break;
+    }
+  }
+  REQUIRE(allInRange);
+}
+
+/**
+ * Optimize for the Schaffer N.1 function using NSGA-II optimizer.
+ */
+TEST_CASE("NSGA2SchafferN1TestVectorBounds", "[NSGA2Test]")
+{
+  SchafferFunctionN1<arma::mat> SCH;
+  const arma::vec lowerBound = {-1000};
+  const arma::vec upperBound = {1000};
+
+  NSGA2 opt(20, 5000, 0.5, 0.5, 1e-3, 1e-6, lowerBound, upperBound);
+
+  typedef decltype(SCH.objectiveA) ObjectiveTypeA;
+  typedef decltype(SCH.objectiveB) ObjectiveTypeB;
+
+  arma::mat coords = SCH.GetInitialPoint();
+  std::tuple<ObjectiveTypeA, ObjectiveTypeB> objectives = SCH.GetObjectives();
+
+  opt.Optimize(objectives, coords);
+  std::vector<arma::mat> bestFront = opt.Front();
+
+  bool allInRange = true;
+
+  for (arma::mat solution: bestFront)
+  {
+    double val = arma::as_scalar(solution);
+
+    if (val < 0.0 || val > 2.0)
+    {
+      allInRange = false;
+      break;
+    }
+  }
+  REQUIRE(allInRange);
+}
+
+/**
+ * Optimize for the Fonseca Fleming function using NSGA-II optimizer.
+ */
+TEST_CASE("NSGA2FonsecaFlemingTest", "[NSGA2Test]")
+{
+  FonsecaFlemingFunction<arma::mat> FON;
+  const double lowerBound = -4;
+  const double upperBound = 4;
+  const double tolerance = 1e-6;
+  const double strength = 1e-4;
+  const double expectedLowerBound = -1.0 / sqrt(3);
+  const double expectedUpperBound = 1.0 / sqrt(3);
+
+  NSGA2 opt(20, 4000, 0.6, 0.3, strength, tolerance, lowerBound, upperBound);
+
+  typedef decltype(FON.objectiveA) ObjectiveTypeA;
+  typedef decltype(FON.objectiveB) ObjectiveTypeB;
+
+  arma::mat coords = FON.GetInitialPoint();
+  std::tuple<ObjectiveTypeA, ObjectiveTypeB> objectives = FON.GetObjectives();
+
+  opt.Optimize(objectives, coords);
+  std::vector<arma::mat> bestFront = opt.Front();
+
+  bool allInRange = true;
+
+  for (size_t i = 0; i < bestFront.size(); i++)
+  {
+    const arma::mat solution = bestFront[i];
+    double valX = arma::as_scalar(solution(0));
+    double valY = arma::as_scalar(solution(1));
+    double valZ = arma::as_scalar(solution(2));
+
+    if (!IsInBounds(valX, expectedLowerBound, expectedUpperBound) ||
+        !IsInBounds(valY, expectedLowerBound, expectedUpperBound) ||
+        !IsInBounds(valZ, expectedLowerBound, expectedUpperBound))
+    {
+      allInRange = false;
+      break;
+    }
+  }
+  REQUIRE(allInRange);
+}
+
+/**
+ * Optimize for the Fonseca Fleming function using NSGA-II optimizer.
+ */
+TEST_CASE("NSGA2FonsecaFlemingTestVectorBounds", "[NSGA2Test]")
+{
+  FonsecaFlemingFunction<arma::mat> FON;
+  const arma::vec lowerBound = {-4, -4, -4};
+  const arma::vec upperBound = {4, 4, 4};
+  const double tolerance = 1e-6;
+  const double strength = 1e-4;
+  const double expectedLowerBound = -1.0 / sqrt(3);
+  const double expectedUpperBound = 1.0 / sqrt(3);
+
+  NSGA2 opt(20, 4000, 0.6, 0.3, strength, tolerance, lowerBound, upperBound);
+
+  typedef decltype(FON.objectiveA) ObjectiveTypeA;
+  typedef decltype(FON.objectiveB) ObjectiveTypeB;
+
+  arma::mat coords = FON.GetInitialPoint();
+  std::tuple<ObjectiveTypeA, ObjectiveTypeB> objectives = FON.GetObjectives();
+
+  opt.Optimize(objectives, coords);
+  std::vector<arma::mat> bestFront = opt.Front();
+
+  bool allInRange = true;
+
+  for (size_t i = 0; i < bestFront.size(); i++)
+  {
+    const arma::mat solution = bestFront[i];
+    double valX = arma::as_scalar(solution(0));
+    double valY = arma::as_scalar(solution(1));
+    double valZ = arma::as_scalar(solution(2));
+
+    if (!IsInBounds(valX, expectedLowerBound, expectedUpperBound) ||
+        !IsInBounds(valY, expectedLowerBound, expectedUpperBound) ||
+        !IsInBounds(valZ, expectedLowerBound, expectedUpperBound))
+    {
+      allInRange = false;
+      break;
+    }
+  }
+  REQUIRE(allInRange);
+}

tests/sdp_primal_dual_test.cpp

@@ -243,7 +243,7 @@ static void SolveMaxCutFeasibleSDP(const SDP<arma::sp_mat>& sdp)
   ysparse = -1.1 * arma::vec(arma::sum(arma::abs(sdp.C()), 0).t());
   Z = -arma::diagmat(ysparse) + sdp.C();
 
-  PrimalDualSolver<> solver;
+  PrimalDualSolver solver;
 
   solver.Optimize(sdp, X, ysparse, ydense, Z);
   CheckKKT(sdp, X, ysparse, ydense, Z);
@@ -260,7 +260,7 @@ static void SolveMaxCutPositiveSDP(const SDP<arma::sp_mat>& sdp)
   ysparse = arma::randu<arma::vec>(sdp.NumSparseConstraints());
   Z.eye(sdp.N(), sdp.N());
 
-  PrimalDualSolver<> solver;
+  PrimalDualSolver solver;
   solver.Optimize(sdp, X, ysparse, ydense, Z);
   CheckKKT(sdp, X, ysparse, ydense, Z);
 }
@@ -289,7 +289,7 @@ TEST_CASE("DeprecatedSmallLovaszThetaSdp", "[SdpPrimalDualTest]")
   UndirectedGraph::LoadFromEdges(g, "data/johnson8-4-4.csv", true);
   auto sdp = ConstructLovaszThetaSDPFromGraph(g);
 
-  PrimalDualSolver<> solver;
+  PrimalDualSolver solver;
 
   arma::mat X, Z;
   arma::mat ysparse, ydense;
@@ -304,7 +304,7 @@ TEST_CASE("SmallLovaszThetaSdp", "[SdpPrimalDualTest]")
   UndirectedGraph::LoadFromEdges(g, "data/johnson8-4-4.csv", true);
   auto sdp = ConstructLovaszThetaSDPFromGraph(g);
 
-  PrimalDualSolver<> solver;
+  PrimalDualSolver solver;
 
   arma::mat X, Z, ysparse, ydense;
   sdp.GetInitialPoints(X, ysparse, ydense, Z);
@@ -439,7 +439,7 @@ TEST_CASE("LogChebychevApproxSdp","[SdpPrimalDualTest]")
     const arma::mat A0 = RandomFullRowRankMatrix(p0, k0);
     const arma::vec b0 = arma::randu<arma::vec>(p0);
     const auto sdp0 = ConstructLogChebychevApproxSdp(A0, b0);
-    PrimalDualSolver<> solver0;
+    PrimalDualSolver solver0;
     arma::mat X0, Z0;
     arma::mat ysparse0, ydense0;
     sdp0.GetInitialPoints(X0, ysparse0, ydense0, Z0);
@@ -459,7 +459,7 @@ TEST_CASE("LogChebychevApproxSdp","[SdpPrimalDualTest]")
     const arma::mat A1 = RandomFullRowRankMatrix(p1, k1);
     const arma::vec b1 = arma::randu<arma::vec>(p1);
     const auto sdp1 = ConstructLogChebychevApproxSdp(A1, b1);
-    PrimalDualSolver<> solver1;
+    PrimalDualSolver solver1;
     arma::mat X1, Z1;
     arma::mat ysparse1, ydense1;
     sdp1.GetInitialPoints(X1, ysparse1, ydense1, Z1);
@@ -573,7 +573,7 @@ TEST_CASE("CorrelationCoeffToySdp","[SdpPrimalDualTest]")
   sdp.C().zeros();
   sdp.C()(0, 2) = sdp.C()(2, 0) = 1.;
 
-  PrimalDualSolver<> solver;
+  PrimalDualSolver solver;
   arma::mat X, Z;
   arma::mat ysparse, ydense;
   sdp.GetInitialPoints(X, ysparse, ydense, Z);

tests/spsa_test.cpp

@@ -45,8 +45,8 @@ TEST_CASE("SPSASphereFunctionFMatTest", "[SPSATest]")
   arma::fmat coordinates = f.GetInitialPoint<arma::fmat>();
   optimizer.Optimize(f, coordinates);
 
-  REQUIRE(coordinates(0) == Approx(0.0).margin(0.1));
-  REQUIRE(coordinates(1) == Approx(0.0).margin(0.1));
+  REQUIRE(coordinates(0) == Approx(0.0f).margin(0.1));
+  REQUIRE(coordinates(1) == Approx(0.0f).margin(0.1));
 }
 
 /**

tests/swats_test.cpp

@@ -83,8 +83,8 @@ TEST_CASE("SWATSStyblinskiTangFunctionFMatTest", "[SWATSTest]")
   arma::fmat coordinates = f.GetInitialPoint<arma::fmat>();
   optimizer.Optimize(f, coordinates);
 
-  REQUIRE(coordinates(0) == Approx(-2.9).epsilon(0.01));
-  REQUIRE(coordinates(1) == Approx(-2.9).epsilon(0.01));
+  REQUIRE(coordinates(0) == Approx(-2.9).epsilon(0.1));
+  REQUIRE(coordinates(1) == Approx(-2.9).epsilon(0.1));
 }
 
 #if ARMA_VERSION_MAJOR > 9 ||\