anchor.diff:              kappa% eQ = -1.0
           end if
       end if
  -    kappa% total = kappa_pre/nu 
  -    
  -    ! neutrons
  -    nu_n = 0.0  
  -    if (ionic% Yn > 0.0) then
  -    nn = rho*ionic% Yn/(1.0-chi)/Mneutron / density_n
  -    kFn = (threepisquare*nn)**onethird
  -    xFn = hbar*kFn/Mneutron/clight
  -    eFn = sqrt(1.0+xFn**2)
  -    nion = (1.0-ionic%Yn)*rho/Mneutron/ionic% A /density_n    
  -    kappa_n_pre = onethird*pi**2*boltzmann**2*T*nn/Mneutron/eFn      
  -    if (K_components(icond_sf)) then
  +    kappa% electron_total = kappa_pre/nu
  +    if (K_components(icond_sf) .and. ionic% Yn > 0.0) then
  +        nn = rho*ionic% Yn/(1.0-chi)/Mneutron / density_n
  +        nion = (1.0-ionic%Yn)*rho/Mneutron/ionic% A /density_n
           kappa% sf =  sPh(nn,nion,T,ionic)
       end if
  -<<<<<<< HEAD
  -	if (K_components(icond_nQ)) then
  -    	nu_c = n_imp(nn,nion,T,ionic)
  -    	kappa% nQ = kappa_n_pre/nu_c
  -        nu_n = nu_n + nu_c
  -    end if    
  - 	if (K_components(icond_np)) then
  -        nu_c = n_phonon(nn,nion,T,ionic)
  -        kappa% np = kappa_n_pre/nu_c
  -        nu_n = nu_n + nu_c
  -    end if      
  -		kappa% total = kappa% total + kappa_n_pre/nu_n + kappa% sf
  -	end if
  -=======
       if (K_components(icond_kap)) then
           kappa% kap = Rosseland_kappa(rho,T,mu_e,ionic)
           K_opacity = 4.0*onethird*arad*clight*T**3/rho/kappa% kap
       end if
       kappa% total = kappa% sf + K_opacity
       if (nu > 0.0) kappa% total = kappa% total + kappa_pre/nu
  ->>>>>>> nworbde/master
  -    
  + 
  +     ! neutrons
  +    nu_n = 0.0  
  +    nu_c = 0.0
  +    kappa_n_pre = 0.0
  +    
  +    if (ionic% Yn > 0.0) then
  +    nn = rho/amu*ionic% Yn/(1.0-chi)
  +    kFn = (threepisquare*nn)**onethird
  +    nion = rho/amu*(1.0-ionic%Yn)/ionic%A   
  +    mnstar = Mneutron
  +    kappa_n_pre = onethird*pi**2*boltzmann**2*T*nn/mnstar    
  +
  +	call sf_get_results(0.0_dp,kfn/cm_to_fm,Tc)
  +  	! Now set the superfluid reduction factor
  +	Tns = Tc(neutron_1S0)
  +	if (T < Tns) then
  +		tau = T/Tns
  +		v = sqrt(1.0-tau)*(1.456-0.157/sqrt(tau) + 1.764/tau)
  +		R = (0.4186+sqrt(1.007**2+(0.5010*v)**2))**2.5 * exp(1.456-sqrt(1.456**2+v**2))
  +		if (R < 1.0E-10) then
  +		R = 0.0
  +		end if
  +	else
  +		R = 1.0
  +	end if
  +
  +	if (K_components(icond_nQ)) then
  +    	nu_c = n_imp(nn,nion,T,ionic,ierr)
  +    	kappa% nQ = kappa_n_pre/nu_c*R
  +        nu_n = nu_n + nu_c
  +    end if    
  + 	if (K_components(icond_np)) then
  +        nu_c = n_phonon(nn,nion,T,ionic, ierr)        
  +        kappa% np = kappa_n_pre/nu_c*R
  +        nu_n = nu_n + nu_c
  +   	end if   
  +
  +	kappa % neutron_total = kappa_n_pre/nu_n*R
  +
  +   	end if	   
  +
  + 	kappa % total = kappa% total + kappa % neutron_total
  +
       contains
       subroutine clear_kappa()
           kappa% total = 0.0
           kappa% ee  = 0.0
           kappa% ei = 0.0
           kappa% eQ = 0.0
  -<<<<<<< HEAD
  -        kappa% sf = 0.0     
  +        kappa% sf = 0.0    
  +        kappa% kap = 0.0 
           kappa% nQ = 0.0
           kappa% np = 0.0
  -=======
  -        kappa% sf = 0.0    
  -        kappa% kap = 0.0 
  ->>>>>>> nworbde/master
  +        kappa% electron_total = 0.0
  +        kappa% neutron_total = 0.0
       end subroutine clear_kappa
   end subroutine conductivity
   
  
