//
// russian : cp1251 encoding
// http://mech.math.msu.su/~iliagri/ext_value.htm
// version 25.10.2023
//
#include <iostream>
#include <cmath>
using namespace std;

class ext_value3
{
 int n_var;
 double val, // значение
        *der_1, // массив первых производных
        *der_2, // массив вторых производных
        *der_3; // массив третьих производных
 public:
// конструктор по-умолчанию для создания массива объектов
 ext_value3() { n_var=0; val=0; der_1=der_2=der_3=NULL; } 
// деструктор: нужно удалять динамически выделенную память (иначе утечка памяти)
 ~ext_value3() { if(der_1)delete [] der_1; if(der_2) delete [] der_2; if(der_3) delete [] der_3; n_var=0; }
// конструктор копирования
 ext_value3(const ext_value3 &ob);
// конструктор для создания независимой переменной или константы
 ext_value3(int n, double a=0.0, int m=-1);
// установка размерности (для использования конструктора по-умолчанию)
 void set_dim(int n); 
// установка значения 
 void set(double a){ val=a; }
// установка значения i-й компоненты первой производной
 void set(double a,int i);
// установка значения i,j-й компоненты второй производной
 void set(double a,int i,int j);
// установка значения i,j,k-й компоненты третьей производной
 void set(double a,int i,int j,int k);
// получить значение
 double get() { return val; }
// получить значение i-й компоненты первой производной 
 double get(int i); 
// получить значение i,j-й компоненты второй производной
 double get(int i,int j);
// получить значение i,j,k-й компоненты третьей производной
 double get(int i,int j,int k);
 ext_value3& operator =(const ext_value3 &a);
 ext_value3& operator +=(const ext_value3 &a);
 ext_value3& operator +=(double a);
 ext_value3 operator +(const ext_value3 &a);
 ext_value3 operator +(double a);
 friend ext_value3 operator +(double b,const ext_value3 &a);
 ext_value3& operator *=(const ext_value3 &a);
 ext_value3& operator *=(double a);
 ext_value3 operator *(const ext_value3 &a);
 ext_value3 operator *(double a);
 friend ext_value3 operator *(double b,const ext_value3 &a);
 ext_value3& operator -=(const ext_value3 &a);
 ext_value3& operator -=(double a);
 ext_value3 operator -(const ext_value3 &a);
 ext_value3 operator -(double a);
 friend ext_value3 operator -(double b,const ext_value3 &a);
// унарный минус:
 ext_value3 operator -();
// унарный плюс не реализован. Нужен?
 ext_value3 operator /(const ext_value3 &a);
 ext_value3 operator /(double a);
 ext_value3& operator /=(const ext_value3 &a);
 ext_value3& operator /=(double a);
 friend ext_value3 operator /(double b,const ext_value3 &a);
 friend ext_value3 sin(const ext_value3 &a);
 friend ext_value3 cos(const ext_value3 &a);
 friend ext_value3 sqrt(const ext_value3 &a);
 friend ext_value3 exp(const ext_value3 &a);
 friend ext_value3 atan2(const ext_value3 &a,const ext_value3 &b);
 friend ext_value3 asin(const ext_value3 &a);
 friend ext_value3 acos(const ext_value3 &a);
};

class ext_value2
{
 int n_var;
 double val, // значение
        *der_1, // массив первых производных
        *der_2; // массив вторых производных
 public:
// конструктор по-умолчанию для создания массива объектов
 ext_value2() { n_var=0; val=0; der_1=der_2=NULL; } 
// деструктор: нужно удалять динамически выделенную память (иначе утечка памяти)
 ~ext_value2() { if(der_1)delete [] der_1; if(der_2) delete [] der_2; n_var=0; }
// конструктор копирования
 ext_value2(const ext_value2 &ob);
// конструктор для создания независимой переменной или константы
 ext_value2(int n, double a=0.0, int m=-1);
// установка размерности (для использования конструктора по-умолчанию)
 void set_dim(int n);
// установка значения 
 void set(double a){ val=a; }
// установка значения i-й компоненты первой производной
 void set(double a,int i);
// установка значения i,j-й компоненты второй производной
 void set(double a,int i,int j);
// получить значение
 double get() { return val; }
// получить значение i-й компоненты первой производной 
 double get(int i); 
// получить значение i,j-й компоненты второй производной
 double get(int i,int j);
 ext_value2& operator =(const ext_value2 &a);
 ext_value2& operator +=(const ext_value2 &a);
 ext_value2& operator +=(double a);
 ext_value2 operator +(const ext_value2 &a);
 ext_value2 operator +(double a);
 friend ext_value2 operator +(double b,const ext_value2 &a);
 ext_value2& operator *=(const ext_value2 &a);
 ext_value2& operator *=(double a);
 ext_value2 operator *(const ext_value2 &a);
 ext_value2 operator *(double a);
 friend ext_value2 operator *(double b,const ext_value2 &a);
 ext_value2& operator -=(const ext_value2 &a);
 ext_value2& operator -=(double a);
 ext_value2 operator -(const ext_value2 &a);
 ext_value2 operator -(double a);
 friend ext_value2 operator -(double b,const ext_value2 &a);
// унарный минус:
 ext_value2 operator -();
 ext_value2 operator /(const ext_value2 &a);
 ext_value2 operator /(double a);
 ext_value2& operator /=(const ext_value2 &a);
 ext_value2& operator /=(double a);
 friend ext_value2 operator /(double b,const ext_value2 &a);
 friend ext_value2 sin(const ext_value2 &a);
 friend ext_value2 cos(const ext_value2 &a);
 friend ext_value2 sqrt(const ext_value2 &a);
 friend ext_value2 exp(const ext_value2 &a);
 friend ext_value2 atan2(const ext_value2 &a,const ext_value2 &b);
 friend ext_value2 asin(const ext_value2 &a);
 friend ext_value2 acos(const ext_value2 &a);
};

class ext_value1
{
 int n_var;
 double val, // значение
        *der_1; // массив первых производных
 public:
// конструктор по-умолчанию для создания массива объектов
 ext_value1() { n_var=0; val=0; der_1=NULL; } 
// деструктор: нужно удалять динамически выделенную память (иначе утечка памяти)
 ~ext_value1() { if(der_1)delete [] der_1; n_var=0; }
// конструктор копирования
 ext_value1(const ext_value1 &ob);
// конструктор для создания независимой переменной или константы
 ext_value1(int n, double a=0.0, int m=-1);
// установка размерности (для использования конструктора по-умолчанию)
 void set_dim(int n);
// установка значения 
 void set(double a){ val=a; }
// установка значения i-й компоненты первой производной
 void set(double a,int i);
// получить значение
 double get() { return val; }
// получить значение i-й компоненты первой производной 
 double get(int i); 
 ext_value1& operator =(const ext_value1 &a);
 ext_value1& operator +=(const ext_value1 &a);
 ext_value1& operator +=(double a);
 ext_value1 operator +(const ext_value1 &a);
 ext_value1 operator +(double a);
 friend ext_value1 operator +(double b,const ext_value1 &a);
 ext_value1& operator *=(const ext_value1 &a);
 ext_value1& operator *=(double a);
 ext_value1 operator *(const ext_value1 &a);
 ext_value1 operator *(double a);
 friend ext_value1 operator *(double b,const ext_value1 &a);
 ext_value1& operator -=(const ext_value1 &a);
 ext_value1& operator -=(double a);
 ext_value1 operator -(const ext_value1 &a);
 ext_value1 operator -(double a);
 friend ext_value1 operator -(double b,const ext_value1 &a);
// унарный минус:
 ext_value1 operator -();
 ext_value1 operator /(const ext_value1 &a);
 ext_value1 operator /(double a);
 ext_value1& operator /=(const ext_value1 &a);
 ext_value1& operator /=(double a);
 friend ext_value1 operator /(double b,const ext_value1 &a);
 friend ext_value1 sin(const ext_value1 &a);
 friend ext_value1 cos(const ext_value1 &a);
 friend ext_value1 sqrt(const ext_value1 &a);
 friend ext_value1 exp(const ext_value1 &a);
 friend ext_value1 atan2(const ext_value1 &a,const ext_value1 &b);
 friend ext_value1 asin(const ext_value1 &a);
 friend ext_value1 acos(const ext_value1 &a);
};


// конструктор копирования
 ext_value3::ext_value3(const ext_value3 &ob)
 { 
  int i,t;
  val=ob.val; 
  n_var=ob.n_var;
  if(n_var)
	{
	  der_1=new double[n_var];
	  for(i=0;i<n_var;i++) der_1[i]=ob.der_1[i];
	  t=(n_var*(n_var+1))/2;
	  der_2=new double[t];
	  for(i=0;i<t;i++) der_2[i]=ob.der_2[i];
	  t=(n_var*(n_var+1))/2+((n_var-1)*(n_var*n_var-n_var+2))/4;
	  der_3=new double[t];
	  for(i=0;i<t;i++) der_3[i]=ob.der_3[i];
	}
  else	{ der_1=der_2=der_3=NULL; }
 }
 ext_value2::ext_value2(const ext_value2 &ob)
 { 
  int i,t;
  val=ob.val; 
  n_var=ob.n_var;
  if(n_var)
	{
	  der_1=new double[n_var];
	  for(i=0;i<n_var;i++) der_1[i]=ob.der_1[i];
	  t=(n_var*(n_var+1))/2;
	  der_2=new double[t];
	  for(i=0;i<t;i++) der_2[i]=ob.der_2[i];
	}
  else	{ der_1=der_2=NULL; }
 }
 ext_value1::ext_value1(const ext_value1 &ob)
 { 
  int i;
  val=ob.val; 
  n_var=ob.n_var;
  if(n_var)
	{
	  der_1=new double[n_var];
	  for(i=0;i<n_var;i++) der_1[i]=ob.der_1[i];
	}
  else	{ der_1=NULL; }
 }


// конструктор для создания независимой переменной или константы
 ext_value3::ext_value3(int n, double a, int m)
 {
  int i,t;
  n_var=n;
  der_1=new double[n_var];
  for(i=0;i<n_var;i++) der_1[i]=0;
  t=(n_var*(n_var+1))/2;
  der_2=new double[t];
  for(i=0;i<t;i++) der_2[i]=0;
  t=(n_var*(n_var+1))/2+((n_var-1)*(n_var*n_var-n_var+2))/4;
  der_3=new double[t];
  for(i=0;i<t;i++) der_3[i]=0;
  val=a; // задание значения переменной (по-умолчанию 0)
  if(m>=0 && m<n_var) der_1[m]=1;   // является ли переменная независимой (иначе константой)
 }
 ext_value2::ext_value2(int n, double a, int m)
 {
  int i,t;
  n_var=n;
  der_1=new double[n_var];
  for(i=0;i<n_var;i++) der_1[i]=0;
  t=(n_var*(n_var+1))/2;
  der_2=new double[t];
  for(i=0;i<t;i++) der_2[i]=0;
  val=a; // задание значения переменной (по-умолчанию 0)
  if(m>=0 && m<n_var) der_1[m]=1;   // является ли переменная независимой (иначе константой)
 }
 ext_value1::ext_value1(int n, double a, int m)
 {
  int i;
  n_var=n;
  der_1=new double[n_var];
  for(i=0;i<n_var;i++) der_1[i]=0;
  val=a; // задание значения переменной (по-умолчанию 0)
  if(m>=0 && m<n_var) der_1[m]=1;   // является ли переменная независимой (иначе константой)
 }


// установка размерности
 void ext_value3::set_dim(int n)
{ int i,t;
  if(n<0) throw n; // исключительная ситуация ошибка размерности
  else if (n==0) {
	 if(der_1) delete [] der_1; 
	 if(der_2) delete [] der_2;
	 if(der_3) delete [] der_3;
	 n_var=0;
	}
  else if(n != n_var)
	{ 
	 if(der_1) delete [] der_1; 
	 if(der_2) delete [] der_2;
	 if(der_3) delete [] der_3;
         n_var=n;
	 der_1=new double[n_var];
	 for(i=0;i<n_var;i++) der_1[i]=0;
	 t=(n_var*(n_var+1))/2;
	 der_2=new double[t];
	 for(i=0;i<t;i++) der_2[i]=0;
	 t=(n_var*(n_var+1))/2+((n_var-1)*(n_var*n_var-n_var+2))/4;
	 der_3=new double[t];
	 for(i=0;i<t;i++) der_3[i]=0;
	}
// иначе ничего не происходит : может быть следовало бы обнулить все значения ?
}
 void ext_value2::set_dim(int n)
{ int i,t;
  if(n<0) throw n; // исключительная ситуация ошибка размерности
  else if (n==0) {
	 if(der_1) delete [] der_1; 
	 if(der_2) delete [] der_2;
	 n_var=0;
	}
  else if(n != n_var)
	{ 
	 if(der_1) delete [] der_1; 
	 if(der_2) delete [] der_2;
         n_var=n;
	 der_1=new double[n_var];
	 for(i=0;i<n_var;i++) der_1[i]=0;
	 t=(n_var*(n_var+1))/2;
	 der_2=new double[t];
	 for(i=0;i<t;i++) der_2[i]=0;
	}
// иначе ничего не происходит : может быть следовало бы обнулить все значения ?
}
 void ext_value1::set_dim(int n)
{ int i;
  if(n<0) throw n; // исключительная ситуация ошибка размерности
  else if (n==0) {
	 if(der_1) delete [] der_1; 
	 n_var=0;
	}
  else if(n != n_var)
	{ 
	 if(der_1) delete [] der_1; 
         n_var=n;
	 der_1=new double[n_var];
	 for(i=0;i<n_var;i++) der_1[i]=0;
	}
// иначе ничего не происходит : может быть следовало бы обнулить все значения ?
}


// установка значения i-й компоненты первой производной
 void ext_value3::set(double a,int i)
{
  if(n_var<=0) throw n_var; // исключительная ситуация ошибка размерности
  if(i>=0 && i<n_var) der_1[i]=a; // при ошибке индекса ничего не происходит   
}
 void ext_value2::set(double a,int i)
{
  if(n_var<=0) throw n_var; // исключительная ситуация ошибка размерности
  if(i>=0 && i<n_var) der_1[i]=a; // при ошибке индекса ничего не происходит   
}
 void ext_value1::set(double a,int i)
{
  if(n_var<=0) throw n_var; // исключительная ситуация ошибка размерности
  if(i>=0 && i<n_var) der_1[i]=a; // при ошибке индекса ничего не происходит   
}


// установка значения i,j-й компоненты второй производной
 void ext_value3::set(double a,int i,int j)
{
  if(n_var<=0) throw n_var; // исключительная ситуация ошибка размерности
  if(i<0 || i>=n_var) throw i; 
  if(j<0 || j>=n_var) throw j; 
  der_2[i+(j*(j+1))/2]=a;
}
 void ext_value2::set(double a,int i,int j)
{
  if(n_var<=0) throw n_var; // исключительная ситуация ошибка размерности
  if(i<0 || i>=n_var) throw i; 
  if(j<0 || j>=n_var) throw j; 
  der_2[i+(j*(j+1))/2]=a;
}


// установка значения i,j,k-й компоненты третьей производной
 void ext_value3::set(double a,int i,int j,int k)
{
  if(n_var<=0) throw n_var; // исключительная ситуация ошибка размерности
  if(i<0 || i>=n_var) throw i; 
  if(j<0 || j>=n_var) throw j; 
  if(k<0 || k>=n_var) throw j; 
  der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4]=a; 
}


// получить значение i-й компоненты первой производной 
 double ext_value3::get(int i) 
{ 
 if(n_var<=0) {cerr<<"ext_value3 n_var throw "<<n_var<<endl; throw n_var; }
 if(i<0 || i>=n_var) {cerr<<"ext_value3 i throw "<<i<<endl; throw i; }
 return der_1[i]; 
}
 double ext_value2::get(int i) 
{ 
 if(n_var<=0) {cerr<<"ext_value3 n_var throw "<<n_var<<endl; throw n_var; }
 if(i<0 || i>=n_var) {cerr<<"ext_value3 i throw "<<i<<endl; throw i; }
 return der_1[i]; 
}
 double ext_value1::get(int i) 
{ 
 if(n_var<=0) {cerr<<"ext_value3 n_var throw "<<n_var<<endl; throw n_var; }
 if(i<0 || i>=n_var) {cerr<<"ext_value3 i throw "<<i<<endl; throw i; }
 return der_1[i]; 
}


// получить значение i,j-й компоненты второй производной
 double ext_value3::get(int i,int j)
  {
 if(n_var<=0) throw n_var; 
 if(i<0 || i>=n_var) throw i;
 if(j<0 || j>=n_var) throw j;
  int t;
  if(i>j) { t=i; i=j; j=t;}  
  return der_2[i+(j*(j+1))/2];
 }
 double ext_value2::get(int i,int j)
  {
 if(n_var<=0) throw n_var; 
 if(i<0 || i>=n_var) throw i;
 if(j<0 || j>=n_var) throw j;
  int t;
  if(i>j) { t=i; i=j; j=t;}  
  return der_2[i+(j*(j+1))/2];
 }


// получить значение i,j,k-й компоненты третьей производной
 double ext_value3::get(int i,int j,int k)
 {
 if(n_var<=0) throw n_var; 
 if(i<0 || i>=n_var) throw i;
 if(j<0 || j>=n_var) throw j;
 if(k<0 || k>=n_var) throw j;
  int t;
  if(i>j) { t=i; i=j; j=t;} 
  if(i>k) { t=i; i=k; k=t;}
  if(j>k) { t=j; j=k; k=t;}
  return der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4];
 }


// оператор =
 ext_value3& ext_value3::operator =(const ext_value3 &a)
 {
  int i,t;
  val=a.val;
  if (a.n_var==0) {
	 if(der_1) delete [] der_1; 
	 if(der_2) delete [] der_2;
	 if(der_3) delete [] der_3;
	 n_var=0;
	 return *this;
	}
  else if(n_var != a.n_var)
	{
	 if(der_1) delete [] der_1; 
	 if(der_2) delete [] der_2;
	 if(der_3) delete [] der_3;
	 n_var=a.n_var;
	 t=(n_var*(n_var+1))/2;
	 der_1=new double[n_var];
	 t=(n_var*(n_var+1))/2;
	 der_2=new double[t];
	 t=(n_var*(n_var+1))/2+((n_var-1)*(n_var*n_var-n_var+2))/4;
	 der_3=new double[t];
	}
  for(i=0;i<n_var;i++) der_1[i]=a.der_1[i];
  t=(n_var*(n_var+1))/2;
  for(i=0;i<t;i++) der_2[i]=a.der_2[i];
  t=(n_var*(n_var+1))/2+((n_var-1)*(n_var*n_var-n_var+2))/4;
  for(i=0;i<t;i++) der_3[i]=a.der_3[i];
  return *this;
 }
 ext_value2& ext_value2::operator =(const ext_value2 &a)
 {
  int i,t=(n_var*(n_var+1))/2;
  val=a.val;
  if (a.n_var==0) {
	 if(der_1) delete [] der_1; 
	 if(der_2) delete [] der_2;
	 n_var=0;
	 return *this;
	}
  else if(n_var != a.n_var)
	{
	 if(der_1) delete [] der_1; 
	 if(der_2) delete [] der_2;
	 n_var=a.n_var;
	 t=(n_var*(n_var+1))/2;
	 der_1=new double[n_var];
	 der_2=new double[t];
	}
  for(i=0;i<n_var;i++) der_1[i]=a.der_1[i];
  for(i=0;i<t;i++) der_2[i]=a.der_2[i];
  return *this;
 }
 ext_value1& ext_value1::operator =(const ext_value1 &a)
 {
  int i;
  val=a.val;
  if (a.n_var==0) {
	 if(der_1) delete [] der_1; 
	 n_var=0;
	 return *this;
	}
  else if(n_var != a.n_var)
	{
	 if(der_1) delete [] der_1; 
	 n_var=a.n_var;
	 der_1=new double[n_var];
	}
  for(i=0;i<n_var;i++) der_1[i]=a.der_1[i];
  return *this;
 }


// оператор +=
// допустимо сложение с однотипным или с double (константой)
 ext_value3& ext_value3::operator +=(const ext_value3 &a)
 {
  int i,t=(n_var*(n_var+1))/2+((n_var-1)*(n_var*n_var-n_var+2))/4;
  if(n_var!=a.n_var) throw n_var-a.n_var;

  val+=a.val;
  for(i=0;i<n_var;i++)
      der_1[i]+=a.der_1[i];
  for(i=0;i<(n_var*(n_var+1))/2;i++)der_2[i]+=a.der_2[i];
  for(i=0;i<t;i++) der_3[i]+=a.der_3[i];
  return *this;
 }
 ext_value2& ext_value2::operator +=(const ext_value2 &a)
 {
  int i;
  if(n_var!=a.n_var) throw n_var-a.n_var;

  val+=a.val;
  for(i=0;i<n_var;i++)
      der_1[i]+=a.der_1[i];
  for(i=0;i<(n_var*(n_var+1))/2;i++)der_2[i]+=a.der_2[i];
  return *this;
 }
 ext_value1& ext_value1::operator +=(const ext_value1 &a)
 {
  int i;
  if(n_var!=a.n_var) throw n_var-a.n_var;

  val+=a.val;
  for(i=0;i<n_var;i++)
      der_1[i]+=a.der_1[i];
  return *this;
 }
 ext_value3& ext_value3::operator +=(double a)
 {
  val+=a;
  return *this;
 }
 ext_value2& ext_value2::operator +=(double a)
 {
  val+=a;
  return *this;
 }
 ext_value1& ext_value1::operator +=(double a)
 {
  val+=a;
  return *this;
 }


// оператор +
// допустимо сложение с однотипным или с double (константой)
 ext_value3 ext_value3::operator +(const ext_value3 &a)
 {
  int i,t=(n_var*(n_var+1))/2+((n_var-1)*(n_var*n_var-n_var+2))/4;
  ext_value3 tmp(n_var);
  if(n_var!=a.n_var) throw n_var-a.n_var;
  
  tmp.val=val+a.val;
  
  for(i=0;i<n_var;i++)
      tmp.der_1[i]=der_1[i]+a.der_1[i];
	  
  for(i=0;i<(n_var*(n_var+1))/2;i++)
      tmp.der_2[i]=der_2[i]+a.der_2[i];

  for(i=0;i<t;i++)
      tmp.der_3[i]=der_3[i]+a.der_3[i];
  
  return tmp;
 }
 ext_value2 ext_value2::operator +(const ext_value2 &a)
 {
  int i;
  ext_value2 tmp(n_var);
  if(n_var!=a.n_var) throw n_var-a.n_var;
  
  tmp.val=val+a.val;
  for(i=0;i<n_var;i++) tmp.der_1[i]=der_1[i]+a.der_1[i];
  for(i=0;i<(n_var*(n_var+1))/2;i++) tmp.der_2[i]=der_2[i]+a.der_2[i];
  
  return tmp;
 }
 ext_value1 ext_value1::operator +(const ext_value1 &a)
 {
  int i;
  ext_value1 tmp(n_var);
  if(n_var!=a.n_var) throw n_var-a.n_var;
  
  tmp.val=val+a.val;
  for(i=0;i<n_var;i++) tmp.der_1[i]=der_1[i]+a.der_1[i];
  
  return tmp;
 }
 ext_value3 ext_value3::operator +(double a)
 {
  ext_value3 tmp=*this;
  tmp.val=val+a;
  return tmp;
 }
 ext_value2 ext_value2::operator +(double a)
 {
  ext_value2 tmp=*this;
  tmp.val=val+a;
  return tmp;
 }
 ext_value1 ext_value1::operator +(double a)
 {
  ext_value1 tmp=*this;
  tmp.val=val+a;
  return tmp;
 }
 ext_value3 operator +(double b,const ext_value3 &a)
 {
  ext_value3 tmp=a;
  tmp.val=a.val+b;
  return tmp;
 }
 ext_value2 operator +(double b,const ext_value2 &a)
 {
  ext_value2 tmp=a;
  tmp.val=a.val+b;
  return tmp;
 }
 ext_value1 operator +(double b,const ext_value1 &a)
 {
  ext_value1 tmp=a;
  tmp.val=a.val+b;
  return tmp;
 }


// оператор *=
// допустимо умножение на однотипное или double (константу)
 ext_value3& ext_value3::operator *=(const ext_value3 &a)
 {
  int i,j,k;
  if(n_var!=a.n_var) throw n_var-a.n_var;
  for(i=0;i<n_var;i++)
      for(j=i;j<n_var;j++)
          for(k=j;k<n_var;k++) 
              der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4]=
              der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4]*a.val+
              der_2[i+(j*(j+1))/2]*a.der_1[k]+
              der_2[i+(k*(k+1))/2]*a.der_1[j]+
              der_1[i]*a.der_2[j+(k*(k+1))/2]+
              der_2[j+(k*(k+1))/2]*a.der_1[i]+
              der_1[j]*a.der_2[i+(k*(k+1))/2]+
              der_1[k]*a.der_2[i+(j*(j+1))/2]+
              val*a.der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4];
  for(i=0;i<n_var;i++)
      for(j=i;j<n_var;j++)
              der_2[i+(j*(j+1))/2]=der_2[i+(j*(j+1))/2]*a.val+
                                   der_1[i]*a.der_1[j]+
                                   der_1[j]*a.der_1[i]+
                                   val*a.der_2[i+(j*(j+1))/2];
  for(i=0;i<n_var;i++)
      der_1[i]=der_1[i]*a.val+val*a.der_1[i];
  val=val*a.val;
  return *this;
 }
 ext_value2& ext_value2::operator *=(const ext_value2 &a)
 {
  int i,j;
  if(n_var!=a.n_var) throw n_var-a.n_var;

  for(i=0;i<n_var;i++)
      for(j=i;j<n_var;j++)
              der_2[i+(j*(j+1))/2]=der_2[i+(j*(j+1))/2]*a.val+
                                   der_1[i]*a.der_1[j]+
                                   der_1[j]*a.der_1[i]+
                                   val*a.der_2[i+(j*(j+1))/2];
  for(i=0;i<n_var;i++)
      der_1[i]=der_1[i]*a.val+val*a.der_1[i];
  val=val*a.val;
  return *this;
 }
 ext_value1& ext_value1::operator *=(const ext_value1 &a)
 {
  int i;
  if(n_var!=a.n_var) throw n_var-a.n_var;
  for(i=0;i<n_var;i++) der_1[i]=der_1[i]*a.val+val*a.der_1[i];
  val=val*a.val;
  return *this;
 }
 ext_value3& ext_value3::operator *=(double a)
 {
  int i,t=(n_var*(n_var+1))/2+((n_var-1)*(n_var*n_var-n_var+2))/4;
  val*=a;
  for(i=0;i<n_var;i++) der_1[i]*=a;
  for(i=0;i<(n_var*(n_var+1))/2;i++) der_2[i]*=a;
  for(i=0;i<t;i++) der_3[i]*=a;
  return *this;
 }
 ext_value2& ext_value2::operator *=(double a)
 {
  int i;
  val*=a;
  for(i=0;i<n_var;i++) der_1[i]*=a;
  for(i=0;i<(n_var*(n_var+1))/2;i++) der_2[i]*=a;
  return *this;
 }
 ext_value1& ext_value1::operator *=(double a)
 {
  int i;
  val*=a;
  for(i=0;i<n_var;i++) der_1[i]*=a;
  return *this;
 }


// оператор *
// допустимо умножение на однотипное или double (константу)
 ext_value3 ext_value3::operator*(const ext_value3 &a)
 {
  int i,j,k;
  ext_value3 tmp(n_var);
  if(n_var!=a.n_var) throw n_var-a.n_var;
  tmp.val=val*a.val;
  for(i=0;i<n_var;i++)
      tmp.der_1[i]=der_1[i]*a.val+val*a.der_1[i];
  for(i=0;i<n_var;i++)
      for(j=i;j<n_var;j++)
          tmp.der_2[i+(j*(j+1))/2]=der_2[i+(j*(j+1))/2]*a.val+
                                   der_1[i]*a.der_1[j]+
                                   der_1[j]*a.der_1[i]+
                                   val*a.der_2[i+(j*(j+1))/2];
  for(i=0;i<n_var;i++)
      for(j=i;j<n_var;j++)
          for(k=j;k<n_var;k++) 
              tmp.der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4]=
              der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4]*a.val+
              der_2[i+(j*(j+1))/2]*a.der_1[k]+
              der_2[i+(k*(k+1))/2]*a.der_1[j]+
              der_1[i]*a.der_2[j+(k*(k+1))/2]+
              der_2[j+(k*(k+1))/2]*a.der_1[i]+
              der_1[j]*a.der_2[i+(k*(k+1))/2]+
              der_1[k]*a.der_2[i+(j*(j+1))/2]+
              val*a.der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4];
  return tmp;
 }
 ext_value2 ext_value2::operator*(const ext_value2 &a)
 {
  int i,j;
  ext_value2 tmp(n_var);
  if(n_var!=a.n_var) throw n_var-a.n_var;
  tmp.val=val*a.val;
  for(i=0;i<n_var;i++)
      tmp.der_1[i]=der_1[i]*a.val+val*a.der_1[i];
  for(i=0;i<n_var;i++)
      for(j=i;j<n_var;j++)
          tmp.der_2[i+(j*(j+1))/2]=der_2[i+(j*(j+1))/2]*a.val+
                                   der_1[i]*a.der_1[j]+
                                   der_1[j]*a.der_1[i]+
                                   val*a.der_2[i+(j*(j+1))/2];
  return tmp;
 }
 ext_value1 ext_value1::operator*(const ext_value1 &a)
 {
  int i;
  ext_value1 tmp(n_var);
  if(n_var!=a.n_var) throw n_var-a.n_var;
  tmp.val=val*a.val;
  for(i=0;i<n_var;i++)
      tmp.der_1[i]=der_1[i]*a.val+val*a.der_1[i];
  return tmp;
 }
 ext_value3 ext_value3::operator *(double a)
 {
  int i,t=(n_var*(n_var+1))/2+((n_var-1)*(n_var*n_var-n_var+2))/4;
  ext_value3 tmp=*this;
  tmp.val*=a;
  for(i=0;i<n_var;i++) tmp.der_1[i]*=a;
  for(i=0;i<(n_var*(n_var+1))/2;i++) tmp.der_2[i]*=a;
  for(i=0;i<t;i++) tmp.der_3[i]*=a;
  return tmp;
 }
 ext_value2 ext_value2::operator *(double a)
 {
  int i;
  ext_value2 tmp=*this;
  tmp.val*=a;
  for(i=0;i<n_var;i++) tmp.der_1[i]*=a;
  for(i=0;i<(n_var*(n_var+1))/2;i++) tmp.der_2[i]*=a;
  return tmp;
 }
 ext_value1 ext_value1::operator *(double a)
 {
  int i;
  ext_value1 tmp=*this;
  tmp.val*=a;
  for(i=0;i<n_var;i++) tmp.der_1[i]*=a;
  return tmp;
 }
 ext_value3 operator *(double b,const ext_value3 &a)
 {
  int i,t=(a.n_var*(a.n_var+1))/2+((a.n_var-1)*(a.n_var*a.n_var-a.n_var+2))/4;
  ext_value3 tmp=a;
  tmp.val*=b;
  for(i=0;i<a.n_var;i++) tmp.der_1[i]*=b;
  for(i=0;i<(a.n_var*(a.n_var+1))/2;i++) tmp.der_2[i]*=b;
  for(i=0;i<t;i++) tmp.der_3[i]*=b;
  return tmp;
 }
 ext_value2 operator *(double b,const ext_value2 &a)
 {
  int i;
  ext_value2 tmp=a;
  tmp.val*=b;
  for(i=0;i<a.n_var;i++) tmp.der_1[i]*=b;
  for(i=0;i<(a.n_var*(a.n_var+1))/2;i++) tmp.der_2[i]*=b;
  return tmp;
 }
 ext_value1 operator *(double b,const ext_value1 &a)
 {
  int i;
  ext_value1 tmp=a;
  tmp.val*=b;
  for(i=0;i<a.n_var;i++) tmp.der_1[i]*=b;
  return tmp;
 }

// оператор -=
// допустимо вычитание однотипного или double (константы)
 ext_value3& ext_value3::operator -=(const ext_value3 &a)
 {
  int i,t=(n_var*(n_var+1))/2+((n_var-1)*(n_var*n_var-n_var+2))/4;
  if(n_var!=a.n_var) throw n_var-a.n_var;

  val-=a.val;
  for(i=0;i<n_var;i++)
      der_1[i]-=a.der_1[i];
  for(i=0;i<(n_var*(n_var+1))/2;i++)der_2[i]-=a.der_2[i];
  for(i=0;i<t;i++) der_3[i]-=a.der_3[i];
  return *this;
 }
 ext_value2& ext_value2::operator -=(const ext_value2 &a)
 {
  int i;
  if(n_var!=a.n_var) throw n_var-a.n_var;

  val-=a.val;
  for(i=0;i<n_var;i++)
      der_1[i]-=a.der_1[i];
  for(i=0;i<(n_var*(n_var+1))/2;i++)der_2[i]-=a.der_2[i];
  return *this;
 }
 ext_value1& ext_value1::operator -=(const ext_value1 &a)
 {
  int i;
  if(n_var!=a.n_var) throw n_var-a.n_var;

  val-=a.val;
  for(i=0;i<n_var;i++)
      der_1[i]-=a.der_1[i];
  return *this;
 }
 ext_value3& ext_value3::operator -=(double a)
 {
  val-=a;
  return *this;
 }
 ext_value2& ext_value2::operator -=(double a)
 {
  val-=a;
  return *this;
 }
 ext_value1& ext_value1::operator -=(double a)
 {
  val-=a;
  return *this;
 }


// оператор -
// допустимо вычитание с однотипным или с double (константой)
 ext_value3 ext_value3::operator -(const ext_value3 &a)
 {
  int i,t=(n_var*(n_var+1))/2+((n_var-1)*(n_var*n_var-n_var+2))/4;
  ext_value3 tmp(n_var);
  if(n_var!=a.n_var) throw n_var-a.n_var;
  tmp.val=val-a.val;
  for(i=0;i<n_var;i++) tmp.der_1[i]=der_1[i]-a.der_1[i];
  for(i=0;i<(n_var*(n_var+1))/2;i++) tmp.der_2[i]=der_2[i]-a.der_2[i];
  for(i=0;i<t;i++) tmp.der_3[i]=der_3[i]-a.der_3[i];
  return tmp;
 }
 ext_value2 ext_value2::operator -(const ext_value2 &a)
 {
  int i;
  ext_value2 tmp(n_var);
  if(n_var!=a.n_var) throw n_var-a.n_var;
  tmp.val=val-a.val;
  for(i=0;i<n_var;i++) tmp.der_1[i]=der_1[i]-a.der_1[i];
  for(i=0;i<(n_var*(n_var+1))/2;i++) tmp.der_2[i]=der_2[i]-a.der_2[i];
  return tmp;
 }
 ext_value1 ext_value1::operator -(const ext_value1 &a)
 {
  int i;
  ext_value1 tmp(n_var);
  if(n_var!=a.n_var) throw n_var-a.n_var;
  tmp.val=val-a.val;
  for(i=0;i<n_var;i++) tmp.der_1[i]=der_1[i]-a.der_1[i];
  return tmp;
 }
 ext_value3 ext_value3::operator -(double a)
 {
  ext_value3 tmp=*this;
  tmp.val=val-a;
  return tmp;
 }
 ext_value2 ext_value2::operator -(double a)
 {
  ext_value2 tmp=*this;
  tmp.val=val-a;
  return tmp;
 }
 ext_value1 ext_value1::operator -(double a)
 {
  ext_value1 tmp=*this;
  tmp.val=val-a;
  return tmp;
 }
 ext_value3 operator -(double b,const ext_value3 &a)
 {
  int i,n,t;
  ext_value3 tmp(a.n_var);
  n=a.n_var;
  t=(n*(n+1))/2+((n-1)*(n*n-n+2))/4;
  tmp.val=b-a.val;
  for(i=0;i<n;i++) tmp.der_1[i]=-a.der_1[i];
  for(i=0;i<(n*(n+1))/2;i++) tmp.der_2[i]=-a.der_2[i];
  for(i=0;i<t;i++) tmp.der_3[i]=-a.der_3[i];
  return tmp;
 }
 ext_value2 operator -(double b,const ext_value2 &a)
 {
  int i,n;
  ext_value2 tmp(a.n_var);
  n=a.n_var;
  tmp.val=b-a.val;
  for(i=0;i<n;i++) tmp.der_1[i]=-a.der_1[i];
  for(i=0;i<(n*(n+1))/2;i++) tmp.der_2[i]=-a.der_2[i];
  return tmp;
 }
 ext_value1 operator -(double b,const ext_value1 &a)
 {
  int i,n;
  ext_value1 tmp(a.n_var);
  n=a.n_var;
  tmp.val=b-a.val;
  for(i=0;i<n;i++) tmp.der_1[i]=-a.der_1[i];
  return tmp;
 }


// унарный минус
 ext_value3 ext_value3::operator -()
 {
	 int i,t=(n_var*(n_var+1))/2+((n_var-1)*(n_var*n_var-n_var+2))/4;
  ext_value3 tmp=*this;
  tmp.val=-val;
  for(i=0;i<n_var;i++) tmp.der_1[i]=-der_1[i];
  for(i=0;i<(n_var*(n_var+1))/2;i++) tmp.der_2[i]=-der_2[i];
  for(i=0;i<t;i++) tmp.der_3[i]=-der_3[i];
  return tmp;
 }
 ext_value2 ext_value2::operator -()
 {
	 int i;
  ext_value2 tmp=*this;
  tmp.val=-val;
  for(i=0;i<n_var;i++) tmp.der_1[i]=-der_1[i];
  for(i=0;i<(n_var*(n_var+1))/2;i++) tmp.der_2[i]=-der_2[i];
  return tmp;
 }
 ext_value1 ext_value1::operator -()
 {
	 int i;
  ext_value1 tmp=*this;
  tmp.val=-val;
  for(i=0;i<n_var;i++) tmp.der_1[i]=-der_1[i];
  return tmp;
 }

// оператор /
 ext_value3 ext_value3::operator/(const ext_value3 &a)
 {
  int i,j,k;
  double w2,w3,w4;
  ext_value3 tmp(n_var);
  w2=a.val*a.val;
  w3=w2*a.val;
  w4=w2*w2;
  if(n_var!=a.n_var) throw n_var-a.n_var;
  tmp.val=val/a.val;
  for(i=0;i<n_var;i++)
      tmp.der_1[i]=der_1[i]/a.val-val*a.der_1[i]/w2;
  for(i=0;i<n_var;i++)
      for(j=i;j<n_var;j++)
	tmp.der_2[i+(j*(j+1))/2]=
	der_2[i+(j*(j+1))/2]/a.val-
	(der_1[i]*a.der_1[j]+
	 der_1[j]*a.der_1[i]+
	 val*a.der_2[i+(j*(j+1))/2]
	)/w2+
	2*val*a.der_1[i]*a.der_1[j]/w3;
  for(i=0;i<n_var;i++)
      for(j=i;j<n_var;j++)
          for(k=j;k<n_var;k++) 
tmp.der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4]=
der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4]/a.val
-
(
 der_2[i+(j*(j+1))/2]*a.der_1[k]+
 der_2[i+(k*(k+1))/2]*a.der_1[j]+
 der_2[j+(k*(k+1))/2]*a.der_1[i]+
 a.der_2[i+(j*(j+1))/2]*der_1[k]+
 a.der_2[i+(k*(k+1))/2]*der_1[j]+
 a.der_2[j+(k*(k+1))/2]*der_1[i]+
 val*a.der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4]
)/w2
+
(
 2*der_1[i]*a.der_1[j]*a.der_1[k]+
 2*a.der_1[i]*der_1[j]*a.der_1[k]+
 2*a.der_1[i]*a.der_1[j]*der_1[k]+
 2*val*a.der_2[i+(j*(j+1))/2]*a.der_1[k]+
 2*val*a.der_2[j+(k*(k+1))/2]*a.der_1[i]+
 2*val*a.der_2[i+(k*(k+1))/2]*a.der_1[j]
)/w3
-6*val*a.der_1[i]*a.der_1[j]*a.der_1[k]/w4;
  return tmp;
 }
 ext_value2 ext_value2::operator/(const ext_value2 &a)
 {
  int i,j;
  double w2,w3;
  ext_value2 tmp(n_var);
  w2=a.val*a.val;
  w3=w2*a.val;
  if(n_var!=a.n_var) throw n_var-a.n_var;
  tmp.val=val/a.val;
  for(i=0;i<n_var;i++)
      tmp.der_1[i]=der_1[i]/a.val-val*a.der_1[i]/w2;
  for(i=0;i<n_var;i++)
      for(j=i;j<n_var;j++)
          tmp.der_2[i+(j*(j+1))/2]=der_2[i+(j*(j+1))/2]/a.val-
  (der_1[i]*a.der_1[j]+der_1[j]*a.der_1[i]+val*a.der_2[i+(j*(j+1))/2])/w2+
  2*val*a.der_1[i]*a.der_1[j]/w3;
  return tmp;
 }
 ext_value1 ext_value1::operator/(const ext_value1 &a)
 {
  int i;
  ext_value1 tmp(n_var);
  if(n_var!=a.n_var) throw n_var-a.n_var;
  tmp.val=val/a.val;
  for(i=0;i<n_var;i++)
      tmp.der_1[i]=der_1[i]/a.val-val*a.der_1[i]/(a.val*a.val);
  return tmp;
 }
 ext_value3 operator /(double b,const ext_value3 &a)
 {
  int i,j,k,n=a.n_var;
  double w=a.val,w2=a.val*a.val;
  ext_value3 tmp(a.n_var);
  tmp.val=b/a.val;
  for(i=0;i<n;i++)
      tmp.der_1[i]=-b*a.der_1[i]/w2;
  for(i=0;i<n;i++)
      for(j=i;j<n;j++)
	tmp.der_2[i+(j*(j+1))/2]=
	(2*a.der_1[i]*a.der_1[j]/w-a.der_2[i+(j*(j+1))/2])*b/w2;
  for(i=0;i<n;i++)
      for(j=i;j<n;j++)
          for(k=j;k<n;k++) 
		tmp.der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4]=
		(
		- a.der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4]
		+2*
			(
			 a.der_2[i+(j*(j+1))/2]*a.der_1[k]+
			 a.der_2[j+(k*(k+1))/2]*a.der_1[i]+
			 a.der_2[i+(k*(k+1))/2]*a.der_1[j]
			)/w
		-6*a.der_1[i]*a.der_1[j]*a.der_1[k]/w2
		)*b/w2;
  return tmp;
}
 ext_value2 operator /(double b,const ext_value2 &a)
 {
  int i,j,n=a.n_var;
  double w=a.val,w2=a.val*a.val;
  ext_value2 tmp(a.n_var);
  tmp.val=b/a.val;
  for(i=0;i<n;i++)
      tmp.der_1[i]=-b*a.der_1[i]/w2;
  for(i=0;i<n;i++)
      for(j=i;j<n;j++)
	tmp.der_2[i+(j*(j+1))/2]=
	(2*a.der_1[i]*a.der_1[j]/w-a.der_2[i+(j*(j+1))/2])*b/w2;
  return tmp;
}
 ext_value1 operator /(double b,const ext_value1 &a)
 {
  int i,n=a.n_var;
  double w2=a.val*a.val;
  ext_value1 tmp(a.n_var);
  tmp.val=b/a.val;
  for(i=0;i<n;i++)
      tmp.der_1[i]=-b*a.der_1[i]/w2;
  return tmp;
}
 ext_value3 ext_value3::operator /(double a)
 {
  int i,t=(n_var*(n_var+1))/2+((n_var-1)*(n_var*n_var-n_var+2))/4;
  ext_value3 tmp=*this;
  tmp.val/=a;
  for(i=0;i<n_var;i++) tmp.der_1[i]/=a;
  for(i=0;i<(n_var*(n_var+1))/2;i++) tmp.der_2[i]/=a;
  for(i=0;i<t;i++) tmp.der_3[i]/=a;
  return tmp;
 }
 ext_value2 ext_value2::operator /(double a)
 {
  int i;
  ext_value2 tmp=*this;
  tmp.val/=a;
  for(i=0;i<n_var;i++) tmp.der_1[i]/=a;
  for(i=0;i<(n_var*(n_var+1))/2;i++) tmp.der_2[i]/=a;
  return tmp;
 }
 ext_value3& ext_value3::operator /=(const ext_value3 &a)
 {
  int i,j,k;
  double w2,w3,w4;
  w2=a.val*a.val;
  w3=w2*a.val;
  w4=w2*w2;
  if(n_var!=a.n_var) throw n_var-a.n_var;
  for(i=0;i<n_var;i++)
      for(j=i;j<n_var;j++)
          for(k=j;k<n_var;k++) 
der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4]=
der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4]/a.val
-
(
 der_2[i+(j*(j+1))/2]*a.der_1[k]+
 der_2[i+(k*(k+1))/2]*a.der_1[j]+
 der_2[j+(k*(k+1))/2]*a.der_1[i]+
 a.der_2[i+(j*(j+1))/2]*der_1[k]+
 a.der_2[i+(k*(k+1))/2]*der_1[j]+
 a.der_2[j+(k*(k+1))/2]*der_1[i]+
 val*a.der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4]
)/w2
+
(
 2*der_1[i]*a.der_1[j]*a.der_1[k]+
 2*a.der_1[i]*der_1[j]*a.der_1[k]+
 2*a.der_1[i]*a.der_1[j]*der_1[k]+
 2*val*a.der_2[i+(j*(j+1))/2]*a.der_1[k]+
 2*val*a.der_2[j+(k*(k+1))/2]*a.der_1[i]+
 2*val*a.der_2[i+(k*(k+1))/2]*a.der_1[j]
)/w3
-6*val*a.der_1[i]*a.der_1[j]*a.der_1[k]/w4;
  for(i=0;i<n_var;i++)
      for(j=i;j<n_var;j++)
	der_2[i+(j*(j+1))/2]=
	der_2[i+(j*(j+1))/2]/a.val-
	(der_1[i]*a.der_1[j]+
	 der_1[j]*a.der_1[i]+
	 val*a.der_2[i+(j*(j+1))/2]
	)/w2+
	2*val*a.der_1[i]*a.der_1[j]/w3;
  for(i=0;i<n_var;i++)
      der_1[i]=der_1[i]/a.val-val*a.der_1[i]/w2;
  val/=a.val;
  return *this;
 }
 ext_value2& ext_value2::operator /=(const ext_value2 &a)
 {
  int i,j;
  double w2,w3;
  w2=a.val*a.val;
  w3=w2*a.val;
  if(n_var!=a.n_var) throw n_var-a.n_var;
  for(i=0;i<n_var;i++)
      for(j=i;j<n_var;j++)
	der_2[i+(j*(j+1))/2]=
	der_2[i+(j*(j+1))/2]/a.val-
	(der_1[i]*a.der_1[j]+
	 der_1[j]*a.der_1[i]+
	 val*a.der_2[i+(j*(j+1))/2]
	)/w2+
	2*val*a.der_1[i]*a.der_1[j]/w3;
  for(i=0;i<n_var;i++)
      der_1[i]=der_1[i]/a.val-val*a.der_1[i]/w2;
  val/=a.val;
  return *this;
 }
 ext_value1& ext_value1::operator /=(const ext_value1 &a)
 {
  int i;
  double w2=a.val*a.val;
  if(n_var!=a.n_var) throw n_var-a.n_var;
  for(i=0;i<n_var;i++)
      der_1[i]=der_1[i]/a.val-val*a.der_1[i]/w2;
  val/=a.val;
  return *this;
 }
 ext_value1 ext_value1::operator /(double a)
 {
  int i;
  ext_value1 tmp=*this;
  tmp.val/=a;
  for(i=0;i<n_var;i++) tmp.der_1[i]/=a;
  return tmp;
 }
 ext_value3& ext_value3::operator /=(double a)
 {
  int i,t=(n_var*(n_var+1))/2+((n_var-1)*(n_var*n_var-n_var+2))/4;
  val/=a;
  for(i=0;i<n_var;i++) der_1[i]/=a;
  for(i=0;i<(n_var*(n_var+1))/2;i++) der_2[i]/=a;
  for(i=0;i<t;i++) der_3[i]/=a;
  return *this;
 }
 ext_value2& ext_value2::operator /=(double a)
 {
  int i;
  val/=a;
  for(i=0;i<n_var;i++) der_1[i]/=a;
  for(i=0;i<(n_var*(n_var+1))/2;i++) der_2[i]/=a;
  return *this;
 }
 ext_value1& ext_value1::operator /=(double a)
 {
  int i;
  val/=a;
  for(i=0;i<n_var;i++) der_1[i]/=a;
  return *this;
 }


// функция sin
 ext_value3 sin(const ext_value3 &a)
 {
  int i,j,k,n=a.n_var;
  ext_value3 tmp(a.n_var);
  
  tmp.val=sin(a.val);
  
  for(i=0;i<n;i++)
      tmp.der_1[i]=a.der_1[i]*cos(a.val);
  
  for(i=0;i<n;i++)
      for(j=i;j<n;j++)
          tmp.der_2[i+(j*(j+1))/2]=a.der_2[i+(j*(j+1))/2]*cos(a.val)-
                                   a.der_1[i]*a.der_1[j]*sin(a.val);
		  
  for(i=0;i<n;i++)
      for(j=i;j<n;j++)
          for(k=j;k<n;k++) 
              tmp.der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4]=
              a.der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4]*cos(a.val)-
              a.der_1[i]*a.der_1[j]*a.der_1[k]*cos(a.val)-
              sin(a.val)*(a.der_2[j+(k*(k+1))/2]*a.der_1[i]+
              a.der_2[i+(k*(k+1))/2]*a.der_1[j]+
              a.der_2[i+(j*(j+1))/2]*a.der_1[k]);

  return tmp;
 }
 ext_value2 sin(const ext_value2 &a)
 {
  int i,j,n=a.n_var;
  ext_value2 tmp(a.n_var);
  
  tmp.val=sin(a.val);
  for(i=0;i<n;i++)
      tmp.der_1[i]=a.der_1[i]*cos(a.val);
  for(i=0;i<n;i++)
      for(j=i;j<n;j++)
          tmp.der_2[i+(j*(j+1))/2]=a.der_2[i+(j*(j+1))/2]*cos(a.val)-
                                   a.der_1[i]*a.der_1[j]*sin(a.val);
  return tmp;
 }
 ext_value1 sin(const ext_value1 &a)
 {
  int i,n=a.n_var;
  ext_value1 tmp(a.n_var);
  
  tmp.val=sin(a.val);
  for(i=0;i<n;i++)
      tmp.der_1[i]=a.der_1[i]*cos(a.val);
  return tmp;
 }


