関数の実装にバグがあります。
atan(b*x - 1 / 2)
項1 / 2
は整数除算を行い、0 に評価されますが、これはおそらくあなたが望むものではありません。一般に、double 変数で演算を行う場合は double リテラルを使用します。このpow()
関数は(double, int)
オーバーロードの 1 つとして使用するので、問題はありません。オーバーロードもあり(double, double)
ますが、指数が実際に整数である場合、それは望ましくありません。
これは、最も単純なルート検索方法の簡単な実装です-二分法です(後で、OPが二分タグを使用していることに気付きました、完璧です)。
#include <iostream>
#include <cmath>
#include <random>
double f(const double x, const double a, const double b, const double c)
{
return sin((a*x / (1.0 + pow(x, 2))) + 1.0) * atan(b*x - 1.0 / 2.0) + exp(-c*x) * atan(x);
}
double BisectionMethod(
double f(double, double, double, double),
const double a, const double b, const double c,
const std::random_device::result_type entropy)
{
std::mt19937 gen(entropy);
static const auto lower_bound = -1.0;
static const auto upper_bound = 1.0;
std::uniform_real_distribution<> dis(lower_bound, upper_bound);
auto pos_pt = dis(gen);
auto neg_pt = dis(gen);
while (f(pos_pt, a, b, c) < 0.0)
pos_pt = dis(gen);
while (f(neg_pt, a, b, c) > 0.0)
neg_pt = dis(gen);
static const auto about_zero_mag = 1E-8;
for (;;)
{
const auto mid_pt = (pos_pt + neg_pt)/2.0;
const auto f_mid_pt = f(mid_pt, a, b, c);
if (fabs(f_mid_pt) < about_zero_mag)
return mid_pt;
if (f_mid_pt >= 0.0)
pos_pt = mid_pt;
else
neg_pt = mid_pt;
}
}
int main()
{
double a, b, c;
std::random_device rd;
static const auto entropy = rd();
a =10, b = 2.0, c = 0.0;
const auto root1 = BisectionMethod(f, a, b, c, entropy);
std::cout << "a = " << a << ", b = " << b << ", c = " << c << std::endl;
std::cout << "Found root: (" << root1 << ", " << f(root1, a, b, c) << ")" << std::endl;
a =4.5, b = 2.8, c = 1.0;
const auto root2 = BisectionMethod(f, a, b, c, entropy);
std::cout << "a = " << a << ", b = " << b << ", c = " << c << std::endl;
std::cout << "Found root: (" << root2 << ", " << f(root2, a, b, c) << ")" << std::endl;
}
出力:
g++ -O3 -std=c++11 -Wall -Wextra -pedantic main.cpp -o root && ./root
a = 10, b = 2, c = 0
Found root: (0.143042, -2.12425e-09)
a = 4.5, b = 2.8, c = 1
Found root: (0.136172, 5.81247e-09)
これは RNG を使用するため、実行ごとに出力が変わります。視覚的に出力は正しく見えます。
コードは、ルートが -1.0 と 1.0 で囲まれていると想定していますが、これはあなたの場合に当てはまります。より一般的なものにしたい場合は、オーバーフローを処理し、ナンをチェックするためのロジックを追加する必要があります。ルートが -1.0 と 1.0 の間にない場合、これは永久にループします。それにもかかわらず、それはこの質問の特定の問題を解決し、より一般的なものの始まりです。
また、関数には複数のルートがあり、指定されたコードは単一のルートを見つけるだけであることに注意してください。
編集:コードをクリーンアップしました。数値メソッドについて話している場合に望ましいと思われる再現性があるようにentropy
、 の引数として追加されました。BisectionMethod()