positive.code:     sPh = onethird*Cv*vs*Lsph * K_n
  end function sPh
  
  <<<<<<< HEAD
  function n_imp(nn, nion, temperature, ionic) result(nu)
      ! neutron-impurity scattering, following S. Reddy's notes
      !
      use nucchem_def, only: composition_info_type
      use constants_def
      real(dp), intent(in) :: nn, nion, temperature
      type(composition_info_type), intent(in) :: ionic
      real(dp) :: ne		! number density of electrons
      real(dp) :: kFn 	! neutron Fermi wavevector 
      real(dp) :: EFn 	! neutron Fermi energy
      real(dp) :: nu		! neutron-impurity scattering frequency
      real(dp) :: V_0 	! effective neutron-impurity potential at V(R)
      real(dp) :: R_a		! radius of scattering center
      real(dp) :: Lambda_n_imp 	! Coulomb logarithm (Potekhin et al. 1999)
  	real(dp) :: fac, sf_frac, Tc
  
  	V_0 = 20.0*mev_to_ergs
  	R_a = 10.0*fm_to_cm
  	Tc = 1.0E8	
  	Lambda_n_imp = 1.0
      fac = 4.0/27.0/pi/hbar**2/clight
      sf_frac = exp(sqrt(1.0-temperature/Tc))
     
      kFn = (3.0*pi*nn)**onethird
      EFn = hbar*clight*kFn
      ne = ionic%Z*nion 
      nu = fac*EFn*(ne/nn)*(V_0*R_a)**2*Lambda_n_imp*ionic%Q/(ionic%Z)**3
  end function n_imp
  
  function n_phonon(nn, nion, temperature, ionic) result(nu)
      ! neutron-phonon scattering, following S. Reddy's notes
      !
      use nucchem_def, only: composition_info_type
      use constants_def
      real(dp), intent(in) :: nn, nion, temperature
      type(composition_info_type), intent(in) :: ionic
      real(dp) :: ne		! electron number density
      real(dp) :: kFn 	! neutron Fermi wavevector 
      real(dp) :: EFn 	! neutron Fermi energy
      real(dp) :: nu		! neutron-impurity scattering frequency
      real(dp) :: V_0 	! effective neutron-impurity potential at V(R)
      real(dp) :: R_a		! radius of scattering center
      real(dp) :: Lambda_n_phonon 	! Coulomb logarithm (Potekhin et al. 1999)
  	real(dp) :: fac, sf_frac, Tc
  
  	V_0 = 20.0*mev_to_ergs
  	R_a = 10.0*fm_to_cm
  	Tc = 1.0E8
  	Lambda_n_phonon = 1.0
      fac = 4.0/27.0/pi/hbar**2/clight
      ! suppression near critical temperature
      sf_frac = exp(sqrt(1.0-temperature/Tc))
     
      kFn = (3.0*pi*nn)**onethird
      EFn = hbar*clight*kFn
      ne = ionic%Z*nion 
      nu = fac*EFn*(nion/nn)*(V_0*R_a)**2*Lambda_n_phonon
  end function n_phonon
  =======
  function electron_scattering(eta_e,theta,Ye) result(kTh)
      use constants_def
      real(dp), intent(in) :: eta_e   ! electron chemical pot./kT
  
negative.code:     type(composition_info_type), intent(in) :: ionic
      type(conductivity_components), intent(out) :: kappa
      integer, intent(in) :: which_ee, which_eQ
      integer :: ierr
      logical, dimension(num_conductivity_channels), intent(in) :: K_components
      real(dp) :: nn,nion, nu, nu_c, kappa_pre, ne, kF, xF, eF, Gamma_e, K_opacity
      real(dp) :: kappa_n_pre, nu_n, kfn, mnstar, tau, v, R, Tns
      real(dp) :: Tc(max_number_sf_types)
      
      ne = rho/amu*ionic%Ye
      kF = (threepisquare*ne)**onethird
  