Подпрограмма:  mnbbr_c

Назначение

Решение задачи минимизации функции многих переменных при наличии ограничений методом скользящего допуска.

Математическое описание

Для решения задачи
     min  f(x) ,   x ∈ En .

При ограничениях    h i(x)  =  0 ,   i = 1, ..., m ,
                                    g i(x) ≥ 0 ,    i = m + 1, ..., P , 

используется метод скользящего (нежесткого) допуска.

B соответствии с алгоритмом исходная задача заменяется следующей

     min  f(x) ,    x ∈ En 

при ограничениях    Ф^(k) - T (x) ≥ 0,

где  Ф^(k) ≥ 0 - значение критерия допуска для нарушения ограничений решаемой задачи на  k - ом шаге алгоритма, а  T (x) ≥ 0 - функционал над множеством всех функций, задающих ограничения в исходной задаче.

Hа каждом шаге алгоритма задача безусловной минимизации решается методом деформируемого многогранника (Нельдера и Мида).

 Функция  Ф имеет вид:

      Ф(k)  =  min  { Ф(k-1) ,
                                             r+1
                                         [   ∑   || xi(k) - x(k)r+2 ||  ] (m+1) / (r+1) } ,
                                             i =1
      Ф0  =  2 (m+1) A ,
 где 

A - величина, характеризующая размер исходного многогранника (см. замечания по использованию);

m - число ограничений в виде pавенств;

x_i^(k) - вектоp, задающий положение  i - ой вершины многогранника в пространстве  Eⁿ;

r = (n - m) - число степеней свободы целевой функции  f (x);

x^(k)_r + 2 - вектоp, задающий положение вершины, которая соответствует центру тяжести рассматриваемого многогранника при  n = r;

Ф^(k - 1) - значение  Ф на (k - 1) - ом шаге алгоритма;

k = 0, 1,... - индекс, указывающий число полностью законченных шагов алгоритма.

Функционал  T (x) является мерой степени нарушения ограничений и имеет вид:

                           m                   P
       T(x)  =  + [  ∑  hi2(x)   +  ∑     Ui gi2(x) ]1/2 ,
                          i =1               i =m+1
   где  Ui  =  0   при  gi(x) ≥ 0    и   Ui  =  1   при  gi(x) < 0. 

Работа алгоритма заканчивается, если выполнено хотя бы одно из условий:

1. Ф(k) < 1.E - 8 ,
       n+1
2. (  ∑  ( f(xi(k)) - f(x(k)r+2) )2 / n )1/2  <  ACC ,
       i =1
где 

f (x_i^(k)) - значение целевой функции в  i - ой вершине многогранника на  k - ом шаге алгоритма;

f (x^(k)_r + 2) - значение целевой функции в центре тяжести многогранника на  k - ом шаге алгоритма;

ACC - точность вычисления минимума целевой функции.

Д.Химмельблау, Прикладное нелинейное программирование, Изд - во "Мир", Mосква, 1975.

Использование

    int mnbbr_c (integer *n, integer *nmp1, integer *itmax,
            integer *itmt, real *alfa, real *beta, real *gam, real *acc, real *a, 
            real *xyz, real *xx, real *vec, real *ff, integer *maxk,
            integer *maxt, S_fp fun, S_fp funt, integer *ierr)

   Параметры 

   Версии:  нет 

   Вызываемые подпрограммы 

   Замечания по использованию  

               alfa = 1
      0.4 ≤ beta ≤ 0.6
      2.8 ≤ gam  ≤ 3.0 

          int fun(float *x, float *f, float *fe)

      Параметры
      x   - вещественный вектоp длины  n, задающий
              точку пространства, в которой вычисляется
              значение функции;
      f   - вещественная переменная, содержащая значение
              функции в точке  x;
      fe - вещественная  переменная, задающая точность
              вычисления значения функции в точке  x. 

          int funt(float *x, float *f, float *fe)

      Параметры
       x    - вещественный вектор длины  n, задающий
                точку пространства, в которой вычисляется
                значение функционала;
       f    - вещественная переменная, содержащая значение
                функционала в точке  x;
      fte - вещественная  переменная, задающая точность
                вычисления значения функционала в точке  x. 

Пример использования

int main(void)
{
    /* Initialized data */
    static float alfa = 1.f;
    static int nmp1 = 3;
    static int itmt = 200;
    static int maxt = 500;
    static float beta = .5f;
    static float gam = 2.f;
    static int n = 2;
    static int itmax = 50;
    static int maxk = 400;
    static float acc = 1e-6f;
    static float a = .3f;
    static float xx[2] = { 1.f,1.f };

    /* Local variables */
    static int ierr;
    extern int t_c();
    extern int mnbbr_c(int *, int *, int *, int *, float *, float *,
                       float *, float *, float *, float *, float *, float *,
                       float *, int *, int *, U_fp, U_fp, int *);
    extern int f0_c();
    static float ff, vec[16] /* was [2][8] */,
                     xyz[18] /* was [3][6] */;

    printf("\n %5i \n", n);
    printf("\n %10.4e \n", acc);
    printf("\n %5.2f %5.2f %5.2f \n", alfa, beta, gam);
    printf("\n %16.7e %16.7e \n", xx[0], xx[1]);
    mnbbr_c(&n, &nmp1, &itmax, &itmt, &alfa, &beta, &gam, &acc, &a, xyz, xx,
            vec, &ff, &maxk, &maxt, (U_fp)f0_c, (U_fp)t_c, &ierr);

    printf("\n\n %5i \n", maxt);
    printf("\n %5i \n", ierr);
    printf("\n %5i \n", itmax);
    printf("\n %5i \n", maxk);
    printf("\n %16.7e \n", ff);
    printf("\n %16.7e %16.7e \n", xx[0], xx[1]);
/* l1: */
    return 0;
} /* main */

int f0_c(float *x, float *f, float *fe)
{
    /* System generated locals */
    float r__1;

    /* Parameter adjustments */
    --x;

    /* Function Body */
/* Computing 2nd power */
    r__1 = x[2];
    *f = x[1] * 4.f - r__1 * r__1 - 12.f;
    return 0;
} /* f0_c */

int t_c(float *x, float *f, float *fe)
{
    /* System generated locals */
    float r__1, r__2;

    /* Builtin functions */
    double sqrt(double);

    /* Local variables */
    static float r__;

    /* Parameter adjustments */
    --x;

    /* Function Body */
/* Computing 2nd power */
    r__1 = x[1];
/* Computing 2nd power */
    r__2 = x[2];
    *f = 25.f - r__1 * r__1 - r__2 * r__2;
/* Computing 2nd power */
    r__1 = *f;
    *f = r__1 * r__1;
/* Computing 2nd power */
    r__1 = x[1] - 5.f;
/* Computing 2nd power */
    r__2 = x[2] - 5.f;
    r__ = r__1 * r__1 + r__2 * r__2 - 16.f;
    r__ = -r__;
    if (r__ < 0.f) {
/* Computing 2nd power */
        r__1 = r__;
        *f += r__1 * r__1;
    }
    if (x[1] < 0.f) {
/* Computing 2nd power */
        r__1 = x[1];
        *f += r__1 * r__1;
    }
    if (x[2] < 0.f) {
/* Computing 2nd power */
        r__1 = x[2];
        *f += r__1 * r__1;
    }
    *f = (float)sqrt(*f);
    return 0;
} /* t_c */


Результаты:

      ierr   =  65
      ff       =  -0.3199231 + 02

      xx(1)  =  1.0012830 - 00
      xx(2)  =  4.8987190 - 00