У меня есть рабочий Integrator класс, который будет вычислять определенный интеграл базовых c функций одной переменной. Я протестировал интеграцию некоторых базовых c функций, и, похоже, он работает правильно.
Я сейчас нахожусь в той точке, где я хотел бы расширить этот класс, чтобы иметь возможность выполнять несколько интегралов та же функция ... и здесь я столкнулся с препятствием ...
Вот мой класс интегратора и некоторые основные c примеры использования:
Integrator.h
#pragma once
#include <algorithm>
#include <utility>
#include <functional>
struct Limits {
double lower;
double upper;
Limits() : lower{ 0 }, upper{ 0 } {}
Limits(double a, double b) : lower{ a }, upper{ b } {
if (a > b) std::swap(lower, upper);
}
void applyLimits(double a, double b) {
lower = a;
upper = b;
if (a > b) std::swap(lower, upper);
}
};
class Integrator {
private:
Limits limits_;
std::function<double(double)> integrand_;
double dx_;
double dy_;
double integral_;
int step_size_;
public:
Integrator(Limits limits, int stepSize, std::function<double(double)> integrand, double dy = 0)
: limits_{ limits },
step_size_{ stepSize },
integrand_{ integrand },
dx_{ 0 }, dy_{ 0 }
{}
~Integrator() = default;
constexpr double dx() const { return this->dx_; }
constexpr double dy() const { return this->dy_; }
constexpr double integral() const { return this->integral_; }
Limits limits() const { return limits_; }
std::function<double(double)>* integrand() { return &this->integrand_; }
// This is always a 1st order of integration!
constexpr double evaluate() {
double distance = limits_.upper - limits_.lower; // Distance is defined as X0 to XN. (upperLimit - lowerLimit)
dx_ = distance / step_size_; // Calculate the amount of iterations by dividing
// the x-distance by the dx stepsize
integral_ = 0; // Initialize area to zero
for (auto i = 0; i < step_size_; i++) { // For each dx step or iteration calculate the area at Xi
dy_ = integrand_(limits_.lower + i * dx_);
double area = dy_ * dx_; // Where the width along x is defines as dxStepSize*i
integral_ += area; // and height(dy) is f(x) at Xi. Sum all of the results
}
return integral_;
}
};
main. cpp
#include <iostream>
#include <exception>
#include <cmath>
#include "Integrator.h"
constexpr double PI = 3.14159265358979;
constexpr double funcA(double x) {
return x;
}
constexpr double funcB(double x) {
return (x*x);
}
constexpr double funcC(double x) {
return ((0.5*(x*x)) + (3*x) - (1/x));
}
double funcD(double x) {
return sin(x);
}
int main() {
try {
std::cout << "Integration of f(x) = x from a=3.0 to b=5.0\nwith an expected output of 8\n";
Integrator integratorA(Limits(3.0, 5.0), 10000, &funcA);
std::cout << integratorA.evaluate() << '\n';
std::cout << "\n\nIntegration of f(x) = x^2 from a=2.0 to b=20.0\nwith an expected output of 2664\n";
Integrator integratorB(Limits(2.0, 20.0), 10000, &funcB);
std::cout << integratorB.evaluate() << '\n';
std::cout << "\n\nIntegration of f(x) = (1\\2)x^2 + 3x - (1\\x) from a=1.0 to b=10.0\nwith an expected output of 312.6974\n";
Integrator integratorC(Limits(1.0, 10.0), 10000, &funcC);
std::cout << integratorC.evaluate() << '\n';
std::cout << "\n\nIntegration of f(x) = sin(x) from a=0.0 to b=" <<PI<< "\nwith an expected output of 2\n";
Integrator integratorD(Limits(0.0, PI), 10000, &funcD);
std::cout << integratorD.evaluate() << '\n';
} catch (const std::exception& e) {
std::cerr << e.what() << std::endl;
return EXIT_FAILURE;
}
return EXIT_SUCCESS;
}
Выход
Integration of f(x) = x from a=3.0 to b=5.0
with an expected output of 8
7.9998
Integration of f(x) = x^2 from a=2.0 to b=20.0
with an expected output of 2664
2663.64
Integration of f(x) = (1\2)x^2 + 3x - (1\x) from a=1.0 to b=10.0
with an expected output of 312.6974
312.663
Integration of f(x) = sin(x) from a=0.0 to b=3.14159
with an expected output of 2
2
Я думал добавить в этот класс еще одну функцию, похожую на его функцию evaluate() ... Сейчас она выглядит примерно так:
double integrate(Limits limits, double dy) {
double total = 0;
dy_ = dy;
for (int i = 0; i < step_size_; i++) {
double yi = limits_.lower*i*dy_;
double dx = static_cast<double>(yi - limits.lower) / stepSize;
double innerArea = 0;
for (int j = 0; j < step_size_; j++) {
Integrator inner(limits, step_size_, integrand_, dy_);
innerArea += inner.evaluate();
}
double outerArea = innerArea * dy_;
total += outerArea;
}
integral_ = total;
return integral_;
}
И здесь я запутался или в тупике ... Я не уверен, как правильно реализовать эту функцию, когда дело доходит до пределов интегрирования по внутренним и внешним интегралам.
Возьмем, например, следующий интеграл ниже:
Верхний предел внутреннего интеграла основан на y для каждой итерации вычислений ... Это должно выполняться динамически. Внешний интеграл прямолинейный, поскольку он идет от [3,5], а не от [1,y].
Я думаю, что я на правильном пути, но что-то в приведенном выше подходе полностью не соответствует действительности ... Я ' m получение совершенно неверных значений из ожидаемых или предполагаемых значений ...
Любые и все предложения и / или советы приветствуются!
Редактировать - Примечание. Я предоставил неверное изображение выше, которое было обновлено ...
Ожидаемый результат должен быть : 65.582 с правильно поставленной функцией f(x) = 1/2x^2 + 3x - (1/x). И когда я пытаюсь вычислить двойной интеграл, я получаю это ...
А вот добавленный код к программе драйвера или main.cpp ...
std::cout << "\n\nTesting Double Integration of f(x) = (1\\2)x^2 + 3x - (1\\x) from [3,5] and [1,y]\nwith an expected output of 65.582\n";
Integrator integratorE(Limits(3, 5), 1000, &funcC);
double dy = integratorE.limits().upper - integratorE.limits().lower;
integratorE.integrate(Limits(1, integratorE.integral()), dy);
std::cout << integratorE.integral() << '\n';
Однако , он ничего не выводит на консоль ...
Edit
Я не получал вывода, потому что не ждал достаточно долго. Итерации были определены как 1000 с помощью step_size. Это приведет к генерации 1000^1000 итераций ... Я упустил это из виду при создании объекта Integrator. Я изменил это в своем коде, чтобы step_size был 100. И теперь мое приложение выводит значение 2.68306e+189, что явно неверно! Когда я увеличиваю step_size до 500, это дает мне что-то порядка 6.62804e+190, что все еще неверно.