ocn_f.h

图像处理的压缩算法

		calculates all the eigenvalues and eigenvectors of 
		a real unsymmetric matrix.

Example:
	To calculate all the eigenvalues and eigenvectors of the real matrix.
	.
	..
	1.5 0.1 4.5 -1.5
	-22.5 3.5 1 2 .5 -2.5
	-2.5 0.3 4.5 -2.5
	-2.5 0.1 4.5 2.5
	.
	..
	.

void test_nag_real_eigensystem()
{

	int i, j, n;
	matrix a;
	a.SetSize(4,4);
	matrix<complex>r, v;
	r.SetSize(1,4);
	v.SetSize(4,4);
	int iter[4];
  
	a[0][0] = 1.5;	
	a[0][1] = 0.1;
	a[0][2] = 4.5;
	a[0][3] = -1.5;
		
	a[1][0] = -22.5;
	a[1][1] = 3.5;
	a[1][2] = 12.5;
	a[1][3] = -2.5;
		
	a[2][0] = -2.5;
	a[2][1] = 0.3;
	a[2][2] = 4.5;
	a[2][3] = -2.5;
		
	a[3][0] = -2.5;
	a[3][1] = 0.1; 
	a[3][2] = 4.5;
	a[3][3] = 2.5;
	n = 4;
	nag_real_eigensystem(n, a, 4, r, v, 4, iter);
	printf("Eigenvalues\n");
	for (i=0; i<n; i++)
	{
		complex temp = r[0][i];
		printf("(%7.3f, %7.3f) \n", temp.m_re, temp.m_im);
	
	}
	printf("\nEigenvectors\n");
	for (i=0; i<n; i++)
		for (j=0; j<n; j++)
		{
			complex temp = v[i][j];
			printf("(%7.3f, %7.3f)", temp.m_re, temp.m_im);
			if(j%4==3 || j==n-1)
				printf("\n");
		}

}


	The output is following:
	
	Eigenvalues
	
	( 3.000,  4.000)
	( 3.000, -4.000)
	( 4.000,  0.000)
	( 2.000,  0.000)
	
	Eigenvectors
	
	( 0.113, -0.151) ( 0.113, 0.151) ( -0.033, 0.000) ( 0.063, 0.000)
	( 0.945,  0.000) ( 0.945, 0.000) (  0.988, 0.000) ( 0.996, 0.000)
	( 0.189,  0.000) ( 0.189, 0.000) (  0.011, 0.000) ( 0.006, 0.000)
	( 0.113, -0.151) ( 0.113, 0.151) (  0.154, 0.000) ( 0.063, 0.000)
		

Parameters:

Return:
	This function returns NAG error code, 0 if no error.
	
	362: More than _value_ iterations are required to isolate all the eigenvalues.
	11: On entry, n must not be less than 1: n = _value_.
	17: On entry, tda = _value_ while n = _value_. These parameters must satisfy tda = n.  On entry, tdv = _value_ while n = _value_. These parameters must satisfy tdv = n.
	73: Memory allocation failed.

	successfully call of the nag_real_eigensystem function.


*/
	
	int  nag_real_eigensystem(
	int n, // the order of the matrix A.
	double a[], // Input: the n by n matrix A.  Output: the array is overwritten.
	int tda, // the second dimension of the array a as declared in the function from which nag real eigensystem is called.
	complex r[], // the eigenvalues.
	complex v[], // the eigenvectors, stored by columns. The ith column corresponds to the ith eigenvalue. The eigenvectors are normalised so that the sum of the squares of the moduli of the elements is equal to 1an d the element of largest modulus is real. This ensures that real eigenvalues have real eigenvectors.
	int tdv, // the second dimension of the array v as declared in the function from which nag real eigensystem is called.
	int iter[] // contains the number of iterations used to find the ith eigenvalue. If iter[i - 1] is negative, the ith eigenvalue is the second of a pair found simultaneously.
	);


