"Fossies" - the Fresh Open Source Software Archive

Member "fityk-1.3.1/tests/gradient.cpp" (13 May 2016, 4548 Bytes) of package /linux/misc/fityk-1.3.1.tar.gz:


As a special service "Fossies" has tried to format the requested source page into HTML format using (guessed) C and C++ source code syntax highlighting (style: standard) with prefixed line numbers and code folding option. Alternatively you can here view or download the uninterpreted source code file. For more information about "gradient.cpp" see the Fossies "Dox" file reference documentation.

    1 
    2 #include <stdio.h>
    3 #include <boost/scoped_ptr.hpp>
    4 #include "fityk/logic.h"
    5 #include "fityk/data.h"
    6 #include "fityk/fit.h"
    7 
    8 #include "catch.hpp"
    9 
   10 using namespace std;
   11 using namespace fityk;
   12 
   13 
   14 // Box and Betts exponential quadratic sum, equivalent to least-squares of
   15 // f(x) = exp(-0.1*a0*x) - exp(-0.1*a1*x) - a2 * (exp(-0.1*x) - exp(-x))
   16 // with points (1,0), (2,0), ... (10,0)
   17 // 
   18 // domains of a0, a1, a2 are, respectively, (0.9, 1.2), (9, 11.2), (0.9, 1.2)
   19 // minimum: (1,10,1) -> 0
   20 
   21 // modified boxbetts_f() from nlopt-2.3/test/testfuncs.c
   22 static double boxbetts_f(const double *a, double *grad)
   23 {
   24     double f = 0;
   25     if (grad)
   26         grad[0] = grad[1] = grad[2] = 0;
   27     for (int x = 1; x <= 10; ++x) {
   28         double e0 = exp(-0.1*x*a[0]);
   29         double e1 = exp(-0.1*x*a[1]);
   30         double e2 = exp(-0.1*x) - exp(-x);
   31         double g = e0 - e1 - e2 * a[2];
   32         f += g*g;
   33         if (grad) {
   34             grad[0] += (2 * g) * (-0.1*x*e0);
   35             grad[1] += (2 * g) * (0.1*x*e1);
   36             grad[2] += -(2 * g) * e2;
   37         }
   38     }
   39     return f;
   40 }
   41 
   42 static double boxbetts_in_fityk(const double *a, double *grad)
   43 {
   44     boost::scoped_ptr<Fityk> ftk(new Fityk);
   45     Full* priv = ftk->priv();
   46     ftk->set_option_as_number("verbosity", -1);
   47     for (int i = 1; i <= 10; ++i)
   48         priv->dk.data(0)->add_one_point(i, 0, 1);
   49     ftk->execute("define BoxBetts(a0,a1,a2) = "
   50             "exp(-0.1*a0*x) - exp(-0.1*a1*x) - a2 * (exp(-0.1*x) - exp(-x))");
   51     ftk->execute("F = BoxBetts(~0.9, ~11.8, ~1.08)");
   52 
   53     priv->get_fit()->get_dof(priv->dk.datas()); // to invoke update_par_usage()
   54     vector<realt> avec(a, a+3);
   55     return priv->get_fit()->compute_wssr_gradient(avec, priv->dk.datas(), grad);
   56 }
   57 
   58 #if 0
   59 static int test_gradient()
   60 {
   61     const double a[3] = { 0.9, 11.8, 1.08 };
   62     double grad[3], grad_again[3];
   63     double ssr = boxbetts_f(a, grad);
   64     printf("ssr=%.10g  grad: %.10g %.10g %.10g\n",
   65             ssr, grad[0], grad[1], grad[2]);
   66 
   67     double ssr_again = boxbetts_in_fityk(a, grad_again);
   68     printf("ssr=%.10g  grad: %.10g %.10g %.10g\n",
   69             ssr_again, grad_again[0], grad_again[1], grad_again[2]);
   70     return 0;
   71 }
   72 int main() { return test_gradient(); }
   73 #endif
   74 
   75 
   76 TEST_CASE("gradient", "test Fit::compute_wssr_gradient()") {
   77     const double a[3] = { 0.9, 11.8, 1.08 };
   78     double grad[3], grad_again[3];
   79     double ssr = boxbetts_f(a, grad);
   80     double ssr_again = boxbetts_in_fityk(a, grad_again);
   81     REQUIRE(ssr == Approx(ssr_again));
   82     REQUIRE(grad[0] == Approx(grad_again[0]));
   83     REQUIRE(grad[1] == Approx(grad_again[1]));
   84     REQUIRE(grad[2] == Approx(grad_again[2]));
   85 }
   86 
   87 //----------- + some unrelated random tests
   88 
   89 TEST_CASE("set-throws", "test Fityk::set_throws()") {
   90     boost::scoped_ptr<Fityk> fik(new Fityk);
   91     REQUIRE(fik->get_throws() == true);
   92     REQUIRE_THROWS_AS(fik->execute("hi"), SyntaxError);
   93     REQUIRE_THROWS_AS(fik->execute("fit"), ExecuteError);
   94     fik->set_throws(false);
   95     REQUIRE_NOTHROW(fik->execute("hi"));
   96     fik->execute("hi");
   97     REQUIRE(fik->last_error() != "");
   98     REQUIRE(fik->last_error().substr(0,12) == "SyntaxError:");
   99     fik->clear_last_error();
  100     REQUIRE(fik->last_error() == "");
  101 }
  102 
  103 
  104 TEST_CASE("set-get-option", "test Fityk::[gs]et_option_as_*()") {
  105     boost::scoped_ptr<Fityk> fik(new Fityk);
  106 
  107     REQUIRE(fik->get_option_as_string("numeric_format") == "%g"); // string
  108     REQUIRE(fik->get_option_as_string("on_error") == "stop"); // enum
  109     REQUIRE(fik->get_option_as_number("verbosity") == 0.); // integer
  110     REQUIRE(fik->get_option_as_number("fit_replot") == 0.); // bool
  111     REQUIRE(fik->get_option_as_number("epsilon") == 1e-12); // double
  112 
  113     // no such option
  114     REQUIRE_THROWS_AS(fik->get_option_as_string("hi"), ExecuteError);
  115     // non-numeric option
  116     REQUIRE_THROWS_AS(fik->get_option_as_number("on_error"), ExecuteError);
  117     REQUIRE_THROWS_AS(fik->get_option_as_number("cwd"), ExecuteError);
  118 
  119     fik->set_option_as_string("numeric_format", "%.8E");
  120     fik->set_option_as_string("on_error", "nothing");
  121     fik->set_option_as_number("verbosity", -1);
  122     fik->set_option_as_number("fit_replot", 1);
  123     fik->set_option_as_number("epsilon", 1e-10);
  124 
  125     REQUIRE(fik->get_option_as_string("numeric_format") == "%.8E"); // string
  126     REQUIRE(fik->get_option_as_string("on_error") == "nothing"); // enum
  127     REQUIRE(fik->get_option_as_number("verbosity") == -1); // integer
  128     REQUIRE(fik->get_option_as_number("fit_replot") == 1.); // bool
  129     REQUIRE(fik->get_option_as_number("epsilon") == 1e-10); // double
  130 }
  131