kusano 7d535a
#include "f2c.h"
kusano 7d535a
kusano 7d535a
/* Subroutine */ int zsymv_(char *uplo, integer *n, doublecomplex *alpha, 
kusano 7d535a
	doublecomplex *a, integer *lda, doublecomplex *x, integer *incx, 
kusano 7d535a
	doublecomplex *beta, doublecomplex *y, integer *incy)
kusano 7d535a
{
kusano 7d535a
/*  -- LAPACK auxiliary routine (version 2.0) --   
kusano 7d535a
       Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd.,   
kusano 7d535a
       Courant Institute, Argonne National Lab, and Rice University   
kusano 7d535a
       October 31, 1992   
kusano 7d535a
kusano 7d535a
kusano 7d535a
    Purpose   
kusano 7d535a
    =======   
kusano 7d535a
kusano 7d535a
    ZSYMV  performs the matrix-vector  operation   
kusano 7d535a
kusano 7d535a
       y := alpha*A*x + beta*y,   
kusano 7d535a
kusano 7d535a
    where alpha and beta are scalars, x and y are n element vectors and   
kusano 7d535a
    A is an n by n symmetric matrix.   
kusano 7d535a
kusano 7d535a
    Arguments   
kusano 7d535a
    ==========   
kusano 7d535a
kusano 7d535a
    UPLO   - CHARACTER*1   
kusano 7d535a
             On entry, UPLO specifies whether the upper or lower   
kusano 7d535a
             triangular part of the array A is to be referenced as   
kusano 7d535a
             follows:   
kusano 7d535a
kusano 7d535a
                UPLO = 'U' or 'u'   Only the upper triangular part of A   
kusano 7d535a
                                    is to be referenced.   
kusano 7d535a
kusano 7d535a
                UPLO = 'L' or 'l'   Only the lower triangular part of A   
kusano 7d535a
                                    is to be referenced.   
kusano 7d535a
kusano 7d535a
             Unchanged on exit.   
kusano 7d535a
kusano 7d535a
    N      - INTEGER   
kusano 7d535a
             On entry, N specifies the order of the matrix A.   
kusano 7d535a
             N must be at least zero.   
kusano 7d535a
             Unchanged on exit.   
kusano 7d535a
kusano 7d535a
    ALPHA  - COMPLEX*16   
kusano 7d535a
             On entry, ALPHA specifies the scalar alpha.   
kusano 7d535a
             Unchanged on exit.   
kusano 7d535a
kusano 7d535a
    A      - COMPLEX*16 array, dimension ( LDA, N )   
kusano 7d535a
             Before entry, with  UPLO = 'U' or 'u', the leading n by n   
kusano 7d535a
             upper triangular part of the array A must contain the upper 
kusano 7d535a
  
kusano 7d535a
             triangular part of the symmetric matrix and the strictly   
kusano 7d535a
             lower triangular part of A is not referenced.   
kusano 7d535a
             Before entry, with UPLO = 'L' or 'l', the leading n by n   
kusano 7d535a
             lower triangular part of the array A must contain the lower 
kusano 7d535a
  
kusano 7d535a
             triangular part of the symmetric matrix and the strictly   
kusano 7d535a
             upper triangular part of A is not referenced.   
kusano 7d535a
             Unchanged on exit.   
kusano 7d535a
kusano 7d535a
    LDA    - INTEGER   
kusano 7d535a
             On entry, LDA specifies the first dimension of A as declared 
kusano 7d535a
  
kusano 7d535a
             in the calling (sub) program. LDA must be at least   
kusano 7d535a
             max( 1, N ).   
kusano 7d535a
             Unchanged on exit.   
kusano 7d535a
kusano 7d535a
    X      - COMPLEX*16 array, dimension at least   
kusano 7d535a
             ( 1 + ( N - 1 )*abs( INCX ) ).   
kusano 7d535a
             Before entry, the incremented array X must contain the N-   
kusano 7d535a
             element vector x.   
kusano 7d535a
             Unchanged on exit.   
kusano 7d535a
kusano 7d535a
    INCX   - INTEGER   
kusano 7d535a
             On entry, INCX specifies the increment for the elements of   
kusano 7d535a
             X. INCX must not be zero.   
kusano 7d535a
             Unchanged on exit.   
kusano 7d535a
kusano 7d535a
    BETA   - COMPLEX*16   
kusano 7d535a
             On entry, BETA specifies the scalar beta. When BETA is   
kusano 7d535a
             supplied as zero then Y need not be set on input.   
kusano 7d535a
             Unchanged on exit.   
kusano 7d535a
kusano 7d535a
    Y      - COMPLEX*16 array, dimension at least   
kusano 7d535a
             ( 1 + ( N - 1 )*abs( INCY ) ).   
kusano 7d535a
             Before entry, the incremented array Y must contain the n   
kusano 7d535a
             element vector y. On exit, Y is overwritten by the updated   
kusano 7d535a
             vector y.   
kusano 7d535a
kusano 7d535a
    INCY   - INTEGER   
kusano 7d535a
             On entry, INCY specifies the increment for the elements of   
kusano 7d535a
             Y. INCY must not be zero.   
kusano 7d535a
             Unchanged on exit.   
kusano 7d535a
kusano 7d535a
   ===================================================================== 
kusano 7d535a
  
kusano 7d535a
kusano 7d535a
kusano 7d535a
       Test the input parameters.   
kusano 7d535a
kusano 7d535a
    
kusano 7d535a
   Parameter adjustments   
kusano 7d535a
       Function Body */
kusano 7d535a
    /* System generated locals */
kusano 7d535a
    integer a_dim1, a_offset, i__1, i__2, i__3, i__4, i__5;
kusano 7d535a
    doublecomplex z__1, z__2, z__3, z__4;
kusano 7d535a
    /* Local variables */
kusano 7d535a
    static integer info;
kusano 7d535a
    static doublecomplex temp1, temp2;
kusano 7d535a
    static integer i, j;
kusano 7d535a
    extern logical lsame_(char *, char *);
kusano 7d535a
    static integer ix, iy, jx, jy, kx, ky;
kusano 7d535a
    extern /* Subroutine */ int xerbla_(char *, integer *);
kusano 7d535a
kusano 7d535a
kusano 7d535a
#define X(I) x[(I)-1]
kusano 7d535a
#define Y(I) y[(I)-1]
kusano 7d535a
kusano 7d535a
#define A(I,J) a[(I)-1 + ((J)-1)* ( *lda)]
kusano 7d535a
kusano 7d535a
    info = 0;
kusano 7d535a
    if (! lsame_(uplo, "U") && ! lsame_(uplo, "L")) {
kusano 7d535a
	info = 1;
kusano 7d535a
    } else if (*n < 0) {
kusano 7d535a
	info = 2;
kusano 7d535a
    } else if (*lda < max(1,*n)) {
kusano 7d535a
	info = 5;
kusano 7d535a
    } else if (*incx == 0) {
kusano 7d535a
	info = 7;
kusano 7d535a
    } else if (*incy == 0) {
kusano 7d535a
	info = 10;
kusano 7d535a
    }
kusano 7d535a
    if (info != 0) {
kusano 7d535a
	xerbla_("ZSYMV ", &info);
kusano 7d535a
	return 0;
kusano 7d535a
    }
kusano 7d535a
kusano 7d535a
/*     Quick return if possible. */
kusano 7d535a
kusano 7d535a
    if (*n == 0 || alpha->r == 0. && alpha->i == 0. && (beta->r == 1. && 
kusano 7d535a
	    beta->i == 0.)) {
kusano 7d535a
	return 0;
kusano 7d535a
    }
kusano 7d535a
kusano 7d535a
/*     Set up the start points in  X  and  Y. */
kusano 7d535a
kusano 7d535a
    if (*incx > 0) {
kusano 7d535a
	kx = 1;
kusano 7d535a
    } else {
kusano 7d535a
	kx = 1 - (*n - 1) * *incx;
kusano 7d535a
    }
kusano 7d535a
    if (*incy > 0) {
kusano 7d535a
	ky = 1;
kusano 7d535a
    } else {
kusano 7d535a
	ky = 1 - (*n - 1) * *incy;
kusano 7d535a
    }
kusano 7d535a
kusano 7d535a
/*     Start the operations. In this version the elements of A are   
kusano 7d535a
       accessed sequentially with one pass through the triangular part   
kusano 7d535a
       of A.   
kusano 7d535a
kusano 7d535a
       First form  y := beta*y. */
kusano 7d535a
kusano 7d535a
    if (beta->r != 1. || beta->i != 0.) {
kusano 7d535a
	if (*incy == 1) {
kusano 7d535a
	    if (beta->r == 0. && beta->i == 0.) {
kusano 7d535a
		i__1 = *n;
kusano 7d535a
		for (i = 1; i <= *n; ++i) {
kusano 7d535a
		    i__2 = i;
kusano 7d535a
		    Y(i).r = 0., Y(i).i = 0.;
kusano 7d535a
/* L10: */
kusano 7d535a
		}
kusano 7d535a
	    } else {
kusano 7d535a
		i__1 = *n;
kusano 7d535a
		for (i = 1; i <= *n; ++i) {
kusano 7d535a
		    i__2 = i;
kusano 7d535a
		    i__3 = i;
kusano 7d535a
		    z__1.r = beta->r * Y(i).r - beta->i * Y(i).i, 
kusano 7d535a
			    z__1.i = beta->r * Y(i).i + beta->i * Y(i)
kusano 7d535a
			    .r;
kusano 7d535a
		    Y(i).r = z__1.r, Y(i).i = z__1.i;
kusano 7d535a
/* L20: */
kusano 7d535a
		}
kusano 7d535a
	    }
kusano 7d535a
	} else {
kusano 7d535a
	    iy = ky;
kusano 7d535a
	    if (beta->r == 0. && beta->i == 0.) {
kusano 7d535a
		i__1 = *n;
kusano 7d535a
		for (i = 1; i <= *n; ++i) {
kusano 7d535a
		    i__2 = iy;
kusano 7d535a
		    Y(iy).r = 0., Y(iy).i = 0.;
kusano 7d535a
		    iy += *incy;
kusano 7d535a
/* L30: */
kusano 7d535a
		}
kusano 7d535a
	    } else {
kusano 7d535a
		i__1 = *n;
kusano 7d535a
		for (i = 1; i <= *n; ++i) {
kusano 7d535a
		    i__2 = iy;
kusano 7d535a
		    i__3 = iy;
kusano 7d535a
		    z__1.r = beta->r * Y(iy).r - beta->i * Y(iy).i, 
kusano 7d535a
			    z__1.i = beta->r * Y(iy).i + beta->i * Y(iy)
kusano 7d535a
			    .r;
kusano 7d535a
		    Y(iy).r = z__1.r, Y(iy).i = z__1.i;
kusano 7d535a
		    iy += *incy;
kusano 7d535a
/* L40: */
kusano 7d535a
		}
kusano 7d535a
	    }
kusano 7d535a
	}
kusano 7d535a
    }
kusano 7d535a
    if (alpha->r == 0. && alpha->i == 0.) {
kusano 7d535a
	return 0;
kusano 7d535a
    }
kusano 7d535a
    if (lsame_(uplo, "U")) {
kusano 7d535a
kusano 7d535a
/*        Form  y  when A is stored in upper triangle. */
kusano 7d535a
kusano 7d535a
	if (*incx == 1 && *incy == 1) {
kusano 7d535a
	    i__1 = *n;
kusano 7d535a
	    for (j = 1; j <= *n; ++j) {
kusano 7d535a
		i__2 = j;
kusano 7d535a
		z__1.r = alpha->r * X(j).r - alpha->i * X(j).i, z__1.i =
kusano 7d535a
			 alpha->r * X(j).i + alpha->i * X(j).r;
kusano 7d535a
		temp1.r = z__1.r, temp1.i = z__1.i;
kusano 7d535a
		temp2.r = 0., temp2.i = 0.;
kusano 7d535a
		i__2 = j - 1;
kusano 7d535a
		for (i = 1; i <= j-1; ++i) {
kusano 7d535a
		    i__3 = i;
kusano 7d535a
		    i__4 = i;
kusano 7d535a
		    i__5 = i + j * a_dim1;
kusano 7d535a
		    z__2.r = temp1.r * A(i,j).r - temp1.i * A(i,j).i, 
kusano 7d535a
			    z__2.i = temp1.r * A(i,j).i + temp1.i * A(i,j)
kusano 7d535a
			    .r;
kusano 7d535a
		    z__1.r = Y(i).r + z__2.r, z__1.i = Y(i).i + z__2.i;
kusano 7d535a
		    Y(i).r = z__1.r, Y(i).i = z__1.i;
kusano 7d535a
		    i__3 = i + j * a_dim1;
kusano 7d535a
		    i__4 = i;
kusano 7d535a
		    z__2.r = A(i,j).r * X(i).r - A(i,j).i * X(i).i, 
kusano 7d535a
			    z__2.i = A(i,j).r * X(i).i + A(i,j).i * X(
kusano 7d535a
			    i).r;
kusano 7d535a
		    z__1.r = temp2.r + z__2.r, z__1.i = temp2.i + z__2.i;
kusano 7d535a
		    temp2.r = z__1.r, temp2.i = z__1.i;
kusano 7d535a
/* L50: */
kusano 7d535a
		}
kusano 7d535a
		i__2 = j;
kusano 7d535a
		i__3 = j;
kusano 7d535a
		i__4 = j + j * a_dim1;
kusano 7d535a
		z__3.r = temp1.r * A(j,j).r - temp1.i * A(j,j).i, z__3.i = 
kusano 7d535a
			temp1.r * A(j,j).i + temp1.i * A(j,j).r;
kusano 7d535a
		z__2.r = Y(j).r + z__3.r, z__2.i = Y(j).i + z__3.i;
kusano 7d535a
		z__4.r = alpha->r * temp2.r - alpha->i * temp2.i, z__4.i = 
kusano 7d535a
			alpha->r * temp2.i + alpha->i * temp2.r;
kusano 7d535a
		z__1.r = z__2.r + z__4.r, z__1.i = z__2.i + z__4.i;
kusano 7d535a
		Y(j).r = z__1.r, Y(j).i = z__1.i;
kusano 7d535a
/* L60: */
kusano 7d535a
	    }
kusano 7d535a
	} else {
kusano 7d535a
	    jx = kx;
kusano 7d535a
	    jy = ky;
kusano 7d535a
	    i__1 = *n;
kusano 7d535a
	    for (j = 1; j <= *n; ++j) {
kusano 7d535a
		i__2 = jx;
kusano 7d535a
		z__1.r = alpha->r * X(jx).r - alpha->i * X(jx).i, z__1.i =
kusano 7d535a
			 alpha->r * X(jx).i + alpha->i * X(jx).r;
kusano 7d535a
		temp1.r = z__1.r, temp1.i = z__1.i;
kusano 7d535a
		temp2.r = 0., temp2.i = 0.;
kusano 7d535a
		ix = kx;
kusano 7d535a
		iy = ky;
kusano 7d535a
		i__2 = j - 1;
kusano 7d535a
		for (i = 1; i <= j-1; ++i) {
kusano 7d535a
		    i__3 = iy;
kusano 7d535a
		    i__4 = iy;
kusano 7d535a
		    i__5 = i + j * a_dim1;
kusano 7d535a
		    z__2.r = temp1.r * A(i,j).r - temp1.i * A(i,j).i, 
kusano 7d535a
			    z__2.i = temp1.r * A(i,j).i + temp1.i * A(i,j)
kusano 7d535a
			    .r;
kusano 7d535a
		    z__1.r = Y(iy).r + z__2.r, z__1.i = Y(iy).i + z__2.i;
kusano 7d535a
		    Y(iy).r = z__1.r, Y(iy).i = z__1.i;
kusano 7d535a
		    i__3 = i + j * a_dim1;
kusano 7d535a
		    i__4 = ix;
kusano 7d535a
		    z__2.r = A(i,j).r * X(ix).r - A(i,j).i * X(ix).i, 
kusano 7d535a
			    z__2.i = A(i,j).r * X(ix).i + A(i,j).i * X(
kusano 7d535a
			    ix).r;
kusano 7d535a
		    z__1.r = temp2.r + z__2.r, z__1.i = temp2.i + z__2.i;
kusano 7d535a
		    temp2.r = z__1.r, temp2.i = z__1.i;
kusano 7d535a
		    ix += *incx;
kusano 7d535a
		    iy += *incy;
kusano 7d535a
/* L70: */
kusano 7d535a
		}
kusano 7d535a
		i__2 = jy;
kusano 7d535a
		i__3 = jy;
kusano 7d535a
		i__4 = j + j * a_dim1;
kusano 7d535a
		z__3.r = temp1.r * A(j,j).r - temp1.i * A(j,j).i, z__3.i = 
kusano 7d535a
			temp1.r * A(j,j).i + temp1.i * A(j,j).r;
kusano 7d535a
		z__2.r = Y(jy).r + z__3.r, z__2.i = Y(jy).i + z__3.i;
kusano 7d535a
		z__4.r = alpha->r * temp2.r - alpha->i * temp2.i, z__4.i = 
kusano 7d535a
			alpha->r * temp2.i + alpha->i * temp2.r;
kusano 7d535a
		z__1.r = z__2.r + z__4.r, z__1.i = z__2.i + z__4.i;
kusano 7d535a
		Y(jy).r = z__1.r, Y(jy).i = z__1.i;
kusano 7d535a
		jx += *incx;
kusano 7d535a
		jy += *incy;
kusano 7d535a
/* L80: */
kusano 7d535a
	    }
kusano 7d535a
	}
kusano 7d535a
    } else {
kusano 7d535a
kusano 7d535a
/*        Form  y  when A is stored in lower triangle. */
kusano 7d535a
kusano 7d535a
	if (*incx == 1 && *incy == 1) {
kusano 7d535a
	    i__1 = *n;
kusano 7d535a
	    for (j = 1; j <= *n; ++j) {
kusano 7d535a
		i__2 = j;
kusano 7d535a
		z__1.r = alpha->r * X(j).r - alpha->i * X(j).i, z__1.i =
kusano 7d535a
			 alpha->r * X(j).i + alpha->i * X(j).r;
kusano 7d535a
		temp1.r = z__1.r, temp1.i = z__1.i;
kusano 7d535a
		temp2.r = 0., temp2.i = 0.;
kusano 7d535a
		i__2 = j;
kusano 7d535a
		i__3 = j;
kusano 7d535a
		i__4 = j + j * a_dim1;
kusano 7d535a
		z__2.r = temp1.r * A(j,j).r - temp1.i * A(j,j).i, z__2.i = 
kusano 7d535a
			temp1.r * A(j,j).i + temp1.i * A(j,j).r;
kusano 7d535a
		z__1.r = Y(j).r + z__2.r, z__1.i = Y(j).i + z__2.i;
kusano 7d535a
		Y(j).r = z__1.r, Y(j).i = z__1.i;
kusano 7d535a
		i__2 = *n;
kusano 7d535a
		for (i = j + 1; i <= *n; ++i) {
kusano 7d535a
		    i__3 = i;
kusano 7d535a
		    i__4 = i;
kusano 7d535a
		    i__5 = i + j * a_dim1;
kusano 7d535a
		    z__2.r = temp1.r * A(i,j).r - temp1.i * A(i,j).i, 
kusano 7d535a
			    z__2.i = temp1.r * A(i,j).i + temp1.i * A(i,j)
kusano 7d535a
			    .r;
kusano 7d535a
		    z__1.r = Y(i).r + z__2.r, z__1.i = Y(i).i + z__2.i;
kusano 7d535a
		    Y(i).r = z__1.r, Y(i).i = z__1.i;
kusano 7d535a
		    i__3 = i + j * a_dim1;
kusano 7d535a
		    i__4 = i;
kusano 7d535a
		    z__2.r = A(i,j).r * X(i).r - A(i,j).i * X(i).i, 
kusano 7d535a
			    z__2.i = A(i,j).r * X(i).i + A(i,j).i * X(
kusano 7d535a
			    i).r;
kusano 7d535a
		    z__1.r = temp2.r + z__2.r, z__1.i = temp2.i + z__2.i;
kusano 7d535a
		    temp2.r = z__1.r, temp2.i = z__1.i;
kusano 7d535a
/* L90: */
kusano 7d535a
		}
kusano 7d535a
		i__2 = j;
kusano 7d535a
		i__3 = j;
kusano 7d535a
		z__2.r = alpha->r * temp2.r - alpha->i * temp2.i, z__2.i = 
kusano 7d535a
			alpha->r * temp2.i + alpha->i * temp2.r;
kusano 7d535a
		z__1.r = Y(j).r + z__2.r, z__1.i = Y(j).i + z__2.i;
kusano 7d535a
		Y(j).r = z__1.r, Y(j).i = z__1.i;
kusano 7d535a
/* L100: */
kusano 7d535a
	    }
kusano 7d535a
	} else {
kusano 7d535a
	    jx = kx;
kusano 7d535a
	    jy = ky;
kusano 7d535a
	    i__1 = *n;
kusano 7d535a
	    for (j = 1; j <= *n; ++j) {
kusano 7d535a
		i__2 = jx;
kusano 7d535a
		z__1.r = alpha->r * X(jx).r - alpha->i * X(jx).i, z__1.i =
kusano 7d535a
			 alpha->r * X(jx).i + alpha->i * X(jx).r;
kusano 7d535a
		temp1.r = z__1.r, temp1.i = z__1.i;
kusano 7d535a
		temp2.r = 0., temp2.i = 0.;
kusano 7d535a
		i__2 = jy;
kusano 7d535a
		i__3 = jy;
kusano 7d535a
		i__4 = j + j * a_dim1;
kusano 7d535a
		z__2.r = temp1.r * A(j,j).r - temp1.i * A(j,j).i, z__2.i = 
kusano 7d535a
			temp1.r * A(j,j).i + temp1.i * A(j,j).r;
kusano 7d535a
		z__1.r = Y(jy).r + z__2.r, z__1.i = Y(jy).i + z__2.i;
kusano 7d535a
		Y(jy).r = z__1.r, Y(jy).i = z__1.i;
kusano 7d535a
		ix = jx;
kusano 7d535a
		iy = jy;
kusano 7d535a
		i__2 = *n;
kusano 7d535a
		for (i = j + 1; i <= *n; ++i) {
kusano 7d535a
		    ix += *incx;
kusano 7d535a
		    iy += *incy;
kusano 7d535a
		    i__3 = iy;
kusano 7d535a
		    i__4 = iy;
kusano 7d535a
		    i__5 = i + j * a_dim1;
kusano 7d535a
		    z__2.r = temp1.r * A(i,j).r - temp1.i * A(i,j).i, 
kusano 7d535a
			    z__2.i = temp1.r * A(i,j).i + temp1.i * A(i,j)
kusano 7d535a
			    .r;
kusano 7d535a
		    z__1.r = Y(iy).r + z__2.r, z__1.i = Y(iy).i + z__2.i;
kusano 7d535a
		    Y(iy).r = z__1.r, Y(iy).i = z__1.i;
kusano 7d535a
		    i__3 = i + j * a_dim1;
kusano 7d535a
		    i__4 = ix;
kusano 7d535a
		    z__2.r = A(i,j).r * X(ix).r - A(i,j).i * X(ix).i, 
kusano 7d535a
			    z__2.i = A(i,j).r * X(ix).i + A(i,j).i * X(
kusano 7d535a
			    ix).r;
kusano 7d535a
		    z__1.r = temp2.r + z__2.r, z__1.i = temp2.i + z__2.i;
kusano 7d535a
		    temp2.r = z__1.r, temp2.i = z__1.i;
kusano 7d535a
/* L110: */
kusano 7d535a
		}
kusano 7d535a
		i__2 = jy;
kusano 7d535a
		i__3 = jy;
kusano 7d535a
		z__2.r = alpha->r * temp2.r - alpha->i * temp2.i, z__2.i = 
kusano 7d535a
			alpha->r * temp2.i + alpha->i * temp2.r;
kusano 7d535a
		z__1.r = Y(jy).r + z__2.r, z__1.i = Y(jy).i + z__2.i;
kusano 7d535a
		Y(jy).r = z__1.r, Y(jy).i = z__1.i;
kusano 7d535a
		jx += *incx;
kusano 7d535a
		jy += *incy;
kusano 7d535a
/* L120: */
kusano 7d535a
	    }
kusano 7d535a
	}
kusano 7d535a
    }
kusano 7d535a
kusano 7d535a
    return 0;
kusano 7d535a
kusano 7d535a
/*     End of ZSYMV */
kusano 7d535a
kusano 7d535a
} /* zsymv_ */
kusano 7d535a