fail by default when negative eigenvalues are present

LICENSE.md

 Copyright holders:
 ------------------
 
-| Year         | Name          |
-|--------------|---------------|
-| 2010 -- 2016 | Ole Klein     |
-| 2010         | Adrian Ngo    |
-| 2010         |Ronnie Schwede |
+| Year         | Name           |
+|--------------|----------------|
+| 2010 -- 2019 | Ole Klein      |
+| 2010         | Adrian Ngo     |
+| 2017         | Lukas Riedel   |
+| 2010         | Ronnie Schwede |
 
 The module dune-randomfield is free software and documentation.
 

dune/randomfield/fieldtraits.hh

@@ -135,6 +135,7 @@ namespace Dune {
         const RF           variance;
         const std::string  covariance;
         const bool         periodic;
+        const bool         approximate;
         const bool         verbose;
         const unsigned int cgIterations;
 
@@ -173,6 +174,7 @@ namespace Dune {
           variance       (config.get<RF>                           ("stochastic.variance")),
           covariance     (config.get<std::string>                  ("stochastic.covariance")),
           periodic       (config.get<bool>                         ("randomField.periodic",false)),
+          approximate    (config.get<bool>                         ("randomField.approximate",false)),
           verbose        (config.get<bool>                         ("randomField.verbose",false)),
           cgIterations   (config.get<unsigned int>                 ("randomField.cgIterations",100)),
           embeddingFactor(config.get<unsigned int>                 ("randomField.embeddingFactor",2)),

dune/randomfield/matrix.hh

@@ -24,12 +24,12 @@ namespace Dune {
           const std::shared_ptr<Traits> traits;
 
           int rank, commSize;
-          std::array<RF,dim>        extensions;
-          unsigned int              level;
-          std::array<RF,dim>        meshsize;
-          RF                        variance;
-          std::string               covariance;
-          unsigned int              cgIterations;
+          std::array<RF,dim> extensions;
+          unsigned int       level;
+          std::array<RF,dim> meshsize;
+          RF                 variance;
+          std::string        covariance;
+          unsigned int       cgIterations;
 
           ptrdiff_t allocLocal, localN0, local0Start;
 
@@ -264,6 +264,15 @@ namespace Dune {
               std::cout << small << " small, " << smallNegative << " small negative and "
                 << negative << " large negative eigenvalues in covariance matrix, smallest "
                 << smallest << std::endl;
+
+            if (negative > 0 && !(*traits).approximate)
+            {
+              if (rank == 0)
+                std::cerr << "negative eigenvalues in covariance matrix, "
+                  << "consider increasing embeddingFactor, or alternatively "
+                  << "allow generation of approximate samples" << std::endl;
+              DUNE_THROW(Dune::Exception,"negative eigenvalues in covariance matrix");
+            }
           }
 
           template<typename Covariance>

src/fieldgenerator.cc

@@ -377,6 +377,11 @@ void printFullExample()
     << "# sets embeddingFactor = 1, i.e., behavior can't be controlled per\n"
     << "# boundary segment and correlation length must be small enough\n"
     << "periodic = 0\n"
+    << "# accept approximate results (1) or not (0)\n"
+    << "# simply sets negative eigenvalues to zero if they occur\n"
+    << "approximate = 0\n"
+    << "# whether information should be printed (1) or not (0)\n"
+    << "verbose = 0\n"
     << "# Conjugate Gradients iterations for matrix inverse multiplication\n"
     << "cgIterations = 100\n"
     << "\n"
@@ -386,9 +391,8 @@ void printFullExample()
     << "# possible values:\n"
     << "#     exponential, gaussian, spherical,\n"
     << "#     separapleExponential, matern32, matern52,\n"
-    << "#     dampedOscillation, cauchy, invCauchy,\n"
-    << "#     yaglom12, yaglom32, cubic,\n"
-    << "#     holeEffect, whiteNoise\n"
+    << "#     dampedOscillation, cauchy, cubic,\n"
+    << "#     whiteNoise\n"
     << "covariance = exponential\n"
     << "# variance of random field\n"
     << "variance = 1.\n"