Skip to content

Commit

Permalink
expand F functions around 0 in one/both arguments
Browse files Browse the repository at this point in the history
  • Loading branch information
Alexander Voigt authored and Alexander Voigt committed Jan 31, 2017
1 parent 115807d commit e573649
Show file tree
Hide file tree
Showing 2 changed files with 136 additions and 4 deletions.
55 changes: 51 additions & 4 deletions src/threshold_loop_functions.cpp
Expand Up @@ -66,6 +66,9 @@ double F1(double x)
{
const double x2 = sqr(x);

if (is_equal(x, 0.))
return 0.;

if (is_equal(x, 1., 0.01))
return (16 + 41*x - 44*x2 + 21*cube(x) - 4*quad(x))/30.;

Expand All @@ -76,6 +79,9 @@ double F2(double x)
{
const double x2 = sqr(x);

if (is_equal(x, 0.))
return 0.;

if (is_equal(x, 1., 0.01))
return (-5 + 216*x - 226*x2 + 104*cube(x) - 19*quad(x))/70.;

Expand All @@ -86,6 +92,9 @@ double F3(double x)
{
const double x2 = sqr(x);

if (is_equal(x, 0.))
return 0.;

if (is_equal(x, 1., 0.01))
return (-27 + 218*x - 142*x2 + 48*cube(x) - 7*quad(x))/90.;

Expand All @@ -96,6 +105,9 @@ double F4(double x)
{
const double x2 = sqr(x);

if (is_equal(x, 0.))
return 0.;

if (is_equal(x, 1., 0.01))
return (31 + 22*x - 42*x2 + 24*cube(x) - 5*quad(x))/30.;

Expand All @@ -107,6 +119,9 @@ double F5(double x)
const double x2 = sqr(x);
const double x4 = quad(x);

if (is_equal(x, 0.))
return 0.;

if (is_equal(x, 1., 0.01))
return (13 + 165*x - 174*x2 + 81*cube(x) - 15*x4)/70.;

Expand All @@ -117,6 +132,9 @@ double F6(double x)
{
const double x2 = sqr(x);

if (is_equal(x, 0.))
return -0.75;

if (is_equal(x, 1., 0.01))
return (-103 + 128*x - 26*x2 + quad(x))/120.;

Expand All @@ -128,6 +146,9 @@ double F7(double x)
const double x2 = sqr(x);
const double x4 = quad(x);

if (is_equal(x, 0.))
return -1.5;

if (is_equal(x, 1., 0.01))
return -1.8642857142857139 + 2.057142857142856*x
+ 1.7142857142857153*x2 - 1.1428571428571432*cube(x)
Expand Down Expand Up @@ -197,14 +218,32 @@ static double F8_x1_x2(double x1, double x2)

double F8(double x1, double x2)
{
if (is_equal(x1, 0.) && is_equal(x2, 0.))
return -2.;

if (is_equal(x1, 1., 0.01) && is_equal(x2, 1., 0.01))
return F8_1_1(x1, x2);

if (is_equal(x1, 1., 0.01))
if (is_equal(x1, 1., 0.01)) {
if (is_equal(x2, 0., 0.01)) {
return -2.333333333333332 + 5.66666666666667*sqr(x2) +
x1*(2.6666666666666643 - 5.333333333333339*sqr(x2)) +
sqr(x1)*(-0.33333333333333215 + 1.6666666666666696*sqr(x2));
}

return F8_1_x2(x1, x2);
}

if (is_equal(x2, 1., 0.01)) {
if (is_equal(x1, 0., 0.01)) {
return -2.3333333333333335 + 2.6666666666666665*x2 -
0.3333333333333333*sqr(x2) +
sqr(x1)*(5.666666666666667 - 5.333333333333334*x2 +
1.6666666666666667*sqr(x2));
}

if (is_equal(x2, 1., 0.01))
return F8_1_x2(x2, x1);
}

if (is_equal(x1, x2, 0.00001))
return F8_x1_x2(x1, x2);
Expand Down Expand Up @@ -292,11 +331,19 @@ double F9(double x1, double x2)
if (is_equal(x1, 1., 0.01) && is_equal(x2, 1., 0.01))
return F9_1_1(x1, x2);

if (is_equal(x1, 1., 0.01))
if (is_equal(x1, 1., 0.01)) {
if (is_equal(x2, 0., 0.01))
return 2. - 2.*(-1 + x1) + 5./3.*sqr(-1 + x1);

return F9_1_x2(x1, x2);
}

if (is_equal(x2, 1., 0.01)) {
if (is_equal(x1, 0., 0.01))
return 2. - 2.*(-1 + x2) + 5./3.*sqr(-1 + x2);

if (is_equal(x2, 1., 0.01))
return F9_1_x2(x2, x1);
}

if (is_equal(x1, x2, 0.00001))
return F9_x1_x2(x1, x2);
Expand Down
85 changes: 85 additions & 0 deletions test/test_threshold_loop_functions.cpp
Expand Up @@ -237,6 +237,14 @@ BOOST_AUTO_TEST_CASE(test_F1)
x = 1.001; BOOST_CHECK_CLOSE(F1(x), F1_bare(x), 1e-5);
x = 1.0001; BOOST_CHECK_CLOSE(F1(x), F1_bare(x), 1e-5);

x = 0.1; BOOST_CHECK_CLOSE(F1(x), F1_bare(x), 1e-5);
x = 0.02; BOOST_CHECK_CLOSE(F1(x), F1_bare(x), 1e-5);
x = 0.01; BOOST_CHECK_CLOSE(F1(x), F1_bare(x), 1e-5);
x = 0.001; BOOST_CHECK_CLOSE(F1(x), F1_bare(x), 1e-5);
x = 0.0001; BOOST_CHECK_CLOSE(F1(x), F1_bare(x), 1e-5);
x = 0.00001; BOOST_CHECK_CLOSE(F1(x), F1_bare(x), 1e-5);

BOOST_CHECK(!std::isnan(F1(0)));
BOOST_CHECK(!std::isnan(F1(1)));
}

Expand All @@ -252,6 +260,14 @@ BOOST_AUTO_TEST_CASE(test_F2)
x = 1.001; BOOST_CHECK_CLOSE(F2(x), F2_bare(x), 1e-5);
x = 1.0001; BOOST_CHECK_CLOSE(F2(x), F2_bare(x), 1e-5);

x = 0.1; BOOST_CHECK_CLOSE(F2(x), F2_bare(x), 1e-5);
x = 0.02; BOOST_CHECK_CLOSE(F2(x), F2_bare(x), 1e-5);
x = 0.01; BOOST_CHECK_CLOSE(F2(x), F2_bare(x), 1e-5);
x = 0.001; BOOST_CHECK_CLOSE(F2(x), F2_bare(x), 1e-5);
x = 0.0001; BOOST_CHECK_CLOSE(F2(x), F2_bare(x), 1e-5);
x = 0.00001; BOOST_CHECK_CLOSE(F2(x), F2_bare(x), 1e-5);

BOOST_CHECK(!std::isnan(F2(0)));
BOOST_CHECK(!std::isnan(F2(1)));
}

Expand All @@ -267,6 +283,14 @@ BOOST_AUTO_TEST_CASE(test_F3)
x = 1.001; BOOST_CHECK_CLOSE(F3(x), F3_bare(x), 1e-5);
x = 1.0001; BOOST_CHECK_CLOSE(F3(x), F3_bare(x), 1e-5);

x = 0.1; BOOST_CHECK_CLOSE(F3(x), F3_bare(x), 1e-5);
x = 0.02; BOOST_CHECK_CLOSE(F3(x), F3_bare(x), 1e-5);
x = 0.01; BOOST_CHECK_CLOSE(F3(x), F3_bare(x), 1e-5);
x = 0.001; BOOST_CHECK_CLOSE(F3(x), F3_bare(x), 1e-5);
x = 0.0001; BOOST_CHECK_CLOSE(F3(x), F3_bare(x), 1e-5);
x = 0.00001; BOOST_CHECK_CLOSE(F3(x), F3_bare(x), 1e-5);

BOOST_CHECK(!std::isnan(F3(0)));
BOOST_CHECK(!std::isnan(F3(1)));
}

Expand All @@ -282,6 +306,14 @@ BOOST_AUTO_TEST_CASE(test_F4)
x = 1.001; BOOST_CHECK_CLOSE(F4(x), F4_bare(x), 1e-5);
x = 1.0001; BOOST_CHECK_CLOSE(F4(x), F4_bare(x), 1e-5);

x = 0.1; BOOST_CHECK_CLOSE(F4(x), F4_bare(x), 1e-5);
x = 0.02; BOOST_CHECK_CLOSE(F4(x), F4_bare(x), 1e-5);
x = 0.01; BOOST_CHECK_CLOSE(F4(x), F4_bare(x), 1e-5);
x = 0.001; BOOST_CHECK_CLOSE(F4(x), F4_bare(x), 1e-5);
x = 0.0001; BOOST_CHECK_CLOSE(F4(x), F4_bare(x), 1e-5);
x = 0.00001; BOOST_CHECK_CLOSE(F4(x), F4_bare(x), 1e-5);

BOOST_CHECK(!std::isnan(F4(0)));
BOOST_CHECK(!std::isnan(F4(1)));
}

Expand All @@ -297,6 +329,14 @@ BOOST_AUTO_TEST_CASE(test_F5)
x = 1.001; BOOST_CHECK_CLOSE(F5(x), F5_bare(x), 1e-5);
x = 1.0001; BOOST_CHECK_CLOSE(F5(x), F5_bare(x), 5e-4);

x = 0.1; BOOST_CHECK_CLOSE(F5(x), F5_bare(x), 1e-5);
x = 0.02; BOOST_CHECK_CLOSE(F5(x), F5_bare(x), 1e-5);
x = 0.01; BOOST_CHECK_CLOSE(F5(x), F5_bare(x), 1e-5);
x = 0.001; BOOST_CHECK_CLOSE(F5(x), F5_bare(x), 1e-5);
x = 0.0001; BOOST_CHECK_CLOSE(F5(x), F5_bare(x), 1e-5);
x = 0.00001; BOOST_CHECK_CLOSE(F5(x), F5_bare(x), 1e-5);

BOOST_CHECK(!std::isnan(F5(0)));
BOOST_CHECK(!std::isnan(F5(1)));
}

Expand All @@ -312,6 +352,14 @@ BOOST_AUTO_TEST_CASE(test_F6)
x = 1.001; BOOST_CHECK_CLOSE(F6(x), F6_bare(x), 1e-5);
x = 1.0001; BOOST_CHECK_CLOSE(F6(x), F6_bare(x), 1e-5);

x = 0.1; BOOST_CHECK_CLOSE(F6(x), F6_bare(x), 1e-5);
x = 0.02; BOOST_CHECK_CLOSE(F6(x), F6_bare(x), 1e-5);
x = 0.01; BOOST_CHECK_CLOSE(F6(x), F6_bare(x), 1e-5);
x = 0.001; BOOST_CHECK_CLOSE(F6(x), F6_bare(x), 1e-5);
x = 0.0001; BOOST_CHECK_CLOSE(F6(x), F6_bare(x), 1e-5);
x = 0.00001; BOOST_CHECK_CLOSE(F6(x), F6_bare(x), 1e-5);

BOOST_CHECK(!std::isnan(F6(0)));
BOOST_CHECK(!std::isnan(F6(1)));
}