/**	f02awc
		calculates all the eigenvalues of a complex Hermitian matrix.

Example:
	To calculate all the eigenvalues of the complex Hermitian matrix:
	.
	..
	0.50 0.00 1.84 + 1.38i 2.08 - 1.56i
	0.00 0.50 1.1 2 +0.84i -0.56 + 0.42i
	1.84 - 1.38i 1.12 - 0.84i 0.50 0.00
	2.08 + 1.56i -0.56 - 0.42i 0.00 0.50
	.
	..
	.

void test_nag_hermitian_eigenvalues()
{	
	int i, j, n;
	matrix<complex>a;
	a.SetSize(4,4);
	double r[4];
	n = 4;

	a[0][0] = 0.50 + 0.00i;
	a[0][1] =  0.00 + 0.00i;
	a[0][2] = 1.84 + 1.38i;
	a[0][3] = 2.08 - 1.56i;
	
	a[1][0] = 0.00 + 0.00i;
	a[1][1] = 0.50 + 0.00i;
	a[1][2] = 1.12 + 0.84i;
	a[1][3] = -0.56 + 0.42i;
	
	a[2][0] = 1.84 - 1.38i;
	a[2][1] = 1.12 - 0.84i;
	a[2][2] = 0.50 + 0.00i;
	a[2][3] = 0.00 + 0.00i;
	
	a[3][0] = 2.08 + 1.56i;
	a[3][1] = -0.56 - 0.42i;
	a[3][2] = 0.00 + 0.00i;
	a[3][3] = 0.50 + 0.00i;	

	nag_hermitian_eigenvalues(n, a, 4, r);
	printf("Eigenvalues\n");
	for (i=0; i<n; i++)
		printf("%9.4f", r[i]);
	printf("\n");
}

	The output is following:
	
	Eigenvalues
	
	-3.0000 -1.0000 2.0000 4.0000
	

Parameters:

Return:
	This function returns NAG error code, 0 if no error.

	11: On entry, n must not be less than 1: n = _value_.
	17: On entry tda = _value_ while n = _value_. These parameters must satisfy tda = n.
	73: Memory allocation failed.
	362: More than _value_ iterations are required to isolate all the eigenvalues.

	successfully call of the nag_hermitian_eigenvalues function.

*/

	int  nag_hermitian_eigenvalues(
	int n, // the order of the matrix A.
	complex a[], // Input: the elements of the lower triangle of the n by n complex Hermitian matrix A. Elements of the array above the diagonal need not be set.  Output: the array is overwritten.
	int tda, // the last dimension of the array a as declared in the function from which nag hermitian eigenvalues is called.
	double r[] // the eigenvalues in ascending order.
	);

