1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
|
// Licensed to the .NET Foundation under one or more agreements.
// The .NET Foundation licenses this file to you under the MIT license.
// See the LICENSE file in the project root for more information.
//
// The Adams-Moulton Predictor Corrector Method adapted from Conte and de Boor
// original source: adams_d.c
using System;
using System.Runtime.CompilerServices;
using System.Diagnostics;
#if XUNIT_PERF
using Xunit;
using Microsoft.Xunit.Performance;
#endif // XUNIT_PERF
#if XUNIT_PERF
[assembly: OptimizeForBenchmarks]
[assembly: MeasureInstructionsRetired]
#endif // XUNIT_PERF
namespace Benchstone.BenchF
{
public static class Adams
{
#if DEBUG
public static int Iterations = 1;
#else
public static int Iterations = 200000;
#endif // DEBUG
static double g_xn, g_yn, g_dn, g_en;
const double g_xn_base = 0.09999999E+01;
const double g_yn_base = 0.71828180E+00;
const double g_dn_base = 0.21287372E-08;
const double g_en_base = 0.74505806E-08;
[MethodImpl(MethodImplOptions.NoInlining)]
private static void Bench()
{
double[] f = new double[5];
double xn, yn, dn, en, yxn, h, fnp, ynp, y0, x0, nz;
int i, k, n, nstep;
#if VERBOSE
Console.WriteLine(" ADAMS-MOULTON METHOD ");
#endif // VERBOSE
n = 4;
h = 1.0 / 32.0;
nstep = 32;
y0 = 0.0;
x0 = 0.0;
xn = 0.0;
yn = 0.0;
dn = 0.0;
en = 0.0;
nz = 0;
f[1] = x0 + y0;
#if VERBOSE
Console.WriteLine("{0}, {1}, {2}, {3}, {4}", nz, x0, y0, dn, en);
#endif // VERBOSE
xn = x0;
for (i = 2; i <= 4; i++)
{
k = i - 1;
xn = xn + h;
yn = Soln(xn);
f[i] = xn + yn;
#if VERBOSE
Console.WriteLine("{0}, {1}, {2}, {3}, {4}", k, xn, yn, dn, en);
#endif // VERBOSE
}
for (k = 4; k <= nstep; k++)
{
ynp = yn + (h / 24) * (55 * f[n] - 59 * f[n - 1] + 37 * f[n - 2] - 9 * f[n - 3]);
xn = xn + h;
fnp = xn + ynp;
yn = yn + (h / 24) * (9 * fnp + 19 * f[n] - 5 * f[n - 1] + f[n - 2]);
dn = (yn - ynp) / 14;
f[n - 3] = f[n - 2];
f[n - 2] = f[n - 1];
f[n - 1] = f[n];
f[n] = xn + yn;
yxn = Soln(xn);
en = yn - yxn;
#if VERBOSE
Console.WriteLine("{0}, {1}, {2}, {3}, {4}", k, xn, yn, dn, en);
#endif // VERBOSE
}
// Escape calculated values:
g_xn = xn;
g_yn = yn;
g_dn = dn;
g_en = en;
}
private static double Soln(double x)
{
return (System.Math.Exp(x) - 1.0 - (x));
}
[MethodImpl(MethodImplOptions.NoOptimization | MethodImplOptions.NoInlining)]
private static void TestBench()
{
for (int l = 1; l <= Iterations; l++)
{
Bench();
}
}
#if XUNIT_PERF
[Benchmark]
public static void Test()
{
foreach (var iteration in Benchmark.Iterations)
{
using (iteration.StartMeasurement())
{
TestBench();
}
}
}
#endif // XUNIT_PERF
[MethodImpl(MethodImplOptions.NoOptimization)]
public static int Main(string[] argv)
{
if (argv.Length > 0)
{
Iterations = Int32.Parse(argv[0]);
}
Stopwatch sw = Stopwatch.StartNew();
TestBench();
sw.Stop();
bool result = true;
// Note: we can't check xn or yn better because of the precision
// with which original results are given
result &= System.Math.Abs(g_xn_base - g_xn) <= 1.5e-7;
result &= System.Math.Abs(g_yn_base - g_yn) <= 1.5e-7;
result &= System.Math.Abs(g_dn) <= 2.5e-9;
// Actual error is much bigger than base error;
// this is likely due to the fact that the original program was written in Fortran
// and was running on a mainframe with a non-IEEE floating point arithmetic
// (it's still beyond the published precision of yn)
result &= System.Math.Abs(g_en) <= 5.5e-8;
Console.WriteLine(result ? "Passed" : "Failed");
Console.WriteLine(" BASE.....P1 1/4 (ADAMS-MOULTON), XN = {0}", g_xn_base);
Console.WriteLine(" VERIFY...P1 1/4 (ADAMS-MOULTON), XN = {0}\n", g_xn);
Console.WriteLine(" BASE.....P1 2/4 (ADAMS-MOULTON), YN = {0}", g_yn_base);
Console.WriteLine(" VERIFY...P1 2/4 (ADAMS-MOULTON), YN = {0}\n", g_yn);
Console.WriteLine(" BASE.....P1 3/4 (ADAMS-MOULTON), DN = {0}", g_dn_base);
Console.WriteLine(" VERIFY...P1 3/4 (ADAMS-MOULTON), DN = {0}\n", g_dn);
Console.WriteLine(" BASE.....P1 4/4 (ADAMS-MOULTON), EN = {0}", g_en_base);
Console.WriteLine(" VERIFY...P1 4/4 (ADAMS-MOULTON), EN = {0}\n", g_en);
Console.WriteLine("Test iterations: {0}; Total time: {1} sec", Iterations, sw.Elapsed.TotalSeconds);
return (result ? 100 : -1);
}
}
}
|