Expand All @@ -327,6 +375,14 @@ BOOST_AUTO_TEST_CASE(test_F7)
x = 1.001; BOOST_CHECK_CLOSE(F7(x), F7_bare(x), 1e-5);
x = 1.0001; BOOST_CHECK_CLOSE(F7(x), F7_bare(x), 1e-5);

x = 0.1; BOOST_CHECK_CLOSE(F7(x), F7_bare(x), 1e-5);
x = 0.02; BOOST_CHECK_CLOSE(F7(x), F7_bare(x), 1e-5);
x = 0.01; BOOST_CHECK_CLOSE(F7(x), F7_bare(x), 1e-5);
x = 0.001; BOOST_CHECK_CLOSE(F7(x), F7_bare(x), 1e-5);
x = 0.0001; BOOST_CHECK_CLOSE(F7(x), F7_bare(x), 1e-5);
x = 0.00001; BOOST_CHECK_CLOSE(F7(x), F7_bare(x), 1e-5);

BOOST_CHECK(!std::isnan(F7(0)));
BOOST_CHECK(!std::isnan(F7(1)));
}

Expand Down Expand Up @@ -360,6 +416,21 @@ BOOST_AUTO_TEST_CASE(test_F8)
x = 2.001; BOOST_CHECK_CLOSE(F8(x,x), F8_bare(x,x + 0.0001), 1e-2);
x = 2.0001; BOOST_CHECK_CLOSE(F8(x,x), F8_bare(x,x + 0.0001), 1e-2);

