Advertisement
Not a member of Pastebin yet?
Sign Up,
it unlocks many cool features!
- #include <iostream>
- #include <cmath>
- #include <complex>
- #include <vector>
- #include <fstream>
- using namespace std;
- ofstream out(“output.txt”);
- complex<double> a0(0, 1), a1(6, 0), a2(3, 0), a3(20, 0), a4(10, 0); //коэффициенты исходного полинома
- vector<complex<double>> P = { a4, a3, a2, a1, a0 };
- complex<double> li(0, 1); //комплексная единица
- complex<double> Horner(complex<double> z, vector<complex<double>> poly){ //функция, раскладывающая многочлен по схеме Горнера
- complex<double> answer(0, 0);
- for (int i = 0; i < poly.size(); i++) { answer = answer * z + poly[i]; }
- return answer;
- }
- double func(complex<double> x, vector<complex<double>> poly){ // вычисляет значение функции
- complex<double> h = Horner(x, poly);
- return abs(h)*abs(h);
- }
- complex<double> CoordDescent(complex<double> z, vector<complex<double>> poly, double eps){ // координатный спуск
- double a, b = 0.5; // шаг по вещественной и мнимой части
- double current; //служебная переменная для вычисления функции
- double F = func(z, poly);
- while (F > eps){
- if (F > (current = func(z + b, poly))) { z += b; F = current; }
- else if (F > (current = func(z - b, poly))) { z -= b; F = current; }
- else { b /= 2; }
- if (F <= eps) break;
- if (F > (current = func(z + li * a, poly))) { z += li * a; F = current; }
- else if (F > (current = func(z - li * a, poly))) { z -= li * a; F = current; }
- else { a /= 2; }
- }
- return z;
- }
- complex<double> deriv(complex<double> z, vector<complex<double>> poly){
- complex<double> answer(0, 0);
- double n = poly.size() - 1;
- for (int i = 0; i < n; i++) { answer = answer * z + poly[i] * (n - i); }
- return answer;
- }
- double fi;
- complex<double> grad(complex<double> z, vector<complex<double>> poly) { //вычисление градиента
- complex<double> h = Horner(z, poly);
- complex<double> dh = deriv(z, poly);
- complex<double> gr = 2.0*(real(dh*conj(h)) - li * imag(dh*conj(h)));
- fi = real(abs(h)*abs(h));
- return gr;
- }
- complex<double> CrushingStep(complex<double> z, vector<complex<double>> poly, double eps){ //градиентный метод с дроблением шага
- double a, delta = 0.5;
- double current;
- complex<double> zk;
- complex<double> zi = z;
- complex<double> gr = grad(z, poly);
- while (fi > eps){
- while (true){
- zk = zi - a * gr;
- if (((current = func(zk, poly)) - fi) <= -a * delta * (real(gr)*real(gr) + imag(gr)*imag(gr))){
- zi = zk; fi = current; break;
- }
- else {
- a /= 2;
- }
- }
- gr = grad(zi, poly);
- a *= 2;
- }
- return zi;
- }
- complex<double> ConstantStep(complex<double> z, vector<complex<double>> poly, double eps){ //градиентный метод с постоянным шагом
- complex<double> zi = z;
- complex<double> gr = grad(z, poly);
- double a = 0.00005;
- while (fi > eps){
- zi -= a * gr;
- gr = grad(zi, poly);
- }
- return zi;
- }
- complex<double> GivenStep(complex<double> z, vector<complex<double>> poly, double eps){ //градиентный метод с заданным шагом
- complex<double> zi = z;
- complex<double> gr = grad(zi, poly);
- double k = 1001;
- while (fi > eps){
- k++;
- zi -= 1 / k * gr;
- gr = grad(zi, poly);
- }
- return zi;
- }
- vector<complex<double>> polynom = P;
- void division(complex<double> z){ //деление полинома на корень
- vector<complex<double>> temp = polynom;
- polynom.clear();
- polynom.push_back(temp[0]);
- for (int i = 0; i < temp.size() - 2; i++) { polynom.push_back(z*polynom[i] + temp[i + 1]); }
- }
- void SearchRoot(complex<double> x0, vector<complex<double>> poly, double eps, int method){// поиск корней
- complex<double> root, n_root;
- double E = 2 * eps;
- auto alch = GivenStep;
- polynom = P;
- switch (method){
- case 1: alch = CoordDescent; out << “Координатный спуск : “ << endl; break;
- case 2: alch = CrushingStep; out << “Дробление шага : “ << endk; break;
- case 3: alch = ConstantStep; out << “Постоянный шаг : “ << endl; break;
- case 4: alch = GivenStep; out << “Заданный шаг : “ << endl; break;
- default: break;
- }
- fi = 1;
- root = alch(x0, poly, eps);
- division(root);
- for (int i = 1; i < poly.size() - 2; i++){
- n_root = alch(x0, polynom, E);
- root = alch(n_root, poly, E);
- division(root);
- E *= 2;
- }
- n_root = -polynom[1] / polynom[0];
- root = alch(n_root, P, eps);
- }
- int main(){
- complex<double> point(0, 0); //точка, из которой мы вызываем методы
- SearchRoot(point, P, 1e-6, 1);
- SearchRoot(point, P, 1e-6, 2);
- SearchRoot(point, P, 1e-6, 3);
- SearchRoot(point, P, 1e-6, 4);
- out.close();
- return 0;
- }
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement