c     this program computes the fractional energy shifts and splits for
c     the Dirac model to order (Z*fine)^2 in units of 
c     (Z*fine)^(-2) *dE/E

      PROGRAM paschen_dirac

      implicit none

      integer           l,k
      double precision  dEj1,dEj2,n,j1,j2
      character*2       str

c     the table header

      WRITE(6,*) '  l       j1          dEj1      j2          dEj2'

c     set n=4

      n = 4.0d0

c     the l=0 state (j=1/2 only); zero the -1/2 state

      l    = 0
      j1   = 0.0d0
      dEj1 = 0.0d0
      j2   = float(l) + 0.50d0
      dEj2 = - (1.0d0/(j2+0.5d0) - (3.0d0/(4.0d0*n)))/n

      WRITE(6,600) l,'s',j1,dEj1,j2,dEj2

c     the l=1 to l=3 states

      DO 13 l=1,3
        j1 = float(l) - 0.50d0
        dEj1 = - (1.0d0/(j1+0.5d0) - (3.0d0/(4.0d0*n)))/n
        j2 = float(l) + 0.50d0
        dEj2 = - (1.0d0/(j2+0.5d0) - (3.0d0/(4.0d0*n)))/n
        IF (l.eq.1) str = 'p'
        IF (l.eq.2) str = 'd'
        IF (l.eq.3) str = 'f'
        WRITE(6,600) l,str,j1,dEj1,j2,dEj2
 13   CONTINUE

 600  FORMAT(1x,i3,2x,a1,2x,f5.2,2x,f13.5,2x,f5.2,2x,f13.5)

      STOP
      END