x = 0.1; BOOST_CHECK_CLOSE(F8(x,1), F8_bare(x,1.000001), 2e-2);
x = 0.02; BOOST_CHECK_SMALL(F8(x,1) - F8_bare(x,1.00001), 1e-2);
x = 0.01; BOOST_CHECK_SMALL(F8(x,1) - F8_bare(x,1.00001), 1e-2);
x = 0.001; BOOST_CHECK_SMALL(F8(x,1) - F8_bare(x,1.00001), 1e-2);
x = 0.0001; BOOST_CHECK_SMALL(F8(x,1) - F8_bare(x,1.00001), 1e-2);

x = 0.1; BOOST_CHECK_CLOSE(F8(1,x), F8_bare(1.000001,x), 2e-2);
x = 0.02; BOOST_CHECK_SMALL(F8(1,x) - F8_bare(1.00001,x), 1e-2);
x = 0.01; BOOST_CHECK_SMALL(F8(1,x) - F8_bare(1.00001,x), 1e-2);
x = 0.001; BOOST_CHECK_SMALL(F8(1,x) - F8_bare(1.00001,x), 1e-2);
x = 0.0001; BOOST_CHECK_SMALL(F8(1,x) - F8_bare(1.00001,x), 1e-2);

BOOST_CHECK(!std::isnan(F8(0,0)));
BOOST_CHECK(!std::isnan(F8(0,1)));
BOOST_CHECK(!std::isnan(F8(1,0)));
BOOST_CHECK(!std::isnan(F8(1,1)));
BOOST_CHECK(!std::isnan(F8(2,2)));
}
Expand Down Expand Up @@ -394,6 +465,20 @@ BOOST_AUTO_TEST_CASE(test_F9)
x = 2.001; BOOST_CHECK_CLOSE(F9(x,x), F9_bare(x,x + 0.0001), 1e-2);
x = 2.0001; BOOST_CHECK_CLOSE(F9(x,x), F9_bare(x,x + 0.0001), 1e-2);