// функция cos
 ext_value3 cos(const ext_value3 &a)
 {
  int i,j,k,n=a.n_var;
  ext_value3 tmp(a.n_var);
  
  tmp.val=cos(a.val);
  
  for(i=0;i<n;i++)
      tmp.der_1[i]=-a.der_1[i]*sin(a.val);

  for(i=0;i<n;i++)
      for(j=i;j<n;j++)
          tmp.der_2[i+(j*(j+1))/2]=-a.der_1[i]*a.der_1[j]*cos(a.val)-
                                   a.der_2[i+(j*(j+1))/2]*sin(a.val);
		
  for(i=0;i<n;i++)
      for(j=i;j<n;j++)
          for(k=j;k<n;k++) 
              tmp.der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4]=
              a.der_1[i]*a.der_1[j]*a.der_1[k]*sin(a.val)-
              cos(a.val)*(a.der_2[j+(k*(k+1))/2]*a.der_1[i]+
              a.der_2[i+(k*(k+1))/2]*a.der_1[j]+
              a.der_2[i+(j*(j+1))/2]*a.der_1[k])-
              a.der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4]*sin(a.val);

  return tmp;
 }
 ext_value2 cos(const ext_value2 &a)
 {
  int i,j,n=a.n_var;
  ext_value2 tmp(a.n_var);
  
  tmp.val=cos(a.val);
  for(i=0;i<n;i++)
      tmp.der_1[i]=-a.der_1[i]*sin(a.val);
  for(i=0;i<n;i++)
      for(j=i;j<n;j++)
          tmp.der_2[i+(j*(j+1))/2]=-a.der_1[i]*a.der_1[j]*cos(a.val)-
                                   a.der_2[i+(j*(j+1))/2]*sin(a.val);
  return tmp;
 }
 ext_value1 cos(const ext_value1 &a)
 {
  int i,n=a.n_var;
  ext_value1 tmp(a.n_var);
  
  tmp.val=cos(a.val);
  for(i=0;i<n;i++)
      tmp.der_1[i]=-a.der_1[i]*sin(a.val);
  return tmp;
 }

// функция sqrt
 ext_value3 sqrt(const ext_value3 &a)
 {
  int i,j,k,n=a.n_var;
  ext_value3 tmp(a.n_var);
  double w3,w5;
// Это не самый лучший способ отслеживать проблему, но хоть что-то:
  if(a.val<=0){cerr<<"ext_value3 sqrt a.val<=0";throw a.val; }
  tmp.val=sqrt(a.val);
   w3=a.val*tmp.val; 	// это a^3/2
   w5=a.val*w3;		// это a^5/2
  for(i=0;i<n;i++)
      tmp.der_1[i]=a.der_1[i]/(2*tmp.val);
  for(i=0;i<n;i++)
      for(j=i;j<n;j++)
          tmp.der_2[i+(j*(j+1))/2]=
            a.der_2[i+(j*(j+1))/2]/(2*tmp.val)
           -a.der_1[i]*a.der_1[j]/(4*w3);
  for(i=0;i<n;i++)
      for(j=i;j<n;j++)
          for(k=j;k<n;k++) 
              tmp.der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4]
              =a.der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4]/(2*tmp.val)
              -(a.der_2[j+(k*(k+1))/2]*a.der_1[i]+
                a.der_2[i+(k*(k+1))/2]*a.der_1[j]+
                a.der_2[i+(j*(j+1))/2]*a.der_1[k])/(4*w3)
              +3*a.der_1[i]*a.der_1[j]*a.der_1[k]/(8*w5);
  return tmp;
 }
 ext_value2 sqrt(const ext_value2 &a)
 {
  int i,j,n=a.n_var;
  ext_value2 tmp(a.n_var);
  double w3;
// Это не самый лучший способ отслеживать проблему, но хоть что-то:
  if(a.val<=0){cerr<<"ext_value3 sqrt a.val<=0";throw a.val; }
  tmp.val=sqrt(a.val);
   w3=a.val*tmp.val; 	// это a^3/2
  for(i=0;i<n;i++)
      tmp.der_1[i]=a.der_1[i]/(2*tmp.val);
  for(i=0;i<n;i++)
      for(j=i;j<n;j++)
          tmp.der_2[i+(j*(j+1))/2]=
            a.der_2[i+(j*(j+1))/2]/(2*tmp.val)
           -a.der_1[i]*a.der_1[j]/(4*w3);
  return tmp;
 }
 ext_value1 sqrt(const ext_value1 &a)
 {
  int i,n=a.n_var;
  ext_value1 tmp(a.n_var);
// Это не самый лучший способ отслеживать проблему, но хоть что-то:
  if(a.val<=0){cerr<<"ext_value3 sqrt a.val<=0";throw a.val; }
  tmp.val=sqrt(a.val);
  for(i=0;i<n;i++)
      tmp.der_1[i]=a.der_1[i]/(2*tmp.val);
  return tmp;
 }

// функция exp
 ext_value3 exp(const ext_value3 &a)
 {
  int i,j,k,n=a.n_var;
  ext_value3 tmp(a.n_var);
      tmp.val=exp(a.val);
  for(i=0;i<n;i++)
      tmp.der_1[i]=a.der_1[i]*tmp.val;
  for(i=0;i<n;i++)
      for(j=i;j<n;j++)
          tmp.der_2[i+(j*(j+1))/2]=
          ( a.der_1[i]*a.der_1[j]+a.der_2[i+(j*(j+1))/2] )
          *tmp.val;
  for(i=0;i<n;i++)
      for(j=i;j<n;j++)
          for(k=j;k<n;k++) 
              tmp.der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4]=
              (a.der_1[i]*a.der_1[j]*a.der_1[k]+
               a.der_2[j+(k*(k+1))/2]*a.der_1[i]+
               a.der_2[i+(k*(k+1))/2]*a.der_1[j]+
               a.der_2[i+(j*(j+1))/2]*a.der_1[k]+
               a.der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4]
              )*tmp.val;
  return tmp;
 }
 ext_value2 exp(const ext_value2 &a)
 {
  int i,j,n=a.n_var;
  ext_value2 tmp(a.n_var);
      tmp.val=exp(a.val);
  for(i=0;i<n;i++)
      tmp.der_1[i]=a.der_1[i]*tmp.val;
  for(i=0;i<n;i++)
      for(j=i;j<n;j++)
          tmp.der_2[i+(j*(j+1))/2]=
          ( a.der_1[i]*a.der_1[j]+a.der_2[i+(j*(j+1))/2] )
          *tmp.val;
  return tmp;
 }
 ext_value1 exp(const ext_value1 &a)
 {
  int i,n=a.n_var;
  ext_value1 tmp(a.n_var);
      tmp.val=exp(a.val);
  for(i=0;i<n;i++)
      tmp.der_1[i]=a.der_1[i]*tmp.val;
  return tmp;
 }

 ext_value3 atan2(const ext_value3 &a,const ext_value3 &b)
 {
  int i,j,k,n=a.n_var;
  double r2,r4,r6;
  if(b.n_var!=a.n_var) throw b.n_var-a.n_var; // Ошибка размерности
  ext_value3 tmp(a.n_var);
      tmp.val=atan2(a.val,b.val);
      r2=a.val*a.val+b.val*b.val;
      r4=r2*r2;
      r6=r4*r2;
  for(i=0;i<n;i++)
      tmp.der_1[i]=(b.val*a.der_1[i]-a.val*b.der_1[i])/r2;

  for(i=0;i<n;i++)
      for(j=i;j<n;j++)
          tmp.der_2[i+(j*(j+1))/2]=
          ( b.der_1[j]*a.der_1[i]-b.der_1[i]*a.der_1[j]+
            b.val*a.der_2[i+(j*(j+1))/2]-
            a.val*b.der_2[i+(j*(j+1))/2] )/r2
         -2*(b.val*a.der_1[i]-a.val*b.der_1[i])
           *(a.val*a.der_1[j]+b.val*b.der_1[j])/r4;

  for(i=0;i<n;i++)
      for(j=i;j<n;j++)
          for(k=j;k<n;k++) 
              tmp.der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4]=
              ( b.der_2[j+(k*(k+1))/2]*a.der_1[i]
               -a.der_2[j+(k*(k+1))/2]*b.der_1[i]
               +a.der_2[i+(k*(k+1))/2]*b.der_1[j]
               -b.der_2[i+(k*(k+1))/2]*a.der_1[j]
               +a.der_2[i+(j*(j+1))/2]*b.der_1[k]
               -b.der_2[i+(j*(j+1))/2]*a.der_1[k]
               +b.val*a.der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4]
               -a.val*b.der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4] )/r2
             -( b.der_1[j]*a.der_1[i]-b.der_1[i]*a.der_1[j]
               +b.val*a.der_2[i+(j*(j+1))/2]
               -a.val*b.der_2[i+(j*(j+1))/2] )*2*
              ( a.val*a.der_1[k]+b.val*b.der_1[k] )/ r4
             -( b.der_1[k]*a.der_1[i]-b.der_1[i]*a.der_1[k]
               +b.val*a.der_2[i+(k*(k+1))/2]
               -a.val*b.der_2[i+(k*(k+1))/2] )*2*
              ( a.val*a.der_1[j]+b.val*b.der_1[j] )/ r4
             -2*( b.val*a.der_1[i]-a.val*b.der_1[i] )
               *( a.der_1[j]*a.der_1[k]+b.der_1[j]*b.der_1[k]
                 +a.val*a.der_2[j+(k*(k+1))/2]
                 +b.val*b.der_2[j+(k*(k+1))/2] )/r4
             +8*( b.val*a.der_1[i]-a.val*b.der_1[i] )
               *( a.val*a.der_1[j]+b.val*b.der_1[j] )
               *( a.val*a.der_1[k]+b.val*b.der_1[k] )/r6;
  return tmp;
 }
 ext_value2 atan2(const ext_value2 &a,const ext_value2 &b)
 {
  int i,j,n=a.n_var;
  double r2,r4;
  if(b.n_var!=a.n_var) throw b.n_var-a.n_var; // Ошибка размерности
  ext_value2 tmp(a.n_var);
      tmp.val=atan2(a.val,b.val);
      r2=a.val*a.val+b.val*b.val;
      r4=r2*r2;
  for(i=0;i<n;i++)
      tmp.der_1[i]=(b.val*a.der_1[i]-a.val*b.der_1[i])/r2;
  for(i=0;i<n;i++)
      for(j=i;j<n;j++)
          tmp.der_2[i+(j*(j+1))/2]=
          ( b.der_1[j]*a.der_1[i]-b.der_1[i]*a.der_1[j]+
            b.val*a.der_2[i+(j*(j+1))/2]-
            a.val*b.der_2[i+(j*(j+1))/2] )/r2
         -2*(b.val*a.der_1[i]-a.val*b.der_1[i])
           *(a.val*a.der_1[j]+b.val*b.der_1[j])/r4;
  return tmp;
 }
 ext_value1 atan2(const ext_value1 &a,const ext_value1 &b)
 {
  int i,n=a.n_var;
  double r2;
  if(b.n_var!=a.n_var) throw b.n_var-a.n_var; // Ошибка размерности
  ext_value1 tmp(a.n_var);
      tmp.val=atan2(a.val,b.val);
      r2=a.val*a.val+b.val*b.val;
  for(i=0;i<n;i++)
      tmp.der_1[i]=(b.val*a.der_1[i]-a.val*b.der_1[i])/r2;
  return tmp;
 }