/**	f02axc
		calculates all the eigenvalues and eigenvectors of 
		a complex Hermitian matrix.


Example:
	To calculate the eigenvalues and eigenvectors of the complex Hermitian matrix:
	.
	..
	0.50 0.00 1.84 + 1.38i 2.08 - 1.56i
	0.00 0.50 1.1 2 +0.84i -0.56 + 0.42i
	1.84 - 1.38i 1.12 - 0.84i 0.50 0.00
	2.08 + 1.56i -0.56 - 0.42i 0.00 0.50
	.
	..
	.

void test_nag_hermitian_eigensystem()
{	
	int i, j, n;
	matrix<complex> a, v;
	a.SetSize(4,4);
	v.SetSize(4,4);
	double r[4];
	n = 4;

	a[0][0] = 0.50 + 0.00i;
	a[0][1] =  0.00 + 0.00i;
	a[0][2] = 1.84 + 1.38i;
	a[0][3] = 2.08 - 1.56i;
	
	a[1][0] = 0.00 + 0.00i;
	a[1][1] = 0.50 + 0.00i;
	a[1][2] = 1.12 + 0.84i;
	a[1][3] = -0.56 + 0.42i;
	
	a[2][0] = 1.84 - 1.38i;
	a[2][1] = 1.12 - 0.84i;
	a[2][2] = 0.50 + 0.00i;
	a[2][3] = 0.00 + 0.00i;
	
	a[3][0] = 2.08 + 1.56i;
	a[3][1] = -0.56 - 0.42i;
	a[3][2] = 0.00 + 0.00i;
	a[3][3] = 0.50 + 0.00i;
	nag_hermitian_eigensystem(n, a, 4, r, v, 4);
	printf("Eigenvalues\n");
	for (i=0; i<n; i++)
		printf("%9.4f", r[i]);
	printf("\nEigenvectors\n");
	for (i=0; i<n; i++)
		for (j=0; j<n; j++)
		{
			complex temp = v[i][j];
			printf("(%7.3f %7.3f )",temp.m_re, temp.m_im);
			if(j%4==3 || j==n-1)
				printf("\n");
		}
}

	The output is following:
		
	Eigenvalues
	
	-3.0000 -1.0000 2.0000 4.0000
	
	Eigenvectors
	
	(  0.700  0.000 ) ( -0.100 0.000 ) ( -0.100  0.000 ) ( 0.700  0.000 )
	(  0.100  0.000 ) (  0.700 0.000 ) (  0.700  0.000 ) ( 0.100  0.000 )
	( -0.400  0.300 ) ( -0.400 0.300 ) (  0.400 -0.300 ) ( 0.400 -0.300 )
	( -0.400 -0.300 ) (  0.400 0.300 ) ( -0.400 -0.300 ) ( 0.400  0.300 )

Parameters:

Return:
	This function returns NAG error code, 0 if no error.
	
	11: On entry, n must not be less than 1: n = _value_.
	17: On entry tda = _value_ while n = _value_. These parameters must satisfy tda = n.  On entry tdv = _value_ while n = _value_. These parameters must satisfy tdv = n.
	73: Memory allocation failed.

	successfully call of the nag_hermitian_eigensystem function.

*/

	int  nag_hermitian_eigensystem(
	int n, // the order of the matrix A.
	complex a[], // the elements of the lower triangle of the n by n complex Hermitian matrix A. Elements of the array above the diagonal need not be set. 
	int tda, // the last dimension of the array a as declared in the function from which nag hermitian eigensystem is called.
	double r[], // the eigenvalues in ascending order.
	complex v[], // the eigenvectors, stored by columns. The ith column corresponds to the ith eigenvector. The eigenvectors are normalised so that the sum of the squares of the moduli of the elements is equal to 1an d the element of largest modulus is real.
	int tdv // the last dimension of the array v as declared in the function from which nag hermitian eigensystem is called.
	);