x = 0.1; BOOST_CHECK_CLOSE(F9(x,1), F9_bare(x,1.00001), 1e-2);
x = 0.02; BOOST_CHECK_CLOSE(F9(x,1), F9_bare(x,1.00001), 1e-2);
x = 0.01; BOOST_CHECK_CLOSE(F9(x,1), F9_bare(x,1.00001), 1e-2);
x = 0.001; BOOST_CHECK_CLOSE(F9(x,1), F9_bare(x,1.00001), 1e-2);
x = 0.0001; BOOST_CHECK_CLOSE(F9(x,1), F9_bare(x,1.00001), 1e-2);

x = 0.1; BOOST_CHECK_CLOSE(F9(1,x), F9_bare(1.00001,x), 1e-2);
x = 0.02; BOOST_CHECK_CLOSE(F9(1,x), F9_bare(1.00001,x), 1e-2);
x = 0.01; BOOST_CHECK_CLOSE(F9(1,x), F9_bare(1.00001,x), 1e-2);
x = 0.001; BOOST_CHECK_CLOSE(F9(1,x), F9_bare(1.00001,x), 1e-2);
x = 0.0001; BOOST_CHECK_CLOSE(F9(1,x), F9_bare(1.00001,x), 1e-2);

BOOST_CHECK(!std::isnan(F9(0,1)));
BOOST_CHECK(!std::isnan(F9(1,0)));
BOOST_CHECK(!std::isnan(F9(1,1)));
BOOST_CHECK(!std::isnan(F9(2,2)));
}
Expand Down

0 comments on commit e573649

Please sign in to comment.