SUBROUTINE ROTGEN( X, Y, CS, SN )
C
      DOUBLE PRECISION   X, Y, CS, SN
C
      DOUBLE PRECISION   WMACH
      COMMON    /SOLMCH/ WMACH(15)
C
C  *********************************************************************
C  ROTGEN  GENERATES A PLANE ROTATION THAT REDUCES THE VECTOR (X, Y)
C  TO THE VECTOR (A, 0),  WHERE  A  IS DEFINED AS FOLLOWS...
C
C     IF BOTH X AND Y ARE NEGLIGIBLY SMALL, OR
C     IF Y IS NEGLIGIBLE RELATIVE TO X,
C     THEN  A = X,  AND THE IDENTITY ROTATION IS RETURNED.
C
C     IF X IS NEGLIGIBLE RELATIVE TO Y,
C     THEN  A = Y,  AND THE SWAP ROTATION IS RETURNED.
C
C     OTHERWISE,  A = SIGN(X) * SQRT( X**2 + Y**2 ).
C
C  IN ALL CASES,  X  AND  Y  ARE OVERWRITTEN BY  A  AND  0,
C  AND  CS  WILL LIE IN THE CLOSED INTERVAL  (0, 1).  ALSO,
C  THE ABSOLUTE VALUE OF  CS  AND  SN  (IF NONZERO) WILL BE NO LESS
C  THAN THE MACHINE PRECISION,  EPS.
C
C  ROTGEN  GUARDS AGAINST OVERFLOW AND UNDERFLOW.
C  IT IS ASSUMED THAT  FLMIN .LT. EPS**2  (I.E.  RTMIN .LT. EPS).
C
C  SYSTEMS OPTIMIZATION LABORATORY, STANFORD UNIVERSITY.
C  VERSION OF JANUARY 1982.
C  *********************************************************************
C
      DOUBLE PRECISION   A, B, EPS, ONE, RTMIN, ZERO
      DOUBLE PRECISION   DSQRT
      DOUBLE PRECISION   DABS, DMAX1
      DATA               ONE/1.0D+0/, ZERO/0.0D+0/
C
      IF (Y .EQ. ZERO) GO TO 100
      IF (X .EQ. ZERO) GO TO 200
C
      EPS    = WMACH(3)
      RTMIN  = WMACH(6)
      A      = DABS(X)
      B      = DABS(Y)
      IF (DMAX1(A,B) .LE. RTMIN) GO TO 100
      IF (A .LT. B) GO TO 50
      IF (B .LE. EPS*A) GO TO 100
      A  = A * DSQRT( ONE + (B/A)**2 )
      GO TO 60
C
   50 IF (A .LE. EPS*B) GO TO 200
      A  = B * DSQRT( ONE + (A/B)**2 )
C
   60 IF (X .LT. ZERO) A = - A
      CS = X/A
      SN = Y/A
      X  = A
      GO TO 300
C
  100 CS = ONE
      SN = ZERO
      GO TO 300
C
  200 CS = ZERO
      SN = ONE
      X  = Y
  300 Y  = ZERO
      RETURN
C
C  END OF ROTGEN
      END