Контрольная работа по Информатике

КОНТРОЛЬНАЯ РАБОТА №2

Задание 1. Методы нахождения безусловного минимума

 

Составить программу нахождения минимума функции n переменных методом Флетчера-Ривза. Для одномерной минимизации по текущему направлению уточнения минимума производить методом деления пополам (MDP). Проверить программу функцией Пауэлла (П):

Анализ  функции:

Как видно, данная функция принимает минимальное  значение в точке 

X=(0, 0, 0, 0), при этом F(X) = 0

Проверим  это при помощи программы (листинг 1, рис.1).

Алгоритм  метода Флетчера-Ривса:

1. Задается начальная точка х₀, точность решения задачи минимизации e.
2. Вычисляем градиент целевой функции g₀=Ѧ(x₀).
3. Вычисляем вектор направления поиска p₀=-g₀.
4. Решаем задачу одномерной минимизации целевой функции

 f(х₀+ ap₀) относительно шага a одним из ранее рассмотренных методов.

5. Вычисляем новую точку х₁=х₀+a_minp₀. Для к=1,2,…, n -1 последовательно повторяем п. 6 – 10.
6.  Вычисляем градиент целевой функции g_k=Ñf(x_k).
7.  Вычисляем b_k

8. Вычисляем вектор направления поиска

p_k=-g_k+ b_kp_k-1

9. Решаем задачу одномерной минимизации целевой функции ¦(х_k+ap_k) относительно a.
10. Вычисляем новую точку х_k+1=х_k+a_{k min}p_k
11. Проверяем условие сходимости ||x_n–x₀|| < e  или  ||g_n|| < e.
12. Если заданная точность решения не достигнута, то переходим к п. 2, с заменой x₀ на x_n.

 

 

 

 

 

 

 

 

Результаты  работы программы:

 

Рис.1. Результаты работы программы

 

 

 

 

 

 

Листинг 1.

/*Составить программу  нахождения минимума функции  n переменных (функцию Пауэла)

методом Флетчера-Ривса. Одномерную минимизацию по текущему направлению,

реализовать методом  деления отрезка пополам.*/

 

#include <iostream>

#include <math.h>

 

#define PI 3.141592654

#define N 4

 

using namespace std;

 

double f(const int, double *);

double Grad(const int, double *,int);

double MDP(const int, double *,double *);

void copyVector(int razmer, double *istochnik, double *priemnik);

void copyVectorMinus(int razmer, double *istochnik, double *priemnik);

 

int main(){

    double M, k=0, alfa_k, beta, norm = 0, normpr = 0;

    double grx[N], x[N], xx[N], eps1, eps2;

    double xpr[N];// начальное приближение

    double pk[N];

    int q = 0;

    cout << "Vvedite nachalnie znachenia peremennih: " << endl;

    for(int j = 0; j < N; j++){

        cout << "x[" << j << "] = ";

        cin >> xpr[j];

        }

 

   

    cout<<"Pogreshnosty:   " << endl;

    cout<<"eps1:   ";

    cin>>eps1;

    cout<<"eps2:   ";

    cin>>eps2;

    cout<<"vvedite predelnoe chislo iteracii: ";

    cin>>M;

   

    //зададим начальное значение переменной хi

    copyVector(N, xpr, x);

       

    do{

       //определим градиент функции

       

        norm = 0;

        for(int j=0; j < N; j++){

            grx[j] = Grad(N, x, j);

            norm += grx[j] * grx[j];

            }

        if (pow(norm, 0.5) < eps1){

            copyVector(N, x, xx);

            cout<<"||f'(x)|| < eps1"<<endl;

            q=1;

            }

        if( k >= M){

            copyVector(N, x, xx);

            cout<<"k = M"<<endl;

            q = 1;

            }

        if(q==0){

            if (k==0)

                copyVectorMinus(N, grx, pk);

            else{

                beta = norm / normpr;

                for(int j=0; j < N; j++)

                    pk[j]=- grx[j] + beta * pk[j];

                           

                }

            if (pow(norm, 0.5)<eps1){

                copyVector(N, x, xx);

                cout<<"||f'(x)|| < eps1"<<endl;

                q=1;

                }

            if(k>=M){

                copyVector(N, x, xx);

                cout << "k = M" << endl;

                q = 1;

                }

            if(q == 0){

                if (k==0) copyVectorMinus(N, grx, pk);

                else{

                    beta = norm/normpr;

                    for(int j=0; j < N; j++)

                        pk[j]=-grx[j]+beta*pk[j];

                    }

                normpr = 0;

                for(int j=0; j < N; j++) {

                    normpr += grx[j] * grx[j];

                    }

                alfa_k = MDP(N, pk, x);

                for(int j=0; j < N; j++) 

                    x[j] -= alfa_k * pk[j];

                double delta_x = 0;

                for(int i = 0; i<N; i++)

                    delta_x += (x[i] - xpr[i])*(x[i] - xpr[i]);

                          

                if(((pow(delta_x, 0.5))<eps2)&&(abs(f(N, x)-f(N, xpr)))<eps2){

                    copyVector(N, x, xx);

                    cout<<"\n||xk-x|| < eps2 && |f(xk)-f(x)| < eps2 "<<endl<<"____________________________________"<<endl;

                    q = 1;

                    }

                copyVector(N, x, xpr);

                }

            }

        for(int j = 0; j <N; j++)

        ++k;

        }

    while(q == 0);

    for(int j = 0; j <N; j++)

        cout << "xx{" << j << "} = " << xx[j] << endl;

    cout << "f(x) = " << f(N, xx) << endl;

    system("pause");

    return 0;

    }

 

