//    ##################################################################
//    ##           К У Р С О В О Й  П Р О Е К Т                       ##
//    ##  по предмету  :    Исследование      операций                ##
//    ##  по теме      :    Метод деформируемого многогранника        ##
//    ##  студент      :    Ратаев  М.Ю.                              ##
//    ##  преподаватель:    Лемешко Б.Ю.                              ##
//    ##################################################################

#include "ini.h"  //    инициализация всех библиотек
   Matrix X(7,2);            // указатель на массив координат

   float  f( Vector x )
      { return ff( x[0],x[1] );}

   void	 Simplex();          // инициализация
			     // построение первоначального правильного симплекса
			     //  с одной из вершин в начале координат

    void Podgotovka( void ); // вычисление центра тяжести и
			     // MIN и MAX значения функции
    void reflection( void ); // отражение
    void extension( void );  // растяжение
    void compression( void );// сжатие
    void reduction( void );  // редукция
   float KritOstan( void );  // критерий останова

    void draw( void );       // графический вывод результатов вычислений



// описания  функций    Simplex

void Simplex()
{
  int i,j;
  float d1,d2;
// выделение памяти
//  X=new float [(N+5)*N];
//  Nkoor=N;
// построение первоначального правильного симплекса
//  с одной из вершин в начале координат
  for( i=0; i<Nkoor; i++ )
     X[0][i]=0;
  d1=t*(sqrt(Nkoor+1)+Nkoor-1)/(Nkoor*sqrt(2.));
  d2=t*(sqrt(Nkoor+1)-1)/(Nkoor*sqrt(2.));
  for( i=1; i<Nkoor+1; i++ )
    for( j=0; j<Nkoor; j++ )
       if(i-1==j) X[i][j]=d1;
	   else   X[i][j]=d2;

  return;
}
void Podgotovka()            // вычисление центра тяжести и
					     // MIN и MAX значения функции
{
   int i,j;
   float fmmm,sum=0;
   fXh=f(X[0]);
   h=0;
   fXl=fXh;
   l=0;
   for( j=1; j<Nkoor+1; j++ )
     { fmmm=f(X[j]);
       if( fmmm <  fXl ) {  fXl=fmmm; l=j; }
       if( fmmm >= fXh ) {  fXh=fmmm; h=j; }
     }
//подсчёт координат центра тяжести симплекса
    for( i=0; i<Nkoor; i++ )
       {   sum=0;
	   for( j=0; j<Nkoor+1; j++ )
		       sum+=X[j][i];
		 X[Nkoor+1][i]=(sum-X[h][i])/Nkoor;
	}
return;
}


void reflection()        // отражение ( проектирование вершины с MAX
			     //             значением функции через центр
			     //             тяжести  )
{ int i;
 // for( i=0; i<Nkoor; i++ )
      X[Nkoor+2]=X[Nkoor+1]+(X[Nkoor+1]-X[h])*alfa;
      // изменение флагов
  flagK=0;
//  K++;
  return;
}


void extension()         // растяжение ( спроектированный через
			     //              центр тяжести вектор
			     //              растягивается )
{ int i;
   if( f(X[Nkoor+2]) <= f(X[l]))
      { //   for( i=0; i<Nkoor; i++ )
	       X[Nkoor+3]=X[Nkoor+1]+(X[Nkoor+2]-X[Nkoor+1])*gamma;

	 if( f(X[Nkoor+3]) < f(X[l]))
	     {   // for( i=0; i<Nkoor; i++ )
			X[h]=X[Nkoor+3];
	     }
	    else
	     {   // for( i=0; i<Nkoor; i++ )
			X[h]=X[Nkoor+2];
	      }
	   flagK=1;K++;
       }
return;
}

void compression()       // сжатие

{  int i,j;
   int flagE=1;              //флаг сжатия : 0- ни разу не выполнилось
			     //              1- выполнилось хотя бы раз

  if( f(X[Nkoor+2]) < f(X[h])) flagE=0;
//   for( j=0; j<Nkoor+1; j++ )
//    {  if( j!=h )
//	 {
//		  if( f(X[Nkoor+2]) > f(X[j]));
//	     else  flagE=1;
//	 }
//    }
   if( flagE==0 )
     {
	// for( i=0; i<Nkoor; i++ )
		  X[Nkoor+4]=X[Nkoor+1]+(X[h]-X[Nkoor+1])*betta;

	 //for( i=0; i<Nkoor; i++ )
		  X[h]=X[Nkoor+4];
	 flagK=1;  K++;
      }
return;
}



void reduction()         // редукция  ( уменьшение  рёбер симплекса
			     //            в delta раз  )
			     //
{  int i,j;

//  if( f(X[Nkoor+2]) > f(X[h]))
//       {
	       for( j=0; j<Nkoor+1; j++ )
		  //for( i=0; i<Nkoor; i++ )
			 X[j]=X[l]+(X[j]-X[l])*delta;
	 t*=delta;
//	 K--;
//       }
//    else X[h]=X[Nkoor+2];
  K++;
  flagK=1;
 return;
 }



float KritOstan()        // критерий останова
 {
  int i,j;
  float sum=0;
    for( j=0; j<Nkoor+1 ; j++ )
       sum+=(f(X[j])-f(X[Nkoor+1]))*(f(X[j])-f(X[Nkoor+1]));
  return sqrt((sum/Nkoor) );
 }



void draw()
{
  VIDEO(K);

return;
}

FILE *outfile;
//    головная   программа
void main( void )
{ int i=0;
  char ima[20];
  printf("\n\n введите имя выходного файла :");
  scanf("%s",ima);
  outfile=fopen( ima, "a");
  fprintf(outfile,"\n       МЕТОД ДЕФОРМИРУЕМОГО МНОГОГРАННИКА \n ");
  fprintf(outfile,"\n файл выходных данных для функции : ............................ \n");
  fprintf(outfile,"\n .............................................................. \n\n");

  iniGraf();


 while (flagPeresch==1)
 { t=T1;
   i=0;
 Simplex();
 while((i<98)&&((ee=KritOstan()) > e) )
   { i++;
	 Podgotovka();       // вычисление центра тяжести и
					     // MIN и MAX значения функции
	     reflection();       // отражение
	     extension();        // растяжение
    if(flagK!=1)
      {
	     compression();      // сжатие
	   if(flagK!=1)
	      {
		     reduction();        // редукция
	      }
      }
     Xx[K]=X[Nkoor+1][0]; //Nkoor+1
     Yy[K]=X[Nkoor+1][1];
     Zz[K]=fXl;
     XxX[K][0][0]=X[0][0] ;
     XxX[K][0][1]=X[0][1] ;
     XxX[K][1][0]=X[1][0] ;
     XxX[K][1][1]=X[1][1] ;
     XxX[K][2][0]=X[2][0] ;
     XxX[K][2][1]=X[2][1] ;

   }
	 flagPeresch=0;
	 colT=K;
 //вывод результатов в файл

 fprintf(outfile,"alfa=%6.4f betta=%6.4f gamma=%6.4f delta=%6.4f e=%7.5f t=%6.4f \n",alfa,betta,gamma,delta,e,T1);
 fprintf(outfile,"экстремум  : f( %f , %f ) = %f \n",Xx[K],Yy[K],Zz[K]);
 fprintf(outfile,"с точностью:  %f    ",ee);
 fprintf(outfile,"за  %d  шагов   \n\n",K);

	 draw();             // графический вывод результатов вычислений
	 K=1;
 }
   // КОНЕЦ
   DeiniGraf();
 fclose(outfile);
}





