subroutine eldenc
C     =================
C
C     compare the actual electon density to that which follows
C     from the values used in the opacity table, interpolated to
C     the actual temperature and mass density
C
      INCLUDE 'IMPLIC.FOR'
      INCLUDE 'BASICS.FOR'
      INCLUDE 'MODELQ.FOR'
      INCLUDE 'ATOMIC.FOR'
      common/eletab/elecgr(mtabt,mtabr)
      common/eospar/anmol(600,mdepth),
     *              anato(100,mdepth),
     *              anion(100,mdepth)
      common/hmolab/anh2(mdepth),anhm(mdepth)
      dimension elcon(31,mdepth)
c
      if(ipelch.gt.0) then
      write(6,600)
  600 format(/' -------------------------'/
     *        ' CHECK OF ELECTRON DENSITY'/
     *        ' -------------------------'/
     * '  ID     TEMP      ACTUAL     LTE EOS  interpol.op.tab.'/)
c     
      do id=1,nd
         t=temp(id)
         rho=dens(id)
         if(numtemp.eq.nd) then
            opac=elecgr(id,1)
            go to 10
         end if
c
         TTAB1=TEMPVEC(1)
         TTAB2=TEMPVEC(NUMTEMP)
         TL=LOG(T)
         DELTAT=(TL-TTAB1)/(TTAB2-TTAB1)*FLOAT(numtemp-1)
         JT = 1 + IDINT(DELTAT)
         IF(JT.LT.1) JT = 1
         IF(JT.GT.numtemp-1) JT = numtemp-1
         ju = jt+1
         t1i=tempvec(jt)
         t2i=tempvec(jt+1)
         dti=(TL-T1i)/(T2i-T1i)
         if(deltat.lt.0.) dti = 0.d0
C
         if(numrh(1).ne.1) then
c
c        lower temperature
c
            numrho=numrh(jt)
            rtab1=rhomat(jt,1)
            rtab2=rhomat(jt,numrho)
            RL = LOG(RHO)
            DELTAR=(RL-RTAB1)/(RTAB2-RTAB1)*FLOAT(numrho-1)
            JR = 1 + IDINT(DELTAR)
            IF(JR.LT.1) JR = 1
            IF(JR.GT.(numrho-1)) JR = numrho-1
            r1i=rhomat(jt,jr)
            r2i=rhomat(jt,jr+1)
            dri=(RL-R1i)/(R2i-R1i)
            if(deltar.lt.0.) dri = 0.d0
            opr1=elecgr(jt,jr)+
     *           dri*(elecgr(jt,jr+1)-elecgr(jt,jr))
c
c           higher temperature
c
            ju=jt+1
            numrho=numrh(ju)
            rtab1=rhomat(ju,1)
            rtab2=rhomat(ju,numrho)
            RL = LOG(RHO)
            DELTAR=(RL-RTAB1)/(RTAB2-RTAB1)*FLOAT(numrho-1)
            JR = 1 + IDINT(DELTAR)
            IF(JR.LT.1) JR = 1
            IF(JR.GT.(numrho-1)) JR = numrho-1
            r1i=rhomat(ju,jr)
            r2i=rhomat(ju,jr+1)
            dri=(RL-R1i)/(R2i-R1i)
            if(deltar.lt.0.) dri = 0.d0
            opr2=elecgr(ju,jr)+
     *           dri*(elecgr(ju,jr+1)-elecgr(ju,jr))
            opac=opr1+(opr2-opr1)*dti
          else
            jr=1
            opac=elecgr(jt,jr)+(elecgr(ju,jr)-elecgr(jt,jr))*dti
         end if
   10    continue
         elecg=exp(opac)
         call rhonen(id,t,rho,an,ane)
         write(6,601) id,t,elec(id),ane,elecg
  601    format(i4,f10.1,1p3e12.4)
      end do
      end if
c
c     electron donors
c
      if(idisk.eq.0.or.ipeldo.eq.0) return
      do id=1,nd
         t=temp(id)
         if(ifmol.gt.0.and.t.lt.tmolim) then
            rho=dens(id)
            call rhonen(id,t,rho,an,ane)
            aein=ane
            call moleq(id,t,an,aein,ane,enrg,entt,wm,1)
            do ia=1,30
               elcon(ia,id)=anion(ia,id)/elec(id)
            end do
            elcon(31,id)=-anhm(id)/elec(id)
          else
            call state(2,id,t,elec(id))
            do ia=1,30
               iat=iatex(ia)
               if(iat.gt.0) then
                  qs=0.
                  n1=n0a(iat)
                  if(ia.eq.1) n1=nfirst(ielh)
                  do i=n1,nka(iat)
                     ch=iz(iel(i))-1
                     if(ilk(i).gt.0) ch=iz(ilk(i))
                     qs=qs+ch*popul(i,id)
                  end do
                  elcon(ia,id)=qs/elec(id)
                else 
                  elcon(ia,id)=rr(ia,99)/elec(id)
               end if
            end do
            if(ielhm.gt.0) then
               elcon(31,id)=-popul(nfirst(ielhm),id)/elec(id)
             else
               aref=dens(id)/wmm(id)/ytot(id)
               elcon(31,id)=-qm*rr(1,1)*aref
            end if
         end if
      end do
c
      write(6,611)
      do id=1,nd
      write(6,612) id,temp(id),elcon(31,id),elcon(1,id),elcon(2,id),
     *             elcon(6,id),elcon(7,id),elcon(8,id),
     *             (elcon(i,id),i=11,15),elcon(20,id),elcon(26,id)
      end do
  611 format(/' RELATIVE CONTRIBUTIONS OF INDIVIDUAL ELECTRON DONORS'//
     * '  ID   TEMP      H-       H      He     C      N      O',
     * '     Na     Mg     Al',
     * '     Si     S      Ca     Fe')
  612 format(i3,f9.1,1pe10.2,0p12f7.3)  
c
      return
      end