2.4 Программирование

Ниже представлен листинг программы, реализующей метод наискорейшего спуска. Так как программа содержит очень большое число мелких функций, то все они были вынесены в отдельный файл «optim.cpp», претендующий на звание «библиотеки». Данный файл вставляется в листинг основной программы через свой заголовочный файл «optim.h». Исходный код основной программы хранится в «minimize.cpp».

// Файл: «minimize.cpp»

//------------------------------------------------------

#include <vcl.h>

#include <stdio.h>

#include <conio.h>

#pragma hdrstop

// prototypes

double func(int k, ...);

double f(double *x);

#include "optim.h"

//------------------------------------------------------

double func(int k, ...)

{

va_list par;

double *pp;

int i;

if(_var==NULL)

{

_inp=-2;

return 0.;

}

va_start(par,k);

pp=va_arg(par,double*);

for(i=0;i<k;i++)

{

_var[i]=*pp;

pp++;

}

va_end(par);

return f(_var);

}

double f(double *x)

{

return 4*x[0]*x[1]+5+(x[0]-3)*(x[0]-3)+(3*x[1]+2)*(3*x[1]+2);

}

#pragma argsused

int main(int argc, char* argv[])

{

int m=2,i;

double *dot=NULL, eps;

if((dot=set_array(m))==NULL)

{

printf("Failed to allocate memory for operational array\n");

return -1;

}

printf("Input the initital estimate\n Size of the estimate is %d\n",m);

for(i=0;i<m;i++)

{

printf("Variable %d: ",i+1);

scanf("%lf",&dot[i]);

}

printf("Input the precision: ");

scanf("%lf",&eps);

gradientLowering(func,m,dot,eps);

printf("Minimum of the function searched by gradient lowering method\n");

printf("[");

for(i=0;i<m;i++)

{

printf("%8.4lf",dot[i]);

if(i!=m-1)

printf(";");

}

printf("]\n");

free(dot);

getch();

return 0;

}

//------------------------------------------------------

// Файл-заголовок «optim.h»

//------------------------------------------------------

#include <math.h>

#include <stdlib.h>

#include <stdarg.h>

double *_lt=NULL;

double *_rt=NULL;

double *_var=NULL;

int _inp=0;

int setTempMas(int m);

void unsetTempMas();

double * set_array(int m);

int shEqu(int k,double *dot1,double *dot2);

int getMin(double (*func)(int k, ...),int m, double *left,double *right,double *res,double eps);

double derOfFunc(double (*func)(int k, ...),int m,double *dot,int p,double h);

int walk(double (*func)(int k, ...),int m,double *sDot, double *eDot, double *tDot, double *h, double eps);

void gradientLowering(double (*func)(int k, ...),int m, double *sDot, double eps);

void alarm();

#include "optim.cpp"

//------------------------------------------------------

// Файл «optim.cpp»

//------------------------------------------------------

void alarm()

{

switch(_inp)

{

case -1:printf("\nError 1: failed to allocate memory for operational arrays\n");

break;

case -2:printf("\nError 2: there was an error reading the arguments of the function f(x1,x2...xn)\n");

break;

case -4:printf("\nError 3: failed to allocate memory in gradientLowering()\n");

break;

case -6:printf("\nError 4: happened breakpoint in derOfFunc()\n");

break;

default:break;

}

unsetTempMas();

}

void gradientLowering(double (*func)(int k, ...),int m, double *sDot, double eps)

{

double *chkPoint=NULL, *eDot=NULL, *tDot=NULL, *grad=NULL;

int i,check;

chkPoint=set_array(m);

eDot=set_array(m);

tDot=set_array(m);

grad=set_array(m);

if(chkPoint==NULL || eDot==NULL || tDot==NULL || grad==NULL)

_inp=-4;

if(setTempMas(m)!=-1 && _inp>0)

{

for(i=0;i<m;i++)

chkPoint[i]=sDot[i];

while(1)

{

check=0;

for(i=0;i<m;i++)

{

grad[i]=derOfFunc(func,m,sDot,i+1,0.9);

grad[i]=-grad[i];

}

walk(func,m,sDot,eDot,tDot,grad,eps);

for(i=0;i<m;i++)

{

sDot[i]=tDot[i];

if(fabs(sDot[i]-chkPoint[i])<eps)

check=1;

else

check=0;

}

if(check)

break;

for(i=0;i<m;i++)

chkPoint[i]=tDot[i];

}

}

free(chkPoint);

free(eDot);

free(tDot);

free(grad);

alarm();

}

