det_inv.c

很好的一个约束遗传算法优化程序

#include "genocop.h"

/********************************************************************************/
/*                                                                              */
/*           FUNCTION NAME     :   inverse()                                    */
/*                                                                              */
/*           SYNOPSIS          :   void inverse(a,a_inverse,n)                  */
/*                                                                              */
/*           DESCRIPTION       :   This function finds the inverse matrix of a  */
/*                                  given float matrix                          */
/*                                                                              */
/*           FUNCTIONS CALLED  :   det(),                                       */
/*                                 matrix().                                    */
/*                                                                              */
/*           CALLING FUNCITONS :   main()                                       */
/*                                                                              */
/*           AUTHOR            :   Swarnalatha Swaminathan                      */
/*                                                                              */
/*           DATE              :   1/17/92                                      */
/*                                                                              */
/*                                                                              */
/*           REV            DATE            BY           DESCRIPTION            */
/*           ---            ----            --           -----------            */
/*                                                                              */
/*                                                                              */
/********************************************************************************/




void 
inverse(a, a_inverse, n)
    MATRIX          a, a_inverse;	    /* the original and the inverse
					     * matrices */
    int             n;			    /* the dimension of the matrix */
{
    int             i, j, k, k1, l, l1, n1; /* counter variables */
    float           d1,			    /* determinat of the matrix */
                    v, w;
    MATRIX          b;			    /* a temporary float matrix */

    b = matrix(1, n, 1, n);
    for (i = 1; i <= n; i++)
	for (j = 1; j <= n; j++)
	    b[i][j] = 0.0;

    d1 = det(a, n);

    if (d1 == 0.0)
    {
	printf("Variables are dependant : No inverse exists.\n\n");
	exit(1);
    } else
    {
	v = -1.0;
	for (i = 1; i <= n; i++)
	{
	    v = -v;
	    w = -1.0;
	    for (j = 1; j <= n; j++)
	    {
		w = -w;
		n1 = n - 1;
		for (k = 1; k <= n1; k++)
		{
		    k1 = k;
		    if (k >= i)
			k1 = k + 1;
		    for (l = 1; l <= n1; l++)
		    {
			l1 = l;
			if (l >= j)
			    l1 = l + 1;
			b[k][l] = a[k1][l1];
		    }
		}
		a_inverse[j][i] = det(b, n1) / d1 * v * w;
	    }
	}
    }
    free_matrix(b, 1, n, 1);
}

/********************************************************************************/
/*                                                                              */
/*           FUNCTION NAME     :   det()                                        */
/*                                                                              */
/*           SYNOPSIS          :   float det(input_matrix,nl)                   */
/*                                                                              */
/*           DESCRIPTION       :   This function returns the determinant of a   */
/*                                  given float matrix                          */
/*                                                                              */
/*           FUNCTIONS CALLED  :   matrix().                                    */
/*                                                                              */
/*           CALLING FUNCITONS :   inverse(),                                   */
/*                                 p_equalities().                              */
/*                                                                              */
/*           AUTHOR            :   Swarnalatha Swaminathan                      */
/*                                                                              */
/*           DATE              :   1/17/92                                      */
/*                                                                              */
/*                                                                              */
/*           REV            DATE            BY           DESCRIPTION            */
/*           ---            ----            --           -----------            */
/*                                                                              */
/*                                                                              */
/********************************************************************************/




float 
det(input_matrix, nl)
    MATRIX          input_matrix;	    /* input float matrix */
    int             nl;			    /* order of the matrix */
{
    int             ia, ib, i1, j1, n2, nz, n3, l, l1, m, m1;
    float           d, x;
    MATRIX          temp_input_matrix, b1;

    temp_input_matrix = matrix(1, nl, 1, nl);
    b1 = matrix(1, nl, 1, nl);

    for (i1 = 1; i1 <= nl; i1++)
	for (j1 = 1; j1 <= nl; j1++)
	    temp_input_matrix[i1][j1] = input_matrix[i1][j1];

    d = 1.0;
    for (ia = 1; ia <= nl; ia++)
    {
	x = -1.0;
	n2 = nl - ia + 1;
	n3 = n2 - 1;
	nz = 0;
	for (ib = 1; ib <= n2; ib++)
	    if (temp_input_matrix[ib][1] != 0.0)
		nz = ib;
	if (nz == 0.0)
	    d = 0.0;
	else
	{
	    for (i1 = 1; i1 <= nz; i1++)
		x = -x;
	    d *= temp_input_matrix[nz][1] * x;

	    for (l = 1; l <= n3; l++)
	    {
		l1 = l;
		if (l >= nz)
		    l1 = l + 1;
		for (m = 1; m <= n3; m++)
		{
		    m1 = m + 1;
		    b1[l][m] = temp_input_matrix[l1][m1] - temp_input_matrix[l1][1]
			/ temp_input_matrix[nz][1] * temp_input_matrix[nz][m1];
		}
	    }
	}
	for (i1 = 1; i1 <= n3; i1++)
	    for (j1 = 1; j1 <= n3; j1++)
		temp_input_matrix[i1][j1] = b1[i1][j1];

    }
    free_matrix(b1, 1, nl, 1);
    free_matrix(temp_input_matrix, 1, nl, 1);

    return (d);
}