//сюда записываем  функцию, минимум которой необходимо  найти

//Функция Пауэла

double f(const int n, double *x){

    double g = pow(x[0]+10*x[1], 2);

    g += 5*pow(x[2]-x[3],2);

    g += pow(x[1]-2*x[2],4);

    g += 10*pow(x[0]+x[3],4); 

    return g;

}

 

void copyVector(int razmer, double *istochnik, double *priemnik){

    for (int i = 0; i < razmer; i++)

        priemnik[i] = istochnik[i];

    }

void copyVectorMinus(int razmer, double *istochnik, double *priemnik){

    for (int i = 0; i < razmer; i++)

        priemnik[i] = -istochnik[i];

    }

//производная по i-ой переменной

double Grad(const int n, double *x, int i){

    double h = 0.1e-7, ly, *lx;

    lx = new double [n];

    for(int j = 0; j<n; j++)

        lx[j] = x[j];

    ly = f(n, x);

    lx[i] = lx[i] - h;

    ly = (ly - f(n, lx)) / h;

    return ly;

    }

 

//находим оптимальное alfa_k

double MDP(const int n, double *pk, double *x){

    double a, b, ll, delta;

    double e=0.00001, l=0.0001, t=1.0, kk=1;

    double alfa_k=0, alfa_p=0, alfa_xk=0, alfa_yk=0;

    int p=0;

    //алгоритмом Свена находим отрезок локализации

    double *x0, *x1, *x2, *xp, *xk, *yk;

    x0 = new double [n];

    x1 = new double [n];

    x2 = new double [n];

    xp = new double [n];

    xk = new double [n];

    yk = new double [n];

 

    for(int j = 0; j < n; j++){

      // cout << pk[j] << "_______________x[j]" << x[j] << endl;

        x0[j] = x[j]-alfa_k*pk[j];

        x1[j] = x[j]-(alfa_k-t)*pk[j];

        x2[j] = x[j]-(alfa_k+t)*pk[j];

        }

 

    double f0 = f(n, x0);

    double f1 = f(n, x1);

    double f2 = f(n, x2);

   

    if (f1 >= f0 && f2 >= f0){

        a = alfa_k-t;

        b = alfa_k+t;

        p = 1;

        }

 

    if (f1 <= f0 && f2 <= f0){

        cout << "Funkziya ne unimodal'na " << endl ;

        p=1;

        }

 

    alfa_p = alfa_k;

 

    if (f1 > f0 && f0 > f2){

        delta=t;

        a=alfa_k ;

        alfa_k += t;

        }

   

    if (f1 < f0 && f0 < f2){

        delta=-t;

        b=alfa_k ;

        alfa_k=alfa_k-t;

        }

   

    double fp ;

    while (p != 1){        

        for(int j = 0; j < n; j++){

            x0[j] = x[j]-alfa_k*pk[j];

            xp[j] = x[j]-alfa_p*pk[j];

            }

        f0 = f(n, x0);

        fp = f(n, xp);

 

        if (f0 < fp )

            if(delta*t >0)

                a=alfa_p;

            else

                b=alfa_p;

        if (f0 > fp){

            if(delta*t >0)

                b=alfa_k;

            else if(delta*t<0)

                a=alfa_k;

            p=1;

            }

        kk++;

        alfa_p = alfa_k;

        alfa_k = alfa_p + pow(2.0, kk) * delta;

        }

    do{

        // методом деления пополам (методом дихотомии) находим alfa_k

        alfa_xk = (a + b)/2 - e;

        alfa_yk = (a + b)/2 + e;

        for(int j = 0; j < n; j++){

            xk[j] = x[j]-(alfa_xk)*pk[j];

            yk[j] = x[j]-(alfa_yk)*pk[j];

            }

        double fxk = f(n, xk);

        double fyk = f(n, yk);

 

        if (fyk >= fxk)

            b = alfa_yk;

        if (fyk < fxk)

            a = alfa_xk;

        ll = b-a;

        }

    while(l < ll);

    return (a + b) / 2.0;

    }

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Задание 2. Методы нахождения условного минимума

 

Составить программу нахождения минимума функции Пауэлла при наличии ограничения x₁²+x₂=1. Для методов, в которых на каждой итерации производится одномерная минимизация по текущему направлению метод уточнения минимума  - деление пополам. Сформулировав составную целевую функцию по методу штрафных функций (МШФ).