int walk(double (*func)(int k, ...),int m,double *sDot, double *eDot, double *tDot, double *h, double eps)

{

int i,t;

if(sDot==NULL || eDot==NULL || tDot==NULL || h==NULL)

return -1;

while(1)

{

for(i=0;i<m;i++)

eDot[i]=sDot[i]+h[i];

getMin(func,m,sDot,eDot,tDot,eps);

if((t=shEqu(m,sDot,tDot)) && t!=-1)

for(i=0;i<m;i++)

{

sDot[i]=tDot[i];

h[i]=-h[i];

}

else

if((t=shEqu(m,eDot,tDot)) && t!=-1)

for(i=0;i<m;i++)

sDot[i]=eDot[i];

else

break;

for(i=0;i<m;i++)

h[i]*=2;

}

if(t==-1)

return -1;

else

return 0;

}

int setTempMas(int m)

{

_lt=set_array(m);

_rt=set_array(m);

_var=set_array(m);

if(_lt==NULL || _rt==NULL || _var==NULL)

{

_inp=-1;

unsetTempMas();

return -1;

}

_inp=1;

return 0;

}

void unsetTempMas()

{

if(_lt!=NULL)

free(_lt);

if(_rt!=NULL)

free(_rt);

if(_var!=NULL)

free(_var);

_lt=_rt=_var=NULL;

_inp=0;

}

double *set_array(int m)

{

double *res=NULL;

res=(double *)calloc(m,sizeof(double));

if(res==NULL)

return NULL;

return res;

}

int shEqu(int k,double *dot1,double *dot2)

{

int i;

if(dot1==NULL || dot2==NULL)

return -1;

for(i=0;i<k;i++)

if(dot1[i]!=dot2[i])

return 0;

return 1;

}

int getMin(double (*func)(int k, ...),int m, double *left,double *right,double *res,double eps)

{

double length=0;

int i;

for(i=0;i<m;i++)

length+=(right[i]-left[i])*(right[i]-left[i]);

length=sqrt(fabs(length));

if(length<eps)

{

for(i=0;i<m;i++)

res[i]=(left[i]+right[i])/2;

return 0;

}

for(i=0;i<m;i++)

{

_lt[i]=(left[i]*2+right[i])/3;

_rt[i]=(left[i]+2*right[i])/3;

}

if(func(m,_lt)<func(m,_rt))

return getMin(func,m,left,_rt,res,eps);

else

return getMin(func,m,_lt,right,res,eps);

}

double derOfFunc(double (*func)(int k, ...),int m,double *dot,int p,double h)

{

int i;

if(_inp==0)

setTempMas(m);

if(dot==NULL)

_inp=-6;

if(p<=0 || p>m)

_inp=-6;

if(_inp<0)

return -1.;

if(h<=0.)

h=0.0001;

for(i=0;i<m;i++)

{

_lt[i]=dot[i];

_rt[i]=dot[i];

}

_lt[p-1]=_lt[p-1]-h;

_rt[p-1]=_rt[p-1]+h;

return (func(m,_rt)-func(m,_lt))/2*h;

}

//------------------------------------------------------

Основная программа представляет собой вызов функции gradientLowering, которая выполняет наискорейший спуск. Тело функции хранится в «optim.cpp». Функции необходимо знать математическую функцию, для которой нужно найти минимум, число аргументов в ней, массив с начальным приближением и точность вычисления. Размер массива с приближением должен соответствовать (или, по крайней мере, быть не меньше) числа переменных в функции, так как весь упор делается на это число. В случае удачи в массив с приближением запишется минимум. Функция ничего не возвращает, однако, в структуре используется аппарат по отлову всевозможных ошибок, кроме ошибок, связанных с вычислением функции.

Для создания массива с приближением используется функция set_array, тело которой хранится в «optim.cpp». Функции необходимо передавать желаемый размер массива. Функция возвращает указатель на массив в случае удачи и NULL значение в противном случае. Данная функция используется в некоторых функциях файла «optim.cpp» и может использоваться сама по себе.

