**********************************************************************
C
C Copyright (C) 1992 Roland W. Freund and Noel M. Nachtigal
C All rights reserved.
C
C This code is part of a copyrighted package. For details, see the
C file `cpyrit.doc' in the top-level directory.
C
C *****************************************************************
C ANY USE OF THIS CODE CONSTITUTES ACCEPTANCE OF THE TERMS OF THE
C COPYRIGHT NOTICE
C *****************************************************************
C
C**********************************************************************
C
SUBROUTINE DAXPBY (N,DZ,DA,DX,DB,DY)
C
C Purpose:
C This subroutine computes DZ = DA * DX + DB * DY. Several special
C cases are handled separately:
C DA = 0.0, DB = 0.0 => DZ = 0.0
C DA = 0.0, DB = 1.0 => DZ = DY (this is COPY)
C DA = 0.0, DB = -1.0 => DZ = -DY
C DA = 0.0, DB = DB => DZ = DB * DY (this is SCAL)
C DA = 1.0, DB = 0.0 => DZ = DX (this is COPY)
C DA = 1.0, DB = 1.0 => DZ = DX + DY
C DA = 1.0, DB = -1.0 => DZ = DX - DY
C DA = 1.0, DB = DB => DZ = DX + DB * DY (this is AXPY)
C DA = -1.0, DB = 0.0 => DZ = -DX
C DA = -1.0, DB = 1.0 => DZ = -DX + DY
C DA = -1.0, DB = -1.0 => DZ = -DX - DY
C DA = -1.0, DB = DB => DZ = -DX + DB * DY
C DA = DA, DB = 0.0 => DZ = DA * DX (this is SCAL)
C DA = DA, DB = 1.0 => DZ = DA * DX + DY (this is AXPY)
C DA = DA, DB = -1.0 => DZ = DA * DX - DY
C DA = DA, DB = DB => DZ = DA * DX + DB * DY
C DZ may be the same as DX or DY.
C
C Parameters:
C N = the dimension of the vectors (input).
C DZ = the vector result (output).
C DA = scalar multiplier for DX (input).
C DX = one of the vectors (input).
C DB = scalar multiplier for DY (input).
C DY = the other vector (input).
C
C Noel M. Nachtigal
C March 23, 1993
C
C**********************************************************************
C
INTEGER N
DOUBLE PRECISION DA, DB, DX(N), DY(N), DZ(N)
C
C Local variables.
C
INTEGER I
C
IF (N.LE.0) RETURN
C
IF (DA.EQ.0.0D0) THEN
IF (DB.EQ.0.0D0) THEN
C DA = 0.0, DB = 0.0 => DZ = 0.0.
DO 10 I = 1, N
DZ(I) = 0.0D0
10 CONTINUE
ELSE IF (DB.EQ.1.0D0) THEN
C DA = 0.0, DB = 1.0 => DZ = DY (this is COPY).
DO 20 I = 1, N
DZ(I) = DY(I)
20 CONTINUE
ELSE IF (DB.EQ.-1.0D0) THEN
C DA = 0.0, DB = -1.0 => DZ = -DY.
DO 30 I = 1, N
DZ(I) = -DY(I)
30 CONTINUE
ELSE
C DA = 0.0, DB = DB => DZ = DB * DY (this is SCAL).
DO 40 I = 1, N
DZ(I) = DB * DY(I)
40 CONTINUE
END IF
ELSE IF (DA.EQ.1.0D0) THEN
IF (DB.EQ.0.0D0) THEN
C DA = 1.0, DB = 0.0 => DZ = DX (this is COPY).
DO 50 I = 1, N
DZ(I) = DX(I)
50 CONTINUE
ELSE IF (DB.EQ.1.0D0) THEN
C DA = 1.0, DB = 1.0 => DZ = DX + DY.
DO 60 I = 1, N
DZ(I) = DX(I) + DY(I)
60 CONTINUE
ELSE IF (DB.EQ.-1.0D0) THEN
C DA = 1.0, DB = -1.0 => DZ = DX - DY.
DO 70 I = 1, N
DZ(I) = DX(I) - DY(I)
70 CONTINUE
ELSE
C DA = 1.0, DB = DB => DZ = DX + DB * DY (this is AXPY).
DO 80 I = 1, N
DZ(I) = DX(I) + DB * DY(I)
80 CONTINUE
END IF
ELSE IF (DA.EQ.-1.0D0) THEN
IF (DB.EQ.0.0D0) THEN
C DA = -1.0, DB = 0.0 => DZ = -DX
DO 90 I = 1, N
DZ(I) = -DX(I)
90 CONTINUE
ELSE IF (DB.EQ.1.0D0) THEN
C DA = -1.0, DB = 1.0 => DZ = -DX + DY
DO 100 I = 1, N
DZ(I) = -DX(I) + DY(I)
100 CONTINUE
ELSE IF (DB.EQ.-1.0D0) THEN
C DA = -1.0, DB = -1.0 => DZ = -DX - DY.
DO 110 I = 1, N
DZ(I) = -DX(I) - DY(I)
110 CONTINUE
ELSE
C DA = -1.0, DB = DB => DZ = -DX + DB * DY
DO 120 I = 1, N
DZ(I) = -DX(I) + DB * DY(I)
120 CONTINUE
END IF
ELSE
IF (DB.EQ.0.0D0) THEN
C DA = DA, DB = 0.0 => DZ = DA * DX (this is SCAL).
DO 130 I = 1, N
DZ(I) = DA * DX(I)
130 CONTINUE
ELSE IF (DB.EQ.1.0D0) THEN
C DA = DA, DB = 1.0 => DZ = DA * DX + DY (this is AXPY)
DO 140 I = 1, N
DZ(I) = DA * DX(I) + DY(I)
140 CONTINUE
ELSE IF (DB.EQ.-1.0D0) THEN
C DA = DA, DB = -1.0 => DZ = DA * DX - DY.
DO 150 I = 1, N
DZ(I) = DA * DX(I) - DY(I)
150 CONTINUE
ELSE
C DA = DA, DB = DB => DZ = DA * DX + DB * DY.
DO 160 I = 1, N
DZ(I) = DA * DX(I) + DB * DY(I)
160 CONTINUE
END IF
END IF
C
RETURN
END
C
C**********************************************************************