      FUNCTION QFADD(X,Y,IX,IY)
      DOUBLE PRECISION QFADD,X,Y,FEXPC,DTEN(-3/3),FIEXP,XX,YY,DABS
      COMMON/CIEXP/IEXP,FEXPC
      DATA DTEN/0.0D0,1.0D-24,1.0D-12,1.0D0,1.0D12,1.0D24,1.0D36/
      XX=FIEXP(X)
      IXX=IX+IEXP
      IS=-IEXP
      YY=FIEXP(Y)
      IYY=IY+IEXP
      IF(DABS(XX)-1.0D0)11,12,12
   11 IXX=IXX-12
      XX=XX*1.0D12
      IS=IS+12
   12 I=(IXX-IYY)/12
      IF(I)25,10,15
   10 QFADD=FIEXP(XX+YY)
      ISS=IEXP
      IEXP=IXX+IEXP
      FEXPC=DTEN((IS-ISS)/12)
      GO TO 35
   15 IF(I-2)16,17,18
   16 QFADD=FIEXP(XX+YY*DTEN(-I))
      ISS=IEXP
      IEXP=IXX+IEXP
      FEXPC=DTEN((IS-ISS)/12)
      GO TO 35
   17 IF(DABS(YY)-1.0D0)18,18,16
   18 QFADD=XX
      IEXP=IXX
      FEXPC=DTEN(IS/12)
      GO TO 35
   25 IF(I+3)28,27,16
   27 IF(DABS(YY)-1.0D0)16,28,28
   28 QFADD=YY
      IEXP=IYY
      FEXPC=0.0D0
   35 RETURN
      END
      FUNCTION QFSQRT(X,I)
      DOUBLE PRECISION QFSQRT,DSQRT,X,FIEXP,FEXPC,XX
      COMMON/CIEXP/IEXP,FEXPC
      XX=FIEXP(X)
      IXX=(I+IEXP)/12
      II=IXX/2
      IF(IXX-2*II)11,15,10
   10 QFSQRT=DSQRT(XX*1.0D12)
      GO TO 20
   11 QFSQRT=DSQRT(XX*1.0D-12)
      GO TO 20
   15 QFSQRT=DSQRT(XX)
   20 IEXP=II*12
      RETURN
      END
      FUNCTION QCUBRT(X,I)
      IMPLICIT DOUBLE PRECISION (A-H,O-R,T-Z)
      COMMON/CIEXP/IEXP,FEXPC
      DIMENSION DTEN3(-2/2)
      DATA DTEN3/1.0D-8,1.0D-4,1.0D0,1.0D4,1.0D8/
      XX=X
      II=I/12
      I3=3*(II/3)
      IF(XX)15,10,15
   10 QCUBRT=0.0D0
      IEXP=0
      GO TO 20
   15 QCUBRT=DSIGN(DEXP(DLOG(DABS(XX))/3.0D0),XX)
      QCUBRT=FIEXP(QCUBRT*DTEN3(II-I3))
      IEXP=IEXP+I3*12
   20 RETURN
      END
      FUNCTION QATANB(Y,X,IY,IX)
      IMPLICIT DOUBLE PRECISION (A-H,O-R,T-Z)
      COMMON/CIEXP/IEXP,FEXPC
      REAL FLOAT,ALOG10,ABS
      PI=3.14159265358979324D0
      IF(X)4,15,4
    4 IF(Y)10,5,10
    5 IF(X)7,12,12
    7 QATANB=PI
      GO TO 30
   10 XX=FIEXP(X)
      IXX=IEXP+IX
      YY=FIEXP(Y)
      IYY=IEXP+IY
      R=YY/DABS(XX)
      IYX=IYY-IXX
      ST=ALOG10(ABS(SNGL(R)))+FLOAT(IYX)
      IF(ST+20.0E0)11,11,14
   11 IF(XX)13,12,12
   12 QATANB=0.0D0
      GO TO 30
   13 QATANB=DSIGN(PI,Y)
      GO TO 30
   14 IF(ST-20.0E0)16,15,15
   15 QATANB=DSIGN(PI/2.0D0,Y)
      GO TO 30
   16 QATANB=QATAN2(QFLOAT(R,IYX),DSIGN(1.0D0,XX))
   30 RETURN
      END
      FUNCTION QATAN2(Y,X)
      DOUBLE PRECISION QATAN2,Y,X,QATAN,DABS,PI,PI2,R
      PI=3.14159265358979324D0
      IF(DABS(Y)-DABS(X))10,10,15
   10 R=Y/X
      IF(X)22,20,20
   20 QATAN2=QATAN(R)
      GO TO 30
   22 IF(Y)23,24,24
   23 QATAN2=QATAN(R)-PI
      GO TO 30
   24 QATAN2=PI+QATAN(R)
      GO TO 30
   15 R=X/Y
      PI2=PI/2.0D0
      IF(Y)27,25,25
   25 QATAN2=PI2-QATAN(R)
      GO TO 30
   27 QATAN2=-PI2-QATAN(R)
   30 RETURN
      END
      DOUBLE PRECISION FUNCTION QATAN(X)
      DOUBLE PRECISION X,PI,DATAN,DSQRT,Z,DS2,A
      Z=X
      PI=3.14159265358979324D0
      IF(X)2,1,3
    1 QATAN=0
      RETURN
    2 M1=-1
      X=-X
      GO TO 4
    3 M1=1
    4 IF(X-1.D0)6,6,5
    5 M2=-1
      X=1.D0/X
      GO TO 7
    6 M2=1
    7 DS2=DSQRT(2.0D0)
      A=-1.0D0+DS2
      IF(X-A)9,9,8
    8 M3=-1
      X=(1.D0-X)/(1.D0+X)
      GO TO 10
    9 M3=1
   10 IF(X+((1.0D0+DS2)*(1.0D0-DSQRT(4.0D0-2.0D0*DS2))))12,12,11
   11 M4=-1
      X=(A-X)/(1.0D0+A*X)
      GO TO 13
   12 M4=1
   13 QATAN=DATAN(X)
      IF(M4.LT.0)QATAN=(PI/8.D0)-QATAN
      IF(M3.LT.0)QATAN=(PI/4.D0)-QATAN
      IF(M2.LT.0)QATAN=(PI/2.D0)-QATAN
      IF(M1.LT.0)QATAN=-QATAN
      X=Z
      RETURN
      END
      FUNCTION QEXP(X)
      IMPLICIT DOUBLE PRECISION (A-H,O-R,T-Z)
      REAL FLOAT
      COMMON/CIEXP/IEXP,FEXPC
      A=27.6310211159285482D0
      XX=DABS(X)
      IX=IDINT(XX/A)
      IEXP=12*IX
      XX=XX-DBLE(FLOAT(IX))*A
      IF(X)20,15,10
   15 QEXP=1.0D0
      GO TO 25
   20 XX=-XX
      IEXP=-IEXP
   10 QEXP=DEXP(XX)
   25 RETURN
      END
      FUNCTION FIEXP(X)
      DOUBLE PRECISION FIEXP,X,FEXPC,DABS,XX
      COMMON/CIEXP/IEXP,FEXPC
      XX=DABS(X)
      IF(XX-1.0D12)20,10,10
   10 IF(XX-1.0D24)16,15,15
   15 IEXP=24
      FIEXP=X*1.0D-24
      FEXPC=1.0D-24
      GO TO 25
   16 IEXP=12
      FIEXP=X*1.0D-12
      FEXPC=1.0D-12
      GO TO 25
   20 IF(XX-1.0D-12)22,22,21
   21 IEXP=0
      FIEXP=X
      FEXPC=1.0D0
      GO TO 25
   22 IF(XX-1.0D-24)24,24,23
   23 IEXP=-12
      FIEXP=X*1.0D+12
      FEXPC=1.0D+12
      GO TO 25
   24 IEXP=-24
      FIEXP=X*1.0D+24
      FEXPC=1.0D+24
   25 RETURN
      END
      FUNCTION QFLOAT(X,I)
      DOUBLE PRECISION X,QFLOAT,DTEN(-1/1),XX
      DATA DTEN/1.0D-12,1.0D0,1.0D12/
      II=I/12
      XX=X
   15 IF(II-2)10,11,12
   12 XX=XX*1.0D24
      II=II-2
      GO TO 15
   11 QFLOAT=XX*1.0D24
      GO TO 35
   10 IF(II+2)13,14,16
   13 XX=XX*1.0D-24
      II=II+2
      GO TO 10
   14 QFLOAT=XX*1.0D-24
      GO TO 35
   16 QFLOAT=XX*DTEN(II)
   35 RETURN
      END