Для объявления арифметической функции используется функция func. Для простоты, в данной реализации все арифметические операции функции (2.1) помещены в отдельную функцию f, который нужно передавать указатель на массив с аргументами. Функция func в данной реализации фактически генерирует массив, который передаёт потом f. Благодаря особенности func, исходная арифметическая формула может не иметь ограничений на число аргументов, главное знать их число. В основе func лежит механизм макросредств, обеспечивающих доступ к переменному списку параметров [2, стр. 249]. Данной функции необходимо передавать число аргументов (это обязательный параметр) и указатель на вектор с аргументами, который должен соответствовать числу аргументов (или, по крайней мере, быть не меньше этого числа). Функция перехватывает ошибку, когда ей передают NULL значение в качестве указателя.

Перейдем к непосредственному описанию функций файла «optim.cpp». В данном файле хранятся следующие функции:

• int setTempMas(int m);

• void unsetTempMas();

• double * set_array(int m);

• int shEqu(int k,double *dot1,double *dot2);

• int getMin(double (*func)(int k, ...),int m, double *left,double *right,double *res,double eps);

• double derOfFunc(double (*func)(int k, ...),int m,double *dot,int p,double h);

• int walk(double (*func)(int k, ...),int m,double *sDot, double *eDot, double *tDot, double *h, double eps);

• void gradientLowering(double (*func)(int k, ...),int m, double *sDot, double eps);

• void alarm().

Функция setTempMas инициализирует внутренние временные массивы, которыми пользуются некоторые функции файла «optim.cpp». Имена этих массивов зарезервированы и не могут использоваться, если подключен данный файл. Данные массивы имеют имена _lt, _rt, _var. Массивами _lt и _rt пользуются в своих целях функции getMin и derOfFunc; массивом _var пользуется func.

Функция unsetTempMas освобождает память от временных массивов.

Функция shEqu сравнивает две точки и возвращает 1, если это одна и та же точка, в противном случае 0. Используется в связке с функцией walk.

Функция getMin выполняет троичный поиск по алгоритму рисунка 2.4. В данной функции реализована рекурсия. Возвращает 0, если минимум найден, при этом результат записывается в массив res. Используется в связке с функцией walk.

Функция derOfFunc возвращает значение производной в точке. Функции необходимо передавать дифференцируемую функцию, число переменных, точку, в которой рассчитывается производная, номер переменной, по которой нужно дифференцировать, шаг. Здесь подразумевается, что переменные нумеруются от единицы, поэтому необходимо иметь это в виду при передаче номера. Если функции передан отрицательный или нулевой шаг, то для вычисления производной используется шаг 0,0001, рекомендуемый в некоторых источниках [4].

Функция walk фактически реализует поиск минимума на антиградиентном направлении. Функции необходимо передавать арифметическую функцию, число переменных, начальную точку, два временных массива (во второй в случае удачи записывается результат), шаговый массив (здесь это градиент) и точность.

Функция alarm служит для «отлова» ошибок, которые могут произойти в процессе работы программы, кроме ошибок, связанных с вычислением функции. Для оценки ситуации она пользуется глобальной переменной _inp, чьё имя зарезервировано и не может быть использовано, если подключен «optim.h». Все ошибки имеют свои коды, по которым alarm их вычисляет и выдаёт предупреждающие сообщения. В задачи alarm также входит очистка временных массивов функции setTempMas.

В данной реализации могут быть обнаружены ошибки указанные в таблице 2.

Таблица 2

Код ошибки	Сообщение для пользователя	Возможные причины
1	Не удалось выделить память под временные массивы.	Передан ненатуральный размер при вызове функции setTempMas, либо недостаточно памяти.
2	Не удалось прочитать аргументы функции.	Функции func передан NULL указатель.
3	Не удалось выделить память в функции gradientLowering().	Передан ненатуральный размер числа переменных, либо недостаточно памяти.
4	Аварийный останов при вычислении производной.	Функции derOfFunc передан NULL указатель вместо реальной точки, либо передан номер несуществующей переменной.

Математическая функция может иметь любую сложность, главное, чтобы она определялась через func.