Минимизации методом случайного поиска

Функция

[pic 1]

Листинг кода минимизации методом случайного поиска для функции с ограничениями

#include "stdafx.h"

#include 

#include 

#include 

#include

#include 

using namespace std;

#define K 2// количество переменных

double x_k[K];

double x_k_1[K];

double vec[K];

double S;//шаг

double E;//точность

double Y;//коэф уменьшения шага

int N;//число неудачных попыток

int j;//счетчик

double n_vec;

int C;// увеличение штрафа

double P;//штраф

double y(double *x, int n)

{

        return (pow(x[0], 2) + pow(x[1], 2) - 20 * x[0] - 30 * x[1]); //функция

}

const double goldenRatio = (1 + sqrt(5)) / 2;

double G1(double *x, int n)//функция ограничения 1

{

        double limit;

        limit = 2 * x[0] + 3 * x[1] - 13;

        if (limit <= 0) limit = 0;

        return limit;

};

double G2(double *x, int n)//функция ограничения 2

{

        double limit;

        limit = 2*x[0] +  x[1] - 10;

        if (limit <= 0) limit = 0;

        return limit;

};

//double G3(double *x, int n)//функция ограничения 3

//{

//        double limit;

//        limit =  x[1];

//        if (limit <= 0) limit = 0;

//        return limit;

//};

double function(double *x, int n, double s)// функция для рачета шага

{

        double m = (pow((x[0] + s*vec[0] / n_vec), 2) + pow((x[1] + s*vec[1] / n_vec), 2) - 20 * (x[0] + (s*vec[0] / n_vec)) - 30 * (x[1] + (s*vec[1] / n_vec)));

        if (G1(x_k, K) > 0) m = m + (P / 2)*(pow((2 * (x[0] + (s*vec[0] / n_vec)) + 3 * (x[1] + (s*vec[1] / n_vec)) - 13), 2));

        if (G2(x_k, K) > 0) m = m + (P / 2)*(pow((2 * (x[0] + (s*vec[0] / n_vec)) +  (x[1] + (s*vec[1] / n_vec)) - 10), 2));

        //if (G3(x_k, K) > 0) m = m + (P / 2)*pow((x[1] + (s*vec[1] / n_vec)), 2);

        return m;//функция

}

double Step()// одномерная минимизация

{

        int i = 0;

        double a, b;

        double accuracy = 0.00001;

        double x1, x2; // Точки, делящие текущий отрезок в отношении золотого сечения

        a = -10;b = 20;

        while (abs(b - a) > accuracy)

        {

                i++;

                x1 = b - (b - a) / goldenRatio;

                x2 = a + (b - a) / goldenRatio;

                if (function(x_k, K, x1) >= function(x_k, K, x2))

                        a = x1;

                else

                        b = x2;

        }

        return (a + b) / 2;

}

double GPenalty(double *x, int n)//штрафная функция

{

        double g;

        g = (P / 2)*(pow(G1(x, K), 2) + pow(G2(x, K), 2));

        return g;

};

double curr_func(double *x, int n)//штрафная функция

{

        double f;

        f = y(x, K) + GPenalty(x, K);

        return f;

};

int main(void)

{

        setlocale(LC_CTYPE, "Russian");

        P = 0.1;

        C = 2;

        E = 0.00001;

        Y = 5;