/**	f02bjc
		calculates all the eigenvalues and, if required, all the 
		eigenvectors of the generalized eigenproblem Ax = lamdaBx 
		where A and B are real, square matrices, using the QZ algorithm.

Example:
	To find all the eigenvalues and eigenvectors of Ax = ?Bx where
	A =
	.
	..
	3.9 4.3 4.3 4.4
	12.5 21.5 21.5 26.0
	-34.5 -47.5 -43.5 -46.0
	-0.5 7.5 3.5 6.0
	.
	..
	and
	B =
	.
	..
	1 1 1 1
	2 3 3 3
	-3 -5 -4 -4
	1 4 3 4
	.
	..
	.

void test_nag_real_general_eigensystem()
{
	int i, j, k, n, ip, iter[8];
	matrix<complex> alfa;
	alfa.SetSize(1,8);
	double beta[8], eps1;
	matrix a, b, z;
	a.SetSize(8,8);
	b.SetSize(8,8);
	z.SetSize(8,8);
	BOOL matz;	
	n = 4;
	
	a[0][0] = 3.9;
	a[0][1] = 12.5;
	a[0][2] = -34.5;
	a[0][3] = -0.5;
	
	a[1][0] = 4.3;  
	a[1][1] = 21.5;
	a[1][2] = -47.5;
	a[1][3] = 7.5;
	
	a[2][0] = 4.3;
	a[2][1] = 21.5;
	a[2][2] = -43.5;
	a[2][3] = 3.5;
	
	a[3][0] = 4.4;
	a[3][1] = 26.0;
	a[3][2] = -46.0;
	a[3][3] = 6.0;
	
	b[0][0] = 1.0;
	b[0][1] = 2.0; 
	b[0][2] = -3.0; 
	b[0][3] = 1.0;
	
	b[1][0] = 1.0; 
	b[1][1] = 3.0; 
	b[1][2] = -5.0;
	b[1][3] = 4.0;
	
	b[2][0] = 1.0;
	b[2][1] = 3.0;
	b[2][2] = -4.0;
	b[2][3] = 3.0;
	
	b[3][0] = 1.0;
	b[3][1] = 3.0;
	b[3][2] = -4.0;
	b[3][3] = 4.0;
	  	
	matz = TRUE;
	eps1 = 1e-8;
	nag_real_general_eigensystem(n, a, 8, b, 8, eps1, alfa, beta, matz, z, 8, iter);
	ip = 0;
	for (i=0; i<n; ++i)
	{
		complex tmp = alfa[0][i];
		printf("Eigensolution %4ld\n",i+1);
		printf("alfa[%ld].re %7.3f",i,tmp.m_re);
		printf(" alfa[%ld].im %7.3f",i,tmp.m_im);
		printf(" beta[%ld] %7.3f\n",i,beta[i]);
		if (beta[i] == 0.0)
			printf("lambda is infinite");
		else if (tmp.m_im == 0.0)
		{
			printf("lambda %7.3f\n",tmp.m_re/beta[i]);
			printf("Eigenvector\n");
			for (j=0; j<n; ++j)
				printf("%7.3f\n", z[j][i]);
		}
		else
		{
			complex temp = alfa[0][i];
			printf("lambda %7.3f %7.3f\n",temp.m_re/beta[i], temp.m_im/beta[i] );
			printf("Eigenvector\n");
			k = pow(-1, ip+2);
			for (j=0; j<n; ++j)
			{
				printf("%7.3f",z[j][i-ip]);
				printf("%7.3f\n",k*z[j][i-ip+1]);
			}
			ip = 1-ip;
		}
	}
	printf("Number of iterations (machine-dependent)\n");
	for (i=0; i<n; ++i)
		printf("%2ld",iter[i]);
	printf("\n");	
}

	The output is following:
	
	Eigensolution 1
	
	alfa[0].re 3.801 alfa[0].im 0.000 beta[0] 1.900
	
	lambda 2.000
	
	Eigenvector
	
	0.996
	0.006
	0.063
	0.063
	
	Eigensolution 2
	
	alfa[1].re 1.563 alfa[1].im 2.084 beta[1] 0.521
	
	lambda 3.000 4.000
	
	Eigenvector
	
	0.945 0.000
	0.189 0.000
	0.113 -0.151
	0.113 -0.151
	
	Eigensolution 3
	
	alfa[2].re 3.030 alfa[2].im -4.040 beta[2] 1.010
	
	lambda 3.000 -4.000
	
	Eigenvector
	
	0.945 0.000
	0.189 0.000
	0.113 0.151
	0.113 0.151
	
	Eigensolution 4
	
	alfa[3].re 4.000 alfa[3].im 0.000 beta[3] 1.000
	
	lambda 4.000
	
	Eigenvector
	
	0.988
	0.011
	-0.033
	0.154
	
	Number of iterations (machine-dependent)
	0 0 5 0


Parameters:

Return:
	This function returns NAG error code, 0 if no error.
	
	11: On entry, n must not be less than 1: n = _value_.
	17: On entry tda = _value_ while n = _value_. These parameters must satisfy tda = n.  On entry tdb = _value_ while n = _value_. These parameters must satisfy tdb = n.  On entry tdv = _value_ while n = _value_. These parameters must satisfy tdv = n.
	374: More than n