//---------------------------------------------------------------------------
//   Klasa obsługująca macierze   NxN

#include <except.h>
#include <iostream.h>
#include <math.h>
//#include <bcd.h>
#include <vcl.h>


#pragma hdrstop

#include "Mattrix.h"

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

#pragma package(smart_init)


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


//----------------------------------------------------------
__fastcall HMattrix::HMattrix( int dim)
{
try {                                 // TEST FOR EXCEPTIONS.
	aM = new long double* [dim];        // STEP 1: SET UP THE ROWS.
	for (int j = 0; j < dim; j++)
		aM[j] = new long double[dim];  // STEP 2: SET UP THE COLUMNS

} catch (...) {  // ENTER THIS BLOCK ONLY IF xalloc IS THROWN.
		// YOU COULD REQUEST OTHER ACTIONS BEFORE TERMINATING
		MessageBox(0, L"Brak pamięci (alokacja dla macierzy). Bye ...",
        	L"Błąd!", MB_OK);
		exit(-1);
}

	for (int i = 0; i < dim; i++)
		for (int j = 0; j < dim; j++)
			aM[i][j] = (long double)0;

	nDim=dim;
}
//----------------------------------------------------------
__fastcall HMattrix::~HMattrix()
{
	if (aM){
		for (int i = 0; i < nDim;  i++)
			delete[] aM[i];     // STEP 1: DELETE THE COLUMNS

		delete[] aM;           // STEP 2: DELETE THE ROWS

	}
}

//----------------------------------------------------------
void __fastcall
HMattrix::StoreColumn(long double *aCol, int col)
{
	for (int j = 0; j < nDim; j++)
		aM[col][j] = aCol[j];
}

//----------------------------------------------------------
void __fastcall
HMattrix::StoreRow(long double *aRow, int row)
{
	for (int i = 0; i < nDim; i++)
		aM[i][row] = aRow[i];
}
//----------------------------------------------------------
void __fastcall
HMattrix::StoreElement(long double val, int i, int j)    // i=row, j=col
{
	aM[i][j]= val;
}

//----------------------------------------------------------
void __fastcall
HMattrix::StoreElement(long double *val, int i, int j)    // i=row, j=col
{
	aM[i][j]= *val;
}
//----------------------------------------------------------
void __fastcall
HMattrix::StoreMattrix(HMattrix *M)
{
	int dim=nDim;
	if (M->nDim<dim)
		dim=M->nDim;
	for(int i=0; i<dim; i++)
		for (int j = 0; j < dim; j++)
			aM[i][j] = M->M(i,j);
}

//----------------------------------------------------------
long double ** __fastcall HMattrix::GetM()
{
	return aM;
}
//----------------------------------------------------------

long double __fastcall HMattrix::M(int i, int j)
{
	long double val=0;
	if (i<nDim && j<nDim && aM)
		val=aM[i][j];
	return val;
}

//----------------------------------------------------------
//int
//HMattrix::DeleteRowCol(bcd **aMx, int dim, int row, int col, bcd ***aMNew)
long double ** __fastcall HMattrix::DeleteRowCol(long double **aMx, int dim, int row, int col)
{
	int nI, nJ; //, retval=0;
	long double **aMtmp=NULL;

try {                                 // TEST FOR EXCEPTIONS.
	if(aMx && col<dim && row<dim){
		aMtmp = new long double* [dim-1];        // STEP 1: SET UP THE ROWS.
		for (int j = 0; j < dim-1; j++)
			aMtmp[j] = new long double[dim-1];  // STEP 2: SET UP THE COLUMNS

		for (int i = 0; i < dim-1; i++)
			for (int j = 0; j < dim-1; j++){
				nI=(i<row ? i : i+1);
				nJ=(j<col ? j : j+1);
				aMtmp[i][j] = aMx[nI][nJ];
			}

//		*aMNew=aMtmp;

		//**** retval=1;
	}
} catch (xalloc) {  // ENTER THIS BLOCK ONLY IF xalloc IS THROWN.
		// YOU COULD REQUEST OTHER ACTIONS BEFORE TERMINATING
		MessageBox(0, L"Brak pamięci (alokacja dla macierzy). Bye ...",
							L"Błąd!", MB_OK);
//		cout << "Brak pamięci (alokacja dla macierzy).";
		//*********retval=0;
}

//	return retval;
	return aMtmp;
}

//----------------------------------------------------------
void __fastcall
HMattrix::XchCol(int col, long double *aCol)
{
	StoreColumn(aCol, col);
}
//----------------------------------------------------------
void __fastcall
HMattrix::deallocate(long double **aX, int dim)
{
	if (aX){
		for (int i = 0; i < dim;  i++)
			if( aX[i])
				delete[] aX[i];     // STEP 1: DELETE THE COLUMNS

		delete[] aX	;           // STEP 2: DELETE THE ROWS
	}
}
//----------------------------------------------------------

long double __fastcall HMattrix::Det()
{
	long double nDet=(long double)0;

	if( nDim==1 )
		nDet=aM[0][0];
	else if (nDim==2)
		nDet=aM[0][0]*aM[1][1]-aM[0][1]*aM[1][0];
	else
		for (int i = 0; i<nDim; i++)
			nDet+=(pow((double)-1, 1+i))*aM[0][i]*Minor(aM,0,i, nDim);

	return(nDet);
}

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

long double  __fastcall
HMattrix::Det(long double **aMn, int dim)
{
	long double nDet=(long double)0;

	if( dim==1 )
		nDet=aMn[0][0];
	else if (dim==2)
		nDet=aMn[0][0]*aMn[1][1]-aMn[0][1]*aMn[1][0];
	else
		for (int i = 0; i<dim; i++)
			nDet+=(pow((double)-1, 1+i))*aMn[0][i]*Minor(aMn, 0,i, dim);

	return(nDet);
}
//----------------------------------------------------------