Метод штрафных функций основан на преобразовании исходной задачи с ограничениями  в одну задачу безусловной оптимизации  или в их последовательность. С  помощью функций-ограничений строят штрафную функцию, которая прибавляется к целевой функции исходной задачи, так, чтобы нарушение какого-либо из ограничений исходной задачи было невыгодным с точки зрения полученной задачи безусловной оптимизации.

В качестве штрафной будем использовать функцию:

 

В результате для  оптимизации получим функцию:

 

Результаты работы программы представлены на рис. 2

Алгоритм работы программы:

Начальный этап. Выбрать  . Выбрать начальную точку x1, параметр штрафа и число . Положить и перейти к основному этапу.

Основной этап. Первый шаг. При начальной точке  xk и параметре штрафа решить следующую задачу:

минимизировать

 

,(2.19)

 

где - целое.

Положить xk+1 равным оптимальному решению этой задачи и перейти ко второму шагу.

Второй шаг. Если , то остановиться. В противном случае положить . Заменить k на k+1 и перейти к первому шагу.

 

Рис. 2. Результаты работы программы

 

Листинг 2.

/*Составить программу  нахождения минимума функции  n переменных 

методом Флетчера-Ривза. Одномерную минимизацию по текущему направлению,

реализовать методом  деления отрезка пополам.*/

 

 

 

#include <iostream>

#include <math.h>

 

#define PI 3.141592654

#define N 4

 

using namespace std;

 

double FShtraf(int n, double *xin);             //штрафная функция

double f(const int n, double *X, double r);/целевая функция оптимизации n-переменных, Х - вектор входных параметров

double f1(const int n, double *xin);

void optFR(int n, double *xin, double *xout, double r);    //метод Флетчера-Ривса функции n переменных. xin - вектор входных параметров, xout -

double Grad(const int, double *,int, double r);           //поиск градиента направлений

double MDP(const int, double *,double *, double r);       //метод деления пополам

void copyVector(int razmer, double *istochnik, double *priemnik);       //копирование вектора istochnik в вектор priemnik

void copyVectorMinus(int razmer, double *istochnik, double *priemnik);  //копирование вектора -istochnik в вектор priemnik

 

int main(){

    double xpr[N];// начальное приближение

    double x_out[N]; //вектор выходных значений

    double r_shtraf = 0.01;

    double beta = 5;

    double eps = 0.00001;

 

    cout << "Vvedite nachalnie znachenia peremennih: " << endl;

    for(int j = 0; j < N; j++){

        cout << "x[" << j << "] = ";

        cin >> xpr[j];

        }

    cout << "Nachalnoe znachenie strafa (0.1): "   ;

    cin >> r_shtraf;

    cout << "Nachalnoe znachenie beta (10): "   ;

    cin >> beta;

   

    while(r_shtraf*FShtraf(N, x_out) >= eps){ 

        optFR(N, xpr, x_out, r_shtraf);

        cout << "r_shtraf=" << r_shtraf << " | x_out[0] = " << x_out[0] << "; x_out[1] = " << x_out[1] ;

        cout << " | f(x) = " << f1(N,x_out) << " | FShtraf = "<< FShtraf(N, x_out) << endl;

        if(r_shtraf*FShtraf(N, x_out) >= eps)

            r_shtraf *= beta;

        else

            break; 

        copyVector(N, x_out, xpr);

    }                         

    for(int j = 0; j <N; j++)

        cout << "xout{" << j << "} = " << x_out[j] << endl;

    cout << "f(x) = " << f(N, x_out, r_shtraf) << endl; 

    system("pause");

    return 0;

    }

//Штрафная функция

double FShtraf(int n, double *xin){

    return pow(pow(xin[0],2)-xin[1],2);

    }

//сюда записываем  функцию, минимум которой необходимо  найти

double f1(const int n, double *x){

    double g = pow(x[0]+10*x[1], 2);

    g += 5*pow(x[2]-x[3],2);

    g += pow(x[1]-2*x[2],4);

    g += 10*pow(x[0]+x[3],4);

    return g;

}

//Функция Пауэла

 

double f(const int n, double *x, double r){

    double g = pow(x[0]+10*x[1], 2);

    g += 5*pow(x[2]-x[3],2);

    g += pow(x[1]-2*x[2],4);

    g += 10*pow(x[0]+x[3],4); 

    g += pow(pow(x[0],2)+x[1]-1,2);

 

    return g;

}

 

//Метод Флетчера-Ривса

void optFR(int n, double *xin, double *xout, double r_shtraf){

    double k=0, alfa_k, beta, norm = 0, normpr = 0;

    double grx[N], x[N];

    double pk[N];

    double eps = 0.001;

    int M = 7;

    int q = 0;

   

  

        //зададим начальное значение переменной хi

    copyVector(N, xin, x);

       

    do{

       //определим градиент функции

       

        norm = 0;

        for(int j=0; j < N; j++){

            grx[j] = Grad(N, x, j, r_shtraf);

            norm += grx[j] * grx[j];

            }

        if (pow(norm, 0.5) < eps){

            copyVector(N, x, xout);

        //    cout<<"||f'(x)|| < eps"<<endl;

            q=1;

            }

        if( k >= M){