// функция asin
 ext_value3 asin(const ext_value3 &a)
 {
  int i,j,k,n=a.n_var;
  ext_value3 tmp(a.n_var);
  double cos_a,cos_a2,cos_a3,cos_a5;
      tmp.val=asin(a.val);
      cos_a2=1.0-a.val*a.val;
      cos_a=sqrt(cos_a2);
      cos_a3=cos_a*cos_a2;
      cos_a5=cos_a3*cos_a2;
  for(i=0;i<n;i++)
      tmp.der_1[i]=a.der_1[i]/cos_a;
  for(i=0;i<n;i++)
      for(j=i;j<n;j++)
          tmp.der_2[i+(j*(j+1))/2]=
          a.val*a.der_1[i]*a.der_1[j]/cos_a3+a.der_2[i+(j*(j+1))/2]/cos_a;
  for(i=0;i<n;i++)
      for(j=i;j<n;j++)
          for(k=j;k<n;k++) 
              tmp.der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4]=
              a.der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4]/cos_a+
              a.val*(a.der_2[j+(k*(k+1))/2]*a.der_1[i]+
                     a.der_2[i+(k*(k+1))/2]*a.der_1[j]+
                     a.der_2[i+(j*(j+1))/2]*a.der_1[k])/cos_a3+
              a.der_1[i]*a.der_1[j]*a.der_1[k]/cos_a3+
              3.0*a.val*a.val*a.der_1[i]*a.der_1[j]*a.der_1[k]/cos_a5;
  return tmp;
 }
 ext_value2 asin(const ext_value2 &a)
 {
  int i,j,n=a.n_var;
  ext_value2 tmp(a.n_var);
  double cos_a,cos_a2,cos_a3;
      tmp.val=asin(a.val);
      cos_a2=1.0-a.val*a.val;
      cos_a=sqrt(cos_a2);
      cos_a3=cos_a*cos_a2;
  for(i=0;i<n;i++)
      tmp.der_1[i]=a.der_1[i]/cos_a;
  for(i=0;i<n;i++)
      for(j=i;j<n;j++)
          tmp.der_2[i+(j*(j+1))/2]=
          a.val*a.der_1[i]*a.der_1[j]/cos_a3+a.der_2[i+(j*(j+1))/2]/cos_a;
  return tmp;
 }
 ext_value1 asin(const ext_value1 &a)
 {
  int i,n=a.n_var;
  ext_value1 tmp(a.n_var);
  double cos_a,cos_a2;
      tmp.val=asin(a.val);
      cos_a2=1.0-a.val*a.val;
      cos_a=sqrt(cos_a2);
  for(i=0;i<n;i++)
      tmp.der_1[i]=a.der_1[i]/cos_a;
  return tmp;
 }

// функция acos
 ext_value3 acos(const ext_value3 &a)
 {
  int i,j,k,n=a.n_var;
  ext_value3 tmp(a.n_var);
  double sin_a,sin_a2,sin_a3,sin_a5;
      tmp.val=acos(a.val);
      sin_a2=1.0-a.val*a.val;
      sin_a=sqrt(sin_a2);
      sin_a3=sin_a*sin_a2;
      sin_a5=sin_a3*sin_a2;
  for(i=0;i<n;i++)
      tmp.der_1[i]=-a.der_1[i]/sin_a;
  for(i=0;i<n;i++)
      for(j=i;j<n;j++)
          tmp.der_2[i+(j*(j+1))/2]=
          -a.val*a.der_1[i]*a.der_1[j]/sin_a3
          -a.der_2[i+(j*(j+1))/2]/sin_a;
  for(i=0;i<n;i++)
      for(j=i;j<n;j++)
          for(k=j;k<n;k++) 
              tmp.der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4]=
              -a.der_3[i+(j*(j+1))/2+(k*(k*k+k+2))/4]/sin_a
              -a.val*(a.der_2[j+(k*(k+1))/2]*a.der_1[i]+
                      a.der_2[i+(k*(k+1))/2]*a.der_1[j]+
                      a.der_2[i+(j*(j+1))/2]*a.der_1[k])/sin_a3
              -a.der_1[i]*a.der_1[j]*a.der_1[k]/sin_a3
              -3.0*a.val*a.val*a.der_1[i]*a.der_1[j]*a.der_1[k]/sin_a5;
  return tmp;
 }
 ext_value2 acos(const ext_value2 &a)
 {
  int i,j,n=a.n_var;
  ext_value2 tmp(a.n_var);
  double sin_a,sin_a2,sin_a3;
      tmp.val=acos(a.val);
      sin_a2=1.0-a.val*a.val;
      sin_a=sqrt(sin_a2);
      sin_a3=sin_a*sin_a2;
  for(i=0;i<n;i++)
      tmp.der_1[i]=-a.der_1[i]/sin_a;
  for(i=0;i<n;i++)
      for(j=i;j<n;j++)
          tmp.der_2[i+(j*(j+1))/2]=
          -a.val*a.der_1[i]*a.der_1[j]/sin_a3
          -a.der_2[i+(j*(j+1))/2]/sin_a;
  return tmp;
 }
 ext_value1 acos(const ext_value1 &a)
 {
  int i,n=a.n_var;
  ext_value1 tmp(a.n_var);
  double sin_a,sin_a2;
      tmp.val=acos(a.val);
      sin_a2=1.0-a.val*a.val;
      sin_a=sqrt(sin_a2);
  for(i=0;i<n;i++)
      tmp.der_1[i]=-a.der_1[i]/sin_a;
  return tmp;
 }