long double __fastcall
HMattrix::Minor(long double **aMx, int row, int col, int dim)
{
	long double **aMn, mDet=0;

//	if (DeleteRowCol(aMx, dim, row, col, &aMn))
	aMn=DeleteRowCol(aMx, dim, row, col);

		if(aMn){
			mDet=Det(aMn, dim-1);
			deallocate(aMn, dim-1);
		}

	return( mDet );
}

//----------------------------------------------------------
//	rozwiązanie równania stopnia nDim
int __fastcall
equation(HMattrix *M, long double *aB, int nDim, long double *result)
{
//	HMattrix M(nDim, aM);
	long double nDet=M->Det();
	long double nDetN;

//	return 1;

	if (nDet==(long double)0)
		return 0;

	HMattrix *aMn;
	for (int i=0; i<nDim; i++) {
		aMn = new HMattrix(nDim);
		aMn->StoreMattrix(M);
		aMn->XchCol(i, aB);
		nDetN=aMn->Det();
		result[nDim-i-1]=nDetN/nDet;
		delete aMn;
	}
	return 1;
}
//----------------------------------------------------------
//------------------------------------------------------
//    wartość wielomianu 3-stopnia dla punktu x
double __fastcall Polynom3Value(double x, long double *Cn)
{
	return(x*x*x*Cn[3]+x*x*Cn[2]+x*Cn[1]+Cn[0]);
}
//------------------------------------------------------
// wielomian stopnia trzeciego
bool __fastcall Polynom3(ePoint *aData, int count, long double *result)
{
	long double nX6 = eSup(aData,count,6);
	long double nX5 = eSup(aData,count,5);
	long double nX4 = eSup(aData,count,4);
	long double nX3 = eSup(aData,count,3);
	long double nX2 = eSup(aData,count,2);
	long double nX1 = eSup(aData,count,1);
	long double nX3Y= eSup(aData,count,3,1);
	long double nX2Y= eSup(aData,count,2,1);
	long double nXY = eSup(aData,count,1,1);
	long double nY1 = eSup(aData,count,0,1);

	// macierz dla m.n.kw.

	HMattrix *aM = new HMattrix(4);
	long double r0[] =  {nX6, nX5, nX4, nX3};
	long double r1[] =	{nX5, nX4, nX3, nX2};
	long double r2[] =	{nX4, nX3, nX2, nX1};
	long double r3[] =	{nX3, nX2, nX1, count};
	aM->StoreRow(r0,0);
	aM->StoreRow(r1,1);
	aM->StoreRow(r2,2);
	aM->StoreRow(r3,3);

	long double aB[]={nX3Y, nX2Y, nXY, nY1};
//	bcd aCN[4]={bcd(0),bcd(0),bcd(0),bcd(0)};==result

	bool retval=equation(aM, aB,4, result );

   delete aM;
   return(retval);
}
//----------------------------------------------------------

// ax^2+b
bool __fastcall Polynom2ab(ePoint *aData, int count, long double *result)
{
	long double nX3 = eSup(aData,count,3);
	long double nX2 = eSup(aData,count,2);
	long double nX1 = eSup(aData,count,1);
	long double nXY = eSup(aData,count,1,1);
	long double nY1 = eSup(aData,count,0,1);

	// macierz dla m.n.kw.

	HMattrix *aM = new HMattrix(2);
	long double r0[] = 	{ nX3,  nX1};
	long double r1[] =   { nX2,  count};
//	long double r2[] =
//	long double r3[] =
	aM->StoreRow(r0,0);
	aM->StoreRow(r1,1);
//	aM->StoreRow(r2,2);
//	aM->StoreRow(r3,3);

//	long double aB[]={nX3Y, nX2Y, nXY, nY1};
	long double aB[]={ nXY, nY1};

	bool retval=equation(aM, aB,2, result );

   delete aM;
   return(retval);
}
// ax^2+bx+c
bool __fastcall Polynom2abc(ePoint *aData, int count, long double *result)
{
	long double nX4 = eSup(aData,count,4);
	long double nX3 = eSup(aData,count,3);
	long double nX2 = eSup(aData,count,2);
	long double nX1 = eSup(aData,count,1);
	long double nXY = eSup(aData,count,1,1);
	long double nY1 = eSup(aData,count,0,1);
	long double nX2Y= eSup(aData,count,2,1);

	// macierz dla m.n.kw.

	HMattrix *aM = new HMattrix(3);
	long double r0[] =   {nX4, nX3, nX2};
	long double r1[] =   {nX3, nX2, nX1};
	long double r2[] =   {nX2, nX1, count};
//	long double r3[] =
	aM->StoreRow(r0,0);
	aM->StoreRow(r1,1);
	aM->StoreRow(r2,2);
//	aM->StoreRow(r3,3);

//	long double aB[]={nX3Y, nX2Y, nXY, nY1};
	long double aB[]={ nX2Y, nXY, nY1};

	bool retval=equation(aM, aB,3, result );

   delete aM;
   return(retval);
}
//----------------------------------------------------------
//------------------------------------------------------
//		trochę matematyki

long double __fastcall eSup(ePoint *aData, int count, int nI, int nJ)
{
	long double nSuma=0;
	for (int i = 0; i<count; i++)
		if( aData[i].x != (long double)0 && aData[i].y != (long double)0){
			nSuma += powl(aData[i].x, nI) * powl(aData[i].y, nJ);
		}
	return nSuma;
}
//------------------------------------------------------
//		przypisywanie wartosci do ePoint
void __fastcall  AssignPointValue(ePoint *dt, double x, double y)
{
	dt->x=(long double) x;
	dt->y=(long double) y;
}

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